Annotation of rpl/lapack/lapack/zlags2.f, revision 1.18
1.9 bertrand 1: *> \brief \b ZLAGS2
2: *
3: * =========== DOCUMENTATION ===========
4: *
1.15 bertrand 5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
1.9 bertrand 7: *
8: *> \htmlonly
1.15 bertrand 9: *> Download ZLAGS2 + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlags2.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlags2.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlags2.f">
1.9 bertrand 15: *> [TXT]</a>
1.15 bertrand 16: *> \endhtmlonly
1.9 bertrand 17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE ZLAGS2( UPPER, A1, A2, A3, B1, B2, B3, CSU, SNU, CSV,
22: * SNV, CSQ, SNQ )
1.15 bertrand 23: *
1.9 bertrand 24: * .. Scalar Arguments ..
25: * LOGICAL UPPER
26: * DOUBLE PRECISION A1, A3, B1, B3, CSQ, CSU, CSV
27: * COMPLEX*16 A2, B2, SNQ, SNU, SNV
28: * ..
1.15 bertrand 29: *
1.9 bertrand 30: *
31: *> \par Purpose:
32: * =============
33: *>
34: *> \verbatim
35: *>
36: *> ZLAGS2 computes 2-by-2 unitary matrices U, V and Q, such
37: *> that if ( UPPER ) then
38: *>
39: *> U**H *A*Q = U**H *( A1 A2 )*Q = ( x 0 )
40: *> ( 0 A3 ) ( x x )
41: *> and
42: *> V**H*B*Q = V**H *( B1 B2 )*Q = ( x 0 )
43: *> ( 0 B3 ) ( x x )
44: *>
45: *> or if ( .NOT.UPPER ) then
46: *>
47: *> U**H *A*Q = U**H *( A1 0 )*Q = ( x x )
48: *> ( A2 A3 ) ( 0 x )
49: *> and
50: *> V**H *B*Q = V**H *( B1 0 )*Q = ( x x )
51: *> ( B2 B3 ) ( 0 x )
52: *> where
53: *>
54: *> U = ( CSU SNU ), V = ( CSV SNV ),
55: *> ( -SNU**H CSU ) ( -SNV**H CSV )
56: *>
57: *> Q = ( CSQ SNQ )
58: *> ( -SNQ**H CSQ )
59: *>
60: *> The rows of the transformed A and B are parallel. Moreover, if the
61: *> input 2-by-2 matrix A is not zero, then the transformed (1,1) entry
62: *> of A is not zero. If the input matrices A and B are both not zero,
63: *> then the transformed (2,2) element of B is not zero, except when the
64: *> first rows of input A and B are parallel and the second rows are
65: *> zero.
66: *> \endverbatim
67: *
68: * Arguments:
69: * ==========
70: *
71: *> \param[in] UPPER
72: *> \verbatim
73: *> UPPER is LOGICAL
74: *> = .TRUE.: the input matrices A and B are upper triangular.
75: *> = .FALSE.: the input matrices A and B are lower triangular.
76: *> \endverbatim
77: *>
78: *> \param[in] A1
79: *> \verbatim
80: *> A1 is DOUBLE PRECISION
81: *> \endverbatim
82: *>
83: *> \param[in] A2
84: *> \verbatim
85: *> A2 is COMPLEX*16
86: *> \endverbatim
87: *>
88: *> \param[in] A3
89: *> \verbatim
90: *> A3 is DOUBLE PRECISION
91: *> On entry, A1, A2 and A3 are elements of the input 2-by-2
92: *> upper (lower) triangular matrix A.
93: *> \endverbatim
94: *>
95: *> \param[in] B1
96: *> \verbatim
97: *> B1 is DOUBLE PRECISION
98: *> \endverbatim
99: *>
100: *> \param[in] B2
101: *> \verbatim
102: *> B2 is COMPLEX*16
103: *> \endverbatim
104: *>
105: *> \param[in] B3
106: *> \verbatim
107: *> B3 is DOUBLE PRECISION
108: *> On entry, B1, B2 and B3 are elements of the input 2-by-2
109: *> upper (lower) triangular matrix B.
110: *> \endverbatim
111: *>
112: *> \param[out] CSU
113: *> \verbatim
114: *> CSU is DOUBLE PRECISION
115: *> \endverbatim
116: *>
117: *> \param[out] SNU
118: *> \verbatim
119: *> SNU is COMPLEX*16
120: *> The desired unitary matrix U.
121: *> \endverbatim
122: *>
123: *> \param[out] CSV
124: *> \verbatim
125: *> CSV is DOUBLE PRECISION
126: *> \endverbatim
127: *>
128: *> \param[out] SNV
129: *> \verbatim
130: *> SNV is COMPLEX*16
131: *> The desired unitary matrix V.
132: *> \endverbatim
133: *>
134: *> \param[out] CSQ
135: *> \verbatim
136: *> CSQ is DOUBLE PRECISION
137: *> \endverbatim
138: *>
139: *> \param[out] SNQ
140: *> \verbatim
141: *> SNQ is COMPLEX*16
142: *> The desired unitary matrix Q.
143: *> \endverbatim
144: *
145: * Authors:
146: * ========
147: *
1.15 bertrand 148: *> \author Univ. of Tennessee
149: *> \author Univ. of California Berkeley
150: *> \author Univ. of Colorado Denver
151: *> \author NAG Ltd.
1.9 bertrand 152: *
153: *> \ingroup complex16OTHERauxiliary
154: *
155: * =====================================================================
1.1 bertrand 156: SUBROUTINE ZLAGS2( UPPER, A1, A2, A3, B1, B2, B3, CSU, SNU, CSV,
157: $ SNV, CSQ, SNQ )
158: *
1.18 ! bertrand 159: * -- LAPACK auxiliary routine --
1.1 bertrand 160: * -- LAPACK is a software package provided by Univ. of Tennessee, --
161: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
162: *
163: * .. Scalar Arguments ..
164: LOGICAL UPPER
165: DOUBLE PRECISION A1, A3, B1, B3, CSQ, CSU, CSV
166: COMPLEX*16 A2, B2, SNQ, SNU, SNV
167: * ..
168: *
169: * =====================================================================
170: *
171: * .. Parameters ..
172: DOUBLE PRECISION ZERO, ONE
173: PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
174: * ..
175: * .. Local Scalars ..
1.15 bertrand 176: DOUBLE PRECISION A, AUA11, AUA12, AUA21, AUA22, AVB12, AVB11,
177: $ AVB21, AVB22, CSL, CSR, D, FB, FC, S1, S2,
1.1 bertrand 178: $ SNL, SNR, UA11R, UA22R, VB11R, VB22R
179: COMPLEX*16 B, C, D1, R, T, UA11, UA12, UA21, UA22, VB11,
180: $ VB12, VB21, VB22
181: * ..
182: * .. External Subroutines ..
183: EXTERNAL DLASV2, ZLARTG
184: * ..
185: * .. Intrinsic Functions ..
186: INTRINSIC ABS, DBLE, DCMPLX, DCONJG, DIMAG
187: * ..
188: * .. Statement Functions ..
189: DOUBLE PRECISION ABS1
190: * ..
191: * .. Statement Function definitions ..
192: ABS1( T ) = ABS( DBLE( T ) ) + ABS( DIMAG( T ) )
193: * ..
194: * .. Executable Statements ..
195: *
196: IF( UPPER ) THEN
197: *
198: * Input matrices A and B are upper triangular matrices
199: *
200: * Form matrix C = A*adj(B) = ( a b )
201: * ( 0 d )
202: *
203: A = A1*B3
204: D = A3*B1
205: B = A2*B1 - A1*B2
206: FB = ABS( B )
207: *
208: * Transform complex 2-by-2 matrix C to real matrix by unitary
209: * diagonal matrix diag(1,D1).
210: *
211: D1 = ONE
212: IF( FB.NE.ZERO )
213: $ D1 = B / FB
214: *
215: * The SVD of real 2 by 2 triangular C
216: *
217: * ( CSL -SNL )*( A B )*( CSR SNR ) = ( R 0 )
218: * ( SNL CSL ) ( 0 D ) ( -SNR CSR ) ( 0 T )
219: *
220: CALL DLASV2( A, FB, D, S1, S2, SNR, CSR, SNL, CSL )
221: *
222: IF( ABS( CSL ).GE.ABS( SNL ) .OR. ABS( CSR ).GE.ABS( SNR ) )
223: $ THEN
224: *
1.8 bertrand 225: * Compute the (1,1) and (1,2) elements of U**H *A and V**H *B,
226: * and (1,2) element of |U|**H *|A| and |V|**H *|B|.
1.1 bertrand 227: *
228: UA11R = CSL*A1
229: UA12 = CSL*A2 + D1*SNL*A3
230: *
231: VB11R = CSR*B1
232: VB12 = CSR*B2 + D1*SNR*B3
233: *
234: AUA12 = ABS( CSL )*ABS1( A2 ) + ABS( SNL )*ABS( A3 )
235: AVB12 = ABS( CSR )*ABS1( B2 ) + ABS( SNR )*ABS( B3 )
236: *
1.8 bertrand 237: * zero (1,2) elements of U**H *A and V**H *B
1.1 bertrand 238: *
239: IF( ( ABS( UA11R )+ABS1( UA12 ) ).EQ.ZERO ) THEN
240: CALL ZLARTG( -DCMPLX( VB11R ), DCONJG( VB12 ), CSQ, SNQ,
241: $ R )
242: ELSE IF( ( ABS( VB11R )+ABS1( VB12 ) ).EQ.ZERO ) THEN
243: CALL ZLARTG( -DCMPLX( UA11R ), DCONJG( UA12 ), CSQ, SNQ,
244: $ R )
245: ELSE IF( AUA12 / ( ABS( UA11R )+ABS1( UA12 ) ).LE.AVB12 /
246: $ ( ABS( VB11R )+ABS1( VB12 ) ) ) THEN
247: CALL ZLARTG( -DCMPLX( UA11R ), DCONJG( UA12 ), CSQ, SNQ,
248: $ R )
249: ELSE
250: CALL ZLARTG( -DCMPLX( VB11R ), DCONJG( VB12 ), CSQ, SNQ,
251: $ R )
252: END IF
253: *
254: CSU = CSL
255: SNU = -D1*SNL
256: CSV = CSR
257: SNV = -D1*SNR
258: *
259: ELSE
260: *
1.8 bertrand 261: * Compute the (2,1) and (2,2) elements of U**H *A and V**H *B,
262: * and (2,2) element of |U|**H *|A| and |V|**H *|B|.
1.1 bertrand 263: *
264: UA21 = -DCONJG( D1 )*SNL*A1
265: UA22 = -DCONJG( D1 )*SNL*A2 + CSL*A3
266: *
267: VB21 = -DCONJG( D1 )*SNR*B1
268: VB22 = -DCONJG( D1 )*SNR*B2 + CSR*B3
269: *
270: AUA22 = ABS( SNL )*ABS1( A2 ) + ABS( CSL )*ABS( A3 )
271: AVB22 = ABS( SNR )*ABS1( B2 ) + ABS( CSR )*ABS( B3 )
272: *
1.8 bertrand 273: * zero (2,2) elements of U**H *A and V**H *B, and then swap.
1.1 bertrand 274: *
275: IF( ( ABS1( UA21 )+ABS1( UA22 ) ).EQ.ZERO ) THEN
276: CALL ZLARTG( -DCONJG( VB21 ), DCONJG( VB22 ), CSQ, SNQ,
277: $ R )
278: ELSE IF( ( ABS1( VB21 )+ABS( VB22 ) ).EQ.ZERO ) THEN
279: CALL ZLARTG( -DCONJG( UA21 ), DCONJG( UA22 ), CSQ, SNQ,
280: $ R )
281: ELSE IF( AUA22 / ( ABS1( UA21 )+ABS1( UA22 ) ).LE.AVB22 /
282: $ ( ABS1( VB21 )+ABS1( VB22 ) ) ) THEN
283: CALL ZLARTG( -DCONJG( UA21 ), DCONJG( UA22 ), CSQ, SNQ,
284: $ R )
285: ELSE
286: CALL ZLARTG( -DCONJG( VB21 ), DCONJG( VB22 ), CSQ, SNQ,
287: $ R )
288: END IF
289: *
290: CSU = SNL
291: SNU = D1*CSL
292: CSV = SNR
293: SNV = D1*CSR
294: *
295: END IF
296: *
297: ELSE
298: *
299: * Input matrices A and B are lower triangular matrices
300: *
301: * Form matrix C = A*adj(B) = ( a 0 )
302: * ( c d )
303: *
304: A = A1*B3
305: D = A3*B1
306: C = A2*B3 - A3*B2
307: FC = ABS( C )
308: *
309: * Transform complex 2-by-2 matrix C to real matrix by unitary
310: * diagonal matrix diag(d1,1).
311: *
312: D1 = ONE
313: IF( FC.NE.ZERO )
314: $ D1 = C / FC
315: *
316: * The SVD of real 2 by 2 triangular C
317: *
318: * ( CSL -SNL )*( A 0 )*( CSR SNR ) = ( R 0 )
319: * ( SNL CSL ) ( C D ) ( -SNR CSR ) ( 0 T )
320: *
321: CALL DLASV2( A, FC, D, S1, S2, SNR, CSR, SNL, CSL )
322: *
323: IF( ABS( CSR ).GE.ABS( SNR ) .OR. ABS( CSL ).GE.ABS( SNL ) )
324: $ THEN
325: *
1.8 bertrand 326: * Compute the (2,1) and (2,2) elements of U**H *A and V**H *B,
327: * and (2,1) element of |U|**H *|A| and |V|**H *|B|.
1.1 bertrand 328: *
329: UA21 = -D1*SNR*A1 + CSR*A2
330: UA22R = CSR*A3
331: *
332: VB21 = -D1*SNL*B1 + CSL*B2
333: VB22R = CSL*B3
334: *
335: AUA21 = ABS( SNR )*ABS( A1 ) + ABS( CSR )*ABS1( A2 )
336: AVB21 = ABS( SNL )*ABS( B1 ) + ABS( CSL )*ABS1( B2 )
337: *
1.8 bertrand 338: * zero (2,1) elements of U**H *A and V**H *B.
1.1 bertrand 339: *
340: IF( ( ABS1( UA21 )+ABS( UA22R ) ).EQ.ZERO ) THEN
341: CALL ZLARTG( DCMPLX( VB22R ), VB21, CSQ, SNQ, R )
342: ELSE IF( ( ABS1( VB21 )+ABS( VB22R ) ).EQ.ZERO ) THEN
343: CALL ZLARTG( DCMPLX( UA22R ), UA21, CSQ, SNQ, R )
344: ELSE IF( AUA21 / ( ABS1( UA21 )+ABS( UA22R ) ).LE.AVB21 /
345: $ ( ABS1( VB21 )+ABS( VB22R ) ) ) THEN
346: CALL ZLARTG( DCMPLX( UA22R ), UA21, CSQ, SNQ, R )
347: ELSE
348: CALL ZLARTG( DCMPLX( VB22R ), VB21, CSQ, SNQ, R )
349: END IF
350: *
351: CSU = CSR
352: SNU = -DCONJG( D1 )*SNR
353: CSV = CSL
354: SNV = -DCONJG( D1 )*SNL
355: *
356: ELSE
357: *
1.8 bertrand 358: * Compute the (1,1) and (1,2) elements of U**H *A and V**H *B,
359: * and (1,1) element of |U|**H *|A| and |V|**H *|B|.
1.1 bertrand 360: *
361: UA11 = CSR*A1 + DCONJG( D1 )*SNR*A2
362: UA12 = DCONJG( D1 )*SNR*A3
363: *
364: VB11 = CSL*B1 + DCONJG( D1 )*SNL*B2
365: VB12 = DCONJG( D1 )*SNL*B3
366: *
367: AUA11 = ABS( CSR )*ABS( A1 ) + ABS( SNR )*ABS1( A2 )
368: AVB11 = ABS( CSL )*ABS( B1 ) + ABS( SNL )*ABS1( B2 )
369: *
1.8 bertrand 370: * zero (1,1) elements of U**H *A and V**H *B, and then swap.
1.1 bertrand 371: *
372: IF( ( ABS1( UA11 )+ABS1( UA12 ) ).EQ.ZERO ) THEN
373: CALL ZLARTG( VB12, VB11, CSQ, SNQ, R )
374: ELSE IF( ( ABS1( VB11 )+ABS1( VB12 ) ).EQ.ZERO ) THEN
375: CALL ZLARTG( UA12, UA11, CSQ, SNQ, R )
376: ELSE IF( AUA11 / ( ABS1( UA11 )+ABS1( UA12 ) ).LE.AVB11 /
377: $ ( ABS1( VB11 )+ABS1( VB12 ) ) ) THEN
378: CALL ZLARTG( UA12, UA11, CSQ, SNQ, R )
379: ELSE
380: CALL ZLARTG( VB12, VB11, CSQ, SNQ, R )
381: END IF
382: *
383: CSU = SNR
384: SNU = DCONJG( D1 )*CSR
385: CSV = SNL
386: SNV = DCONJG( D1 )*CSL
387: *
388: END IF
389: *
390: END IF
391: *
392: RETURN
393: *
394: * End of ZLAGS2
395: *
396: END
CVSweb interface <joel.bertrand@systella.fr>