![]() ![]() | ![]() |
Mise à jour de lapack vers la version 3.3.0.
1: SUBROUTINE DTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A, 2: + B, LDB ) 3: * 4: * -- LAPACK routine (version 3.3.0) -- 5: * 6: * -- Contributed by Fred Gustavson of the IBM Watson Research Center -- 7: * November 2010 8: * 9: * -- LAPACK is a software package provided by Univ. of Tennessee, -- 10: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 11: * 12: * .. 13: * .. Scalar Arguments .. 14: CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO 15: INTEGER LDB, M, N 16: DOUBLE PRECISION ALPHA 17: * .. 18: * .. Array Arguments .. 19: DOUBLE PRECISION A( 0: * ), B( 0: LDB-1, 0: * ) 20: * .. 21: * 22: * Purpose 23: * ======= 24: * 25: * Level 3 BLAS like routine for A in RFP Format. 26: * 27: * DTFSM solves the matrix equation 28: * 29: * op( A )*X = alpha*B or X*op( A ) = alpha*B 30: * 31: * where alpha is a scalar, X and B are m by n matrices, A is a unit, or 32: * non-unit, upper or lower triangular matrix and op( A ) is one of 33: * 34: * op( A ) = A or op( A ) = A'. 35: * 36: * A is in Rectangular Full Packed (RFP) Format. 37: * 38: * The matrix X is overwritten on B. 39: * 40: * Arguments 41: * ========== 42: * 43: * TRANSR (input) CHARACTER*1 44: * = 'N': The Normal Form of RFP A is stored; 45: * = 'T': The Transpose Form of RFP A is stored. 46: * 47: * SIDE (input) CHARACTER*1 48: * On entry, SIDE specifies whether op( A ) appears on the left 49: * or right of X as follows: 50: * 51: * SIDE = 'L' or 'l' op( A )*X = alpha*B. 52: * 53: * SIDE = 'R' or 'r' X*op( A ) = alpha*B. 54: * 55: * Unchanged on exit. 56: * 57: * UPLO (input) CHARACTER*1 58: * On entry, UPLO specifies whether the RFP matrix A came from 59: * an upper or lower triangular matrix as follows: 60: * UPLO = 'U' or 'u' RFP A came from an upper triangular matrix 61: * UPLO = 'L' or 'l' RFP A came from a lower triangular matrix 62: * 63: * Unchanged on exit. 64: * 65: * TRANS (input) CHARACTER*1 66: * On entry, TRANS specifies the form of op( A ) to be used 67: * in the matrix multiplication as follows: 68: * 69: * TRANS = 'N' or 'n' op( A ) = A. 70: * 71: * TRANS = 'T' or 't' op( A ) = A'. 72: * 73: * Unchanged on exit. 74: * 75: * DIAG (input) CHARACTER*1 76: * On entry, DIAG specifies whether or not RFP A is unit 77: * triangular as follows: 78: * 79: * DIAG = 'U' or 'u' A is assumed to be unit triangular. 80: * 81: * DIAG = 'N' or 'n' A is not assumed to be unit 82: * triangular. 83: * 84: * Unchanged on exit. 85: * 86: * M (input) INTEGER 87: * On entry, M specifies the number of rows of B. M must be at 88: * least zero. 89: * Unchanged on exit. 90: * 91: * N (input) INTEGER 92: * On entry, N specifies the number of columns of B. N must be 93: * at least zero. 94: * Unchanged on exit. 95: * 96: * ALPHA (input) DOUBLE PRECISION 97: * On entry, ALPHA specifies the scalar alpha. When alpha is 98: * zero then A is not referenced and B need not be set before 99: * entry. 100: * Unchanged on exit. 101: * 102: * A (input) DOUBLE PRECISION array, dimension (NT) 103: * NT = N*(N+1)/2. On entry, the matrix A in RFP Format. 104: * RFP Format is described by TRANSR, UPLO and N as follows: 105: * If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even; 106: * K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If 107: * TRANSR = 'T' then RFP is the transpose of RFP A as 108: * defined when TRANSR = 'N'. The contents of RFP A are defined 109: * by UPLO as follows: If UPLO = 'U' the RFP A contains the NT 110: * elements of upper packed A either in normal or 111: * transpose Format. If UPLO = 'L' the RFP A contains 112: * the NT elements of lower packed A either in normal or 113: * transpose Format. The LDA of RFP A is (N+1)/2 when 114: * TRANSR = 'T'. When TRANSR is 'N' the LDA is N+1 when N is 115: * even and is N when is odd. 116: * See the Note below for more details. Unchanged on exit. 117: * 118: * B (input/output) DOUBLE PRECISION array, dimension (LDB,N) 119: * Before entry, the leading m by n part of the array B must 120: * contain the right-hand side matrix B, and on exit is 121: * overwritten by the solution matrix X. 122: * 123: * LDB (input) INTEGER 124: * On entry, LDB specifies the first dimension of B as declared 125: * in the calling (sub) program. LDB must be at least 126: * max( 1, m ). 127: * Unchanged on exit. 128: * 129: * Further Details 130: * =============== 131: * 132: * We first consider Rectangular Full Packed (RFP) Format when N is 133: * even. We give an example where N = 6. 134: * 135: * AP is Upper AP is Lower 136: * 137: * 00 01 02 03 04 05 00 138: * 11 12 13 14 15 10 11 139: * 22 23 24 25 20 21 22 140: * 33 34 35 30 31 32 33 141: * 44 45 40 41 42 43 44 142: * 55 50 51 52 53 54 55 143: * 144: * 145: * Let TRANSR = 'N'. RFP holds AP as follows: 146: * For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last 147: * three columns of AP upper. The lower triangle A(4:6,0:2) consists of 148: * the transpose of the first three columns of AP upper. 149: * For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first 150: * three columns of AP lower. The upper triangle A(0:2,0:2) consists of 151: * the transpose of the last three columns of AP lower. 152: * This covers the case N even and TRANSR = 'N'. 153: * 154: * RFP A RFP A 155: * 156: * 03 04 05 33 43 53 157: * 13 14 15 00 44 54 158: * 23 24 25 10 11 55 159: * 33 34 35 20 21 22 160: * 00 44 45 30 31 32 161: * 01 11 55 40 41 42 162: * 02 12 22 50 51 52 163: * 164: * Now let TRANSR = 'T'. RFP A in both UPLO cases is just the 165: * transpose of RFP A above. One therefore gets: 166: * 167: * 168: * RFP A RFP A 169: * 170: * 03 13 23 33 00 01 02 33 00 10 20 30 40 50 171: * 04 14 24 34 44 11 12 43 44 11 21 31 41 51 172: * 05 15 25 35 45 55 22 53 54 55 22 32 42 52 173: * 174: * 175: * We then consider Rectangular Full Packed (RFP) Format when N is 176: * odd. We give an example where N = 5. 177: * 178: * AP is Upper AP is Lower 179: * 180: * 00 01 02 03 04 00 181: * 11 12 13 14 10 11 182: * 22 23 24 20 21 22 183: * 33 34 30 31 32 33 184: * 44 40 41 42 43 44 185: * 186: * 187: * Let TRANSR = 'N'. RFP holds AP as follows: 188: * For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last 189: * three columns of AP upper. The lower triangle A(3:4,0:1) consists of 190: * the transpose of the first two columns of AP upper. 191: * For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first 192: * three columns of AP lower. The upper triangle A(0:1,1:2) consists of 193: * the transpose of the last two columns of AP lower. 194: * This covers the case N odd and TRANSR = 'N'. 195: * 196: * RFP A RFP A 197: * 198: * 02 03 04 00 33 43 199: * 12 13 14 10 11 44 200: * 22 23 24 20 21 22 201: * 00 33 34 30 31 32 202: * 01 11 44 40 41 42 203: * 204: * Now let TRANSR = 'T'. RFP A in both UPLO cases is just the 205: * transpose of RFP A above. One therefore gets: 206: * 207: * RFP A RFP A 208: * 209: * 02 12 22 00 01 00 10 20 30 40 50 210: * 03 13 23 33 11 33 11 21 31 41 51 211: * 04 14 24 34 44 43 44 22 32 42 52 212: * 213: * Reference 214: * ========= 215: * 216: * ===================================================================== 217: * 218: * .. 219: * .. Parameters .. 220: DOUBLE PRECISION ONE, ZERO 221: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 222: * .. 223: * .. Local Scalars .. 224: LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR, 225: + NOTRANS 226: INTEGER M1, M2, N1, N2, K, INFO, I, J 227: * .. 228: * .. External Functions .. 229: LOGICAL LSAME 230: EXTERNAL LSAME 231: * .. 232: * .. External Subroutines .. 233: EXTERNAL XERBLA, DGEMM, DTRSM 234: * .. 235: * .. Intrinsic Functions .. 236: INTRINSIC MAX, MOD 237: * .. 238: * .. Executable Statements .. 239: * 240: * Test the input parameters. 241: * 242: INFO = 0 243: NORMALTRANSR = LSAME( TRANSR, 'N' ) 244: LSIDE = LSAME( SIDE, 'L' ) 245: LOWER = LSAME( UPLO, 'L' ) 246: NOTRANS = LSAME( TRANS, 'N' ) 247: IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'T' ) ) THEN 248: INFO = -1 249: ELSE IF( .NOT.LSIDE .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN 250: INFO = -2 251: ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN 252: INFO = -3 253: ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'T' ) ) THEN 254: INFO = -4 255: ELSE IF( .NOT.LSAME( DIAG, 'N' ) .AND. .NOT.LSAME( DIAG, 'U' ) ) 256: + THEN 257: INFO = -5 258: ELSE IF( M.LT.0 ) THEN 259: INFO = -6 260: ELSE IF( N.LT.0 ) THEN 261: INFO = -7 262: ELSE IF( LDB.LT.MAX( 1, M ) ) THEN 263: INFO = -11 264: END IF 265: IF( INFO.NE.0 ) THEN 266: CALL XERBLA( 'DTFSM ', -INFO ) 267: RETURN 268: END IF 269: * 270: * Quick return when ( (N.EQ.0).OR.(M.EQ.0) ) 271: * 272: IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) ) 273: + RETURN 274: * 275: * Quick return when ALPHA.EQ.(0D+0) 276: * 277: IF( ALPHA.EQ.ZERO ) THEN 278: DO 20 J = 0, N - 1 279: DO 10 I = 0, M - 1 280: B( I, J ) = ZERO 281: 10 CONTINUE 282: 20 CONTINUE 283: RETURN 284: END IF 285: * 286: IF( LSIDE ) THEN 287: * 288: * SIDE = 'L' 289: * 290: * A is M-by-M. 291: * If M is odd, set NISODD = .TRUE., and M1 and M2. 292: * If M is even, NISODD = .FALSE., and M. 293: * 294: IF( MOD( M, 2 ).EQ.0 ) THEN 295: MISODD = .FALSE. 296: K = M / 2 297: ELSE 298: MISODD = .TRUE. 299: IF( LOWER ) THEN 300: M2 = M / 2 301: M1 = M - M2 302: ELSE 303: M1 = M / 2 304: M2 = M - M1 305: END IF 306: END IF 307: * 308: * 309: IF( MISODD ) THEN 310: * 311: * SIDE = 'L' and N is odd 312: * 313: IF( NORMALTRANSR ) THEN 314: * 315: * SIDE = 'L', N is odd, and TRANSR = 'N' 316: * 317: IF( LOWER ) THEN 318: * 319: * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'L' 320: * 321: IF( NOTRANS ) THEN 322: * 323: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and 324: * TRANS = 'N' 325: * 326: IF( M.EQ.1 ) THEN 327: CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA, 328: + A, M, B, LDB ) 329: ELSE 330: CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA, 331: + A( 0 ), M, B, LDB ) 332: CALL DGEMM( 'N', 'N', M2, N, M1, -ONE, A( M1 ), 333: + M, B, LDB, ALPHA, B( M1, 0 ), LDB ) 334: CALL DTRSM( 'L', 'U', 'T', DIAG, M2, N, ONE, 335: + A( M ), M, B( M1, 0 ), LDB ) 336: END IF 337: * 338: ELSE 339: * 340: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and 341: * TRANS = 'T' 342: * 343: IF( M.EQ.1 ) THEN 344: CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ALPHA, 345: + A( 0 ), M, B, LDB ) 346: ELSE 347: CALL DTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA, 348: + A( M ), M, B( M1, 0 ), LDB ) 349: CALL DGEMM( 'T', 'N', M1, N, M2, -ONE, A( M1 ), 350: + M, B( M1, 0 ), LDB, ALPHA, B, LDB ) 351: CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ONE, 352: + A( 0 ), M, B, LDB ) 353: END IF 354: * 355: END IF 356: * 357: ELSE 358: * 359: * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'U' 360: * 361: IF( .NOT.NOTRANS ) THEN 362: * 363: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and 364: * TRANS = 'N' 365: * 366: CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA, 367: + A( M2 ), M, B, LDB ) 368: CALL DGEMM( 'T', 'N', M2, N, M1, -ONE, A( 0 ), M, 369: + B, LDB, ALPHA, B( M1, 0 ), LDB ) 370: CALL DTRSM( 'L', 'U', 'T', DIAG, M2, N, ONE, 371: + A( M1 ), M, B( M1, 0 ), LDB ) 372: * 373: ELSE 374: * 375: * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and 376: * TRANS = 'T' 377: * 378: CALL DTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA, 379: + A( M1 ), M, B( M1, 0 ), LDB ) 380: CALL DGEMM( 'N', 'N', M1, N, M2, -ONE, A( 0 ), M, 381: + B( M1, 0 ), LDB, ALPHA, B, LDB ) 382: CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ONE, 383: + A( M2 ), M, B, LDB ) 384: * 385: END IF 386: * 387: END IF 388: * 389: ELSE 390: * 391: * SIDE = 'L', N is odd, and TRANSR = 'T' 392: * 393: IF( LOWER ) THEN 394: * 395: * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'L' 396: * 397: IF( NOTRANS ) THEN 398: * 399: * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and 400: * TRANS = 'N' 401: * 402: IF( M.EQ.1 ) THEN 403: CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA, 404: + A( 0 ), M1, B, LDB ) 405: ELSE 406: CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA, 407: + A( 0 ), M1, B, LDB ) 408: CALL DGEMM( 'T', 'N', M2, N, M1, -ONE, 409: + A( M1*M1 ), M1, B, LDB, ALPHA, 410: + B( M1, 0 ), LDB ) 411: CALL DTRSM( 'L', 'L', 'N', DIAG, M2, N, ONE, 412: + A( 1 ), M1, B( M1, 0 ), LDB ) 413: END IF 414: * 415: ELSE 416: * 417: * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and 418: * TRANS = 'T' 419: * 420: IF( M.EQ.1 ) THEN 421: CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ALPHA, 422: + A( 0 ), M1, B, LDB ) 423: ELSE 424: CALL DTRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA, 425: + A( 1 ), M1, B( M1, 0 ), LDB ) 426: CALL DGEMM( 'N', 'N', M1, N, M2, -ONE, 427: + A( M1*M1 ), M1, B( M1, 0 ), LDB, 428: + ALPHA, B, LDB ) 429: CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ONE, 430: + A( 0 ), M1, B, LDB ) 431: END IF 432: * 433: END IF 434: * 435: ELSE 436: * 437: * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'U' 438: * 439: IF( .NOT.NOTRANS ) THEN 440: * 441: * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and 442: * TRANS = 'N' 443: * 444: CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA, 445: + A( M2*M2 ), M2, B, LDB ) 446: CALL DGEMM( 'N', 'N', M2, N, M1, -ONE, A( 0 ), M2, 447: + B, LDB, ALPHA, B( M1, 0 ), LDB ) 448: CALL DTRSM( 'L', 'L', 'N', DIAG, M2, N, ONE, 449: + A( M1*M2 ), M2, B( M1, 0 ), LDB ) 450: * 451: ELSE 452: * 453: * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and 454: * TRANS = 'T' 455: * 456: CALL DTRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA, 457: + A( M1*M2 ), M2, B( M1, 0 ), LDB ) 458: CALL DGEMM( 'T', 'N', M1, N, M2, -ONE, A( 0 ), M2, 459: + B( M1, 0 ), LDB, ALPHA, B, LDB ) 460: CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ONE, 461: + A( M2*M2 ), M2, B, LDB ) 462: * 463: END IF 464: * 465: END IF 466: * 467: END IF 468: * 469: ELSE 470: * 471: * SIDE = 'L' and N is even 472: * 473: IF( NORMALTRANSR ) THEN 474: * 475: * SIDE = 'L', N is even, and TRANSR = 'N' 476: * 477: IF( LOWER ) THEN 478: * 479: * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'L' 480: * 481: IF( NOTRANS ) THEN 482: * 483: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L', 484: * and TRANS = 'N' 485: * 486: CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA, 487: + A( 1 ), M+1, B, LDB ) 488: CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( K+1 ), 489: + M+1, B, LDB, ALPHA, B( K, 0 ), LDB ) 490: CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ONE, 491: + A( 0 ), M+1, B( K, 0 ), LDB ) 492: * 493: ELSE 494: * 495: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L', 496: * and TRANS = 'T' 497: * 498: CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA, 499: + A( 0 ), M+1, B( K, 0 ), LDB ) 500: CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( K+1 ), 501: + M+1, B( K, 0 ), LDB, ALPHA, B, LDB ) 502: CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ONE, 503: + A( 1 ), M+1, B, LDB ) 504: * 505: END IF 506: * 507: ELSE 508: * 509: * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'U' 510: * 511: IF( .NOT.NOTRANS ) THEN 512: * 513: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U', 514: * and TRANS = 'N' 515: * 516: CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA, 517: + A( K+1 ), M+1, B, LDB ) 518: CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), M+1, 519: + B, LDB, ALPHA, B( K, 0 ), LDB ) 520: CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ONE, 521: + A( K ), M+1, B( K, 0 ), LDB ) 522: * 523: ELSE 524: * 525: * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U', 526: * and TRANS = 'T' 527: CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA, 528: + A( K ), M+1, B( K, 0 ), LDB ) 529: CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), M+1, 530: + B( K, 0 ), LDB, ALPHA, B, LDB ) 531: CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ONE, 532: + A( K+1 ), M+1, B, LDB ) 533: * 534: END IF 535: * 536: END IF 537: * 538: ELSE 539: * 540: * SIDE = 'L', N is even, and TRANSR = 'T' 541: * 542: IF( LOWER ) THEN 543: * 544: * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'L' 545: * 546: IF( NOTRANS ) THEN 547: * 548: * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L', 549: * and TRANS = 'N' 550: * 551: CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA, 552: + A( K ), K, B, LDB ) 553: CALL DGEMM( 'T', 'N', K, N, K, -ONE, 554: + A( K*( K+1 ) ), K, B, LDB, ALPHA, 555: + B( K, 0 ), LDB ) 556: CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ONE, 557: + A( 0 ), K, B( K, 0 ), LDB ) 558: * 559: ELSE 560: * 561: * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L', 562: * and TRANS = 'T' 563: * 564: CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA, 565: + A( 0 ), K, B( K, 0 ), LDB ) 566: CALL DGEMM( 'N', 'N', K, N, K, -ONE, 567: + A( K*( K+1 ) ), K, B( K, 0 ), LDB, 568: + ALPHA, B, LDB ) 569: CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ONE, 570: + A( K ), K, B, LDB ) 571: * 572: END IF 573: * 574: ELSE 575: * 576: * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'U' 577: * 578: IF( .NOT.NOTRANS ) THEN 579: * 580: * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U', 581: * and TRANS = 'N' 582: * 583: CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA, 584: + A( K*( K+1 ) ), K, B, LDB ) 585: CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), K, B, 586: + LDB, ALPHA, B( K, 0 ), LDB ) 587: CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ONE, 588: + A( K*K ), K, B( K, 0 ), LDB ) 589: * 590: ELSE 591: * 592: * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U', 593: * and TRANS = 'T' 594: * 595: CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA, 596: + A( K*K ), K, B( K, 0 ), LDB ) 597: CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), K, 598: + B( K, 0 ), LDB, ALPHA, B, LDB ) 599: CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ONE, 600: + A( K*( K+1 ) ), K, B, LDB ) 601: * 602: END IF 603: * 604: END IF 605: * 606: END IF 607: * 608: END IF 609: * 610: ELSE 611: * 612: * SIDE = 'R' 613: * 614: * A is N-by-N. 615: * If N is odd, set NISODD = .TRUE., and N1 and N2. 616: * If N is even, NISODD = .FALSE., and K. 617: * 618: IF( MOD( N, 2 ).EQ.0 ) THEN 619: NISODD = .FALSE. 620: K = N / 2 621: ELSE 622: NISODD = .TRUE. 623: IF( LOWER ) THEN 624: N2 = N / 2 625: N1 = N - N2 626: ELSE 627: N1 = N / 2 628: N2 = N - N1 629: END IF 630: END IF 631: * 632: IF( NISODD ) THEN 633: * 634: * SIDE = 'R' and N is odd 635: * 636: IF( NORMALTRANSR ) THEN 637: * 638: * SIDE = 'R', N is odd, and TRANSR = 'N' 639: * 640: IF( LOWER ) THEN 641: * 642: * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'L' 643: * 644: IF( NOTRANS ) THEN 645: * 646: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and 647: * TRANS = 'N' 648: * 649: CALL DTRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA, 650: + A( N ), N, B( 0, N1 ), LDB ) 651: CALL DGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ), 652: + LDB, A( N1 ), N, ALPHA, B( 0, 0 ), 653: + LDB ) 654: CALL DTRSM( 'R', 'L', 'N', DIAG, M, N1, ONE, 655: + A( 0 ), N, B( 0, 0 ), LDB ) 656: * 657: ELSE 658: * 659: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and 660: * TRANS = 'T' 661: * 662: CALL DTRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA, 663: + A( 0 ), N, B( 0, 0 ), LDB ) 664: CALL DGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ), 665: + LDB, A( N1 ), N, ALPHA, B( 0, N1 ), 666: + LDB ) 667: CALL DTRSM( 'R', 'U', 'N', DIAG, M, N2, ONE, 668: + A( N ), N, B( 0, N1 ), LDB ) 669: * 670: END IF 671: * 672: ELSE 673: * 674: * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'U' 675: * 676: IF( NOTRANS ) THEN 677: * 678: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and 679: * TRANS = 'N' 680: * 681: CALL DTRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA, 682: + A( N2 ), N, B( 0, 0 ), LDB ) 683: CALL DGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ), 684: + LDB, A( 0 ), N, ALPHA, B( 0, N1 ), 685: + LDB ) 686: CALL DTRSM( 'R', 'U', 'N', DIAG, M, N2, ONE, 687: + A( N1 ), N, B( 0, N1 ), LDB ) 688: * 689: ELSE 690: * 691: * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and 692: * TRANS = 'T' 693: * 694: CALL DTRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA, 695: + A( N1 ), N, B( 0, N1 ), LDB ) 696: CALL DGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ), 697: + LDB, A( 0 ), N, ALPHA, B( 0, 0 ), LDB ) 698: CALL DTRSM( 'R', 'L', 'N', DIAG, M, N1, ONE, 699: + A( N2 ), N, B( 0, 0 ), LDB ) 700: * 701: END IF 702: * 703: END IF 704: * 705: ELSE 706: * 707: * SIDE = 'R', N is odd, and TRANSR = 'T' 708: * 709: IF( LOWER ) THEN 710: * 711: * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'L' 712: * 713: IF( NOTRANS ) THEN 714: * 715: * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and 716: * TRANS = 'N' 717: * 718: CALL DTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA, 719: + A( 1 ), N1, B( 0, N1 ), LDB ) 720: CALL DGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ), 721: + LDB, A( N1*N1 ), N1, ALPHA, B( 0, 0 ), 722: + LDB ) 723: CALL DTRSM( 'R', 'U', 'T', DIAG, M, N1, ONE, 724: + A( 0 ), N1, B( 0, 0 ), LDB ) 725: * 726: ELSE 727: * 728: * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and 729: * TRANS = 'T' 730: * 731: CALL DTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA, 732: + A( 0 ), N1, B( 0, 0 ), LDB ) 733: CALL DGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ), 734: + LDB, A( N1*N1 ), N1, ALPHA, B( 0, N1 ), 735: + LDB ) 736: CALL DTRSM( 'R', 'L', 'T', DIAG, M, N2, ONE, 737: + A( 1 ), N1, B( 0, N1 ), LDB ) 738: * 739: END IF 740: * 741: ELSE 742: * 743: * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'U' 744: * 745: IF( NOTRANS ) THEN 746: * 747: * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and 748: * TRANS = 'N' 749: * 750: CALL DTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA, 751: + A( N2*N2 ), N2, B( 0, 0 ), LDB ) 752: CALL DGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ), 753: + LDB, A( 0 ), N2, ALPHA, B( 0, N1 ), 754: + LDB ) 755: CALL DTRSM( 'R', 'L', 'T', DIAG, M, N2, ONE, 756: + A( N1*N2 ), N2, B( 0, N1 ), LDB ) 757: * 758: ELSE 759: * 760: * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and 761: * TRANS = 'T' 762: * 763: CALL DTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA, 764: + A( N1*N2 ), N2, B( 0, N1 ), LDB ) 765: CALL DGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ), 766: + LDB, A( 0 ), N2, ALPHA, B( 0, 0 ), 767: + LDB ) 768: CALL DTRSM( 'R', 'U', 'T', DIAG, M, N1, ONE, 769: + A( N2*N2 ), N2, B( 0, 0 ), LDB ) 770: * 771: END IF 772: * 773: END IF 774: * 775: END IF 776: * 777: ELSE 778: * 779: * SIDE = 'R' and N is even 780: * 781: IF( NORMALTRANSR ) THEN 782: * 783: * SIDE = 'R', N is even, and TRANSR = 'N' 784: * 785: IF( LOWER ) THEN 786: * 787: * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'L' 788: * 789: IF( NOTRANS ) THEN 790: * 791: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L', 792: * and TRANS = 'N' 793: * 794: CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA, 795: + A( 0 ), N+1, B( 0, K ), LDB ) 796: CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ), 797: + LDB, A( K+1 ), N+1, ALPHA, B( 0, 0 ), 798: + LDB ) 799: CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ONE, 800: + A( 1 ), N+1, B( 0, 0 ), LDB ) 801: * 802: ELSE 803: * 804: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L', 805: * and TRANS = 'T' 806: * 807: CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA, 808: + A( 1 ), N+1, B( 0, 0 ), LDB ) 809: CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ), 810: + LDB, A( K+1 ), N+1, ALPHA, B( 0, K ), 811: + LDB ) 812: CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ONE, 813: + A( 0 ), N+1, B( 0, K ), LDB ) 814: * 815: END IF 816: * 817: ELSE 818: * 819: * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'U' 820: * 821: IF( NOTRANS ) THEN 822: * 823: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U', 824: * and TRANS = 'N' 825: * 826: CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA, 827: + A( K+1 ), N+1, B( 0, 0 ), LDB ) 828: CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ), 829: + LDB, A( 0 ), N+1, ALPHA, B( 0, K ), 830: + LDB ) 831: CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ONE, 832: + A( K ), N+1, B( 0, K ), LDB ) 833: * 834: ELSE 835: * 836: * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U', 837: * and TRANS = 'T' 838: * 839: CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA, 840: + A( K ), N+1, B( 0, K ), LDB ) 841: CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ), 842: + LDB, A( 0 ), N+1, ALPHA, B( 0, 0 ), 843: + LDB ) 844: CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ONE, 845: + A( K+1 ), N+1, B( 0, 0 ), LDB ) 846: * 847: END IF 848: * 849: END IF 850: * 851: ELSE 852: * 853: * SIDE = 'R', N is even, and TRANSR = 'T' 854: * 855: IF( LOWER ) THEN 856: * 857: * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'L' 858: * 859: IF( NOTRANS ) THEN 860: * 861: * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L', 862: * and TRANS = 'N' 863: * 864: CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA, 865: + A( 0 ), K, B( 0, K ), LDB ) 866: CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ), 867: + LDB, A( ( K+1 )*K ), K, ALPHA, 868: + B( 0, 0 ), LDB ) 869: CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ONE, 870: + A( K ), K, B( 0, 0 ), LDB ) 871: * 872: ELSE 873: * 874: * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L', 875: * and TRANS = 'T' 876: * 877: CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA, 878: + A( K ), K, B( 0, 0 ), LDB ) 879: CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ), 880: + LDB, A( ( K+1 )*K ), K, ALPHA, 881: + B( 0, K ), LDB ) 882: CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ONE, 883: + A( 0 ), K, B( 0, K ), LDB ) 884: * 885: END IF 886: * 887: ELSE 888: * 889: * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'U' 890: * 891: IF( NOTRANS ) THEN 892: * 893: * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U', 894: * and TRANS = 'N' 895: * 896: CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA, 897: + A( ( K+1 )*K ), K, B( 0, 0 ), LDB ) 898: CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ), 899: + LDB, A( 0 ), K, ALPHA, B( 0, K ), LDB ) 900: CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ONE, 901: + A( K*K ), K, B( 0, K ), LDB ) 902: * 903: ELSE 904: * 905: * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U', 906: * and TRANS = 'T' 907: * 908: CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA, 909: + A( K*K ), K, B( 0, K ), LDB ) 910: CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ), 911: + LDB, A( 0 ), K, ALPHA, B( 0, 0 ), LDB ) 912: CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ONE, 913: + A( ( K+1 )*K ), K, B( 0, 0 ), LDB ) 914: * 915: END IF 916: * 917: END IF 918: * 919: END IF 920: * 921: END IF 922: END IF 923: * 924: RETURN 925: * 926: * End of DTFSM 927: * 928: END