![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b STPTTF 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download STPTTF + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/stpttf.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/stpttf.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/stpttf.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE STPTTF( TRANSR, UPLO, N, AP, ARF, INFO ) 00022 * 00023 * .. Scalar Arguments .. 00024 * CHARACTER TRANSR, UPLO 00025 * INTEGER INFO, N 00026 * .. 00027 * .. Array Arguments .. 00028 * REAL AP( 0: * ), ARF( 0: * ) 00029 * 00030 * 00031 *> \par Purpose: 00032 * ============= 00033 *> 00034 *> \verbatim 00035 *> 00036 *> STPTTF copies a triangular matrix A from standard packed format (TP) 00037 *> to rectangular full packed format (TF). 00038 *> \endverbatim 00039 * 00040 * Arguments: 00041 * ========== 00042 * 00043 *> \param[in] TRANSR 00044 *> \verbatim 00045 *> TRANSR is CHARACTER*1 00046 *> = 'N': ARF in Normal format is wanted; 00047 *> = 'T': ARF in Conjugate-transpose format is wanted. 00048 *> \endverbatim 00049 *> 00050 *> \param[in] UPLO 00051 *> \verbatim 00052 *> UPLO is CHARACTER*1 00053 *> = 'U': A is upper triangular; 00054 *> = 'L': A is lower triangular. 00055 *> \endverbatim 00056 *> 00057 *> \param[in] N 00058 *> \verbatim 00059 *> N is INTEGER 00060 *> The order of the matrix A. N >= 0. 00061 *> \endverbatim 00062 *> 00063 *> \param[in] AP 00064 *> \verbatim 00065 *> AP is REAL array, dimension ( N*(N+1)/2 ), 00066 *> On entry, the upper or lower triangular matrix A, packed 00067 *> columnwise in a linear array. The j-th column of A is stored 00068 *> in the array AP as follows: 00069 *> if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; 00070 *> if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n. 00071 *> \endverbatim 00072 *> 00073 *> \param[out] ARF 00074 *> \verbatim 00075 *> ARF is REAL array, dimension ( N*(N+1)/2 ), 00076 *> On exit, the upper or lower triangular matrix A stored in 00077 *> RFP format. For a further discussion see Notes below. 00078 *> \endverbatim 00079 *> 00080 *> \param[out] INFO 00081 *> \verbatim 00082 *> INFO is INTEGER 00083 *> = 0: successful exit 00084 *> < 0: if INFO = -i, the i-th argument had an illegal value 00085 *> \endverbatim 00086 * 00087 * Authors: 00088 * ======== 00089 * 00090 *> \author Univ. of Tennessee 00091 *> \author Univ. of California Berkeley 00092 *> \author Univ. of Colorado Denver 00093 *> \author NAG Ltd. 00094 * 00095 *> \date November 2011 00096 * 00097 *> \ingroup realOTHERcomputational 00098 * 00099 *> \par Further Details: 00100 * ===================== 00101 *> 00102 *> \verbatim 00103 *> 00104 *> We first consider Rectangular Full Packed (RFP) Format when N is 00105 *> even. We give an example where N = 6. 00106 *> 00107 *> AP is Upper AP is Lower 00108 *> 00109 *> 00 01 02 03 04 05 00 00110 *> 11 12 13 14 15 10 11 00111 *> 22 23 24 25 20 21 22 00112 *> 33 34 35 30 31 32 33 00113 *> 44 45 40 41 42 43 44 00114 *> 55 50 51 52 53 54 55 00115 *> 00116 *> 00117 *> Let TRANSR = 'N'. RFP holds AP as follows: 00118 *> For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last 00119 *> three columns of AP upper. The lower triangle A(4:6,0:2) consists of 00120 *> the transpose of the first three columns of AP upper. 00121 *> For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first 00122 *> three columns of AP lower. The upper triangle A(0:2,0:2) consists of 00123 *> the transpose of the last three columns of AP lower. 00124 *> This covers the case N even and TRANSR = 'N'. 00125 *> 00126 *> RFP A RFP A 00127 *> 00128 *> 03 04 05 33 43 53 00129 *> 13 14 15 00 44 54 00130 *> 23 24 25 10 11 55 00131 *> 33 34 35 20 21 22 00132 *> 00 44 45 30 31 32 00133 *> 01 11 55 40 41 42 00134 *> 02 12 22 50 51 52 00135 *> 00136 *> Now let TRANSR = 'T'. RFP A in both UPLO cases is just the 00137 *> transpose of RFP A above. One therefore gets: 00138 *> 00139 *> 00140 *> RFP A RFP A 00141 *> 00142 *> 03 13 23 33 00 01 02 33 00 10 20 30 40 50 00143 *> 04 14 24 34 44 11 12 43 44 11 21 31 41 51 00144 *> 05 15 25 35 45 55 22 53 54 55 22 32 42 52 00145 *> 00146 *> 00147 *> We then consider Rectangular Full Packed (RFP) Format when N is 00148 *> odd. We give an example where N = 5. 00149 *> 00150 *> AP is Upper AP is Lower 00151 *> 00152 *> 00 01 02 03 04 00 00153 *> 11 12 13 14 10 11 00154 *> 22 23 24 20 21 22 00155 *> 33 34 30 31 32 33 00156 *> 44 40 41 42 43 44 00157 *> 00158 *> 00159 *> Let TRANSR = 'N'. RFP holds AP as follows: 00160 *> For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last 00161 *> three columns of AP upper. The lower triangle A(3:4,0:1) consists of 00162 *> the transpose of the first two columns of AP upper. 00163 *> For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first 00164 *> three columns of AP lower. The upper triangle A(0:1,1:2) consists of 00165 *> the transpose of the last two columns of AP lower. 00166 *> This covers the case N odd and TRANSR = 'N'. 00167 *> 00168 *> RFP A RFP A 00169 *> 00170 *> 02 03 04 00 33 43 00171 *> 12 13 14 10 11 44 00172 *> 22 23 24 20 21 22 00173 *> 00 33 34 30 31 32 00174 *> 01 11 44 40 41 42 00175 *> 00176 *> Now let TRANSR = 'T'. RFP A in both UPLO cases is just the 00177 *> transpose of RFP A above. One therefore gets: 00178 *> 00179 *> RFP A RFP A 00180 *> 00181 *> 02 12 22 00 01 00 10 20 30 40 50 00182 *> 03 13 23 33 11 33 11 21 31 41 51 00183 *> 04 14 24 34 44 43 44 22 32 42 52 00184 *> \endverbatim 00185 *> 00186 * ===================================================================== 00187 SUBROUTINE STPTTF( TRANSR, UPLO, N, AP, ARF, INFO ) 00188 * 00189 * -- LAPACK computational routine (version 3.4.0) -- 00190 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00191 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00192 * November 2011 00193 * 00194 * .. Scalar Arguments .. 00195 CHARACTER TRANSR, UPLO 00196 INTEGER INFO, N 00197 * .. 00198 * .. Array Arguments .. 00199 REAL AP( 0: * ), ARF( 0: * ) 00200 * 00201 * ===================================================================== 00202 * 00203 * .. Parameters .. 00204 * .. 00205 * .. Local Scalars .. 00206 LOGICAL LOWER, NISODD, NORMALTRANSR 00207 INTEGER N1, N2, K, NT 00208 INTEGER I, J, IJ 00209 INTEGER IJP, JP, LDA, JS 00210 * .. 00211 * .. External Functions .. 00212 LOGICAL LSAME 00213 EXTERNAL LSAME 00214 * .. 00215 * .. External Subroutines .. 00216 EXTERNAL XERBLA 00217 * .. 00218 * .. Intrinsic Functions .. 00219 INTRINSIC MOD 00220 * .. 00221 * .. Executable Statements .. 00222 * 00223 * Test the input parameters. 00224 * 00225 INFO = 0 00226 NORMALTRANSR = LSAME( TRANSR, 'N' ) 00227 LOWER = LSAME( UPLO, 'L' ) 00228 IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'T' ) ) THEN 00229 INFO = -1 00230 ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN 00231 INFO = -2 00232 ELSE IF( N.LT.0 ) THEN 00233 INFO = -3 00234 END IF 00235 IF( INFO.NE.0 ) THEN 00236 CALL XERBLA( 'STPTTF', -INFO ) 00237 RETURN 00238 END IF 00239 * 00240 * Quick return if possible 00241 * 00242 IF( N.EQ.0 ) 00243 $ RETURN 00244 * 00245 IF( N.EQ.1 ) THEN 00246 IF( NORMALTRANSR ) THEN 00247 ARF( 0 ) = AP( 0 ) 00248 ELSE 00249 ARF( 0 ) = AP( 0 ) 00250 END IF 00251 RETURN 00252 END IF 00253 * 00254 * Size of array ARF(0:NT-1) 00255 * 00256 NT = N*( N+1 ) / 2 00257 * 00258 * Set N1 and N2 depending on LOWER 00259 * 00260 IF( LOWER ) THEN 00261 N2 = N / 2 00262 N1 = N - N2 00263 ELSE 00264 N1 = N / 2 00265 N2 = N - N1 00266 END IF 00267 * 00268 * If N is odd, set NISODD = .TRUE. 00269 * If N is even, set K = N/2 and NISODD = .FALSE. 00270 * 00271 * set lda of ARF^C; ARF^C is (0:(N+1)/2-1,0:N-noe) 00272 * where noe = 0 if n is even, noe = 1 if n is odd 00273 * 00274 IF( MOD( N, 2 ).EQ.0 ) THEN 00275 K = N / 2 00276 NISODD = .FALSE. 00277 LDA = N + 1 00278 ELSE 00279 NISODD = .TRUE. 00280 LDA = N 00281 END IF 00282 * 00283 * ARF^C has lda rows and n+1-noe cols 00284 * 00285 IF( .NOT.NORMALTRANSR ) 00286 $ LDA = ( N+1 ) / 2 00287 * 00288 * start execution: there are eight cases 00289 * 00290 IF( NISODD ) THEN 00291 * 00292 * N is odd 00293 * 00294 IF( NORMALTRANSR ) THEN 00295 * 00296 * N is odd and TRANSR = 'N' 00297 * 00298 IF( LOWER ) THEN 00299 * 00300 * N is odd, TRANSR = 'N', and UPLO = 'L' 00301 * 00302 IJP = 0 00303 JP = 0 00304 DO J = 0, N2 00305 DO I = J, N - 1 00306 IJ = I + JP 00307 ARF( IJ ) = AP( IJP ) 00308 IJP = IJP + 1 00309 END DO 00310 JP = JP + LDA 00311 END DO 00312 DO I = 0, N2 - 1 00313 DO J = 1 + I, N2 00314 IJ = I + J*LDA 00315 ARF( IJ ) = AP( IJP ) 00316 IJP = IJP + 1 00317 END DO 00318 END DO 00319 * 00320 ELSE 00321 * 00322 * N is odd, TRANSR = 'N', and UPLO = 'U' 00323 * 00324 IJP = 0 00325 DO J = 0, N1 - 1 00326 IJ = N2 + J 00327 DO I = 0, J 00328 ARF( IJ ) = AP( IJP ) 00329 IJP = IJP + 1 00330 IJ = IJ + LDA 00331 END DO 00332 END DO 00333 JS = 0 00334 DO J = N1, N - 1 00335 IJ = JS 00336 DO IJ = JS, JS + J 00337 ARF( IJ ) = AP( IJP ) 00338 IJP = IJP + 1 00339 END DO 00340 JS = JS + LDA 00341 END DO 00342 * 00343 END IF 00344 * 00345 ELSE 00346 * 00347 * N is odd and TRANSR = 'T' 00348 * 00349 IF( LOWER ) THEN 00350 * 00351 * N is odd, TRANSR = 'T', and UPLO = 'L' 00352 * 00353 IJP = 0 00354 DO I = 0, N2 00355 DO IJ = I*( LDA+1 ), N*LDA - 1, LDA 00356 ARF( IJ ) = AP( IJP ) 00357 IJP = IJP + 1 00358 END DO 00359 END DO 00360 JS = 1 00361 DO J = 0, N2 - 1 00362 DO IJ = JS, JS + N2 - J - 1 00363 ARF( IJ ) = AP( IJP ) 00364 IJP = IJP + 1 00365 END DO 00366 JS = JS + LDA + 1 00367 END DO 00368 * 00369 ELSE 00370 * 00371 * N is odd, TRANSR = 'T', and UPLO = 'U' 00372 * 00373 IJP = 0 00374 JS = N2*LDA 00375 DO J = 0, N1 - 1 00376 DO IJ = JS, JS + J 00377 ARF( IJ ) = AP( IJP ) 00378 IJP = IJP + 1 00379 END DO 00380 JS = JS + LDA 00381 END DO 00382 DO I = 0, N1 00383 DO IJ = I, I + ( N1+I )*LDA, LDA 00384 ARF( IJ ) = AP( IJP ) 00385 IJP = IJP + 1 00386 END DO 00387 END DO 00388 * 00389 END IF 00390 * 00391 END IF 00392 * 00393 ELSE 00394 * 00395 * N is even 00396 * 00397 IF( NORMALTRANSR ) THEN 00398 * 00399 * N is even and TRANSR = 'N' 00400 * 00401 IF( LOWER ) THEN 00402 * 00403 * N is even, TRANSR = 'N', and UPLO = 'L' 00404 * 00405 IJP = 0 00406 JP = 0 00407 DO J = 0, K - 1 00408 DO I = J, N - 1 00409 IJ = 1 + I + JP 00410 ARF( IJ ) = AP( IJP ) 00411 IJP = IJP + 1 00412 END DO 00413 JP = JP + LDA 00414 END DO 00415 DO I = 0, K - 1 00416 DO J = I, K - 1 00417 IJ = I + J*LDA 00418 ARF( IJ ) = AP( IJP ) 00419 IJP = IJP + 1 00420 END DO 00421 END DO 00422 * 00423 ELSE 00424 * 00425 * N is even, TRANSR = 'N', and UPLO = 'U' 00426 * 00427 IJP = 0 00428 DO J = 0, K - 1 00429 IJ = K + 1 + J 00430 DO I = 0, J 00431 ARF( IJ ) = AP( IJP ) 00432 IJP = IJP + 1 00433 IJ = IJ + LDA 00434 END DO 00435 END DO 00436 JS = 0 00437 DO J = K, N - 1 00438 IJ = JS 00439 DO IJ = JS, JS + J 00440 ARF( IJ ) = AP( IJP ) 00441 IJP = IJP + 1 00442 END DO 00443 JS = JS + LDA 00444 END DO 00445 * 00446 END IF 00447 * 00448 ELSE 00449 * 00450 * N is even and TRANSR = 'T' 00451 * 00452 IF( LOWER ) THEN 00453 * 00454 * N is even, TRANSR = 'T', and UPLO = 'L' 00455 * 00456 IJP = 0 00457 DO I = 0, K - 1 00458 DO IJ = I + ( I+1 )*LDA, ( N+1 )*LDA - 1, LDA 00459 ARF( IJ ) = AP( IJP ) 00460 IJP = IJP + 1 00461 END DO 00462 END DO 00463 JS = 0 00464 DO J = 0, K - 1 00465 DO IJ = JS, JS + K - J - 1 00466 ARF( IJ ) = AP( IJP ) 00467 IJP = IJP + 1 00468 END DO 00469 JS = JS + LDA + 1 00470 END DO 00471 * 00472 ELSE 00473 * 00474 * N is even, TRANSR = 'T', and UPLO = 'U' 00475 * 00476 IJP = 0 00477 JS = ( K+1 )*LDA 00478 DO J = 0, K - 1 00479 DO IJ = JS, JS + J 00480 ARF( IJ ) = AP( IJP ) 00481 IJP = IJP + 1 00482 END DO 00483 JS = JS + LDA 00484 END DO 00485 DO I = 0, K - 1 00486 DO IJ = I, I + ( K+I )*LDA, LDA 00487 ARF( IJ ) = AP( IJP ) 00488 IJP = IJP + 1 00489 END DO 00490 END DO 00491 * 00492 END IF 00493 * 00494 END IF 00495 * 00496 END IF 00497 * 00498 RETURN 00499 * 00500 * End of STPTTF 00501 * 00502 END