1*> \brief \b ZLATM6 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8* Definition: 9* =========== 10* 11* SUBROUTINE ZLATM6( TYPE, N, A, LDA, B, X, LDX, Y, LDY, ALPHA, 12* BETA, WX, WY, S, DIF ) 13* 14* .. Scalar Arguments .. 15* INTEGER LDA, LDX, LDY, N, TYPE 16* COMPLEX*16 ALPHA, BETA, WX, WY 17* .. 18* .. Array Arguments .. 19* DOUBLE PRECISION DIF( * ), S( * ) 20* COMPLEX*16 A( LDA, * ), B( LDA, * ), X( LDX, * ), 21* $ Y( LDY, * ) 22* .. 23* 24* 25*> \par Purpose: 26* ============= 27*> 28*> \verbatim 29*> 30*> ZLATM6 generates test matrices for the generalized eigenvalue 31*> problem, their corresponding right and left eigenvector matrices, 32*> and also reciprocal condition numbers for all eigenvalues and 33*> the reciprocal condition numbers of eigenvectors corresponding to 34*> the 1th and 5th eigenvalues. 35*> 36*> Test Matrices 37*> ============= 38*> 39*> Two kinds of test matrix pairs 40*> (A, B) = inverse(YH) * (Da, Db) * inverse(X) 41*> are used in the tests: 42*> 43*> Type 1: 44*> Da = 1+a 0 0 0 0 Db = 1 0 0 0 0 45*> 0 2+a 0 0 0 0 1 0 0 0 46*> 0 0 3+a 0 0 0 0 1 0 0 47*> 0 0 0 4+a 0 0 0 0 1 0 48*> 0 0 0 0 5+a , 0 0 0 0 1 49*> and Type 2: 50*> Da = 1+i 0 0 0 0 Db = 1 0 0 0 0 51*> 0 1-i 0 0 0 0 1 0 0 0 52*> 0 0 1 0 0 0 0 1 0 0 53*> 0 0 0 (1+a)+(1+b)i 0 0 0 0 1 0 54*> 0 0 0 0 (1+a)-(1+b)i, 0 0 0 0 1 . 55*> 56*> In both cases the same inverse(YH) and inverse(X) are used to compute 57*> (A, B), giving the exact eigenvectors to (A,B) as (YH, X): 58*> 59*> YH: = 1 0 -y y -y X = 1 0 -x -x x 60*> 0 1 -y y -y 0 1 x -x -x 61*> 0 0 1 0 0 0 0 1 0 0 62*> 0 0 0 1 0 0 0 0 1 0 63*> 0 0 0 0 1, 0 0 0 0 1 , where 64*> 65*> a, b, x and y will have all values independently of each other. 66*> \endverbatim 67* 68* Arguments: 69* ========== 70* 71*> \param[in] TYPE 72*> \verbatim 73*> TYPE is INTEGER 74*> Specifies the problem type (see futher details). 75*> \endverbatim 76*> 77*> \param[in] N 78*> \verbatim 79*> N is INTEGER 80*> Size of the matrices A and B. 81*> \endverbatim 82*> 83*> \param[out] A 84*> \verbatim 85*> A is COMPLEX*16 array, dimension (LDA, N). 86*> On exit A N-by-N is initialized according to TYPE. 87*> \endverbatim 88*> 89*> \param[in] LDA 90*> \verbatim 91*> LDA is INTEGER 92*> The leading dimension of A and of B. 93*> \endverbatim 94*> 95*> \param[out] B 96*> \verbatim 97*> B is COMPLEX*16 array, dimension (LDA, N). 98*> On exit B N-by-N is initialized according to TYPE. 99*> \endverbatim 100*> 101*> \param[out] X 102*> \verbatim 103*> X is COMPLEX*16 array, dimension (LDX, N). 104*> On exit X is the N-by-N matrix of right eigenvectors. 105*> \endverbatim 106*> 107*> \param[in] LDX 108*> \verbatim 109*> LDX is INTEGER 110*> The leading dimension of X. 111*> \endverbatim 112*> 113*> \param[out] Y 114*> \verbatim 115*> Y is COMPLEX*16 array, dimension (LDY, N). 116*> On exit Y is the N-by-N matrix of left eigenvectors. 117*> \endverbatim 118*> 119*> \param[in] LDY 120*> \verbatim 121*> LDY is INTEGER 122*> The leading dimension of Y. 123*> \endverbatim 124*> 125*> \param[in] ALPHA 126*> \verbatim 127*> ALPHA is COMPLEX*16 128*> \endverbatim 129*> 130*> \param[in] BETA 131*> \verbatim 132*> BETA is COMPLEX*16 133*> \verbatim 134*> Weighting constants for matrix A. 135*> \endverbatim 136*> 137*> \param[in] WX 138*> \verbatim 139*> WX is COMPLEX*16 140*> Constant for right eigenvector matrix. 141*> \endverbatim 142*> 143*> \param[in] WY 144*> \verbatim 145*> WY is COMPLEX*16 146*> Constant for left eigenvector matrix. 147*> \endverbatim 148*> 149*> \param[out] S 150*> \verbatim 151*> S is DOUBLE PRECISION array, dimension (N) 152*> S(i) is the reciprocal condition number for eigenvalue i. 153*> \endverbatim 154*> 155*> \param[out] DIF 156*> \verbatim 157*> DIF is DOUBLE PRECISION array, dimension (N) 158*> DIF(i) is the reciprocal condition number for eigenvector i. 159*> \endverbatim 160* 161* Authors: 162* ======== 163* 164*> \author Univ. of Tennessee 165*> \author Univ. of California Berkeley 166*> \author Univ. of Colorado Denver 167*> \author NAG Ltd. 168* 169*> \date November 2011 170* 171*> \ingroup complex16_matgen 172* 173* ===================================================================== 174 SUBROUTINE ZLATM6( TYPE, N, A, LDA, B, X, LDX, Y, LDY, ALPHA, 175 $ BETA, WX, WY, S, DIF ) 176* 177* -- LAPACK computational routine (version 3.4.0) -- 178* -- LAPACK is a software package provided by Univ. of Tennessee, -- 179* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 180* November 2011 181* 182* .. Scalar Arguments .. 183 INTEGER LDA, LDX, LDY, N, TYPE 184 COMPLEX*16 ALPHA, BETA, WX, WY 185* .. 186* .. Array Arguments .. 187 DOUBLE PRECISION DIF( * ), S( * ) 188 COMPLEX*16 A( LDA, * ), B( LDA, * ), X( LDX, * ), 189 $ Y( LDY, * ) 190* .. 191* 192* ===================================================================== 193* 194* .. Parameters .. 195 DOUBLE PRECISION RONE, TWO, THREE 196 PARAMETER ( RONE = 1.0D+0, TWO = 2.0D+0, THREE = 3.0D+0 ) 197 COMPLEX*16 ZERO, ONE 198 PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ), 199 $ ONE = ( 1.0D+0, 0.0D+0 ) ) 200* .. 201* .. Local Scalars .. 202 INTEGER I, INFO, J 203* .. 204* .. Local Arrays .. 205 DOUBLE PRECISION RWORK( 50 ) 206 COMPLEX*16 WORK( 26 ), Z( 8, 8 ) 207* .. 208* .. Intrinsic Functions .. 209 INTRINSIC CDABS, DBLE, DCMPLX, DCONJG, SQRT 210* .. 211* .. External Subroutines .. 212 EXTERNAL ZGESVD, ZLACPY, ZLAKF2 213* .. 214* .. Executable Statements .. 215* 216* Generate test problem ... 217* (Da, Db) ... 218* 219 DO 20 I = 1, N 220 DO 10 J = 1, N 221* 222 IF( I.EQ.J ) THEN 223 A( I, I ) = DCMPLX( I ) + ALPHA 224 B( I, I ) = ONE 225 ELSE 226 A( I, J ) = ZERO 227 B( I, J ) = ZERO 228 END IF 229* 230 10 CONTINUE 231 20 CONTINUE 232 IF( TYPE.EQ.2 ) THEN 233 A( 1, 1 ) = DCMPLX( RONE, RONE ) 234 A( 2, 2 ) = DCONJG( A( 1, 1 ) ) 235 A( 3, 3 ) = ONE 236 A( 4, 4 ) = DCMPLX( DBLE( ONE+ALPHA ), DBLE( ONE+BETA ) ) 237 A( 5, 5 ) = DCONJG( A( 4, 4 ) ) 238 END IF 239* 240* Form X and Y 241* 242 CALL ZLACPY( 'F', N, N, B, LDA, Y, LDY ) 243 Y( 3, 1 ) = -DCONJG( WY ) 244 Y( 4, 1 ) = DCONJG( WY ) 245 Y( 5, 1 ) = -DCONJG( WY ) 246 Y( 3, 2 ) = -DCONJG( WY ) 247 Y( 4, 2 ) = DCONJG( WY ) 248 Y( 5, 2 ) = -DCONJG( WY ) 249* 250 CALL ZLACPY( 'F', N, N, B, LDA, X, LDX ) 251 X( 1, 3 ) = -WX 252 X( 1, 4 ) = -WX 253 X( 1, 5 ) = WX 254 X( 2, 3 ) = WX 255 X( 2, 4 ) = -WX 256 X( 2, 5 ) = -WX 257* 258* Form (A, B) 259* 260 B( 1, 3 ) = WX + WY 261 B( 2, 3 ) = -WX + WY 262 B( 1, 4 ) = WX - WY 263 B( 2, 4 ) = WX - WY 264 B( 1, 5 ) = -WX + WY 265 B( 2, 5 ) = WX + WY 266 A( 1, 3 ) = WX*A( 1, 1 ) + WY*A( 3, 3 ) 267 A( 2, 3 ) = -WX*A( 2, 2 ) + WY*A( 3, 3 ) 268 A( 1, 4 ) = WX*A( 1, 1 ) - WY*A( 4, 4 ) 269 A( 2, 4 ) = WX*A( 2, 2 ) - WY*A( 4, 4 ) 270 A( 1, 5 ) = -WX*A( 1, 1 ) + WY*A( 5, 5 ) 271 A( 2, 5 ) = WX*A( 2, 2 ) + WY*A( 5, 5 ) 272* 273* Compute condition numbers 274* 275 S( 1 ) = RONE / SQRT( ( RONE+THREE*CDABS( WY )*CDABS( WY ) ) / 276 $ ( RONE+CDABS( A( 1, 1 ) )*CDABS( A( 1, 1 ) ) ) ) 277 S( 2 ) = RONE / SQRT( ( RONE+THREE*CDABS( WY )*CDABS( WY ) ) / 278 $ ( RONE+CDABS( A( 2, 2 ) )*CDABS( A( 2, 2 ) ) ) ) 279 S( 3 ) = RONE / SQRT( ( RONE+TWO*CDABS( WX )*CDABS( WX ) ) / 280 $ ( RONE+CDABS( A( 3, 3 ) )*CDABS( A( 3, 3 ) ) ) ) 281 S( 4 ) = RONE / SQRT( ( RONE+TWO*CDABS( WX )*CDABS( WX ) ) / 282 $ ( RONE+CDABS( A( 4, 4 ) )*CDABS( A( 4, 4 ) ) ) ) 283 S( 5 ) = RONE / SQRT( ( RONE+TWO*CDABS( WX )*CDABS( WX ) ) / 284 $ ( RONE+CDABS( A( 5, 5 ) )*CDABS( A( 5, 5 ) ) ) ) 285* 286 CALL ZLAKF2( 1, 4, A, LDA, A( 2, 2 ), B, B( 2, 2 ), Z, 8 ) 287 CALL ZGESVD( 'N', 'N', 8, 8, Z, 8, RWORK, WORK, 1, WORK( 2 ), 1, 288 $ WORK( 3 ), 24, RWORK( 9 ), INFO ) 289 DIF( 1 ) = RWORK( 8 ) 290* 291 CALL ZLAKF2( 4, 1, A, LDA, A( 5, 5 ), B, B( 5, 5 ), Z, 8 ) 292 CALL ZGESVD( 'N', 'N', 8, 8, Z, 8, RWORK, WORK, 1, WORK( 2 ), 1, 293 $ WORK( 3 ), 24, RWORK( 9 ), INFO ) 294 DIF( 5 ) = RWORK( 8 ) 295* 296 RETURN 297* 298* End of ZLATM6 299* 300 END 301