![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b CLATM6 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * SUBROUTINE CLATM6( TYPE, N, A, LDA, B, X, LDX, Y, LDY, ALPHA, 00012 * BETA, WX, WY, S, DIF ) 00013 * 00014 * .. Scalar Arguments .. 00015 * INTEGER LDA, LDX, LDY, N, TYPE 00016 * COMPLEX ALPHA, BETA, WX, WY 00017 * .. 00018 * .. Array Arguments .. 00019 * REAL DIF( * ), S( * ) 00020 * COMPLEX A( LDA, * ), B( LDA, * ), X( LDX, * ), 00021 * $ Y( LDY, * ) 00022 * .. 00023 * 00024 * 00025 *> \par Purpose: 00026 * ============= 00027 *> 00028 *> \verbatim 00029 *> 00030 *> CLATM6 generates test matrices for the generalized eigenvalue 00031 *> problem, their corresponding right and left eigenvector matrices, 00032 *> and also reciprocal condition numbers for all eigenvalues and 00033 *> the reciprocal condition numbers of eigenvectors corresponding to 00034 *> the 1th and 5th eigenvalues. 00035 *> 00036 *> Test Matrices 00037 *> ============= 00038 *> 00039 *> Two kinds of test matrix pairs 00040 *> (A, B) = inverse(YH) * (Da, Db) * inverse(X) 00041 *> are used in the tests: 00042 *> 00043 *> Type 1: 00044 *> Da = 1+a 0 0 0 0 Db = 1 0 0 0 0 00045 *> 0 2+a 0 0 0 0 1 0 0 0 00046 *> 0 0 3+a 0 0 0 0 1 0 0 00047 *> 0 0 0 4+a 0 0 0 0 1 0 00048 *> 0 0 0 0 5+a , 0 0 0 0 1 00049 *> and Type 2: 00050 *> Da = 1+i 0 0 0 0 Db = 1 0 0 0 0 00051 *> 0 1-i 0 0 0 0 1 0 0 0 00052 *> 0 0 1 0 0 0 0 1 0 0 00053 *> 0 0 0 (1+a)+(1+b)i 0 0 0 0 1 0 00054 *> 0 0 0 0 (1+a)-(1+b)i, 0 0 0 0 1 . 00055 *> 00056 *> In both cases the same inverse(YH) and inverse(X) are used to compute 00057 *> (A, B), giving the exact eigenvectors to (A,B) as (YH, X): 00058 *> 00059 *> YH: = 1 0 -y y -y X = 1 0 -x -x x 00060 *> 0 1 -y y -y 0 1 x -x -x 00061 *> 0 0 1 0 0 0 0 1 0 0 00062 *> 0 0 0 1 0 0 0 0 1 0 00063 *> 0 0 0 0 1, 0 0 0 0 1 , where 00064 *> 00065 *> a, b, x and y will have all values independently of each other. 00066 *> \endverbatim 00067 * 00068 * Arguments: 00069 * ========== 00070 * 00071 *> \param[in] TYPE 00072 *> \verbatim 00073 *> TYPE is INTEGER 00074 *> Specifies the problem type (see futher details). 00075 *> \endverbatim 00076 *> 00077 *> \param[in] N 00078 *> \verbatim 00079 *> N is INTEGER 00080 *> Size of the matrices A and B. 00081 *> \endverbatim 00082 *> 00083 *> \param[out] A 00084 *> \verbatim 00085 *> A is COMPLEX array, dimension (LDA, N). 00086 *> On exit A N-by-N is initialized according to TYPE. 00087 *> \endverbatim 00088 *> 00089 *> \param[in] LDA 00090 *> \verbatim 00091 *> LDA is INTEGER 00092 *> The leading dimension of A and of B. 00093 *> \endverbatim 00094 *> 00095 *> \param[out] B 00096 *> \verbatim 00097 *> B is COMPLEX array, dimension (LDA, N). 00098 *> On exit B N-by-N is initialized according to TYPE. 00099 *> \endverbatim 00100 *> 00101 *> \param[out] X 00102 *> \verbatim 00103 *> X is COMPLEX array, dimension (LDX, N). 00104 *> On exit X is the N-by-N matrix of right eigenvectors. 00105 *> \endverbatim 00106 *> 00107 *> \param[in] LDX 00108 *> \verbatim 00109 *> LDX is INTEGER 00110 *> The leading dimension of X. 00111 *> \endverbatim 00112 *> 00113 *> \param[out] Y 00114 *> \verbatim 00115 *> Y is COMPLEX array, dimension (LDY, N). 00116 *> On exit Y is the N-by-N matrix of left eigenvectors. 00117 *> \endverbatim 00118 *> 00119 *> \param[in] LDY 00120 *> \verbatim 00121 *> LDY is INTEGER 00122 *> The leading dimension of Y. 00123 *> \endverbatim 00124 *> 00125 *> \param[in] ALPHA 00126 *> \verbatim 00127 *> ALPHA is COMPLEX 00128 *> \endverbatim 00129 *> 00130 *> \param[in] BETA 00131 *> \verbatim 00132 *> BETA is COMPLEX 00133 *> 00134 *> Weighting constants for matrix A. 00135 *> \endverbatim 00136 *> 00137 *> \param[in] WX 00138 *> \verbatim 00139 *> WX is COMPLEX 00140 *> Constant for right eigenvector matrix. 00141 *> \endverbatim 00142 *> 00143 *> \param[in] WY 00144 *> \verbatim 00145 *> WY is COMPLEX 00146 *> Constant for left eigenvector matrix. 00147 *> \endverbatim 00148 *> 00149 *> \param[out] S 00150 *> \verbatim 00151 *> S is REAL array, dimension (N) 00152 *> S(i) is the reciprocal condition number for eigenvalue i. 00153 *> \endverbatim 00154 *> 00155 *> \param[out] DIF 00156 *> \verbatim 00157 *> DIF is REAL array, dimension (N) 00158 *> DIF(i) is the reciprocal condition number for eigenvector i. 00159 *> \endverbatim 00160 * 00161 * Authors: 00162 * ======== 00163 * 00164 *> \author Univ. of Tennessee 00165 *> \author Univ. of California Berkeley 00166 *> \author Univ. of Colorado Denver 00167 *> \author NAG Ltd. 00168 * 00169 *> \date November 2011 00170 * 00171 *> \ingroup complex_matgen 00172 * 00173 * ===================================================================== 00174 SUBROUTINE CLATM6( TYPE, N, A, LDA, B, X, LDX, Y, LDY, ALPHA, 00175 $ BETA, WX, WY, S, DIF ) 00176 * 00177 * -- LAPACK computational routine (version 3.4.0) -- 00178 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00179 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00180 * November 2011 00181 * 00182 * .. Scalar Arguments .. 00183 INTEGER LDA, LDX, LDY, N, TYPE 00184 COMPLEX ALPHA, BETA, WX, WY 00185 * .. 00186 * .. Array Arguments .. 00187 REAL DIF( * ), S( * ) 00188 COMPLEX A( LDA, * ), B( LDA, * ), X( LDX, * ), 00189 $ Y( LDY, * ) 00190 * .. 00191 * 00192 * ===================================================================== 00193 * 00194 * .. Parameters .. 00195 REAL RONE, TWO, THREE 00196 PARAMETER ( RONE = 1.0E+0, TWO = 2.0E+0, THREE = 3.0E+0 ) 00197 COMPLEX ZERO, ONE 00198 PARAMETER ( ZERO = ( 0.0E+0, 0.0E+0 ), 00199 $ ONE = ( 1.0E+0, 0.0E+0 ) ) 00200 * .. 00201 * .. Local Scalars .. 00202 INTEGER I, INFO, J 00203 * .. 00204 * .. Local Arrays .. 00205 REAL RWORK( 50 ) 00206 COMPLEX WORK( 26 ), Z( 8, 8 ) 00207 * .. 00208 * .. Intrinsic Functions .. 00209 INTRINSIC CABS, CMPLX, CONJG, REAL, SQRT 00210 * .. 00211 * .. External Subroutines .. 00212 EXTERNAL CGESVD, CLACPY, CLAKF2 00213 * .. 00214 * .. Executable Statements .. 00215 * 00216 * Generate test problem ... 00217 * (Da, Db) ... 00218 * 00219 DO 20 I = 1, N 00220 DO 10 J = 1, N 00221 * 00222 IF( I.EQ.J ) THEN 00223 A( I, I ) = CMPLX( I ) + ALPHA 00224 B( I, I ) = ONE 00225 ELSE 00226 A( I, J ) = ZERO 00227 B( I, J ) = ZERO 00228 END IF 00229 * 00230 10 CONTINUE 00231 20 CONTINUE 00232 IF( TYPE.EQ.2 ) THEN 00233 A( 1, 1 ) = CMPLX( RONE, RONE ) 00234 A( 2, 2 ) = CONJG( A( 1, 1 ) ) 00235 A( 3, 3 ) = ONE 00236 A( 4, 4 ) = CMPLX( REAL( ONE+ALPHA ), REAL( ONE+BETA ) ) 00237 A( 5, 5 ) = CONJG( A( 4, 4 ) ) 00238 END IF 00239 * 00240 * Form X and Y 00241 * 00242 CALL CLACPY( 'F', N, N, B, LDA, Y, LDY ) 00243 Y( 3, 1 ) = -CONJG( WY ) 00244 Y( 4, 1 ) = CONJG( WY ) 00245 Y( 5, 1 ) = -CONJG( WY ) 00246 Y( 3, 2 ) = -CONJG( WY ) 00247 Y( 4, 2 ) = CONJG( WY ) 00248 Y( 5, 2 ) = -CONJG( WY ) 00249 * 00250 CALL CLACPY( 'F', N, N, B, LDA, X, LDX ) 00251 X( 1, 3 ) = -WX 00252 X( 1, 4 ) = -WX 00253 X( 1, 5 ) = WX 00254 X( 2, 3 ) = WX 00255 X( 2, 4 ) = -WX 00256 X( 2, 5 ) = -WX 00257 * 00258 * Form (A, B) 00259 * 00260 B( 1, 3 ) = WX + WY 00261 B( 2, 3 ) = -WX + WY 00262 B( 1, 4 ) = WX - WY 00263 B( 2, 4 ) = WX - WY 00264 B( 1, 5 ) = -WX + WY 00265 B( 2, 5 ) = WX + WY 00266 A( 1, 3 ) = WX*A( 1, 1 ) + WY*A( 3, 3 ) 00267 A( 2, 3 ) = -WX*A( 2, 2 ) + WY*A( 3, 3 ) 00268 A( 1, 4 ) = WX*A( 1, 1 ) - WY*A( 4, 4 ) 00269 A( 2, 4 ) = WX*A( 2, 2 ) - WY*A( 4, 4 ) 00270 A( 1, 5 ) = -WX*A( 1, 1 ) + WY*A( 5, 5 ) 00271 A( 2, 5 ) = WX*A( 2, 2 ) + WY*A( 5, 5 ) 00272 * 00273 * Compute condition numbers 00274 * 00275 S( 1 ) = RONE / SQRT( ( RONE+THREE*CABS( WY )*CABS( WY ) ) / 00276 $ ( RONE+CABS( A( 1, 1 ) )*CABS( A( 1, 1 ) ) ) ) 00277 S( 2 ) = RONE / SQRT( ( RONE+THREE*CABS( WY )*CABS( WY ) ) / 00278 $ ( RONE+CABS( A( 2, 2 ) )*CABS( A( 2, 2 ) ) ) ) 00279 S( 3 ) = RONE / SQRT( ( RONE+TWO*CABS( WX )*CABS( WX ) ) / 00280 $ ( RONE+CABS( A( 3, 3 ) )*CABS( A( 3, 3 ) ) ) ) 00281 S( 4 ) = RONE / SQRT( ( RONE+TWO*CABS( WX )*CABS( WX ) ) / 00282 $ ( RONE+CABS( A( 4, 4 ) )*CABS( A( 4, 4 ) ) ) ) 00283 S( 5 ) = RONE / SQRT( ( RONE+TWO*CABS( WX )*CABS( WX ) ) / 00284 $ ( RONE+CABS( A( 5, 5 ) )*CABS( A( 5, 5 ) ) ) ) 00285 * 00286 CALL CLAKF2( 1, 4, A, LDA, A( 2, 2 ), B, B( 2, 2 ), Z, 8 ) 00287 CALL CGESVD( 'N', 'N', 8, 8, Z, 8, RWORK, WORK, 1, WORK( 2 ), 1, 00288 $ WORK( 3 ), 24, RWORK( 9 ), INFO ) 00289 DIF( 1 ) = RWORK( 8 ) 00290 * 00291 CALL CLAKF2( 4, 1, A, LDA, A( 5, 5 ), B, B( 5, 5 ), Z, 8 ) 00292 CALL CGESVD( 'N', 'N', 8, 8, Z, 8, RWORK, WORK, 1, WORK( 2 ), 1, 00293 $ WORK( 3 ), 24, RWORK( 9 ), INFO ) 00294 DIF( 5 ) = RWORK( 8 ) 00295 * 00296 RETURN 00297 * 00298 * End of CLATM6 00299 * 00300 END