![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b ZLAGS2 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download ZLAGS2 + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlags2.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlags2.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlags2.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE ZLAGS2( UPPER, A1, A2, A3, B1, B2, B3, CSU, SNU, CSV, 00022 * SNV, CSQ, SNQ ) 00023 * 00024 * .. Scalar Arguments .. 00025 * LOGICAL UPPER 00026 * DOUBLE PRECISION A1, A3, B1, B3, CSQ, CSU, CSV 00027 * COMPLEX*16 A2, B2, SNQ, SNU, SNV 00028 * .. 00029 * 00030 * 00031 *> \par Purpose: 00032 * ============= 00033 *> 00034 *> \verbatim 00035 *> 00036 *> ZLAGS2 computes 2-by-2 unitary matrices U, V and Q, such 00037 *> that if ( UPPER ) then 00038 *> 00039 *> U**H *A*Q = U**H *( A1 A2 )*Q = ( x 0 ) 00040 *> ( 0 A3 ) ( x x ) 00041 *> and 00042 *> V**H*B*Q = V**H *( B1 B2 )*Q = ( x 0 ) 00043 *> ( 0 B3 ) ( x x ) 00044 *> 00045 *> or if ( .NOT.UPPER ) then 00046 *> 00047 *> U**H *A*Q = U**H *( A1 0 )*Q = ( x x ) 00048 *> ( A2 A3 ) ( 0 x ) 00049 *> and 00050 *> V**H *B*Q = V**H *( B1 0 )*Q = ( x x ) 00051 *> ( B2 B3 ) ( 0 x ) 00052 *> where 00053 *> 00054 *> U = ( CSU SNU ), V = ( CSV SNV ), 00055 *> ( -SNU**H CSU ) ( -SNV**H CSV ) 00056 *> 00057 *> Q = ( CSQ SNQ ) 00058 *> ( -SNQ**H CSQ ) 00059 *> 00060 *> The rows of the transformed A and B are parallel. Moreover, if the 00061 *> input 2-by-2 matrix A is not zero, then the transformed (1,1) entry 00062 *> of A is not zero. If the input matrices A and B are both not zero, 00063 *> then the transformed (2,2) element of B is not zero, except when the 00064 *> first rows of input A and B are parallel and the second rows are 00065 *> zero. 00066 *> \endverbatim 00067 * 00068 * Arguments: 00069 * ========== 00070 * 00071 *> \param[in] UPPER 00072 *> \verbatim 00073 *> UPPER is LOGICAL 00074 *> = .TRUE.: the input matrices A and B are upper triangular. 00075 *> = .FALSE.: the input matrices A and B are lower triangular. 00076 *> \endverbatim 00077 *> 00078 *> \param[in] A1 00079 *> \verbatim 00080 *> A1 is DOUBLE PRECISION 00081 *> \endverbatim 00082 *> 00083 *> \param[in] A2 00084 *> \verbatim 00085 *> A2 is COMPLEX*16 00086 *> \endverbatim 00087 *> 00088 *> \param[in] A3 00089 *> \verbatim 00090 *> A3 is DOUBLE PRECISION 00091 *> On entry, A1, A2 and A3 are elements of the input 2-by-2 00092 *> upper (lower) triangular matrix A. 00093 *> \endverbatim 00094 *> 00095 *> \param[in] B1 00096 *> \verbatim 00097 *> B1 is DOUBLE PRECISION 00098 *> \endverbatim 00099 *> 00100 *> \param[in] B2 00101 *> \verbatim 00102 *> B2 is COMPLEX*16 00103 *> \endverbatim 00104 *> 00105 *> \param[in] B3 00106 *> \verbatim 00107 *> B3 is DOUBLE PRECISION 00108 *> On entry, B1, B2 and B3 are elements of the input 2-by-2 00109 *> upper (lower) triangular matrix B. 00110 *> \endverbatim 00111 *> 00112 *> \param[out] CSU 00113 *> \verbatim 00114 *> CSU is DOUBLE PRECISION 00115 *> \endverbatim 00116 *> 00117 *> \param[out] SNU 00118 *> \verbatim 00119 *> SNU is COMPLEX*16 00120 *> The desired unitary matrix U. 00121 *> \endverbatim 00122 *> 00123 *> \param[out] CSV 00124 *> \verbatim 00125 *> CSV is DOUBLE PRECISION 00126 *> \endverbatim 00127 *> 00128 *> \param[out] SNV 00129 *> \verbatim 00130 *> SNV is COMPLEX*16 00131 *> The desired unitary matrix V. 00132 *> \endverbatim 00133 *> 00134 *> \param[out] CSQ 00135 *> \verbatim 00136 *> CSQ is DOUBLE PRECISION 00137 *> \endverbatim 00138 *> 00139 *> \param[out] SNQ 00140 *> \verbatim 00141 *> SNQ is COMPLEX*16 00142 *> The desired unitary matrix Q. 00143 *> \endverbatim 00144 * 00145 * Authors: 00146 * ======== 00147 * 00148 *> \author Univ. of Tennessee 00149 *> \author Univ. of California Berkeley 00150 *> \author Univ. of Colorado Denver 00151 *> \author NAG Ltd. 00152 * 00153 *> \date November 2011 00154 * 00155 *> \ingroup complex16OTHERauxiliary 00156 * 00157 * ===================================================================== 00158 SUBROUTINE ZLAGS2( UPPER, A1, A2, A3, B1, B2, B3, CSU, SNU, CSV, 00159 $ SNV, CSQ, SNQ ) 00160 * 00161 * -- LAPACK auxiliary routine (version 3.4.0) -- 00162 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00163 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00164 * November 2011 00165 * 00166 * .. Scalar Arguments .. 00167 LOGICAL UPPER 00168 DOUBLE PRECISION A1, A3, B1, B3, CSQ, CSU, CSV 00169 COMPLEX*16 A2, B2, SNQ, SNU, SNV 00170 * .. 00171 * 00172 * ===================================================================== 00173 * 00174 * .. Parameters .. 00175 DOUBLE PRECISION ZERO, ONE 00176 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 00177 * .. 00178 * .. Local Scalars .. 00179 DOUBLE PRECISION A, AUA11, AUA12, AUA21, AUA22, AVB12, AVB11, 00180 $ AVB21, AVB22, CSL, CSR, D, FB, FC, S1, S2, 00181 $ SNL, SNR, UA11R, UA22R, VB11R, VB22R 00182 COMPLEX*16 B, C, D1, R, T, UA11, UA12, UA21, UA22, VB11, 00183 $ VB12, VB21, VB22 00184 * .. 00185 * .. External Subroutines .. 00186 EXTERNAL DLASV2, ZLARTG 00187 * .. 00188 * .. Intrinsic Functions .. 00189 INTRINSIC ABS, DBLE, DCMPLX, DCONJG, DIMAG 00190 * .. 00191 * .. Statement Functions .. 00192 DOUBLE PRECISION ABS1 00193 * .. 00194 * .. Statement Function definitions .. 00195 ABS1( T ) = ABS( DBLE( T ) ) + ABS( DIMAG( T ) ) 00196 * .. 00197 * .. Executable Statements .. 00198 * 00199 IF( UPPER ) THEN 00200 * 00201 * Input matrices A and B are upper triangular matrices 00202 * 00203 * Form matrix C = A*adj(B) = ( a b ) 00204 * ( 0 d ) 00205 * 00206 A = A1*B3 00207 D = A3*B1 00208 B = A2*B1 - A1*B2 00209 FB = ABS( B ) 00210 * 00211 * Transform complex 2-by-2 matrix C to real matrix by unitary 00212 * diagonal matrix diag(1,D1). 00213 * 00214 D1 = ONE 00215 IF( FB.NE.ZERO ) 00216 $ D1 = B / FB 00217 * 00218 * The SVD of real 2 by 2 triangular C 00219 * 00220 * ( CSL -SNL )*( A B )*( CSR SNR ) = ( R 0 ) 00221 * ( SNL CSL ) ( 0 D ) ( -SNR CSR ) ( 0 T ) 00222 * 00223 CALL DLASV2( A, FB, D, S1, S2, SNR, CSR, SNL, CSL ) 00224 * 00225 IF( ABS( CSL ).GE.ABS( SNL ) .OR. ABS( CSR ).GE.ABS( SNR ) ) 00226 $ THEN 00227 * 00228 * Compute the (1,1) and (1,2) elements of U**H *A and V**H *B, 00229 * and (1,2) element of |U|**H *|A| and |V|**H *|B|. 00230 * 00231 UA11R = CSL*A1 00232 UA12 = CSL*A2 + D1*SNL*A3 00233 * 00234 VB11R = CSR*B1 00235 VB12 = CSR*B2 + D1*SNR*B3 00236 * 00237 AUA12 = ABS( CSL )*ABS1( A2 ) + ABS( SNL )*ABS( A3 ) 00238 AVB12 = ABS( CSR )*ABS1( B2 ) + ABS( SNR )*ABS( B3 ) 00239 * 00240 * zero (1,2) elements of U**H *A and V**H *B 00241 * 00242 IF( ( ABS( UA11R )+ABS1( UA12 ) ).EQ.ZERO ) THEN 00243 CALL ZLARTG( -DCMPLX( VB11R ), DCONJG( VB12 ), CSQ, SNQ, 00244 $ R ) 00245 ELSE IF( ( ABS( VB11R )+ABS1( VB12 ) ).EQ.ZERO ) THEN 00246 CALL ZLARTG( -DCMPLX( UA11R ), DCONJG( UA12 ), CSQ, SNQ, 00247 $ R ) 00248 ELSE IF( AUA12 / ( ABS( UA11R )+ABS1( UA12 ) ).LE.AVB12 / 00249 $ ( ABS( VB11R )+ABS1( VB12 ) ) ) THEN 00250 CALL ZLARTG( -DCMPLX( UA11R ), DCONJG( UA12 ), CSQ, SNQ, 00251 $ R ) 00252 ELSE 00253 CALL ZLARTG( -DCMPLX( VB11R ), DCONJG( VB12 ), CSQ, SNQ, 00254 $ R ) 00255 END IF 00256 * 00257 CSU = CSL 00258 SNU = -D1*SNL 00259 CSV = CSR 00260 SNV = -D1*SNR 00261 * 00262 ELSE 00263 * 00264 * Compute the (2,1) and (2,2) elements of U**H *A and V**H *B, 00265 * and (2,2) element of |U|**H *|A| and |V|**H *|B|. 00266 * 00267 UA21 = -DCONJG( D1 )*SNL*A1 00268 UA22 = -DCONJG( D1 )*SNL*A2 + CSL*A3 00269 * 00270 VB21 = -DCONJG( D1 )*SNR*B1 00271 VB22 = -DCONJG( D1 )*SNR*B2 + CSR*B3 00272 * 00273 AUA22 = ABS( SNL )*ABS1( A2 ) + ABS( CSL )*ABS( A3 ) 00274 AVB22 = ABS( SNR )*ABS1( B2 ) + ABS( CSR )*ABS( B3 ) 00275 * 00276 * zero (2,2) elements of U**H *A and V**H *B, and then swap. 00277 * 00278 IF( ( ABS1( UA21 )+ABS1( UA22 ) ).EQ.ZERO ) THEN 00279 CALL ZLARTG( -DCONJG( VB21 ), DCONJG( VB22 ), CSQ, SNQ, 00280 $ R ) 00281 ELSE IF( ( ABS1( VB21 )+ABS( VB22 ) ).EQ.ZERO ) THEN 00282 CALL ZLARTG( -DCONJG( UA21 ), DCONJG( UA22 ), CSQ, SNQ, 00283 $ R ) 00284 ELSE IF( AUA22 / ( ABS1( UA21 )+ABS1( UA22 ) ).LE.AVB22 / 00285 $ ( ABS1( VB21 )+ABS1( VB22 ) ) ) THEN 00286 CALL ZLARTG( -DCONJG( UA21 ), DCONJG( UA22 ), CSQ, SNQ, 00287 $ R ) 00288 ELSE 00289 CALL ZLARTG( -DCONJG( VB21 ), DCONJG( VB22 ), CSQ, SNQ, 00290 $ R ) 00291 END IF 00292 * 00293 CSU = SNL 00294 SNU = D1*CSL 00295 CSV = SNR 00296 SNV = D1*CSR 00297 * 00298 END IF 00299 * 00300 ELSE 00301 * 00302 * Input matrices A and B are lower triangular matrices 00303 * 00304 * Form matrix C = A*adj(B) = ( a 0 ) 00305 * ( c d ) 00306 * 00307 A = A1*B3 00308 D = A3*B1 00309 C = A2*B3 - A3*B2 00310 FC = ABS( C ) 00311 * 00312 * Transform complex 2-by-2 matrix C to real matrix by unitary 00313 * diagonal matrix diag(d1,1). 00314 * 00315 D1 = ONE 00316 IF( FC.NE.ZERO ) 00317 $ D1 = C / FC 00318 * 00319 * The SVD of real 2 by 2 triangular C 00320 * 00321 * ( CSL -SNL )*( A 0 )*( CSR SNR ) = ( R 0 ) 00322 * ( SNL CSL ) ( C D ) ( -SNR CSR ) ( 0 T ) 00323 * 00324 CALL DLASV2( A, FC, D, S1, S2, SNR, CSR, SNL, CSL ) 00325 * 00326 IF( ABS( CSR ).GE.ABS( SNR ) .OR. ABS( CSL ).GE.ABS( SNL ) ) 00327 $ THEN 00328 * 00329 * Compute the (2,1) and (2,2) elements of U**H *A and V**H *B, 00330 * and (2,1) element of |U|**H *|A| and |V|**H *|B|. 00331 * 00332 UA21 = -D1*SNR*A1 + CSR*A2 00333 UA22R = CSR*A3 00334 * 00335 VB21 = -D1*SNL*B1 + CSL*B2 00336 VB22R = CSL*B3 00337 * 00338 AUA21 = ABS( SNR )*ABS( A1 ) + ABS( CSR )*ABS1( A2 ) 00339 AVB21 = ABS( SNL )*ABS( B1 ) + ABS( CSL )*ABS1( B2 ) 00340 * 00341 * zero (2,1) elements of U**H *A and V**H *B. 00342 * 00343 IF( ( ABS1( UA21 )+ABS( UA22R ) ).EQ.ZERO ) THEN 00344 CALL ZLARTG( DCMPLX( VB22R ), VB21, CSQ, SNQ, R ) 00345 ELSE IF( ( ABS1( VB21 )+ABS( VB22R ) ).EQ.ZERO ) THEN 00346 CALL ZLARTG( DCMPLX( UA22R ), UA21, CSQ, SNQ, R ) 00347 ELSE IF( AUA21 / ( ABS1( UA21 )+ABS( UA22R ) ).LE.AVB21 / 00348 $ ( ABS1( VB21 )+ABS( VB22R ) ) ) THEN 00349 CALL ZLARTG( DCMPLX( UA22R ), UA21, CSQ, SNQ, R ) 00350 ELSE 00351 CALL ZLARTG( DCMPLX( VB22R ), VB21, CSQ, SNQ, R ) 00352 END IF 00353 * 00354 CSU = CSR 00355 SNU = -DCONJG( D1 )*SNR 00356 CSV = CSL 00357 SNV = -DCONJG( D1 )*SNL 00358 * 00359 ELSE 00360 * 00361 * Compute the (1,1) and (1,2) elements of U**H *A and V**H *B, 00362 * and (1,1) element of |U|**H *|A| and |V|**H *|B|. 00363 * 00364 UA11 = CSR*A1 + DCONJG( D1 )*SNR*A2 00365 UA12 = DCONJG( D1 )*SNR*A3 00366 * 00367 VB11 = CSL*B1 + DCONJG( D1 )*SNL*B2 00368 VB12 = DCONJG( D1 )*SNL*B3 00369 * 00370 AUA11 = ABS( CSR )*ABS( A1 ) + ABS( SNR )*ABS1( A2 ) 00371 AVB11 = ABS( CSL )*ABS( B1 ) + ABS( SNL )*ABS1( B2 ) 00372 * 00373 * zero (1,1) elements of U**H *A and V**H *B, and then swap. 00374 * 00375 IF( ( ABS1( UA11 )+ABS1( UA12 ) ).EQ.ZERO ) THEN 00376 CALL ZLARTG( VB12, VB11, CSQ, SNQ, R ) 00377 ELSE IF( ( ABS1( VB11 )+ABS1( VB12 ) ).EQ.ZERO ) THEN 00378 CALL ZLARTG( UA12, UA11, CSQ, SNQ, R ) 00379 ELSE IF( AUA11 / ( ABS1( UA11 )+ABS1( UA12 ) ).LE.AVB11 / 00380 $ ( ABS1( VB11 )+ABS1( VB12 ) ) ) THEN 00381 CALL ZLARTG( UA12, UA11, CSQ, SNQ, R ) 00382 ELSE 00383 CALL ZLARTG( VB12, VB11, CSQ, SNQ, R ) 00384 END IF 00385 * 00386 CSU = SNR 00387 SNU = DCONJG( D1 )*CSR 00388 CSV = SNL 00389 SNV = DCONJG( D1 )*CSL 00390 * 00391 END IF 00392 * 00393 END IF 00394 * 00395 RETURN 00396 * 00397 * End of ZLAGS2 00398 * 00399 END