SRC\clanv2.f
#lines: 145   size: 3 Kb   creation: 18/01/2006 23:36:04   last modification: 08/05/2008 18:37:39   attribute: ARCH    Find   Reload   
1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
27:
28:
29:
30:
31:
32:
33:
34:
35:
36:
37:
38:
39:
40:
41:
42:
43:
44:
45:
46:
47:
48:
49:
50:
51:
52:
53:
54:
55:
56:
57:
58:
59:
60:
61:
62:
63:
64:
65:
66:
67:
68:
69:
70:
71:
72:
73:
74:
75:
76:
77:
78:
79:
80:
81:
82:
83:
84:
85:
86:
87:
88:
89:
90:
91:
92:
93:
94:
95:
96:
97:
98:
99:
100:
101:
102:
103:
104:
105:
106:
107:
108:
109:
110:
111:
112:
113:
114:
115:
116:
117:
118:
119:
120:
121:
122:
123:
124:
125:
126:
127:
128:
129:
130:
131:
132:
133:
134:
135:
136:
137:
138:
139:
140:
141:
142:
143:
144:
145:
 
      SUBROUTINE CLANV2( A, B, C, D, RT1, RT2, CS, SN )
*
*  -- ScaLAPACK routine (version 1.7) --
*     Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
*     Courant Institute, Argonne National Lab, and Rice University
*     May 28, 1999
*
*     .. Scalar Arguments ..
      REAL               CS
      COMPLEX            A, B, C, D, RT1, RT2, SN
*     ..
*
*  Purpose
*  =======
*
*  CLANV2 computes the Schur factorization of a complex 2-by-2
*  nonhermitian matrix in standard form:
*
*       [ A  B ] = [ CS -SN ] [ AA  BB ] [ CS  SN ]
*       [ C  D ]   [ SN  CS ] [  0  DD ] [-SN  CS ]
*
*  Arguments
*  =========
*
*  A       (input/output) COMPLEX
*  B       (input/output) COMPLEX
*  C       (input/output) COMPLEX
*  D       (input/output) COMPLEX
*          On entry, the elements of the input matrix.
*          On exit, they are overwritten by the elements of the
*          standardised Schur form.
*
*  RT1     (output) COMPLEX
*  RT2     (output) COMPLEX
*          The two eigenvalues.
*
*  CS      (output) REAL
*  SN      (output) COMPLEX
*          Parameters of the rotation matrix.
*
*  Further Details
*  ===============
*
*  Implemented by Mark R. Fahey, May 28, 1999
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               RZERO, HALF, RONE
      PARAMETER          ( RZERO = 0.0E+0, HALF = 0.5E+0,
     $                   RONE = 1.0E+0 )
      COMPLEX            ZERO, ONE
      PARAMETER          ( ZERO = ( 0.0E+0, 0.0E+0 ),
     $                   ONE = ( 1.0E+0, 0.0E+0 ) )
*     ..
*     .. Local Scalars ..
      COMPLEX            AA, BB, DD, T, TEMP, TEMP2, U, X, Y
*     ..
*     .. External Functions ..
      COMPLEX            CLADIV
      EXTERNAL           CLADIV
*     ..
*     .. External Subroutines ..
      EXTERNAL           CLARTG
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          REAL, CMPLX, CONJG, AIMAG, SQRT
*     ..
*     .. Executable Statements ..
*
*     Initialize CS and SN
*
      CS = RONE
      SN = ZERO
*
      IF( C.EQ.ZERO ) THEN
         GO TO 10
*
      ELSE IF( B.EQ.ZERO ) THEN
*
*        Swap rows and columns
*
         CS = RZERO
         SN = ONE
         TEMP = D
         D = A
         A = TEMP
         B = -C
         C = ZERO
         GO TO 10
      ELSE IF( ( A-D ).EQ.ZERO ) THEN
         TEMP = SQRT( B*C )
         A = A + TEMP
         D = D - TEMP
         IF( ( B+C ).EQ.ZERO ) THEN
            CS = SQRT( HALF )
            SN = CMPLX( RZERO, RONE )*CS
         ELSE
            TEMP = SQRT( B+C )
            TEMP2 = CLADIV( SQRT( B ), TEMP )
            CS = REAL( TEMP2 )
            SN = CLADIV( SQRT( C ), TEMP )
         END IF
         B = B - C
         C = ZERO
         GO TO 10
      ELSE
*
*        Compute eigenvalue closest to D
*
         T = D
         U = B*C
         X = HALF*( A-T )
         Y = SQRT( X*X+U )
         IF( REAL( X )*REAL( Y )+AIMAG( X )*AIMAG( Y ).LT.RZERO )
     $      Y = -Y
         T = T - CLADIV( U, ( X+Y ) )
*
*        Do one QR step with exact shift T - resulting 2 x 2 in
*        triangular form.
*
         CALL CLARTG( A-T, C, CS, SN, AA )
*
         D = D - T
         BB = CS*B + SN*D
         DD = -CONJG( SN )*B + CS*D
*
         A = AA*CS + BB*CONJG( SN ) + T
         B = -AA*SN + BB*CS
         C = ZERO
         D = T
*
      END IF
*
   10 CONTINUE
*
*     Store eigenvalues in RT1 and RT2.
*
      RT1 = A
      RT2 = D
      RETURN
*
*     End of CLANV2
*
      END