LAPACK: SRC/dsyconv.f Source File

Go to the documentation of this file.
00001       SUBROUTINE DSYCONV( UPLO, WAY, N, A, LDA, IPIV, WORK, INFO )
00002 *
00003 *  -- LAPACK PROTOTYPE routine (version 3.3.0) --
00004 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00005 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00006 *     November 2010
00007 *
00008 *  -- Written by Julie Langou of the Univ. of TN    --
00009 *
00010 *     .. Scalar Arguments ..
00011       CHARACTER          UPLO, WAY
00012       INTEGER            INFO, LDA, N
00013 *     ..
00014 *     .. Array Arguments ..
00015       INTEGER            IPIV( * )
00016       DOUBLE PRECISION   A( LDA, * ), WORK( * )
00017 *     ..
00018 *
00019 *  Purpose
00020 *  =======
00021 *
00022 *  DSYCONV convert A given by TRF into L and D and vice-versa.
00023 *  Get Non-diag elements of D (returned in workspace) and 
00024 *  apply or reverse permutation done in TRF.
00025 *
00026 *  Arguments
00027 *  =========
00028 *
00029 *  UPLO    (input) CHARACTER*1
00030 *          Specifies whether the details of the factorization are stored
00031 *          as an upper or lower triangular matrix.
00032 *          = 'U':  Upper triangular, form is A = U*D*U**T;
00033 *          = 'L':  Lower triangular, form is A = L*D*L**T.
00034 * 
00035 *  WAY     (input) CHARACTER*1
00036 *          = 'C': Convert 
00037 *          = 'R': Revert
00038 *
00039 *  N       (input) INTEGER
00040 *          The order of the matrix A.  N >= 0.
00041 *
00042 *  A       (input) DOUBLE PRECISION array, dimension (LDA,N)
00043 *          The block diagonal matrix D and the multipliers used to
00044 *          obtain the factor U or L as computed by DSYTRF.
00045 *
00046 *  LDA     (input) INTEGER
00047 *          The leading dimension of the array A.  LDA >= max(1,N).
00048 *
00049 *  IPIV    (input) INTEGER array, dimension (N)
00050 *          Details of the interchanges and the block structure of D
00051 *          as determined by DSYTRF.
00052 *
00053 * WORK     (workspace) DOUBLE PRECISION array, dimension (N)
00054 *
00055 * LWORK    (input) INTEGER
00056 *          The length of WORK.  LWORK >=1. 
00057 *          LWORK = N
00058 *
00059 *          If LWORK = -1, then a workspace query is assumed; the routine
00060 *          only calculates the optimal size of the WORK array, returns
00061 *          this value as the first entry of the WORK array, and no error
00062 *          message related to LWORK is issued by XERBLA.
00063 *
00064 *  INFO    (output) INTEGER
00065 *          = 0:  successful exit
00066 *          < 0:  if INFO = -i, the i-th argument had an illegal value
00067 *
00068 *  =====================================================================
00069 *
00070 *     .. Parameters ..
00071       DOUBLE PRECISION   ZERO
00072       PARAMETER          ( ZERO = 0.0D+0 )
00073 *     ..
00074 *     .. External Functions ..
00075       LOGICAL            LSAME
00076       EXTERNAL           LSAME
00077 *
00078 *     .. External Subroutines ..
00079       EXTERNAL           XERBLA
00080 *     .. Local Scalars ..
00081       LOGICAL            UPPER, CONVERT
00082       INTEGER            I, IP, J
00083       DOUBLE PRECISION   TEMP
00084 *     ..
00085 *     .. Executable Statements ..
00086 *
00087       INFO = 0
00088       UPPER = LSAME( UPLO, 'U' )
00089       CONVERT = LSAME( WAY, 'C' )
00090       IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
00091          INFO = -1
00092       ELSE IF( .NOT.CONVERT .AND. .NOT.LSAME( WAY, 'R' ) ) THEN
00093          INFO = -2
00094       ELSE IF( N.LT.0 ) THEN
00095          INFO = -3
00096       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
00097          INFO = -5
00098 
00099       END IF
00100       IF( INFO.NE.0 ) THEN
00101          CALL XERBLA( 'DSYCONV', -INFO )
00102          RETURN
00103       END IF
00104 *
00105 *     Quick return if possible
00106 *
00107       IF( N.EQ.0 )
00108      $   RETURN
00109 *
00110       IF( UPPER ) THEN
00111 *
00112 *      A is UPPER
00113 *
00114 *      Convert A (A is upper)
00115 *
00116 *        Convert VALUE
00117 *
00118          IF ( CONVERT ) THEN
00119             I=N
00120             WORK(1)=ZERO
00121             DO WHILE ( I .GT. 1 )
00122                IF( IPIV(I) .LT. 0 ) THEN
00123                   WORK(I)=A(I-1,I)
00124                   A(I-1,I)=ZERO
00125                   I=I-1
00126                ELSE
00127                   WORK(I)=ZERO
00128                ENDIF
00129                I=I-1
00130             END DO
00131 *
00132 *        Convert PERMUTATIONS
00133 *  
00134          I=N
00135          DO WHILE ( I .GE. 1 )
00136             IF( IPIV(I) .GT. 0) THEN
00137                IP=IPIV(I)
00138                IF( I .LT. N) THEN
00139                   DO 12 J= I+1,N
00140                     TEMP=A(IP,J)
00141                     A(IP,J)=A(I,J)
00142                     A(I,J)=TEMP
00143  12            CONTINUE
00144                ENDIF
00145             ELSE
00146               IP=-IPIV(I)
00147                IF( I .LT. N) THEN
00148              DO 13 J= I+1,N
00149                  TEMP=A(IP,J)
00150                  A(IP,J)=A(I-1,J)
00151                  A(I-1,J)=TEMP
00152  13            CONTINUE
00153                 ENDIF
00154                 I=I-1
00155            ENDIF
00156            I=I-1
00157         END DO
00158 
00159          ELSE
00160 *
00161 *      Revert A (A is upper)
00162 *
00163 *
00164 *        Revert PERMUTATIONS
00165 *  
00166             I=1
00167             DO WHILE ( I .LE. N )
00168                IF( IPIV(I) .GT. 0 ) THEN
00169                   IP=IPIV(I)
00170                   IF( I .LT. N) THEN
00171                   DO J= I+1,N
00172                     TEMP=A(IP,J)
00173                     A(IP,J)=A(I,J)
00174                     A(I,J)=TEMP
00175                   END DO
00176                   ENDIF
00177                ELSE
00178                  IP=-IPIV(I)
00179                  I=I+1
00180                  IF( I .LT. N) THEN
00181                     DO J= I+1,N
00182                        TEMP=A(IP,J)
00183                        A(IP,J)=A(I-1,J)
00184                        A(I-1,J)=TEMP
00185                     END DO
00186                  ENDIF
00187                ENDIF
00188                I=I+1
00189             END DO
00190 *
00191 *        Revert VALUE
00192 *
00193             I=N
00194             DO WHILE ( I .GT. 1 )
00195                IF( IPIV(I) .LT. 0 ) THEN
00196                   A(I-1,I)=WORK(I)
00197                   I=I-1
00198                ENDIF
00199                I=I-1
00200             END DO
00201          END IF
00202       ELSE
00203 *
00204 *      A is LOWER
00205 *
00206          IF ( CONVERT ) THEN
00207 *
00208 *      Convert A (A is lower)
00209 *
00210 *
00211 *        Convert VALUE
00212 *
00213             I=1
00214             WORK(N)=ZERO
00215             DO WHILE ( I .LE. N )
00216                IF( I.LT.N .AND. IPIV(I) .LT. 0 ) THEN
00217                   WORK(I)=A(I+1,I)
00218                   A(I+1,I)=ZERO
00219                   I=I+1
00220                ELSE
00221                   WORK(I)=ZERO
00222                ENDIF
00223                I=I+1
00224             END DO
00225 *
00226 *        Convert PERMUTATIONS
00227 *
00228          I=1
00229          DO WHILE ( I .LE. N )
00230             IF( IPIV(I) .GT. 0 ) THEN
00231                IP=IPIV(I)
00232                IF (I .GT. 1) THEN
00233                DO 22 J= 1,I-1
00234                  TEMP=A(IP,J)
00235                  A(IP,J)=A(I,J)
00236                  A(I,J)=TEMP
00237  22            CONTINUE
00238                ENDIF
00239             ELSE
00240               IP=-IPIV(I)
00241               IF (I .GT. 1) THEN
00242               DO 23 J= 1,I-1
00243                  TEMP=A(IP,J)
00244                  A(IP,J)=A(I+1,J)
00245                  A(I+1,J)=TEMP
00246  23           CONTINUE
00247               ENDIF
00248               I=I+1
00249            ENDIF
00250            I=I+1
00251         END DO
00252          ELSE
00253 *
00254 *      Revert A (A is lower)
00255 *
00256 *
00257 *        Revert PERMUTATIONS
00258 *
00259             I=N
00260             DO WHILE ( I .GE. 1 )
00261                IF( IPIV(I) .GT. 0 ) THEN
00262                   IP=IPIV(I)
00263                   IF (I .GT. 1) THEN
00264                      DO J= 1,I-1
00265                         TEMP=A(I,J)
00266                         A(I,J)=A(IP,J)
00267                         A(IP,J)=TEMP
00268                      END DO
00269                   ENDIF
00270                ELSE
00271                   IP=-IPIV(I)
00272                   I=I-1
00273                   IF (I .GT. 1) THEN
00274                      DO J= 1,I-1
00275                         TEMP=A(I+1,J)
00276                         A(I+1,J)=A(IP,J)
00277                         A(IP,J)=TEMP
00278                      END DO
00279                   ENDIF
00280                ENDIF
00281                I=I-1
00282             END DO
00283 *
00284 *        Revert VALUE
00285 *
00286             I=1
00287             DO WHILE ( I .LE. N-1 )
00288                IF( IPIV(I) .LT. ZERO ) THEN
00289                   A(I+1,I)=WORK(I)
00290                   I=I+1
00291                ENDIF
00292                I=I+1
00293             END DO
00294          END IF
00295       END IF
00296 
00297       RETURN
00298 *
00299 *     End of DSYCONV
00300 *
00301       END