SRC\zlanv2.f

#lines: 145   size: 3 Kb   creation: 18/01/2006 23:36:04   last modification: 08/05/2008 18:38:23   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 ZLANV2( 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 ..
      DOUBLE PRECISION   CS
      COMPLEX*16         A, B, C, D, RT1, RT2, SN
*     ..
*
*  Purpose
*  =======
*
*  ZLANV2 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*16
*  B       (input/output) COMPLEX*16
*  C       (input/output) COMPLEX*16
*  D       (input/output) COMPLEX*16
*          On entry, the elements of the input matrix.
*          On exit, they are overwritten by the elements of the
*          standardised Schur form.
*
*  RT1     (output) COMPLEX*16
*  RT2     (output) COMPLEX*16
*          The two eigenvalues.
*
*  CS      (output) DOUBLE PRECISION
*  SN      (output) COMPLEX*16
*          Parameters of the rotation matrix.
*
*  Further Details
*  ===============
*
*  Implemented by Mark R. Fahey, May 28, 1999
*
*  =====================================================================
*
*     .. Parameters ..
      DOUBLE PRECISION   RZERO, HALF, RONE
      PARAMETER          ( RZERO = 0.0D+0, HALF = 0.5D+0,
     $                   RONE = 1.0D+0 )
      COMPLEX*16         ZERO, ONE
      PARAMETER          ( ZERO = ( 0.0D+0, 0.0D+0 ),
     $                   ONE = ( 1.0D+0, 0.0D+0 ) )
*     ..
*     .. Local Scalars ..
      COMPLEX*16         AA, BB, DD, T, TEMP, TEMP2, U, X, Y
*     ..
*     .. External Functions ..
      COMPLEX*16         ZLADIV
      EXTERNAL           ZLADIV
*     ..
*     .. External Subroutines ..
      EXTERNAL           ZLARTG
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          DBLE, DCMPLX, DCONJG, DIMAG, 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 = DCMPLX( RZERO, RONE )*CS
         ELSE
            TEMP = SQRT( B+C )
            TEMP2 = ZLADIV( SQRT( B ), TEMP )
            CS = DBLE( TEMP2 )
            SN = ZLADIV( 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( DBLE( X )*DBLE( Y )+DIMAG( X )*DIMAG( Y ).LT.RZERO )
     $      Y = -Y
         T = T - ZLADIV( U, ( X+Y ) )
*
*        Do one QR step with exact shift T - resulting 2 x 2 in
*        triangular form.
*
         CALL ZLARTG( A-T, C, CS, SN, AA )
*
         D = D - T
         BB = CS*B + SN*D
         DD = -DCONJG( SN )*B + CS*D
*
         A = AA*CS + BB*DCONJG( 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 ZLANV2
*
      END