![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DTFSM 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download DTFSM + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dtfsm.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dtfsm.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dtfsm.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE DTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A, 00022 * B, LDB ) 00023 * 00024 * .. Scalar Arguments .. 00025 * CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO 00026 * INTEGER LDB, M, N 00027 * DOUBLE PRECISION ALPHA 00028 * .. 00029 * .. Array Arguments .. 00030 * DOUBLE PRECISION A( 0: * ), B( 0: LDB-1, 0: * ) 00031 * .. 00032 * 00033 * 00034 *> \par Purpose: 00035 * ============= 00036 *> 00037 *> \verbatim 00038 *> 00039 *> Level 3 BLAS like routine for A in RFP Format. 00040 *> 00041 *> DTFSM solves the matrix equation 00042 *> 00043 *> op( A )*X = alpha*B or X*op( A ) = alpha*B 00044 *> 00045 *> where alpha is a scalar, X and B are m by n matrices, A is a unit, or 00046 *> non-unit, upper or lower triangular matrix and op( A ) is one of 00047 *> 00048 *> op( A ) = A or op( A ) = A**T. 00049 *> 00050 *> A is in Rectangular Full Packed (RFP) Format. 00051 *> 00052 *> The matrix X is overwritten on B. 00053 *> \endverbatim 00054 * 00055 * Arguments: 00056 * ========== 00057 * 00058 *> \param[in] TRANSR 00059 *> \verbatim 00060 *> TRANSR is CHARACTER*1 00061 *> = 'N': The Normal Form of RFP A is stored; 00062 *> = 'T': The Transpose Form of RFP A is stored. 00063 *> \endverbatim 00064 *> 00065 *> \param[in] SIDE 00066 *> \verbatim 00067 *> SIDE is CHARACTER*1 00068 *> On entry, SIDE specifies whether op( A ) appears on the left 00069 *> or right of X as follows: 00070 *> 00071 *> SIDE = 'L' or 'l' op( A )*X = alpha*B. 00072 *> 00073 *> SIDE = 'R' or 'r' X*op( A ) = alpha*B. 00074 *> 00075 *> Unchanged on exit. 00076 *> \endverbatim 00077 *> 00078 *> \param[in] UPLO 00079 *> \verbatim 00080 *> UPLO is CHARACTER*1 00081 *> On entry, UPLO specifies whether the RFP matrix A came from 00082 *> an upper or lower triangular matrix as follows: 00083 *> UPLO = 'U' or 'u' RFP A came from an upper triangular matrix 00084 *> UPLO = 'L' or 'l' RFP A came from a lower triangular matrix 00085 *> 00086 *> Unchanged on exit. 00087 *> \endverbatim 00088 *> 00089 *> \param[in] TRANS 00090 *> \verbatim 00091 *> TRANS is CHARACTER*1 00092 *> On entry, TRANS specifies the form of op( A ) to be used 00093 *> in the matrix multiplication as follows: 00094 *> 00095 *> TRANS = 'N' or 'n' op( A ) = A. 00096 *> 00097 *> TRANS = 'T' or 't' op( A ) = A'. 00098 *> 00099 *> Unchanged on exit. 00100 *> \endverbatim 00101 *> 00102 *> \param[in] DIAG 00103 *> \verbatim 00104 *> DIAG is CHARACTER*1 00105 *> On entry, DIAG specifies whether or not RFP A is unit 00106 *> triangular as follows: 00107 *> 00108 *> DIAG = 'U' or 'u' A is assumed to be unit triangular. 00109 *> 00110 *> DIAG = 'N' or 'n' A is not assumed to be unit 00111 *> triangular. 00112 *> 00113 *> Unchanged on exit. 00114 *> \endverbatim 00115 *> 00116 *> \param[in] M 00117 *> \verbatim 00118 *> M is INTEGER 00119 *> On entry, M specifies the number of rows of B. M must be at 00120 *> least zero. 00121 *> Unchanged on exit. 00122 *> \endverbatim 00123 *> 00124 *> \param[in] N 00125 *> \verbatim 00126 *> N is INTEGER 00127 *> On entry, N specifies the number of columns of B. N must be 00128 *> at least zero. 00129 *> Unchanged on exit. 00130 *> \endverbatim 00131 *> 00132 *> \param[in] ALPHA 00133 *> \verbatim 00134 *> ALPHA is DOUBLE PRECISION 00135 *> On entry, ALPHA specifies the scalar alpha. When alpha is 00136 *> zero then A is not referenced and B need not be set before 00137 *> entry. 00138 *> Unchanged on exit. 00139 *> \endverbatim 00140 *> 00141 *> \param[in] A 00142 *> \verbatim 00143 *> A is DOUBLE PRECISION array, dimension (NT) 00144 *> NT = N*(N+1)/2. On entry, the matrix A in RFP Format. 00145 *> RFP Format is described by TRANSR, UPLO and N as follows: 00146 *> If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even; 00147 *> K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If 00148 *> TRANSR = 'T' then RFP is the transpose of RFP A as 00149 *> defined when TRANSR = 'N'. The contents of RFP A are defined 00150 *> by UPLO as follows: If UPLO = 'U' the RFP A contains the NT 00151 *> elements of upper packed A either in normal or 00152 *> transpose Format. If UPLO = 'L' the RFP A contains 00153 *> the NT elements of lower packed A either in normal or 00154 *> transpose Format. The LDA of RFP A is (N+1)/2 when 00155 *> TRANSR = 'T'. When TRANSR is 'N' the LDA is N+1 when N is 00156 *> even and is N when is odd. 00157 *> See the Note below for more details. Unchanged on exit. 00158 *> \endverbatim 00159 *> 00160 *> \param[in,out] B 00161 *> \verbatim 00162 *> B is DOUBLE PRECISION array, dimension (LDB,N) 00163 *> Before entry, the leading m by n part of the array B must 00164 *> contain the right-hand side matrix B, and on exit is 00165 *> overwritten by the solution matrix X. 00166 *> \endverbatim 00167 *> 00168 *> \param[in] LDB 00169 *> \verbatim 00170 *> LDB is INTEGER 00171 *> On entry, LDB specifies the first dimension of B as declared 00172 *> in the calling (sub) program. LDB must be at least 00173 *> max( 1, m ). 00174 *> Unchanged on exit. 00175 *> \endverbatim 00176 * 00177 * Authors: 00178 * ======== 00179 * 00180 *> \author Univ. of Tennessee 00181 *> \author Univ. of California Berkeley 00182 *> \author Univ. of Colorado Denver 00183 *> \author NAG Ltd. 00184 * 00185 *> \date November 2011 00186 * 00187 *> \ingroup doubleOTHERcomputational 00188 * 00189 *> \par Further Details: 00190 * ===================== 00191 *> 00192 *> \verbatim 00193 *> 00194 *> We first consider Rectangular Full Packed (RFP) Format when N is 00195 *> even. We give an example where N = 6. 00196 *> 00197 *> AP is Upper AP is Lower 00198 *> 00199 *> 00 01 02 03 04 05 00 00200 *> 11 12 13 14 15 10 11 00201 *> 22 23 24 25 20 21 22 00202 *> 33 34 35 30 31 32 33 00203 *> 44 45 40 41 42 43 44 00204 *> 55 50 51 52 53 54 55 00205 *> 00206 *> 00207 *> Let TRANSR = 'N'. RFP holds AP as follows: 00208 *> For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last 00209 *> three columns of AP upper. The lower triangle A(4:6,0:2) consists of 00210 *> the transpose of the first three columns of AP upper. 00211 *> For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first 00212 *> three columns of AP lower. The upper triangle A(0:2,0:2) consists of 00213 *> the transpose of the last three columns of AP lower. 00214 *> This covers the case N even and TRANSR = 'N'. 00215 *> 00216 *> RFP A RFP A 00217 *> 00218 *> 03 04 05 33 43 53 00219 *> 13 14 15 00 44 54 00220 *> 23 24 25 10 11 55 00221 *> 33 34 35 20 21 22 00222 *> 00 44 45 30 31 32 00223 *> 01 11 55 40 41 42 00224 *> 02 12 22 50 51 52 00225 *> 00226 *> Now let TRANSR = 'T'. RFP A in both UPLO cases is just the 00227 *> transpose of RFP A above. One therefore gets: 00228 *> 00229 *> 00230 *> RFP A RFP A 00231 *> 00232 *> 03 13 23 33 00 01 02 33 00 10 20 30 40 50 00233 *> 04 14 24 34 44 11 12 43 44 11 21 31 41 51 00234 *> 05 15 25 35 45 55 22 53 54 55 22 32 42 52 00235 *> 00236 *> 00237 *> We then consider Rectangular Full Packed (RFP) Format when N is 00238 *> odd. We give an example where N = 5. 00239 *> 00240 *> AP is Upper AP is Lower 00241 *> 00242 *> 00 01 02 03 04 00 00243 *> 11 12 13 14 10 11 00244 *> 22 23 24 20 21 22 00245 *> 33 34 30 31 32 33 00246 *> 44 40 41 42 43 44 00247 *> 00248 *> 00249 *> Let TRANSR = 'N'. RFP holds AP as follows: 00250 *> For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last 00251 *> three columns of AP upper. The lower triangle A(3:4,0:1) consists of 00252 *> the transpose of the first two columns of AP upper. 00253 *> For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first 00254 *> three columns of AP lower. The upper triangle A(0:1,1:2) consists of 00255 *> the transpose of the last two columns of AP lower. 00256 *> This covers the case N odd and TRANSR = 'N'. 00257 *> 00258 *> RFP A RFP A 00259 *> 00260 *> 02 03 04 00 33 43 00261 *> 12 13 14 10 11 44 00262 *> 22 23 24 20 21 22 00263 *> 00 33 34 30 31 32 00264 *> 01 11 44 40 41 42 00265 *> 00266 *> Now let TRANSR = 'T'. RFP A in both UPLO cases is just the 00267 *> transpose of RFP A above. One therefore gets: 00268 *> 00269 *> RFP A RFP A 00270 *> 00271 *> 02 12 22 00 01 00 10 20 30 40 50 00272 *> 03 13 23 33 11 33 11 21 31 41 51 00273 *> 04 14 24 34 44 43 44 22 32 42 52 00274 *> \endverbatim 00275 * 00276 * ===================================================================== 00277 SUBROUTINE DTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A, 00278 $ B, LDB ) 00279 * 00280 * -- LAPACK computational routine (version 3.4.0) -- 00281 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00282 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00283 * November 2011 00284 * 00285 * .. Scalar Arguments .. 00286 CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO 00287 INTEGER LDB, M, N 00288 DOUBLE PRECISION ALPHA 00289 * .. 00290 * .. Array Arguments .. 00291 DOUBLE PRECISION A( 0: * ), B( 0: LDB-1, 0: * ) 00292 * .. 00293 * 00294 * ===================================================================== 00295 * 00296 * .. 00297 * .. Parameters .. 00298 DOUBLE PRECISION ONE, ZERO 00299 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00300 * .. 00301 * .. Local Scalars .. 00302 LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR, 00303 $ NOTRANS 00304 INTEGER M1, M2, N1, N2, K, INFO, I, J 00305 * .. 00306 * .. External Functions .. 00307 LOGICAL LSAME 00308 EXTERNAL LSAME 00309 * .. 00310 * .. External Subroutines .. 00311 EXTERNAL XERBLA, DGEMM, DTRSM 00312 * .. 00313 * .. Intrinsic Functions .. 00314 INTRINSIC MAX, MOD 00315 * .. 00316 * .. Executable Statements .. 00317 * 00318 * Test the input parameters. 00319 * 00320 INFO = 0 00321 NORMALTRANSR = LSAME( TRANSR, 'N' ) 00322 LSIDE = LSAME( SIDE, 'L' ) 00323 LOWER = LSAME( UPLO, 'L' ) 00324 NOTRANS = LSAME( TRANS, 'N' ) 00325 IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'T' ) ) THEN 00326 INFO = -1 00327 ELSE IF( .NOT.LSIDE .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN 00328 INFO = -2 00329 ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN 00330 INFO = -3 00331 ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'T' ) ) THEN 00332 INFO = -4 00333 ELSE IF( .NOT.LSAME( DIAG, 'N' ) .AND. .NOT.LSAME( DIAG, 'U' ) ) 00334 $ THEN 00335 INFO = -5 00336 ELSE IF( M.LT.0 ) THEN 00337 INFO = -6 00338 ELSE IF( N.LT.0 ) THEN 00339 INFO = -7 00340 ELSE IF( LDB.LT.MAX( 1, M ) ) THEN 00341 INFO = -11 00342 END IF 00343 IF( INFO.NE.0 ) THEN 00344 CALL XERBLA( 'DTFSM ', -INFO ) 00345 RETURN 00346 END IF 00347 * 00348 * Quick return when ( (N.EQ.0).OR.(M.EQ.0) ) 00349 * 00350 IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) ) 00351 $ RETURN 00352 * 00353 * Quick return when ALPHA.EQ.(0D+0) 00354 * 00355 IF( ALPHA.EQ.ZERO ) THEN 00356 DO 20 J = 0, N - 1 00357 DO 10 I = 0, M - 1 00358 B( I, J ) = ZERO 00359 10 CONTINUE 00360 20 CONTINUE 00361 RETURN 00362 END IF 00363 * 00364 IF( LSIDE ) THEN 00365 * 00366 * SIDE = 'L' 00367 * 00368 * A is M-by-M. 00369 * If M is odd, set NISODD = .TRUE., and M1 and M2. 00370 * If M is even, NISODD = .FALSE., and M. 00371 * 00372 IF( MOD( M, 2 ).EQ.0 ) THEN 00373 MISODD = .FALSE. 00374 K = M / 2 00375 ELSE 00376 MISODD = .TRUE. 00377 IF( LOWER ) THEN 00378 M2 = M / 2 00379 M1 = M - M2 00380 ELSE 00381 M1 = M / 2 00382 M2 = M - M1 00383 END IF 00384 END IF 00385 * 00386 * 00387 IF( MISODD ) THEN 00388 * 00389 * SIDE = 'L' and N is odd 00390 * 00391 IF( NORMALTRANSR ) THEN 00392 * 00393 * SIDE = 'L', N is odd, and TRANSR = 'N' 00394 * 00395 IF( LOWER ) THEN 00396 * 00397 * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'L' 00398 * 00399 IF( NOTRANS ) THEN 00400 * 00401 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and 00402 * TRANS = 'N' 00403 * 00404 IF( M.EQ.1 ) THEN 00405 CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA, 00406 $ A, M, B, LDB ) 00407 ELSE 00408 CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA, 00409 $ A( 0 ), M, B, LDB ) 00410 CALL DGEMM( 'N', 'N', M2, N, M1, -ONE, A( M1 ), 00411 $ M, B, LDB, ALPHA, B( M1, 0 ), LDB ) 00412 CALL DTRSM( 'L', 'U', 'T', DIAG, M2, N, ONE, 00413 $ A( M ), M, B( M1, 0 ), LDB ) 00414 END IF 00415 * 00416 ELSE 00417 * 00418 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and 00419 * TRANS = 'T' 00420 * 00421 IF( M.EQ.1 ) THEN 00422 CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ALPHA, 00423 $ A( 0 ), M, B, LDB ) 00424 ELSE 00425 CALL DTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA, 00426 $ A( M ), M, B( M1, 0 ), LDB ) 00427 CALL DGEMM( 'T', 'N', M1, N, M2, -ONE, A( M1 ), 00428 $ M, B( M1, 0 ), LDB, ALPHA, B, LDB ) 00429 CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ONE, 00430 $ A( 0 ), M, B, LDB ) 00431 END IF 00432 * 00433 END IF 00434 * 00435 ELSE 00436 * 00437 * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'U' 00438 * 00439 IF( .NOT.NOTRANS ) THEN 00440 * 00441 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and 00442 * TRANS = 'N' 00443 * 00444 CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA, 00445 $ A( M2 ), M, B, LDB ) 00446 CALL DGEMM( 'T', 'N', M2, N, M1, -ONE, A( 0 ), M, 00447 $ B, LDB, ALPHA, B( M1, 0 ), LDB ) 00448 CALL DTRSM( 'L', 'U', 'T', DIAG, M2, N, ONE, 00449 $ A( M1 ), M, B( M1, 0 ), LDB ) 00450 * 00451 ELSE 00452 * 00453 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and 00454 * TRANS = 'T' 00455 * 00456 CALL DTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA, 00457 $ A( M1 ), M, B( M1, 0 ), LDB ) 00458 CALL DGEMM( 'N', 'N', M1, N, M2, -ONE, A( 0 ), M, 00459 $ B( M1, 0 ), LDB, ALPHA, B, LDB ) 00460 CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ONE, 00461 $ A( M2 ), M, B, LDB ) 00462 * 00463 END IF 00464 * 00465 END IF 00466 * 00467 ELSE 00468 * 00469 * SIDE = 'L', N is odd, and TRANSR = 'T' 00470 * 00471 IF( LOWER ) THEN 00472 * 00473 * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'L' 00474 * 00475 IF( NOTRANS ) THEN 00476 * 00477 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and 00478 * TRANS = 'N' 00479 * 00480 IF( M.EQ.1 ) THEN 00481 CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA, 00482 $ A( 0 ), M1, B, LDB ) 00483 ELSE 00484 CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA, 00485 $ A( 0 ), M1, B, LDB ) 00486 CALL DGEMM( 'T', 'N', M2, N, M1, -ONE, 00487 $ A( M1*M1 ), M1, B, LDB, ALPHA, 00488 $ B( M1, 0 ), LDB ) 00489 CALL DTRSM( 'L', 'L', 'N', DIAG, M2, N, ONE, 00490 $ A( 1 ), M1, B( M1, 0 ), LDB ) 00491 END IF 00492 * 00493 ELSE 00494 * 00495 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and 00496 * TRANS = 'T' 00497 * 00498 IF( M.EQ.1 ) THEN 00499 CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ALPHA, 00500 $ A( 0 ), M1, B, LDB ) 00501 ELSE 00502 CALL DTRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA, 00503 $ A( 1 ), M1, B( M1, 0 ), LDB ) 00504 CALL DGEMM( 'N', 'N', M1, N, M2, -ONE, 00505 $ A( M1*M1 ), M1, B( M1, 0 ), LDB, 00506 $ ALPHA, B, LDB ) 00507 CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ONE, 00508 $ A( 0 ), M1, B, LDB ) 00509 END IF 00510 * 00511 END IF 00512 * 00513 ELSE 00514 * 00515 * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'U' 00516 * 00517 IF( .NOT.NOTRANS ) THEN 00518 * 00519 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and 00520 * TRANS = 'N' 00521 * 00522 CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA, 00523 $ A( M2*M2 ), M2, B, LDB ) 00524 CALL DGEMM( 'N', 'N', M2, N, M1, -ONE, A( 0 ), M2, 00525 $ B, LDB, ALPHA, B( M1, 0 ), LDB ) 00526 CALL DTRSM( 'L', 'L', 'N', DIAG, M2, N, ONE, 00527 $ A( M1*M2 ), M2, B( M1, 0 ), LDB ) 00528 * 00529 ELSE 00530 * 00531 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and 00532 * TRANS = 'T' 00533 * 00534 CALL DTRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA, 00535 $ A( M1*M2 ), M2, B( M1, 0 ), LDB ) 00536 CALL DGEMM( 'T', 'N', M1, N, M2, -ONE, A( 0 ), M2, 00537 $ B( M1, 0 ), LDB, ALPHA, B, LDB ) 00538 CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ONE, 00539 $ A( M2*M2 ), M2, B, LDB ) 00540 * 00541 END IF 00542 * 00543 END IF 00544 * 00545 END IF 00546 * 00547 ELSE 00548 * 00549 * SIDE = 'L' and N is even 00550 * 00551 IF( NORMALTRANSR ) THEN 00552 * 00553 * SIDE = 'L', N is even, and TRANSR = 'N' 00554 * 00555 IF( LOWER ) THEN 00556 * 00557 * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'L' 00558 * 00559 IF( NOTRANS ) THEN 00560 * 00561 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L', 00562 * and TRANS = 'N' 00563 * 00564 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA, 00565 $ A( 1 ), M+1, B, LDB ) 00566 CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( K+1 ), 00567 $ M+1, B, LDB, ALPHA, B( K, 0 ), LDB ) 00568 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ONE, 00569 $ A( 0 ), M+1, B( K, 0 ), LDB ) 00570 * 00571 ELSE 00572 * 00573 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L', 00574 * and TRANS = 'T' 00575 * 00576 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA, 00577 $ A( 0 ), M+1, B( K, 0 ), LDB ) 00578 CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( K+1 ), 00579 $ M+1, B( K, 0 ), LDB, ALPHA, B, LDB ) 00580 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ONE, 00581 $ A( 1 ), M+1, B, LDB ) 00582 * 00583 END IF 00584 * 00585 ELSE 00586 * 00587 * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'U' 00588 * 00589 IF( .NOT.NOTRANS ) THEN 00590 * 00591 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U', 00592 * and TRANS = 'N' 00593 * 00594 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA, 00595 $ A( K+1 ), M+1, B, LDB ) 00596 CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), M+1, 00597 $ B, LDB, ALPHA, B( K, 0 ), LDB ) 00598 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ONE, 00599 $ A( K ), M+1, B( K, 0 ), LDB ) 00600 * 00601 ELSE 00602 * 00603 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U', 00604 * and TRANS = 'T' 00605 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA, 00606 $ A( K ), M+1, B( K, 0 ), LDB ) 00607 CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), M+1, 00608 $ B( K, 0 ), LDB, ALPHA, B, LDB ) 00609 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ONE, 00610 $ A( K+1 ), M+1, B, LDB ) 00611 * 00612 END IF 00613 * 00614 END IF 00615 * 00616 ELSE 00617 * 00618 * SIDE = 'L', N is even, and TRANSR = 'T' 00619 * 00620 IF( LOWER ) THEN 00621 * 00622 * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'L' 00623 * 00624 IF( NOTRANS ) THEN 00625 * 00626 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L', 00627 * and TRANS = 'N' 00628 * 00629 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA, 00630 $ A( K ), K, B, LDB ) 00631 CALL DGEMM( 'T', 'N', K, N, K, -ONE, 00632 $ A( K*( K+1 ) ), K, B, LDB, ALPHA, 00633 $ B( K, 0 ), LDB ) 00634 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ONE, 00635 $ A( 0 ), K, B( K, 0 ), LDB ) 00636 * 00637 ELSE 00638 * 00639 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L', 00640 * and TRANS = 'T' 00641 * 00642 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA, 00643 $ A( 0 ), K, B( K, 0 ), LDB ) 00644 CALL DGEMM( 'N', 'N', K, N, K, -ONE, 00645 $ A( K*( K+1 ) ), K, B( K, 0 ), LDB, 00646 $ ALPHA, B, LDB ) 00647 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ONE, 00648 $ A( K ), K, B, LDB ) 00649 * 00650 END IF 00651 * 00652 ELSE 00653 * 00654 * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'U' 00655 * 00656 IF( .NOT.NOTRANS ) THEN 00657 * 00658 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U', 00659 * and TRANS = 'N' 00660 * 00661 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA, 00662 $ A( K*( K+1 ) ), K, B, LDB ) 00663 CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), K, B, 00664 $ LDB, ALPHA, B( K, 0 ), LDB ) 00665 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ONE, 00666 $ A( K*K ), K, B( K, 0 ), LDB ) 00667 * 00668 ELSE 00669 * 00670 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U', 00671 * and TRANS = 'T' 00672 * 00673 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA, 00674 $ A( K*K ), K, B( K, 0 ), LDB ) 00675 CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), K, 00676 $ B( K, 0 ), LDB, ALPHA, B, LDB ) 00677 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ONE, 00678 $ A( K*( K+1 ) ), K, B, LDB ) 00679 * 00680 END IF 00681 * 00682 END IF 00683 * 00684 END IF 00685 * 00686 END IF 00687 * 00688 ELSE 00689 * 00690 * SIDE = 'R' 00691 * 00692 * A is N-by-N. 00693 * If N is odd, set NISODD = .TRUE., and N1 and N2. 00694 * If N is even, NISODD = .FALSE., and K. 00695 * 00696 IF( MOD( N, 2 ).EQ.0 ) THEN 00697 NISODD = .FALSE. 00698 K = N / 2 00699 ELSE 00700 NISODD = .TRUE. 00701 IF( LOWER ) THEN 00702 N2 = N / 2 00703 N1 = N - N2 00704 ELSE 00705 N1 = N / 2 00706 N2 = N - N1 00707 END IF 00708 END IF 00709 * 00710 IF( NISODD ) THEN 00711 * 00712 * SIDE = 'R' and N is odd 00713 * 00714 IF( NORMALTRANSR ) THEN 00715 * 00716 * SIDE = 'R', N is odd, and TRANSR = 'N' 00717 * 00718 IF( LOWER ) THEN 00719 * 00720 * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'L' 00721 * 00722 IF( NOTRANS ) THEN 00723 * 00724 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and 00725 * TRANS = 'N' 00726 * 00727 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA, 00728 $ A( N ), N, B( 0, N1 ), LDB ) 00729 CALL DGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ), 00730 $ LDB, A( N1 ), N, ALPHA, B( 0, 0 ), 00731 $ LDB ) 00732 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N1, ONE, 00733 $ A( 0 ), N, B( 0, 0 ), LDB ) 00734 * 00735 ELSE 00736 * 00737 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and 00738 * TRANS = 'T' 00739 * 00740 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA, 00741 $ A( 0 ), N, B( 0, 0 ), LDB ) 00742 CALL DGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ), 00743 $ LDB, A( N1 ), N, ALPHA, B( 0, N1 ), 00744 $ LDB ) 00745 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N2, ONE, 00746 $ A( N ), N, B( 0, N1 ), LDB ) 00747 * 00748 END IF 00749 * 00750 ELSE 00751 * 00752 * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'U' 00753 * 00754 IF( NOTRANS ) THEN 00755 * 00756 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and 00757 * TRANS = 'N' 00758 * 00759 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA, 00760 $ A( N2 ), N, B( 0, 0 ), LDB ) 00761 CALL DGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ), 00762 $ LDB, A( 0 ), N, ALPHA, B( 0, N1 ), 00763 $ LDB ) 00764 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N2, ONE, 00765 $ A( N1 ), N, B( 0, N1 ), LDB ) 00766 * 00767 ELSE 00768 * 00769 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and 00770 * TRANS = 'T' 00771 * 00772 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA, 00773 $ A( N1 ), N, B( 0, N1 ), LDB ) 00774 CALL DGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ), 00775 $ LDB, A( 0 ), N, ALPHA, B( 0, 0 ), LDB ) 00776 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N1, ONE, 00777 $ A( N2 ), N, B( 0, 0 ), LDB ) 00778 * 00779 END IF 00780 * 00781 END IF 00782 * 00783 ELSE 00784 * 00785 * SIDE = 'R', N is odd, and TRANSR = 'T' 00786 * 00787 IF( LOWER ) THEN 00788 * 00789 * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'L' 00790 * 00791 IF( NOTRANS ) THEN 00792 * 00793 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and 00794 * TRANS = 'N' 00795 * 00796 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA, 00797 $ A( 1 ), N1, B( 0, N1 ), LDB ) 00798 CALL DGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ), 00799 $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, 0 ), 00800 $ LDB ) 00801 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N1, ONE, 00802 $ A( 0 ), N1, B( 0, 0 ), LDB ) 00803 * 00804 ELSE 00805 * 00806 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and 00807 * TRANS = 'T' 00808 * 00809 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA, 00810 $ A( 0 ), N1, B( 0, 0 ), LDB ) 00811 CALL DGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ), 00812 $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, N1 ), 00813 $ LDB ) 00814 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N2, ONE, 00815 $ A( 1 ), N1, B( 0, N1 ), LDB ) 00816 * 00817 END IF 00818 * 00819 ELSE 00820 * 00821 * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'U' 00822 * 00823 IF( NOTRANS ) THEN 00824 * 00825 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and 00826 * TRANS = 'N' 00827 * 00828 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA, 00829 $ A( N2*N2 ), N2, B( 0, 0 ), LDB ) 00830 CALL DGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ), 00831 $ LDB, A( 0 ), N2, ALPHA, B( 0, N1 ), 00832 $ LDB ) 00833 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N2, ONE, 00834 $ A( N1*N2 ), N2, B( 0, N1 ), LDB ) 00835 * 00836 ELSE 00837 * 00838 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and 00839 * TRANS = 'T' 00840 * 00841 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA, 00842 $ A( N1*N2 ), N2, B( 0, N1 ), LDB ) 00843 CALL DGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ), 00844 $ LDB, A( 0 ), N2, ALPHA, B( 0, 0 ), 00845 $ LDB ) 00846 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N1, ONE, 00847 $ A( N2*N2 ), N2, B( 0, 0 ), LDB ) 00848 * 00849 END IF 00850 * 00851 END IF 00852 * 00853 END IF 00854 * 00855 ELSE 00856 * 00857 * SIDE = 'R' and N is even 00858 * 00859 IF( NORMALTRANSR ) THEN 00860 * 00861 * SIDE = 'R', N is even, and TRANSR = 'N' 00862 * 00863 IF( LOWER ) THEN 00864 * 00865 * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'L' 00866 * 00867 IF( NOTRANS ) THEN 00868 * 00869 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L', 00870 * and TRANS = 'N' 00871 * 00872 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA, 00873 $ A( 0 ), N+1, B( 0, K ), LDB ) 00874 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ), 00875 $ LDB, A( K+1 ), N+1, ALPHA, B( 0, 0 ), 00876 $ LDB ) 00877 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ONE, 00878 $ A( 1 ), N+1, B( 0, 0 ), LDB ) 00879 * 00880 ELSE 00881 * 00882 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L', 00883 * and TRANS = 'T' 00884 * 00885 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA, 00886 $ A( 1 ), N+1, B( 0, 0 ), LDB ) 00887 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ), 00888 $ LDB, A( K+1 ), N+1, ALPHA, B( 0, K ), 00889 $ LDB ) 00890 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ONE, 00891 $ A( 0 ), N+1, B( 0, K ), LDB ) 00892 * 00893 END IF 00894 * 00895 ELSE 00896 * 00897 * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'U' 00898 * 00899 IF( NOTRANS ) THEN 00900 * 00901 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U', 00902 * and TRANS = 'N' 00903 * 00904 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA, 00905 $ A( K+1 ), N+1, B( 0, 0 ), LDB ) 00906 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ), 00907 $ LDB, A( 0 ), N+1, ALPHA, B( 0, K ), 00908 $ LDB ) 00909 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ONE, 00910 $ A( K ), N+1, B( 0, K ), LDB ) 00911 * 00912 ELSE 00913 * 00914 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U', 00915 * and TRANS = 'T' 00916 * 00917 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA, 00918 $ A( K ), N+1, B( 0, K ), LDB ) 00919 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ), 00920 $ LDB, A( 0 ), N+1, ALPHA, B( 0, 0 ), 00921 $ LDB ) 00922 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ONE, 00923 $ A( K+1 ), N+1, B( 0, 0 ), LDB ) 00924 * 00925 END IF 00926 * 00927 END IF 00928 * 00929 ELSE 00930 * 00931 * SIDE = 'R', N is even, and TRANSR = 'T' 00932 * 00933 IF( LOWER ) THEN 00934 * 00935 * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'L' 00936 * 00937 IF( NOTRANS ) THEN 00938 * 00939 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L', 00940 * and TRANS = 'N' 00941 * 00942 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA, 00943 $ A( 0 ), K, B( 0, K ), LDB ) 00944 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ), 00945 $ LDB, A( ( K+1 )*K ), K, ALPHA, 00946 $ B( 0, 0 ), LDB ) 00947 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ONE, 00948 $ A( K ), K, B( 0, 0 ), LDB ) 00949 * 00950 ELSE 00951 * 00952 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L', 00953 * and TRANS = 'T' 00954 * 00955 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA, 00956 $ A( K ), K, B( 0, 0 ), LDB ) 00957 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ), 00958 $ LDB, A( ( K+1 )*K ), K, ALPHA, 00959 $ B( 0, K ), LDB ) 00960 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ONE, 00961 $ A( 0 ), K, B( 0, K ), LDB ) 00962 * 00963 END IF 00964 * 00965 ELSE 00966 * 00967 * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'U' 00968 * 00969 IF( NOTRANS ) THEN 00970 * 00971 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U', 00972 * and TRANS = 'N' 00973 * 00974 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA, 00975 $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB ) 00976 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ), 00977 $ LDB, A( 0 ), K, ALPHA, B( 0, K ), LDB ) 00978 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ONE, 00979 $ A( K*K ), K, B( 0, K ), LDB ) 00980 * 00981 ELSE 00982 * 00983 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U', 00984 * and TRANS = 'T' 00985 * 00986 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA, 00987 $ A( K*K ), K, B( 0, K ), LDB ) 00988 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ), 00989 $ LDB, A( 0 ), K, ALPHA, B( 0, 0 ), LDB ) 00990 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ONE, 00991 $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB ) 00992 * 00993 END IF 00994 * 00995 END IF 00996 * 00997 END IF 00998 * 00999 END IF 01000 END IF 01001 * 01002 RETURN 01003 * 01004 * End of DTFSM 01005 * 01006 END