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