1*> \brief \b DLAEDA used by DSTEDC. Computes the Z vector determining the rank-one modification of the diagonal matrix. Used when the original matrix is dense. 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8*> \htmlonly 9*> Download DLAEDA + dependencies 10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaeda.f"> 11*> [TGZ]</a> 12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaeda.f"> 13*> [ZIP]</a> 14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaeda.f"> 15*> [TXT]</a> 16*> \endhtmlonly 17* 18* Definition: 19* =========== 20* 21* SUBROUTINE DLAEDA( N, TLVLS, CURLVL, CURPBM, PRMPTR, PERM, GIVPTR, 22* GIVCOL, GIVNUM, Q, QPTR, Z, ZTEMP, INFO ) 23* 24* .. Scalar Arguments .. 25* INTEGER CURLVL, CURPBM, INFO, N, TLVLS 26* .. 27* .. Array Arguments .. 28* INTEGER GIVCOL( 2, * ), GIVPTR( * ), PERM( * ), 29* $ PRMPTR( * ), QPTR( * ) 30* DOUBLE PRECISION GIVNUM( 2, * ), Q( * ), Z( * ), ZTEMP( * ) 31* .. 32* 33* 34*> \par Purpose: 35* ============= 36*> 37*> \verbatim 38*> 39*> DLAEDA computes the Z vector corresponding to the merge step in the 40*> CURLVLth step of the merge process with TLVLS steps for the CURPBMth 41*> problem. 42*> \endverbatim 43* 44* Arguments: 45* ========== 46* 47*> \param[in] N 48*> \verbatim 49*> N is INTEGER 50*> The dimension of the symmetric tridiagonal matrix. N >= 0. 51*> \endverbatim 52*> 53*> \param[in] TLVLS 54*> \verbatim 55*> TLVLS is INTEGER 56*> The total number of merging levels in the overall divide and 57*> conquer tree. 58*> \endverbatim 59*> 60*> \param[in] CURLVL 61*> \verbatim 62*> CURLVL is INTEGER 63*> The current level in the overall merge routine, 64*> 0 <= curlvl <= tlvls. 65*> \endverbatim 66*> 67*> \param[in] CURPBM 68*> \verbatim 69*> CURPBM is INTEGER 70*> The current problem in the current level in the overall 71*> merge routine (counting from upper left to lower right). 72*> \endverbatim 73*> 74*> \param[in] PRMPTR 75*> \verbatim 76*> PRMPTR is INTEGER array, dimension (N lg N) 77*> Contains a list of pointers which indicate where in PERM a 78*> level's permutation is stored. PRMPTR(i+1) - PRMPTR(i) 79*> indicates the size of the permutation and incidentally the 80*> size of the full, non-deflated problem. 81*> \endverbatim 82*> 83*> \param[in] PERM 84*> \verbatim 85*> PERM is INTEGER array, dimension (N lg N) 86*> Contains the permutations (from deflation and sorting) to be 87*> applied to each eigenblock. 88*> \endverbatim 89*> 90*> \param[in] GIVPTR 91*> \verbatim 92*> GIVPTR is INTEGER array, dimension (N lg N) 93*> Contains a list of pointers which indicate where in GIVCOL a 94*> level's Givens rotations are stored. GIVPTR(i+1) - GIVPTR(i) 95*> indicates the number of Givens rotations. 96*> \endverbatim 97*> 98*> \param[in] GIVCOL 99*> \verbatim 100*> GIVCOL is INTEGER array, dimension (2, N lg N) 101*> Each pair of numbers indicates a pair of columns to take place 102*> in a Givens rotation. 103*> \endverbatim 104*> 105*> \param[in] GIVNUM 106*> \verbatim 107*> GIVNUM is DOUBLE PRECISION array, dimension (2, N lg N) 108*> Each number indicates the S value to be used in the 109*> corresponding Givens rotation. 110*> \endverbatim 111*> 112*> \param[in] Q 113*> \verbatim 114*> Q is DOUBLE PRECISION array, dimension (N**2) 115*> Contains the square eigenblocks from previous levels, the 116*> starting positions for blocks are given by QPTR. 117*> \endverbatim 118*> 119*> \param[in] QPTR 120*> \verbatim 121*> QPTR is INTEGER array, dimension (N+2) 122*> Contains a list of pointers which indicate where in Q an 123*> eigenblock is stored. SQRT( QPTR(i+1) - QPTR(i) ) indicates 124*> the size of the block. 125*> \endverbatim 126*> 127*> \param[out] Z 128*> \verbatim 129*> Z is DOUBLE PRECISION array, dimension (N) 130*> On output this vector contains the updating vector (the last 131*> row of the first sub-eigenvector matrix and the first row of 132*> the second sub-eigenvector matrix). 133*> \endverbatim 134*> 135*> \param[out] ZTEMP 136*> \verbatim 137*> ZTEMP is DOUBLE PRECISION array, dimension (N) 138*> \endverbatim 139*> 140*> \param[out] INFO 141*> \verbatim 142*> INFO is INTEGER 143*> = 0: successful exit. 144*> < 0: if INFO = -i, the i-th argument had an illegal value. 145*> \endverbatim 146* 147* Authors: 148* ======== 149* 150*> \author Univ. of Tennessee 151*> \author Univ. of California Berkeley 152*> \author Univ. of Colorado Denver 153*> \author NAG Ltd. 154* 155*> \ingroup auxOTHERcomputational 156* 157*> \par Contributors: 158* ================== 159*> 160*> Jeff Rutter, Computer Science Division, University of California 161*> at Berkeley, USA 162* 163* ===================================================================== 164 SUBROUTINE DLAEDA( N, TLVLS, CURLVL, CURPBM, PRMPTR, PERM, GIVPTR, 165 $ GIVCOL, GIVNUM, Q, QPTR, Z, ZTEMP, INFO ) 166* 167* -- LAPACK computational routine -- 168* -- LAPACK is a software package provided by Univ. of Tennessee, -- 169* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 170* 171* .. Scalar Arguments .. 172 INTEGER CURLVL, CURPBM, INFO, N, TLVLS 173* .. 174* .. Array Arguments .. 175 INTEGER GIVCOL( 2, * ), GIVPTR( * ), PERM( * ), 176 $ PRMPTR( * ), QPTR( * ) 177 DOUBLE PRECISION GIVNUM( 2, * ), Q( * ), Z( * ), ZTEMP( * ) 178* .. 179* 180* ===================================================================== 181* 182* .. Parameters .. 183 DOUBLE PRECISION ZERO, HALF, ONE 184 PARAMETER ( ZERO = 0.0D0, HALF = 0.5D0, ONE = 1.0D0 ) 185* .. 186* .. Local Scalars .. 187 INTEGER BSIZ1, BSIZ2, CURR, I, K, MID, PSIZ1, PSIZ2, 188 $ PTR, ZPTR1 189* .. 190* .. External Subroutines .. 191 EXTERNAL DCOPY, DGEMV, DROT, XERBLA 192* .. 193* .. Intrinsic Functions .. 194 INTRINSIC DBLE, INT, SQRT 195* .. 196* .. Executable Statements .. 197* 198* Test the input parameters. 199* 200 INFO = 0 201* 202 IF( N.LT.0 ) THEN 203 INFO = -1 204 END IF 205 IF( INFO.NE.0 ) THEN 206 CALL XERBLA( 'DLAEDA', -INFO ) 207 RETURN 208 END IF 209* 210* Quick return if possible 211* 212 IF( N.EQ.0 ) 213 $ RETURN 214* 215* Determine location of first number in second half. 216* 217 MID = N / 2 + 1 218* 219* Gather last/first rows of appropriate eigenblocks into center of Z 220* 221 PTR = 1 222* 223* Determine location of lowest level subproblem in the full storage 224* scheme 225* 226 CURR = PTR + CURPBM*2**CURLVL + 2**( CURLVL-1 ) - 1 227* 228* Determine size of these matrices. We add HALF to the value of 229* the SQRT in case the machine underestimates one of these square 230* roots. 231* 232 BSIZ1 = INT( HALF+SQRT( DBLE( QPTR( CURR+1 )-QPTR( CURR ) ) ) ) 233 BSIZ2 = INT( HALF+SQRT( DBLE( QPTR( CURR+2 )-QPTR( CURR+1 ) ) ) ) 234 DO 10 K = 1, MID - BSIZ1 - 1 235 Z( K ) = ZERO 236 10 CONTINUE 237 CALL DCOPY( BSIZ1, Q( QPTR( CURR )+BSIZ1-1 ), BSIZ1, 238 $ Z( MID-BSIZ1 ), 1 ) 239 CALL DCOPY( BSIZ2, Q( QPTR( CURR+1 ) ), BSIZ2, Z( MID ), 1 ) 240 DO 20 K = MID + BSIZ2, N 241 Z( K ) = ZERO 242 20 CONTINUE 243* 244* Loop through remaining levels 1 -> CURLVL applying the Givens 245* rotations and permutation and then multiplying the center matrices 246* against the current Z. 247* 248 PTR = 2**TLVLS + 1 249 DO 70 K = 1, CURLVL - 1 250 CURR = PTR + CURPBM*2**( CURLVL-K ) + 2**( CURLVL-K-1 ) - 1 251 PSIZ1 = PRMPTR( CURR+1 ) - PRMPTR( CURR ) 252 PSIZ2 = PRMPTR( CURR+2 ) - PRMPTR( CURR+1 ) 253 ZPTR1 = MID - PSIZ1 254* 255* Apply Givens at CURR and CURR+1 256* 257 DO 30 I = GIVPTR( CURR ), GIVPTR( CURR+1 ) - 1 258 CALL DROT( 1, Z( ZPTR1+GIVCOL( 1, I )-1 ), 1, 259 $ Z( ZPTR1+GIVCOL( 2, I )-1 ), 1, GIVNUM( 1, I ), 260 $ GIVNUM( 2, I ) ) 261 30 CONTINUE 262 DO 40 I = GIVPTR( CURR+1 ), GIVPTR( CURR+2 ) - 1 263 CALL DROT( 1, Z( MID-1+GIVCOL( 1, I ) ), 1, 264 $ Z( MID-1+GIVCOL( 2, I ) ), 1, GIVNUM( 1, I ), 265 $ GIVNUM( 2, I ) ) 266 40 CONTINUE 267 PSIZ1 = PRMPTR( CURR+1 ) - PRMPTR( CURR ) 268 PSIZ2 = PRMPTR( CURR+2 ) - PRMPTR( CURR+1 ) 269 DO 50 I = 0, PSIZ1 - 1 270 ZTEMP( I+1 ) = Z( ZPTR1+PERM( PRMPTR( CURR )+I )-1 ) 271 50 CONTINUE 272 DO 60 I = 0, PSIZ2 - 1 273 ZTEMP( PSIZ1+I+1 ) = Z( MID+PERM( PRMPTR( CURR+1 )+I )-1 ) 274 60 CONTINUE 275* 276* Multiply Blocks at CURR and CURR+1 277* 278* Determine size of these matrices. We add HALF to the value of 279* the SQRT in case the machine underestimates one of these 280* square roots. 281* 282 BSIZ1 = INT( HALF+SQRT( DBLE( QPTR( CURR+1 )-QPTR( CURR ) ) ) ) 283 BSIZ2 = INT( HALF+SQRT( DBLE( QPTR( CURR+2 )-QPTR( CURR+ 284 $ 1 ) ) ) ) 285 IF( BSIZ1.GT.0 ) THEN 286 CALL DGEMV( 'T', BSIZ1, BSIZ1, ONE, Q( QPTR( CURR ) ), 287 $ BSIZ1, ZTEMP( 1 ), 1, ZERO, Z( ZPTR1 ), 1 ) 288 END IF 289 CALL DCOPY( PSIZ1-BSIZ1, ZTEMP( BSIZ1+1 ), 1, Z( ZPTR1+BSIZ1 ), 290 $ 1 ) 291 IF( BSIZ2.GT.0 ) THEN 292 CALL DGEMV( 'T', BSIZ2, BSIZ2, ONE, Q( QPTR( CURR+1 ) ), 293 $ BSIZ2, ZTEMP( PSIZ1+1 ), 1, ZERO, Z( MID ), 1 ) 294 END IF 295 CALL DCOPY( PSIZ2-BSIZ2, ZTEMP( PSIZ1+BSIZ2+1 ), 1, 296 $ Z( MID+BSIZ2 ), 1 ) 297* 298 PTR = PTR + 2**( TLVLS-K ) 299 70 CONTINUE 300* 301 RETURN 302* 303* End of DLAEDA 304* 305 END 306