Annotation of rpl/lapack/lapack/zhfrk.f, revision 1.1
1.1 ! bertrand 1: SUBROUTINE ZHFRK( TRANSR, UPLO, TRANS, N, K, ALPHA, A, LDA, BETA,
! 2: + C )
! 3: *
! 4: * -- LAPACK routine (version 3.2.2) --
! 5: *
! 6: * -- Contributed by Julien Langou of the Univ. of Colorado Denver --
! 7: * -- June 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: DOUBLE PRECISION ALPHA, BETA
! 15: INTEGER K, LDA, N
! 16: CHARACTER TRANS, TRANSR, UPLO
! 17: * ..
! 18: * .. Array Arguments ..
! 19: COMPLEX*16 A( LDA, * ), C( * )
! 20: * ..
! 21: *
! 22: * Purpose
! 23: * =======
! 24: *
! 25: * Level 3 BLAS like routine for C in RFP Format.
! 26: *
! 27: * ZHFRK performs one of the Hermitian rank--k operations
! 28: *
! 29: * C := alpha*A*conjg( A' ) + beta*C,
! 30: *
! 31: * or
! 32: *
! 33: * C := alpha*conjg( A' )*A + beta*C,
! 34: *
! 35: * where alpha and beta are real scalars, C is an n--by--n Hermitian
! 36: * matrix and A is an n--by--k matrix in the first case and a k--by--n
! 37: * matrix in the second case.
! 38: *
! 39: * Arguments
! 40: * ==========
! 41: *
! 42: * TRANSR (input) CHARACTER
! 43: * = 'N': The Normal Form of RFP A is stored;
! 44: * = 'C': The Conjugate-transpose Form of RFP A is stored.
! 45: *
! 46: * UPLO (input) CHARACTER
! 47: * On entry, UPLO specifies whether the upper or lower
! 48: * triangular part of the array C is to be referenced as
! 49: * follows:
! 50: *
! 51: * UPLO = 'U' or 'u' Only the upper triangular part of C
! 52: * is to be referenced.
! 53: *
! 54: * UPLO = 'L' or 'l' Only the lower triangular part of C
! 55: * is to be referenced.
! 56: *
! 57: * Unchanged on exit.
! 58: *
! 59: * TRANS (input) CHARACTER
! 60: * On entry, TRANS specifies the operation to be performed as
! 61: * follows:
! 62: *
! 63: * TRANS = 'N' or 'n' C := alpha*A*conjg( A' ) + beta*C.
! 64: *
! 65: * TRANS = 'C' or 'c' C := alpha*conjg( A' )*A + beta*C.
! 66: *
! 67: * Unchanged on exit.
! 68: *
! 69: * N (input) INTEGER
! 70: * On entry, N specifies the order of the matrix C. N must be
! 71: * at least zero.
! 72: * Unchanged on exit.
! 73: *
! 74: * K (input) INTEGER
! 75: * On entry with TRANS = 'N' or 'n', K specifies the number
! 76: * of columns of the matrix A, and on entry with
! 77: * TRANS = 'C' or 'c', K specifies the number of rows of the
! 78: * matrix A. K must be at least zero.
! 79: * Unchanged on exit.
! 80: *
! 81: * ALPHA (input) DOUBLE PRECISION
! 82: * On entry, ALPHA specifies the scalar alpha.
! 83: * Unchanged on exit.
! 84: *
! 85: * A (input) COMPLEX*16 array of DIMENSION (LDA,ka)
! 86: * where KA
! 87: * is K when TRANS = 'N' or 'n', and is N otherwise. Before
! 88: * entry with TRANS = 'N' or 'n', the leading N--by--K part of
! 89: * the array A must contain the matrix A, otherwise the leading
! 90: * K--by--N part of the array A must contain the matrix A.
! 91: * Unchanged on exit.
! 92: *
! 93: * LDA (input) INTEGER
! 94: * On entry, LDA specifies the first dimension of A as declared
! 95: * in the calling (sub) program. When TRANS = 'N' or 'n'
! 96: * then LDA must be at least max( 1, n ), otherwise LDA must
! 97: * be at least max( 1, k ).
! 98: * Unchanged on exit.
! 99: *
! 100: * BETA (input) DOUBLE PRECISION
! 101: * On entry, BETA specifies the scalar beta.
! 102: * Unchanged on exit.
! 103: *
! 104: * C (input/output) COMPLEX*16 array, dimension (N*(N+1)/2)
! 105: * On entry, the matrix A in RFP Format. RFP Format is
! 106: * described by TRANSR, UPLO and N. Note that the imaginary
! 107: * parts of the diagonal elements need not be set, they are
! 108: * assumed to be zero, and on exit they are set to zero.
! 109: *
! 110: * Arguments
! 111: * ==========
! 112: *
! 113: * ..
! 114: * .. Parameters ..
! 115: DOUBLE PRECISION ONE, ZERO
! 116: COMPLEX*16 CZERO
! 117: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
! 118: PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ) )
! 119: * ..
! 120: * .. Local Scalars ..
! 121: LOGICAL LOWER, NORMALTRANSR, NISODD, NOTRANS
! 122: INTEGER INFO, NROWA, J, NK, N1, N2
! 123: COMPLEX*16 CALPHA, CBETA
! 124: * ..
! 125: * .. External Functions ..
! 126: LOGICAL LSAME
! 127: EXTERNAL LSAME
! 128: * ..
! 129: * .. External Subroutines ..
! 130: EXTERNAL XERBLA, ZGEMM, ZHERK
! 131: * ..
! 132: * .. Intrinsic Functions ..
! 133: INTRINSIC MAX, DCMPLX
! 134: * ..
! 135: * .. Executable Statements ..
! 136: *
! 137: *
! 138: * Test the input parameters.
! 139: *
! 140: INFO = 0
! 141: NORMALTRANSR = LSAME( TRANSR, 'N' )
! 142: LOWER = LSAME( UPLO, 'L' )
! 143: NOTRANS = LSAME( TRANS, 'N' )
! 144: *
! 145: IF( NOTRANS ) THEN
! 146: NROWA = N
! 147: ELSE
! 148: NROWA = K
! 149: END IF
! 150: *
! 151: IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'C' ) ) THEN
! 152: INFO = -1
! 153: ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
! 154: INFO = -2
! 155: ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'C' ) ) THEN
! 156: INFO = -3
! 157: ELSE IF( N.LT.0 ) THEN
! 158: INFO = -4
! 159: ELSE IF( K.LT.0 ) THEN
! 160: INFO = -5
! 161: ELSE IF( LDA.LT.MAX( 1, NROWA ) ) THEN
! 162: INFO = -8
! 163: END IF
! 164: IF( INFO.NE.0 ) THEN
! 165: CALL XERBLA( 'ZHFRK ', -INFO )
! 166: RETURN
! 167: END IF
! 168: *
! 169: * Quick return if possible.
! 170: *
! 171: * The quick return case: ((ALPHA.EQ.0).AND.(BETA.NE.ZERO)) is not
! 172: * done (it is in ZHERK for example) and left in the general case.
! 173: *
! 174: IF( ( N.EQ.0 ) .OR. ( ( ( ALPHA.EQ.ZERO ) .OR. ( K.EQ.0 ) ) .AND.
! 175: + ( BETA.EQ.ONE ) ) )RETURN
! 176: *
! 177: IF( ( ALPHA.EQ.ZERO ) .AND. ( BETA.EQ.ZERO ) ) THEN
! 178: DO J = 1, ( ( N*( N+1 ) ) / 2 )
! 179: C( J ) = CZERO
! 180: END DO
! 181: RETURN
! 182: END IF
! 183: *
! 184: CALPHA = DCMPLX( ALPHA, ZERO )
! 185: CBETA = DCMPLX( BETA, ZERO )
! 186: *
! 187: * C is N-by-N.
! 188: * If N is odd, set NISODD = .TRUE., and N1 and N2.
! 189: * If N is even, NISODD = .FALSE., and NK.
! 190: *
! 191: IF( MOD( N, 2 ).EQ.0 ) THEN
! 192: NISODD = .FALSE.
! 193: NK = N / 2
! 194: ELSE
! 195: NISODD = .TRUE.
! 196: IF( LOWER ) THEN
! 197: N2 = N / 2
! 198: N1 = N - N2
! 199: ELSE
! 200: N1 = N / 2
! 201: N2 = N - N1
! 202: END IF
! 203: END IF
! 204: *
! 205: IF( NISODD ) THEN
! 206: *
! 207: * N is odd
! 208: *
! 209: IF( NORMALTRANSR ) THEN
! 210: *
! 211: * N is odd and TRANSR = 'N'
! 212: *
! 213: IF( LOWER ) THEN
! 214: *
! 215: * N is odd, TRANSR = 'N', and UPLO = 'L'
! 216: *
! 217: IF( NOTRANS ) THEN
! 218: *
! 219: * N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
! 220: *
! 221: CALL ZHERK( 'L', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
! 222: + BETA, C( 1 ), N )
! 223: CALL ZHERK( 'U', 'N', N2, K, ALPHA, A( N1+1, 1 ), LDA,
! 224: + BETA, C( N+1 ), N )
! 225: CALL ZGEMM( 'N', 'C', N2, N1, K, CALPHA, A( N1+1, 1 ),
! 226: + LDA, A( 1, 1 ), LDA, CBETA, C( N1+1 ), N )
! 227: *
! 228: ELSE
! 229: *
! 230: * N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'C'
! 231: *
! 232: CALL ZHERK( 'L', 'C', N1, K, ALPHA, A( 1, 1 ), LDA,
! 233: + BETA, C( 1 ), N )
! 234: CALL ZHERK( 'U', 'C', N2, K, ALPHA, A( 1, N1+1 ), LDA,
! 235: + BETA, C( N+1 ), N )
! 236: CALL ZGEMM( 'C', 'N', N2, N1, K, CALPHA, A( 1, N1+1 ),
! 237: + LDA, A( 1, 1 ), LDA, CBETA, C( N1+1 ), N )
! 238: *
! 239: END IF
! 240: *
! 241: ELSE
! 242: *
! 243: * N is odd, TRANSR = 'N', and UPLO = 'U'
! 244: *
! 245: IF( NOTRANS ) THEN
! 246: *
! 247: * N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
! 248: *
! 249: CALL ZHERK( 'L', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
! 250: + BETA, C( N2+1 ), N )
! 251: CALL ZHERK( 'U', 'N', N2, K, ALPHA, A( N2, 1 ), LDA,
! 252: + BETA, C( N1+1 ), N )
! 253: CALL ZGEMM( 'N', 'C', N1, N2, K, CALPHA, A( 1, 1 ),
! 254: + LDA, A( N2, 1 ), LDA, CBETA, C( 1 ), N )
! 255: *
! 256: ELSE
! 257: *
! 258: * N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'C'
! 259: *
! 260: CALL ZHERK( 'L', 'C', N1, K, ALPHA, A( 1, 1 ), LDA,
! 261: + BETA, C( N2+1 ), N )
! 262: CALL ZHERK( 'U', 'C', N2, K, ALPHA, A( 1, N2 ), LDA,
! 263: + BETA, C( N1+1 ), N )
! 264: CALL ZGEMM( 'C', 'N', N1, N2, K, CALPHA, A( 1, 1 ),
! 265: + LDA, A( 1, N2 ), LDA, CBETA, C( 1 ), N )
! 266: *
! 267: END IF
! 268: *
! 269: END IF
! 270: *
! 271: ELSE
! 272: *
! 273: * N is odd, and TRANSR = 'C'
! 274: *
! 275: IF( LOWER ) THEN
! 276: *
! 277: * N is odd, TRANSR = 'C', and UPLO = 'L'
! 278: *
! 279: IF( NOTRANS ) THEN
! 280: *
! 281: * N is odd, TRANSR = 'C', UPLO = 'L', and TRANS = 'N'
! 282: *
! 283: CALL ZHERK( 'U', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
! 284: + BETA, C( 1 ), N1 )
! 285: CALL ZHERK( 'L', 'N', N2, K, ALPHA, A( N1+1, 1 ), LDA,
! 286: + BETA, C( 2 ), N1 )
! 287: CALL ZGEMM( 'N', 'C', N1, N2, K, CALPHA, A( 1, 1 ),
! 288: + LDA, A( N1+1, 1 ), LDA, CBETA,
! 289: + C( N1*N1+1 ), N1 )
! 290: *
! 291: ELSE
! 292: *
! 293: * N is odd, TRANSR = 'C', UPLO = 'L', and TRANS = 'C'
! 294: *
! 295: CALL ZHERK( 'U', 'C', N1, K, ALPHA, A( 1, 1 ), LDA,
! 296: + BETA, C( 1 ), N1 )
! 297: CALL ZHERK( 'L', 'C', N2, K, ALPHA, A( 1, N1+1 ), LDA,
! 298: + BETA, C( 2 ), N1 )
! 299: CALL ZGEMM( 'C', 'N', N1, N2, K, CALPHA, A( 1, 1 ),
! 300: + LDA, A( 1, N1+1 ), LDA, CBETA,
! 301: + C( N1*N1+1 ), N1 )
! 302: *
! 303: END IF
! 304: *
! 305: ELSE
! 306: *
! 307: * N is odd, TRANSR = 'C', and UPLO = 'U'
! 308: *
! 309: IF( NOTRANS ) THEN
! 310: *
! 311: * N is odd, TRANSR = 'C', UPLO = 'U', and TRANS = 'N'
! 312: *
! 313: CALL ZHERK( 'U', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
! 314: + BETA, C( N2*N2+1 ), N2 )
! 315: CALL ZHERK( 'L', 'N', N2, K, ALPHA, A( N1+1, 1 ), LDA,
! 316: + BETA, C( N1*N2+1 ), N2 )
! 317: CALL ZGEMM( 'N', 'C', N2, N1, K, CALPHA, A( N1+1, 1 ),
! 318: + LDA, A( 1, 1 ), LDA, CBETA, C( 1 ), N2 )
! 319: *
! 320: ELSE
! 321: *
! 322: * N is odd, TRANSR = 'C', UPLO = 'U', and TRANS = 'C'
! 323: *
! 324: CALL ZHERK( 'U', 'C', N1, K, ALPHA, A( 1, 1 ), LDA,
! 325: + BETA, C( N2*N2+1 ), N2 )
! 326: CALL ZHERK( 'L', 'C', N2, K, ALPHA, A( 1, N1+1 ), LDA,
! 327: + BETA, C( N1*N2+1 ), N2 )
! 328: CALL ZGEMM( 'C', 'N', N2, N1, K, CALPHA, A( 1, N1+1 ),
! 329: + LDA, A( 1, 1 ), LDA, CBETA, C( 1 ), N2 )
! 330: *
! 331: END IF
! 332: *
! 333: END IF
! 334: *
! 335: END IF
! 336: *
! 337: ELSE
! 338: *
! 339: * N is even
! 340: *
! 341: IF( NORMALTRANSR ) THEN
! 342: *
! 343: * N is even and TRANSR = 'N'
! 344: *
! 345: IF( LOWER ) THEN
! 346: *
! 347: * N is even, TRANSR = 'N', and UPLO = 'L'
! 348: *
! 349: IF( NOTRANS ) THEN
! 350: *
! 351: * N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
! 352: *
! 353: CALL ZHERK( 'L', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
! 354: + BETA, C( 2 ), N+1 )
! 355: CALL ZHERK( 'U', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
! 356: + BETA, C( 1 ), N+1 )
! 357: CALL ZGEMM( 'N', 'C', NK, NK, K, CALPHA, A( NK+1, 1 ),
! 358: + LDA, A( 1, 1 ), LDA, CBETA, C( NK+2 ),
! 359: + N+1 )
! 360: *
! 361: ELSE
! 362: *
! 363: * N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'C'
! 364: *
! 365: CALL ZHERK( 'L', 'C', NK, K, ALPHA, A( 1, 1 ), LDA,
! 366: + BETA, C( 2 ), N+1 )
! 367: CALL ZHERK( 'U', 'C', NK, K, ALPHA, A( 1, NK+1 ), LDA,
! 368: + BETA, C( 1 ), N+1 )
! 369: CALL ZGEMM( 'C', 'N', NK, NK, K, CALPHA, A( 1, NK+1 ),
! 370: + LDA, A( 1, 1 ), LDA, CBETA, C( NK+2 ),
! 371: + N+1 )
! 372: *
! 373: END IF
! 374: *
! 375: ELSE
! 376: *
! 377: * N is even, TRANSR = 'N', and UPLO = 'U'
! 378: *
! 379: IF( NOTRANS ) THEN
! 380: *
! 381: * N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
! 382: *
! 383: CALL ZHERK( 'L', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
! 384: + BETA, C( NK+2 ), N+1 )
! 385: CALL ZHERK( 'U', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
! 386: + BETA, C( NK+1 ), N+1 )
! 387: CALL ZGEMM( 'N', 'C', NK, NK, K, CALPHA, A( 1, 1 ),
! 388: + LDA, A( NK+1, 1 ), LDA, CBETA, C( 1 ),
! 389: + N+1 )
! 390: *
! 391: ELSE
! 392: *
! 393: * N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'C'
! 394: *
! 395: CALL ZHERK( 'L', 'C', NK, K, ALPHA, A( 1, 1 ), LDA,
! 396: + BETA, C( NK+2 ), N+1 )
! 397: CALL ZHERK( 'U', 'C', NK, K, ALPHA, A( 1, NK+1 ), LDA,
! 398: + BETA, C( NK+1 ), N+1 )
! 399: CALL ZGEMM( 'C', 'N', NK, NK, K, CALPHA, A( 1, 1 ),
! 400: + LDA, A( 1, NK+1 ), LDA, CBETA, C( 1 ),
! 401: + N+1 )
! 402: *
! 403: END IF
! 404: *
! 405: END IF
! 406: *
! 407: ELSE
! 408: *
! 409: * N is even, and TRANSR = 'C'
! 410: *
! 411: IF( LOWER ) THEN
! 412: *
! 413: * N is even, TRANSR = 'C', and UPLO = 'L'
! 414: *
! 415: IF( NOTRANS ) THEN
! 416: *
! 417: * N is even, TRANSR = 'C', UPLO = 'L', and TRANS = 'N'
! 418: *
! 419: CALL ZHERK( 'U', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
! 420: + BETA, C( NK+1 ), NK )
! 421: CALL ZHERK( 'L', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
! 422: + BETA, C( 1 ), NK )
! 423: CALL ZGEMM( 'N', 'C', NK, NK, K, CALPHA, A( 1, 1 ),
! 424: + LDA, A( NK+1, 1 ), LDA, CBETA,
! 425: + C( ( ( NK+1 )*NK )+1 ), NK )
! 426: *
! 427: ELSE
! 428: *
! 429: * N is even, TRANSR = 'C', UPLO = 'L', and TRANS = 'C'
! 430: *
! 431: CALL ZHERK( 'U', 'C', NK, K, ALPHA, A( 1, 1 ), LDA,
! 432: + BETA, C( NK+1 ), NK )
! 433: CALL ZHERK( 'L', 'C', NK, K, ALPHA, A( 1, NK+1 ), LDA,
! 434: + BETA, C( 1 ), NK )
! 435: CALL ZGEMM( 'C', 'N', NK, NK, K, CALPHA, A( 1, 1 ),
! 436: + LDA, A( 1, NK+1 ), LDA, CBETA,
! 437: + C( ( ( NK+1 )*NK )+1 ), NK )
! 438: *
! 439: END IF
! 440: *
! 441: ELSE
! 442: *
! 443: * N is even, TRANSR = 'C', and UPLO = 'U'
! 444: *
! 445: IF( NOTRANS ) THEN
! 446: *
! 447: * N is even, TRANSR = 'C', UPLO = 'U', and TRANS = 'N'
! 448: *
! 449: CALL ZHERK( 'U', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
! 450: + BETA, C( NK*( NK+1 )+1 ), NK )
! 451: CALL ZHERK( 'L', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
! 452: + BETA, C( NK*NK+1 ), NK )
! 453: CALL ZGEMM( 'N', 'C', NK, NK, K, CALPHA, A( NK+1, 1 ),
! 454: + LDA, A( 1, 1 ), LDA, CBETA, C( 1 ), NK )
! 455: *
! 456: ELSE
! 457: *
! 458: * N is even, TRANSR = 'C', UPLO = 'U', and TRANS = 'C'
! 459: *
! 460: CALL ZHERK( 'U', 'C', NK, K, ALPHA, A( 1, 1 ), LDA,
! 461: + BETA, C( NK*( NK+1 )+1 ), NK )
! 462: CALL ZHERK( 'L', 'C', NK, K, ALPHA, A( 1, NK+1 ), LDA,
! 463: + BETA, C( NK*NK+1 ), NK )
! 464: CALL ZGEMM( 'C', 'N', NK, NK, K, CALPHA, A( 1, NK+1 ),
! 465: + LDA, A( 1, 1 ), LDA, CBETA, C( 1 ), NK )
! 466: *
! 467: END IF
! 468: *
! 469: END IF
! 470: *
! 471: END IF
! 472: *
! 473: END IF
! 474: *
! 475: RETURN
! 476: *
! 477: * End of ZHFRK
! 478: *
! 479: END
CVSweb interface <joel.bertrand@systella.fr>