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