|
LAPACK 3.3.0
|
|
LAPACK 3.3.0
|
00001 SUBROUTINE DLAGS2( UPPER, A1, A2, A3, B1, B2, B3, CSU, SNU, CSV, 00002 $ SNV, CSQ, SNQ ) 00003 * 00004 * -- LAPACK auxiliary routine (version 3.2) -- 00005 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00006 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00007 * November 2006 00008 * 00009 * .. Scalar Arguments .. 00010 LOGICAL UPPER 00011 DOUBLE PRECISION A1, A2, A3, B1, B2, B3, CSQ, CSU, CSV, SNQ, 00012 $ SNU, SNV 00013 * .. 00014 * 00015 * Purpose 00016 * ======= 00017 * 00018 * DLAGS2 computes 2-by-2 orthogonal matrices U, V and Q, such 00019 * that if ( UPPER ) then 00020 * 00021 * U'*A*Q = U'*( A1 A2 )*Q = ( x 0 ) 00022 * ( 0 A3 ) ( x x ) 00023 * and 00024 * V'*B*Q = V'*( B1 B2 )*Q = ( x 0 ) 00025 * ( 0 B3 ) ( x x ) 00026 * 00027 * or if ( .NOT.UPPER ) then 00028 * 00029 * U'*A*Q = U'*( A1 0 )*Q = ( x x ) 00030 * ( A2 A3 ) ( 0 x ) 00031 * and 00032 * V'*B*Q = V'*( B1 0 )*Q = ( x x ) 00033 * ( B2 B3 ) ( 0 x ) 00034 * 00035 * The rows of the transformed A and B are parallel, where 00036 * 00037 * U = ( CSU SNU ), V = ( CSV SNV ), Q = ( CSQ SNQ ) 00038 * ( -SNU CSU ) ( -SNV CSV ) ( -SNQ CSQ ) 00039 * 00040 * Z' denotes the transpose of Z. 00041 * 00042 * 00043 * Arguments 00044 * ========= 00045 * 00046 * UPPER (input) LOGICAL 00047 * = .TRUE.: the input matrices A and B are upper triangular. 00048 * = .FALSE.: the input matrices A and B are lower triangular. 00049 * 00050 * A1 (input) DOUBLE PRECISION 00051 * A2 (input) DOUBLE PRECISION 00052 * A3 (input) DOUBLE PRECISION 00053 * On entry, A1, A2 and A3 are elements of the input 2-by-2 00054 * upper (lower) triangular matrix A. 00055 * 00056 * B1 (input) DOUBLE PRECISION 00057 * B2 (input) DOUBLE PRECISION 00058 * B3 (input) DOUBLE PRECISION 00059 * On entry, B1, B2 and B3 are elements of the input 2-by-2 00060 * upper (lower) triangular matrix B. 00061 * 00062 * CSU (output) DOUBLE PRECISION 00063 * SNU (output) DOUBLE PRECISION 00064 * The desired orthogonal matrix U. 00065 * 00066 * CSV (output) DOUBLE PRECISION 00067 * SNV (output) DOUBLE PRECISION 00068 * The desired orthogonal matrix V. 00069 * 00070 * CSQ (output) DOUBLE PRECISION 00071 * SNQ (output) DOUBLE PRECISION 00072 * The desired orthogonal matrix Q. 00073 * 00074 * ===================================================================== 00075 * 00076 * .. Parameters .. 00077 DOUBLE PRECISION ZERO 00078 PARAMETER ( ZERO = 0.0D+0 ) 00079 * .. 00080 * .. Local Scalars .. 00081 DOUBLE PRECISION A, AUA11, AUA12, AUA21, AUA22, AVB11, AVB12, 00082 $ AVB21, AVB22, B, C, CSL, CSR, D, R, S1, S2, 00083 $ SNL, SNR, UA11, UA11R, UA12, UA21, UA22, UA22R, 00084 $ VB11, VB11R, VB12, VB21, VB22, VB22R 00085 * .. 00086 * .. External Subroutines .. 00087 EXTERNAL DLARTG, DLASV2 00088 * .. 00089 * .. Intrinsic Functions .. 00090 INTRINSIC ABS 00091 * .. 00092 * .. Executable Statements .. 00093 * 00094 IF( UPPER ) THEN 00095 * 00096 * Input matrices A and B are upper triangular matrices 00097 * 00098 * Form matrix C = A*adj(B) = ( a b ) 00099 * ( 0 d ) 00100 * 00101 A = A1*B3 00102 D = A3*B1 00103 B = A2*B1 - A1*B2 00104 * 00105 * The SVD of real 2-by-2 triangular C 00106 * 00107 * ( CSL -SNL )*( A B )*( CSR SNR ) = ( R 0 ) 00108 * ( SNL CSL ) ( 0 D ) ( -SNR CSR ) ( 0 T ) 00109 * 00110 CALL DLASV2( A, B, D, S1, S2, SNR, CSR, SNL, CSL ) 00111 * 00112 IF( ABS( CSL ).GE.ABS( SNL ) .OR. ABS( CSR ).GE.ABS( SNR ) ) 00113 $ THEN 00114 * 00115 * Compute the (1,1) and (1,2) elements of U'*A and V'*B, 00116 * and (1,2) element of |U|'*|A| and |V|'*|B|. 00117 * 00118 UA11R = CSL*A1 00119 UA12 = CSL*A2 + SNL*A3 00120 * 00121 VB11R = CSR*B1 00122 VB12 = CSR*B2 + SNR*B3 00123 * 00124 AUA12 = ABS( CSL )*ABS( A2 ) + ABS( SNL )*ABS( A3 ) 00125 AVB12 = ABS( CSR )*ABS( B2 ) + ABS( SNR )*ABS( B3 ) 00126 * 00127 * zero (1,2) elements of U'*A and V'*B 00128 * 00129 IF( ( ABS( UA11R )+ABS( UA12 ) ).NE.ZERO ) THEN 00130 IF( AUA12 / ( ABS( UA11R )+ABS( UA12 ) ).LE.AVB12 / 00131 $ ( ABS( VB11R )+ABS( VB12 ) ) ) THEN 00132 CALL DLARTG( -UA11R, UA12, CSQ, SNQ, R ) 00133 ELSE 00134 CALL DLARTG( -VB11R, VB12, CSQ, SNQ, R ) 00135 END IF 00136 ELSE 00137 CALL DLARTG( -VB11R, VB12, CSQ, SNQ, R ) 00138 END IF 00139 * 00140 CSU = CSL 00141 SNU = -SNL 00142 CSV = CSR 00143 SNV = -SNR 00144 * 00145 ELSE 00146 * 00147 * Compute the (2,1) and (2,2) elements of U'*A and V'*B, 00148 * and (2,2) element of |U|'*|A| and |V|'*|B|. 00149 * 00150 UA21 = -SNL*A1 00151 UA22 = -SNL*A2 + CSL*A3 00152 * 00153 VB21 = -SNR*B1 00154 VB22 = -SNR*B2 + CSR*B3 00155 * 00156 AUA22 = ABS( SNL )*ABS( A2 ) + ABS( CSL )*ABS( A3 ) 00157 AVB22 = ABS( SNR )*ABS( B2 ) + ABS( CSR )*ABS( B3 ) 00158 * 00159 * zero (2,2) elements of U'*A and V'*B, and then swap. 00160 * 00161 IF( ( ABS( UA21 )+ABS( UA22 ) ).NE.ZERO ) THEN 00162 IF( AUA22 / ( ABS( UA21 )+ABS( UA22 ) ).LE.AVB22 / 00163 $ ( ABS( VB21 )+ABS( VB22 ) ) ) THEN 00164 CALL DLARTG( -UA21, UA22, CSQ, SNQ, R ) 00165 ELSE 00166 CALL DLARTG( -VB21, VB22, CSQ, SNQ, R ) 00167 END IF 00168 ELSE 00169 CALL DLARTG( -VB21, VB22, CSQ, SNQ, R ) 00170 END IF 00171 * 00172 CSU = SNL 00173 SNU = CSL 00174 CSV = SNR 00175 SNV = CSR 00176 * 00177 END IF 00178 * 00179 ELSE 00180 * 00181 * Input matrices A and B are lower triangular matrices 00182 * 00183 * Form matrix C = A*adj(B) = ( a 0 ) 00184 * ( c d ) 00185 * 00186 A = A1*B3 00187 D = A3*B1 00188 C = A2*B3 - A3*B2 00189 * 00190 * The SVD of real 2-by-2 triangular C 00191 * 00192 * ( CSL -SNL )*( A 0 )*( CSR SNR ) = ( R 0 ) 00193 * ( SNL CSL ) ( C D ) ( -SNR CSR ) ( 0 T ) 00194 * 00195 CALL DLASV2( A, C, D, S1, S2, SNR, CSR, SNL, CSL ) 00196 * 00197 IF( ABS( CSR ).GE.ABS( SNR ) .OR. ABS( CSL ).GE.ABS( SNL ) ) 00198 $ THEN 00199 * 00200 * Compute the (2,1) and (2,2) elements of U'*A and V'*B, 00201 * and (2,1) element of |U|'*|A| and |V|'*|B|. 00202 * 00203 UA21 = -SNR*A1 + CSR*A2 00204 UA22R = CSR*A3 00205 * 00206 VB21 = -SNL*B1 + CSL*B2 00207 VB22R = CSL*B3 00208 * 00209 AUA21 = ABS( SNR )*ABS( A1 ) + ABS( CSR )*ABS( A2 ) 00210 AVB21 = ABS( SNL )*ABS( B1 ) + ABS( CSL )*ABS( B2 ) 00211 * 00212 * zero (2,1) elements of U'*A and V'*B. 00213 * 00214 IF( ( ABS( UA21 )+ABS( UA22R ) ).NE.ZERO ) THEN 00215 IF( AUA21 / ( ABS( UA21 )+ABS( UA22R ) ).LE.AVB21 / 00216 $ ( ABS( VB21 )+ABS( VB22R ) ) ) THEN 00217 CALL DLARTG( UA22R, UA21, CSQ, SNQ, R ) 00218 ELSE 00219 CALL DLARTG( VB22R, VB21, CSQ, SNQ, R ) 00220 END IF 00221 ELSE 00222 CALL DLARTG( VB22R, VB21, CSQ, SNQ, R ) 00223 END IF 00224 * 00225 CSU = CSR 00226 SNU = -SNR 00227 CSV = CSL 00228 SNV = -SNL 00229 * 00230 ELSE 00231 * 00232 * Compute the (1,1) and (1,2) elements of U'*A and V'*B, 00233 * and (1,1) element of |U|'*|A| and |V|'*|B|. 00234 * 00235 UA11 = CSR*A1 + SNR*A2 00236 UA12 = SNR*A3 00237 * 00238 VB11 = CSL*B1 + SNL*B2 00239 VB12 = SNL*B3 00240 * 00241 AUA11 = ABS( CSR )*ABS( A1 ) + ABS( SNR )*ABS( A2 ) 00242 AVB11 = ABS( CSL )*ABS( B1 ) + ABS( SNL )*ABS( B2 ) 00243 * 00244 * zero (1,1) elements of U'*A and V'*B, and then swap. 00245 * 00246 IF( ( ABS( UA11 )+ABS( UA12 ) ).NE.ZERO ) THEN 00247 IF( AUA11 / ( ABS( UA11 )+ABS( UA12 ) ).LE.AVB11 / 00248 $ ( ABS( VB11 )+ABS( VB12 ) ) ) THEN 00249 CALL DLARTG( UA12, UA11, CSQ, SNQ, R ) 00250 ELSE 00251 CALL DLARTG( VB12, VB11, CSQ, SNQ, R ) 00252 END IF 00253 ELSE 00254 CALL DLARTG( VB12, VB11, CSQ, SNQ, R ) 00255 END IF 00256 * 00257 CSU = SNR 00258 SNU = CSR 00259 CSV = SNL 00260 SNV = CSL 00261 * 00262 END IF 00263 * 00264 END IF 00265 * 00266 RETURN 00267 * 00268 * End of DLAGS2 00269 * 00270 END
1.7.3 
