Annotation of rpl/lapack/lapack/dlaeda.f, revision 1.1
1.1 ! bertrand 1: SUBROUTINE DLAEDA( N, TLVLS, CURLVL, CURPBM, PRMPTR, PERM, GIVPTR,
! 2: $ GIVCOL, GIVNUM, Q, QPTR, Z, ZTEMP, INFO )
! 3: *
! 4: * -- LAPACK routine (version 3.2) --
! 5: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 6: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 7: * November 2006
! 8: *
! 9: * .. Scalar Arguments ..
! 10: INTEGER CURLVL, CURPBM, INFO, N, TLVLS
! 11: * ..
! 12: * .. Array Arguments ..
! 13: INTEGER GIVCOL( 2, * ), GIVPTR( * ), PERM( * ),
! 14: $ PRMPTR( * ), QPTR( * )
! 15: DOUBLE PRECISION GIVNUM( 2, * ), Q( * ), Z( * ), ZTEMP( * )
! 16: * ..
! 17: *
! 18: * Purpose
! 19: * =======
! 20: *
! 21: * DLAEDA computes the Z vector corresponding to the merge step in the
! 22: * CURLVLth step of the merge process with TLVLS steps for the CURPBMth
! 23: * problem.
! 24: *
! 25: * Arguments
! 26: * =========
! 27: *
! 28: * N (input) INTEGER
! 29: * The dimension of the symmetric tridiagonal matrix. N >= 0.
! 30: *
! 31: * TLVLS (input) INTEGER
! 32: * The total number of merging levels in the overall divide and
! 33: * conquer tree.
! 34: *
! 35: * CURLVL (input) INTEGER
! 36: * The current level in the overall merge routine,
! 37: * 0 <= curlvl <= tlvls.
! 38: *
! 39: * CURPBM (input) INTEGER
! 40: * The current problem in the current level in the overall
! 41: * merge routine (counting from upper left to lower right).
! 42: *
! 43: * PRMPTR (input) INTEGER array, dimension (N lg N)
! 44: * Contains a list of pointers which indicate where in PERM a
! 45: * level's permutation is stored. PRMPTR(i+1) - PRMPTR(i)
! 46: * indicates the size of the permutation and incidentally the
! 47: * size of the full, non-deflated problem.
! 48: *
! 49: * PERM (input) INTEGER array, dimension (N lg N)
! 50: * Contains the permutations (from deflation and sorting) to be
! 51: * applied to each eigenblock.
! 52: *
! 53: * GIVPTR (input) INTEGER array, dimension (N lg N)
! 54: * Contains a list of pointers which indicate where in GIVCOL a
! 55: * level's Givens rotations are stored. GIVPTR(i+1) - GIVPTR(i)
! 56: * indicates the number of Givens rotations.
! 57: *
! 58: * GIVCOL (input) INTEGER array, dimension (2, N lg N)
! 59: * Each pair of numbers indicates a pair of columns to take place
! 60: * in a Givens rotation.
! 61: *
! 62: * GIVNUM (input) DOUBLE PRECISION array, dimension (2, N lg N)
! 63: * Each number indicates the S value to be used in the
! 64: * corresponding Givens rotation.
! 65: *
! 66: * Q (input) DOUBLE PRECISION array, dimension (N**2)
! 67: * Contains the square eigenblocks from previous levels, the
! 68: * starting positions for blocks are given by QPTR.
! 69: *
! 70: * QPTR (input) INTEGER array, dimension (N+2)
! 71: * Contains a list of pointers which indicate where in Q an
! 72: * eigenblock is stored. SQRT( QPTR(i+1) - QPTR(i) ) indicates
! 73: * the size of the block.
! 74: *
! 75: * Z (output) DOUBLE PRECISION array, dimension (N)
! 76: * On output this vector contains the updating vector (the last
! 77: * row of the first sub-eigenvector matrix and the first row of
! 78: * the second sub-eigenvector matrix).
! 79: *
! 80: * ZTEMP (workspace) DOUBLE PRECISION array, dimension (N)
! 81: *
! 82: * INFO (output) INTEGER
! 83: * = 0: successful exit.
! 84: * < 0: if INFO = -i, the i-th argument had an illegal value.
! 85: *
! 86: * Further Details
! 87: * ===============
! 88: *
! 89: * Based on contributions by
! 90: * Jeff Rutter, Computer Science Division, University of California
! 91: * at Berkeley, USA
! 92: *
! 93: * =====================================================================
! 94: *
! 95: * .. Parameters ..
! 96: DOUBLE PRECISION ZERO, HALF, ONE
! 97: PARAMETER ( ZERO = 0.0D0, HALF = 0.5D0, ONE = 1.0D0 )
! 98: * ..
! 99: * .. Local Scalars ..
! 100: INTEGER BSIZ1, BSIZ2, CURR, I, K, MID, PSIZ1, PSIZ2,
! 101: $ PTR, ZPTR1
! 102: * ..
! 103: * .. External Subroutines ..
! 104: EXTERNAL DCOPY, DGEMV, DROT, XERBLA
! 105: * ..
! 106: * .. Intrinsic Functions ..
! 107: INTRINSIC DBLE, INT, SQRT
! 108: * ..
! 109: * .. Executable Statements ..
! 110: *
! 111: * Test the input parameters.
! 112: *
! 113: INFO = 0
! 114: *
! 115: IF( N.LT.0 ) THEN
! 116: INFO = -1
! 117: END IF
! 118: IF( INFO.NE.0 ) THEN
! 119: CALL XERBLA( 'DLAEDA', -INFO )
! 120: RETURN
! 121: END IF
! 122: *
! 123: * Quick return if possible
! 124: *
! 125: IF( N.EQ.0 )
! 126: $ RETURN
! 127: *
! 128: * Determine location of first number in second half.
! 129: *
! 130: MID = N / 2 + 1
! 131: *
! 132: * Gather last/first rows of appropriate eigenblocks into center of Z
! 133: *
! 134: PTR = 1
! 135: *
! 136: * Determine location of lowest level subproblem in the full storage
! 137: * scheme
! 138: *
! 139: CURR = PTR + CURPBM*2**CURLVL + 2**( CURLVL-1 ) - 1
! 140: *
! 141: * Determine size of these matrices. We add HALF to the value of
! 142: * the SQRT in case the machine underestimates one of these square
! 143: * roots.
! 144: *
! 145: BSIZ1 = INT( HALF+SQRT( DBLE( QPTR( CURR+1 )-QPTR( CURR ) ) ) )
! 146: BSIZ2 = INT( HALF+SQRT( DBLE( QPTR( CURR+2 )-QPTR( CURR+1 ) ) ) )
! 147: DO 10 K = 1, MID - BSIZ1 - 1
! 148: Z( K ) = ZERO
! 149: 10 CONTINUE
! 150: CALL DCOPY( BSIZ1, Q( QPTR( CURR )+BSIZ1-1 ), BSIZ1,
! 151: $ Z( MID-BSIZ1 ), 1 )
! 152: CALL DCOPY( BSIZ2, Q( QPTR( CURR+1 ) ), BSIZ2, Z( MID ), 1 )
! 153: DO 20 K = MID + BSIZ2, N
! 154: Z( K ) = ZERO
! 155: 20 CONTINUE
! 156: *
! 157: * Loop thru remaining levels 1 -> CURLVL applying the Givens
! 158: * rotations and permutation and then multiplying the center matrices
! 159: * against the current Z.
! 160: *
! 161: PTR = 2**TLVLS + 1
! 162: DO 70 K = 1, CURLVL - 1
! 163: CURR = PTR + CURPBM*2**( CURLVL-K ) + 2**( CURLVL-K-1 ) - 1
! 164: PSIZ1 = PRMPTR( CURR+1 ) - PRMPTR( CURR )
! 165: PSIZ2 = PRMPTR( CURR+2 ) - PRMPTR( CURR+1 )
! 166: ZPTR1 = MID - PSIZ1
! 167: *
! 168: * Apply Givens at CURR and CURR+1
! 169: *
! 170: DO 30 I = GIVPTR( CURR ), GIVPTR( CURR+1 ) - 1
! 171: CALL DROT( 1, Z( ZPTR1+GIVCOL( 1, I )-1 ), 1,
! 172: $ Z( ZPTR1+GIVCOL( 2, I )-1 ), 1, GIVNUM( 1, I ),
! 173: $ GIVNUM( 2, I ) )
! 174: 30 CONTINUE
! 175: DO 40 I = GIVPTR( CURR+1 ), GIVPTR( CURR+2 ) - 1
! 176: CALL DROT( 1, Z( MID-1+GIVCOL( 1, I ) ), 1,
! 177: $ Z( MID-1+GIVCOL( 2, I ) ), 1, GIVNUM( 1, I ),
! 178: $ GIVNUM( 2, I ) )
! 179: 40 CONTINUE
! 180: PSIZ1 = PRMPTR( CURR+1 ) - PRMPTR( CURR )
! 181: PSIZ2 = PRMPTR( CURR+2 ) - PRMPTR( CURR+1 )
! 182: DO 50 I = 0, PSIZ1 - 1
! 183: ZTEMP( I+1 ) = Z( ZPTR1+PERM( PRMPTR( CURR )+I )-1 )
! 184: 50 CONTINUE
! 185: DO 60 I = 0, PSIZ2 - 1
! 186: ZTEMP( PSIZ1+I+1 ) = Z( MID+PERM( PRMPTR( CURR+1 )+I )-1 )
! 187: 60 CONTINUE
! 188: *
! 189: * Multiply Blocks at CURR and CURR+1
! 190: *
! 191: * Determine size of these matrices. We add HALF to the value of
! 192: * the SQRT in case the machine underestimates one of these
! 193: * square roots.
! 194: *
! 195: BSIZ1 = INT( HALF+SQRT( DBLE( QPTR( CURR+1 )-QPTR( CURR ) ) ) )
! 196: BSIZ2 = INT( HALF+SQRT( DBLE( QPTR( CURR+2 )-QPTR( CURR+
! 197: $ 1 ) ) ) )
! 198: IF( BSIZ1.GT.0 ) THEN
! 199: CALL DGEMV( 'T', BSIZ1, BSIZ1, ONE, Q( QPTR( CURR ) ),
! 200: $ BSIZ1, ZTEMP( 1 ), 1, ZERO, Z( ZPTR1 ), 1 )
! 201: END IF
! 202: CALL DCOPY( PSIZ1-BSIZ1, ZTEMP( BSIZ1+1 ), 1, Z( ZPTR1+BSIZ1 ),
! 203: $ 1 )
! 204: IF( BSIZ2.GT.0 ) THEN
! 205: CALL DGEMV( 'T', BSIZ2, BSIZ2, ONE, Q( QPTR( CURR+1 ) ),
! 206: $ BSIZ2, ZTEMP( PSIZ1+1 ), 1, ZERO, Z( MID ), 1 )
! 207: END IF
! 208: CALL DCOPY( PSIZ2-BSIZ2, ZTEMP( PSIZ1+BSIZ2+1 ), 1,
! 209: $ Z( MID+BSIZ2 ), 1 )
! 210: *
! 211: PTR = PTR + 2**( TLVLS-K )
! 212: 70 CONTINUE
! 213: *
! 214: RETURN
! 215: *
! 216: * End of DLAEDA
! 217: *
! 218: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dlaeda.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack