Annotation of rpl/lapack/lapack/dlaqz4.f, revision 1.1
1.1 ! bertrand 1: *> \brief \b DLAQZ4
! 2: *
! 3: * =========== DOCUMENTATION ===========
! 4: *
! 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
! 7: *
! 8: *> \htmlonly
! 9: *> Download DLAQZ4 + dependencies
! 10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaqz4.f">
! 11: *> [TGZ]</a>
! 12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaqz4.f">
! 13: *> [ZIP]</a>
! 14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaqz4.f">
! 15: *> [TXT]</a>
! 16: *> \endhtmlonly
! 17: *
! 18: * Definition:
! 19: * ===========
! 20: *
! 21: * SUBROUTINE DLAQZ4( ILSCHUR, ILQ, ILZ, N, ILO, IHI, NSHIFTS,
! 22: * $ NBLOCK_DESIRED, SR, SI, SS, A, LDA, B, LDB, Q, LDQ, Z, LDZ,
! 23: * $ QC, LDQC, ZC, LDZC, WORK, LWORK, INFO )
! 24: * IMPLICIT NONE
! 25: *
! 26: * Function arguments
! 27: * LOGICAL, INTENT( IN ) :: ILSCHUR, ILQ, ILZ
! 28: * INTEGER, INTENT( IN ) :: N, ILO, IHI, LDA, LDB, LDQ, LDZ, LWORK,
! 29: * $ NSHIFTS, NBLOCK_DESIRED, LDQC, LDZC
! 30: *
! 31: * DOUBLE PRECISION, INTENT( INOUT ) :: A( LDA, * ), B( LDB, * ),
! 32: * $ Q( LDQ, * ), Z( LDZ, * ), QC( LDQC, * ), ZC( LDZC, * ),
! 33: * $ WORK( * ), SR( * ), SI( * ), SS( * )
! 34: *
! 35: * INTEGER, INTENT( OUT ) :: INFO
! 36: * ..
! 37: *
! 38: *
! 39: *> \par Purpose:
! 40: * =============
! 41: *>
! 42: *> \verbatim
! 43: *>
! 44: *> DLAQZ4 Executes a single multishift QZ sweep
! 45: *> \endverbatim
! 46: *
! 47: * Arguments:
! 48: * ==========
! 49: *
! 50: *> \param[in] ILSCHUR
! 51: *> \verbatim
! 52: *> ILSCHUR is LOGICAL
! 53: *> Determines whether or not to update the full Schur form
! 54: *> \endverbatim
! 55: *> \param[in] ILQ
! 56: *> \verbatim
! 57: *> ILQ is LOGICAL
! 58: *> Determines whether or not to update the matrix Q
! 59: *> \endverbatim
! 60: *>
! 61: *> \param[in] ILZ
! 62: *> \verbatim
! 63: *> ILZ is LOGICAL
! 64: *> Determines whether or not to update the matrix Z
! 65: *> \endverbatim
! 66: *>
! 67: *> \param[in] N
! 68: *> \verbatim
! 69: *> N is INTEGER
! 70: *> The order of the matrices A, B, Q, and Z. N >= 0.
! 71: *> \endverbatim
! 72: *>
! 73: *> \param[in] ILO
! 74: *> \verbatim
! 75: *> ILO is INTEGER
! 76: *> \endverbatim
! 77: *>
! 78: *> \param[in] IHI
! 79: *> \verbatim
! 80: *> IHI is INTEGER
! 81: *> \endverbatim
! 82: *>
! 83: *> \param[in] NSHIFTS
! 84: *> \verbatim
! 85: *> NSHIFTS is INTEGER
! 86: *> The desired number of shifts to use
! 87: *> \endverbatim
! 88: *>
! 89: *> \param[in] NBLOCK_DESIRED
! 90: *> \verbatim
! 91: *> NBLOCK_DESIRED is INTEGER
! 92: *> The desired size of the computational windows
! 93: *> \endverbatim
! 94: *>
! 95: *> \param[in] SR
! 96: *> \verbatim
! 97: *> SR is DOUBLE PRECISION array. SR contains
! 98: *> the real parts of the shifts to use.
! 99: *> \endverbatim
! 100: *>
! 101: *> \param[in] SI
! 102: *> \verbatim
! 103: *> SI is DOUBLE PRECISION array. SI contains
! 104: *> the imaginary parts of the shifts to use.
! 105: *> \endverbatim
! 106: *>
! 107: *> \param[in] SS
! 108: *> \verbatim
! 109: *> SS is DOUBLE PRECISION array. SS contains
! 110: *> the scale of the shifts to use.
! 111: *> \endverbatim
! 112: *>
! 113: *> \param[in,out] A
! 114: *> \verbatim
! 115: *> A is DOUBLE PRECISION array, dimension (LDA, N)
! 116: *> \endverbatim
! 117: *>
! 118: *> \param[in] LDA
! 119: *> \verbatim
! 120: *> LDA is INTEGER
! 121: *> The leading dimension of the array A. LDA >= max( 1, N ).
! 122: *> \endverbatim
! 123: *>
! 124: *> \param[in,out] B
! 125: *> \verbatim
! 126: *> B is DOUBLE PRECISION array, dimension (LDB, N)
! 127: *> \endverbatim
! 128: *>
! 129: *> \param[in] LDB
! 130: *> \verbatim
! 131: *> LDB is INTEGER
! 132: *> The leading dimension of the array B. LDB >= max( 1, N ).
! 133: *> \endverbatim
! 134: *>
! 135: *> \param[in,out] Q
! 136: *> \verbatim
! 137: *> Q is DOUBLE PRECISION array, dimension (LDQ, N)
! 138: *> \endverbatim
! 139: *>
! 140: *> \param[in] LDQ
! 141: *> \verbatim
! 142: *> LDQ is INTEGER
! 143: *> \endverbatim
! 144: *>
! 145: *> \param[in,out] Z
! 146: *> \verbatim
! 147: *> Z is DOUBLE PRECISION array, dimension (LDZ, N)
! 148: *> \endverbatim
! 149: *>
! 150: *> \param[in] LDZ
! 151: *> \verbatim
! 152: *> LDZ is INTEGER
! 153: *> \endverbatim
! 154: *>
! 155: *> \param[in,out] QC
! 156: *> \verbatim
! 157: *> QC is DOUBLE PRECISION array, dimension (LDQC, NBLOCK_DESIRED)
! 158: *> \endverbatim
! 159: *>
! 160: *> \param[in] LDQC
! 161: *> \verbatim
! 162: *> LDQC is INTEGER
! 163: *> \endverbatim
! 164: *>
! 165: *> \param[in,out] ZC
! 166: *> \verbatim
! 167: *> ZC is DOUBLE PRECISION array, dimension (LDZC, NBLOCK_DESIRED)
! 168: *> \endverbatim
! 169: *>
! 170: *> \param[in] LDZC
! 171: *> \verbatim
! 172: *> LDZ is INTEGER
! 173: *> \endverbatim
! 174: *>
! 175: *> \param[out] WORK
! 176: *> \verbatim
! 177: *> WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
! 178: *> On exit, if INFO >= 0, WORK(1) returns the optimal LWORK.
! 179: *> \endverbatim
! 180: *>
! 181: *> \param[in] LWORK
! 182: *> \verbatim
! 183: *> LWORK is INTEGER
! 184: *> The dimension of the array WORK. LWORK >= max(1,N).
! 185: *>
! 186: *> If LWORK = -1, then a workspace query is assumed; the routine
! 187: *> only calculates the optimal size of the WORK array, returns
! 188: *> this value as the first entry of the WORK array, and no error
! 189: *> message related to LWORK is issued by XERBLA.
! 190: *> \endverbatim
! 191: *>
! 192: *> \param[out] INFO
! 193: *> \verbatim
! 194: *> INFO is INTEGER
! 195: *> = 0: successful exit
! 196: *> < 0: if INFO = -i, the i-th argument had an illegal value
! 197: *> \endverbatim
! 198: *
! 199: * Authors:
! 200: * ========
! 201: *
! 202: *> \author Thijs Steel, KU Leuven
! 203: *
! 204: *> \date May 2020
! 205: *
! 206: *> \ingroup doubleGEcomputational
! 207: *>
! 208: * =====================================================================
! 209: SUBROUTINE DLAQZ4( ILSCHUR, ILQ, ILZ, N, ILO, IHI, NSHIFTS,
! 210: $ NBLOCK_DESIRED, SR, SI, SS, A, LDA, B, LDB, Q,
! 211: $ LDQ, Z, LDZ, QC, LDQC, ZC, LDZC, WORK, LWORK,
! 212: $ INFO )
! 213: IMPLICIT NONE
! 214:
! 215: * Function arguments
! 216: LOGICAL, INTENT( IN ) :: ILSCHUR, ILQ, ILZ
! 217: INTEGER, INTENT( IN ) :: N, ILO, IHI, LDA, LDB, LDQ, LDZ, LWORK,
! 218: $ NSHIFTS, NBLOCK_DESIRED, LDQC, LDZC
! 219:
! 220: DOUBLE PRECISION, INTENT( INOUT ) :: A( LDA, * ), B( LDB, * ),
! 221: $ Q( LDQ, * ), Z( LDZ, * ), QC( LDQC, * ),
! 222: $ ZC( LDZC, * ), WORK( * ), SR( * ), SI( * ),
! 223: $ SS( * )
! 224:
! 225: INTEGER, INTENT( OUT ) :: INFO
! 226:
! 227: * Parameters
! 228: DOUBLE PRECISION :: ZERO, ONE, HALF
! 229: PARAMETER( ZERO = 0.0D0, ONE = 1.0D0, HALF = 0.5D0 )
! 230:
! 231: * Local scalars
! 232: INTEGER :: I, J, NS, ISTARTM, ISTOPM, SHEIGHT, SWIDTH, K, NP,
! 233: $ ISTARTB, ISTOPB, ISHIFT, NBLOCK, NPOS
! 234: DOUBLE PRECISION :: TEMP, V( 3 ), C1, S1, C2, S2, SWAP
! 235: *
! 236: * External functions
! 237: EXTERNAL :: XERBLA, DGEMM, DLAQZ1, DLAQZ2, DLASET, DLARTG, DROT,
! 238: $ DLACPY
! 239:
! 240: INFO = 0
! 241: IF ( NBLOCK_DESIRED .LT. NSHIFTS+1 ) THEN
! 242: INFO = -8
! 243: END IF
! 244: IF ( LWORK .EQ.-1 ) THEN
! 245: * workspace query, quick return
! 246: WORK( 1 ) = N*NBLOCK_DESIRED
! 247: RETURN
! 248: ELSE IF ( LWORK .LT. N*NBLOCK_DESIRED ) THEN
! 249: INFO = -25
! 250: END IF
! 251:
! 252: IF( INFO.NE.0 ) THEN
! 253: CALL XERBLA( 'DLAQZ4', -INFO )
! 254: RETURN
! 255: END IF
! 256:
! 257: * Executable statements
! 258:
! 259: IF ( NSHIFTS .LT. 2 ) THEN
! 260: RETURN
! 261: END IF
! 262:
! 263: IF ( ILO .GE. IHI ) THEN
! 264: RETURN
! 265: END IF
! 266:
! 267: IF ( ILSCHUR ) THEN
! 268: ISTARTM = 1
! 269: ISTOPM = N
! 270: ELSE
! 271: ISTARTM = ILO
! 272: ISTOPM = IHI
! 273: END IF
! 274:
! 275: * Shuffle shifts into pairs of real shifts and pairs
! 276: * of complex conjugate shifts assuming complex
! 277: * conjugate shifts are already adjacent to one
! 278: * another
! 279:
! 280: DO I = 1, NSHIFTS-2, 2
! 281: IF( SI( I ).NE.-SI( I+1 ) ) THEN
! 282: *
! 283: SWAP = SR( I )
! 284: SR( I ) = SR( I+1 )
! 285: SR( I+1 ) = SR( I+2 )
! 286: SR( I+2 ) = SWAP
! 287:
! 288: SWAP = SI( I )
! 289: SI( I ) = SI( I+1 )
! 290: SI( I+1 ) = SI( I+2 )
! 291: SI( I+2 ) = SWAP
! 292:
! 293: SWAP = SS( I )
! 294: SS( I ) = SS( I+1 )
! 295: SS( I+1 ) = SS( I+2 )
! 296: SS( I+2 ) = SWAP
! 297: END IF
! 298: END DO
! 299:
! 300: * NSHFTS is supposed to be even, but if it is odd,
! 301: * then simply reduce it by one. The shuffle above
! 302: * ensures that the dropped shift is real and that
! 303: * the remaining shifts are paired.
! 304:
! 305: NS = NSHIFTS-MOD( NSHIFTS, 2 )
! 306: NPOS = MAX( NBLOCK_DESIRED-NS, 1 )
! 307:
! 308: * The following block introduces the shifts and chases
! 309: * them down one by one just enough to make space for
! 310: * the other shifts. The near-the-diagonal block is
! 311: * of size (ns+1) x ns.
! 312:
! 313: CALL DLASET( 'FULL', NS+1, NS+1, ZERO, ONE, QC, LDQC )
! 314: CALL DLASET( 'FULL', NS, NS, ZERO, ONE, ZC, LDZC )
! 315:
! 316: DO I = 1, NS, 2
! 317: * Introduce the shift
! 318: CALL DLAQZ1( A( ILO, ILO ), LDA, B( ILO, ILO ), LDB, SR( I ),
! 319: $ SR( I+1 ), SI( I ), SS( I ), SS( I+1 ), V )
! 320:
! 321: TEMP = V( 2 )
! 322: CALL DLARTG( TEMP, V( 3 ), C1, S1, V( 2 ) )
! 323: CALL DLARTG( V( 1 ), V( 2 ), C2, S2, TEMP )
! 324:
! 325: CALL DROT( NS, A( ILO+1, ILO ), LDA, A( ILO+2, ILO ), LDA, C1,
! 326: $ S1 )
! 327: CALL DROT( NS, A( ILO, ILO ), LDA, A( ILO+1, ILO ), LDA, C2,
! 328: $ S2 )
! 329: CALL DROT( NS, B( ILO+1, ILO ), LDB, B( ILO+2, ILO ), LDB, C1,
! 330: $ S1 )
! 331: CALL DROT( NS, B( ILO, ILO ), LDB, B( ILO+1, ILO ), LDB, C2,
! 332: $ S2 )
! 333: CALL DROT( NS+1, QC( 1, 2 ), 1, QC( 1, 3 ), 1, C1, S1 )
! 334: CALL DROT( NS+1, QC( 1, 1 ), 1, QC( 1, 2 ), 1, C2, S2 )
! 335:
! 336: * Chase the shift down
! 337: DO J = 1, NS-1-I
! 338:
! 339: CALL DLAQZ2( .TRUE., .TRUE., J, 1, NS, IHI-ILO+1, A( ILO,
! 340: $ ILO ), LDA, B( ILO, ILO ), LDB, NS+1, 1, QC,
! 341: $ LDQC, NS, 1, ZC, LDZC )
! 342:
! 343: END DO
! 344:
! 345: END DO
! 346:
! 347: * Update the rest of the pencil
! 348:
! 349: * Update A(ilo:ilo+ns,ilo+ns:istopm) and B(ilo:ilo+ns,ilo+ns:istopm)
! 350: * from the left with Qc(1:ns+1,1:ns+1)'
! 351: SHEIGHT = NS+1
! 352: SWIDTH = ISTOPM-( ILO+NS )+1
! 353: IF ( SWIDTH > 0 ) THEN
! 354: CALL DGEMM( 'T', 'N', SHEIGHT, SWIDTH, SHEIGHT, ONE, QC, LDQC,
! 355: $ A( ILO, ILO+NS ), LDA, ZERO, WORK, SHEIGHT )
! 356: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT, A( ILO,
! 357: $ ILO+NS ), LDA )
! 358: CALL DGEMM( 'T', 'N', SHEIGHT, SWIDTH, SHEIGHT, ONE, QC, LDQC,
! 359: $ B( ILO, ILO+NS ), LDB, ZERO, WORK, SHEIGHT )
! 360: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT, B( ILO,
! 361: $ ILO+NS ), LDB )
! 362: END IF
! 363: IF ( ILQ ) THEN
! 364: CALL DGEMM( 'N', 'N', N, SHEIGHT, SHEIGHT, ONE, Q( 1, ILO ),
! 365: $ LDQ, QC, LDQC, ZERO, WORK, N )
! 366: CALL DLACPY( 'ALL', N, SHEIGHT, WORK, N, Q( 1, ILO ), LDQ )
! 367: END IF
! 368:
! 369: * Update A(istartm:ilo-1,ilo:ilo+ns-1) and B(istartm:ilo-1,ilo:ilo+ns-1)
! 370: * from the right with Zc(1:ns,1:ns)
! 371: SHEIGHT = ILO-1-ISTARTM+1
! 372: SWIDTH = NS
! 373: IF ( SHEIGHT > 0 ) THEN
! 374: CALL DGEMM( 'N', 'N', SHEIGHT, SWIDTH, SWIDTH, ONE, A( ISTARTM,
! 375: $ ILO ), LDA, ZC, LDZC, ZERO, WORK, SHEIGHT )
! 376: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT, A( ISTARTM,
! 377: $ ILO ), LDA )
! 378: CALL DGEMM( 'N', 'N', SHEIGHT, SWIDTH, SWIDTH, ONE, B( ISTARTM,
! 379: $ ILO ), LDB, ZC, LDZC, ZERO, WORK, SHEIGHT )
! 380: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT, B( ISTARTM,
! 381: $ ILO ), LDB )
! 382: END IF
! 383: IF ( ILZ ) THEN
! 384: CALL DGEMM( 'N', 'N', N, SWIDTH, SWIDTH, ONE, Z( 1, ILO ), LDZ,
! 385: $ ZC, LDZC, ZERO, WORK, N )
! 386: CALL DLACPY( 'ALL', N, SWIDTH, WORK, N, Z( 1, ILO ), LDZ )
! 387: END IF
! 388:
! 389: * The following block chases the shifts down to the bottom
! 390: * right block. If possible, a shift is moved down npos
! 391: * positions at a time
! 392:
! 393: K = ILO
! 394: DO WHILE ( K < IHI-NS )
! 395: NP = MIN( IHI-NS-K, NPOS )
! 396: * Size of the near-the-diagonal block
! 397: NBLOCK = NS+NP
! 398: * istartb points to the first row we will be updating
! 399: ISTARTB = K+1
! 400: * istopb points to the last column we will be updating
! 401: ISTOPB = K+NBLOCK-1
! 402:
! 403: CALL DLASET( 'FULL', NS+NP, NS+NP, ZERO, ONE, QC, LDQC )
! 404: CALL DLASET( 'FULL', NS+NP, NS+NP, ZERO, ONE, ZC, LDZC )
! 405:
! 406: * Near the diagonal shift chase
! 407: DO I = NS-1, 0, -2
! 408: DO J = 0, NP-1
! 409: * Move down the block with index k+i+j-1, updating
! 410: * the (ns+np x ns+np) block:
! 411: * (k:k+ns+np,k:k+ns+np-1)
! 412: CALL DLAQZ2( .TRUE., .TRUE., K+I+J-1, ISTARTB, ISTOPB,
! 413: $ IHI, A, LDA, B, LDB, NBLOCK, K+1, QC, LDQC,
! 414: $ NBLOCK, K, ZC, LDZC )
! 415: END DO
! 416: END DO
! 417:
! 418: * Update rest of the pencil
! 419:
! 420: * Update A(k+1:k+ns+np, k+ns+np:istopm) and
! 421: * B(k+1:k+ns+np, k+ns+np:istopm)
! 422: * from the left with Qc(1:ns+np,1:ns+np)'
! 423: SHEIGHT = NS+NP
! 424: SWIDTH = ISTOPM-( K+NS+NP )+1
! 425: IF ( SWIDTH > 0 ) THEN
! 426: CALL DGEMM( 'T', 'N', SHEIGHT, SWIDTH, SHEIGHT, ONE, QC,
! 427: $ LDQC, A( K+1, K+NS+NP ), LDA, ZERO, WORK,
! 428: $ SHEIGHT )
! 429: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT, A( K+1,
! 430: $ K+NS+NP ), LDA )
! 431: CALL DGEMM( 'T', 'N', SHEIGHT, SWIDTH, SHEIGHT, ONE, QC,
! 432: $ LDQC, B( K+1, K+NS+NP ), LDB, ZERO, WORK,
! 433: $ SHEIGHT )
! 434: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT, B( K+1,
! 435: $ K+NS+NP ), LDB )
! 436: END IF
! 437: IF ( ILQ ) THEN
! 438: CALL DGEMM( 'N', 'N', N, NBLOCK, NBLOCK, ONE, Q( 1, K+1 ),
! 439: $ LDQ, QC, LDQC, ZERO, WORK, N )
! 440: CALL DLACPY( 'ALL', N, NBLOCK, WORK, N, Q( 1, K+1 ), LDQ )
! 441: END IF
! 442:
! 443: * Update A(istartm:k,k:k+ns+npos-1) and B(istartm:k,k:k+ns+npos-1)
! 444: * from the right with Zc(1:ns+np,1:ns+np)
! 445: SHEIGHT = K-ISTARTM+1
! 446: SWIDTH = NBLOCK
! 447: IF ( SHEIGHT > 0 ) THEN
! 448: CALL DGEMM( 'N', 'N', SHEIGHT, SWIDTH, SWIDTH, ONE,
! 449: $ A( ISTARTM, K ), LDA, ZC, LDZC, ZERO, WORK,
! 450: $ SHEIGHT )
! 451: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT,
! 452: $ A( ISTARTM, K ), LDA )
! 453: CALL DGEMM( 'N', 'N', SHEIGHT, SWIDTH, SWIDTH, ONE,
! 454: $ B( ISTARTM, K ), LDB, ZC, LDZC, ZERO, WORK,
! 455: $ SHEIGHT )
! 456: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT,
! 457: $ B( ISTARTM, K ), LDB )
! 458: END IF
! 459: IF ( ILZ ) THEN
! 460: CALL DGEMM( 'N', 'N', N, NBLOCK, NBLOCK, ONE, Z( 1, K ),
! 461: $ LDZ, ZC, LDZC, ZERO, WORK, N )
! 462: CALL DLACPY( 'ALL', N, NBLOCK, WORK, N, Z( 1, K ), LDZ )
! 463: END IF
! 464:
! 465: K = K+NP
! 466:
! 467: END DO
! 468:
! 469: * The following block removes the shifts from the bottom right corner
! 470: * one by one. Updates are initially applied to A(ihi-ns+1:ihi,ihi-ns:ihi).
! 471:
! 472: CALL DLASET( 'FULL', NS, NS, ZERO, ONE, QC, LDQC )
! 473: CALL DLASET( 'FULL', NS+1, NS+1, ZERO, ONE, ZC, LDZC )
! 474:
! 475: * istartb points to the first row we will be updating
! 476: ISTARTB = IHI-NS+1
! 477: * istopb points to the last column we will be updating
! 478: ISTOPB = IHI
! 479:
! 480: DO I = 1, NS, 2
! 481: * Chase the shift down to the bottom right corner
! 482: DO ISHIFT = IHI-I-1, IHI-2
! 483: CALL DLAQZ2( .TRUE., .TRUE., ISHIFT, ISTARTB, ISTOPB, IHI,
! 484: $ A, LDA, B, LDB, NS, IHI-NS+1, QC, LDQC, NS+1,
! 485: $ IHI-NS, ZC, LDZC )
! 486: END DO
! 487:
! 488: END DO
! 489:
! 490: * Update rest of the pencil
! 491:
! 492: * Update A(ihi-ns+1:ihi, ihi+1:istopm)
! 493: * from the left with Qc(1:ns,1:ns)'
! 494: SHEIGHT = NS
! 495: SWIDTH = ISTOPM-( IHI+1 )+1
! 496: IF ( SWIDTH > 0 ) THEN
! 497: CALL DGEMM( 'T', 'N', SHEIGHT, SWIDTH, SHEIGHT, ONE, QC, LDQC,
! 498: $ A( IHI-NS+1, IHI+1 ), LDA, ZERO, WORK, SHEIGHT )
! 499: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT,
! 500: $ A( IHI-NS+1, IHI+1 ), LDA )
! 501: CALL DGEMM( 'T', 'N', SHEIGHT, SWIDTH, SHEIGHT, ONE, QC, LDQC,
! 502: $ B( IHI-NS+1, IHI+1 ), LDB, ZERO, WORK, SHEIGHT )
! 503: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT,
! 504: $ B( IHI-NS+1, IHI+1 ), LDB )
! 505: END IF
! 506: IF ( ILQ ) THEN
! 507: CALL DGEMM( 'N', 'N', N, NS, NS, ONE, Q( 1, IHI-NS+1 ), LDQ,
! 508: $ QC, LDQC, ZERO, WORK, N )
! 509: CALL DLACPY( 'ALL', N, NS, WORK, N, Q( 1, IHI-NS+1 ), LDQ )
! 510: END IF
! 511:
! 512: * Update A(istartm:ihi-ns,ihi-ns:ihi)
! 513: * from the right with Zc(1:ns+1,1:ns+1)
! 514: SHEIGHT = IHI-NS-ISTARTM+1
! 515: SWIDTH = NS+1
! 516: IF ( SHEIGHT > 0 ) THEN
! 517: CALL DGEMM( 'N', 'N', SHEIGHT, SWIDTH, SWIDTH, ONE, A( ISTARTM,
! 518: $ IHI-NS ), LDA, ZC, LDZC, ZERO, WORK, SHEIGHT )
! 519: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT, A( ISTARTM,
! 520: $ IHI-NS ), LDA )
! 521: CALL DGEMM( 'N', 'N', SHEIGHT, SWIDTH, SWIDTH, ONE, B( ISTARTM,
! 522: $ IHI-NS ), LDB, ZC, LDZC, ZERO, WORK, SHEIGHT )
! 523: CALL DLACPY( 'ALL', SHEIGHT, SWIDTH, WORK, SHEIGHT, B( ISTARTM,
! 524: $ IHI-NS ), LDB )
! 525: END IF
! 526: IF ( ILZ ) THEN
! 527: CALL DGEMM( 'N', 'N', N, NS+1, NS+1, ONE, Z( 1, IHI-NS ), LDZ,
! 528: $ ZC, LDZC, ZERO, WORK, N )
! 529: CALL DLACPY( 'ALL', N, NS+1, WORK, N, Z( 1, IHI-NS ), LDZ )
! 530: END IF
! 531:
! 532: END SUBROUTINE
CVSweb interface <joel.bertrand@systella.fr>