LAPACK 3.3.1
Linear Algebra PACKage
|
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