Annotation of rpl/lapack/lapack/zlarfg.f, revision 1.1.1.1
1.1 bertrand 1: SUBROUTINE ZLARFG( N, ALPHA, X, INCX, TAU )
2: *
3: * -- LAPACK auxiliary routine (version 3.2) --
4: * -- LAPACK is a software package provided by Univ. of Tennessee, --
5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
6: * November 2006
7: *
8: * .. Scalar Arguments ..
9: INTEGER INCX, N
10: COMPLEX*16 ALPHA, TAU
11: * ..
12: * .. Array Arguments ..
13: COMPLEX*16 X( * )
14: * ..
15: *
16: * Purpose
17: * =======
18: *
19: * ZLARFG generates a complex elementary reflector H of order n, such
20: * that
21: *
22: * H' * ( alpha ) = ( beta ), H' * H = I.
23: * ( x ) ( 0 )
24: *
25: * where alpha and beta are scalars, with beta real, and x is an
26: * (n-1)-element complex vector. H is represented in the form
27: *
28: * H = I - tau * ( 1 ) * ( 1 v' ) ,
29: * ( v )
30: *
31: * where tau is a complex scalar and v is a complex (n-1)-element
32: * vector. Note that H is not hermitian.
33: *
34: * If the elements of x are all zero and alpha is real, then tau = 0
35: * and H is taken to be the unit matrix.
36: *
37: * Otherwise 1 <= real(tau) <= 2 and abs(tau-1) <= 1 .
38: *
39: * Arguments
40: * =========
41: *
42: * N (input) INTEGER
43: * The order of the elementary reflector.
44: *
45: * ALPHA (input/output) COMPLEX*16
46: * On entry, the value alpha.
47: * On exit, it is overwritten with the value beta.
48: *
49: * X (input/output) COMPLEX*16 array, dimension
50: * (1+(N-2)*abs(INCX))
51: * On entry, the vector x.
52: * On exit, it is overwritten with the vector v.
53: *
54: * INCX (input) INTEGER
55: * The increment between elements of X. INCX > 0.
56: *
57: * TAU (output) COMPLEX*16
58: * The value tau.
59: *
60: * =====================================================================
61: *
62: * .. Parameters ..
63: DOUBLE PRECISION ONE, ZERO
64: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
65: * ..
66: * .. Local Scalars ..
67: INTEGER J, KNT
68: DOUBLE PRECISION ALPHI, ALPHR, BETA, RSAFMN, SAFMIN, XNORM
69: * ..
70: * .. External Functions ..
71: DOUBLE PRECISION DLAMCH, DLAPY3, DZNRM2
72: COMPLEX*16 ZLADIV
73: EXTERNAL DLAMCH, DLAPY3, DZNRM2, ZLADIV
74: * ..
75: * .. Intrinsic Functions ..
76: INTRINSIC ABS, DBLE, DCMPLX, DIMAG, SIGN
77: * ..
78: * .. External Subroutines ..
79: EXTERNAL ZDSCAL, ZSCAL
80: * ..
81: * .. Executable Statements ..
82: *
83: IF( N.LE.0 ) THEN
84: TAU = ZERO
85: RETURN
86: END IF
87: *
88: XNORM = DZNRM2( N-1, X, INCX )
89: ALPHR = DBLE( ALPHA )
90: ALPHI = DIMAG( ALPHA )
91: *
92: IF( XNORM.EQ.ZERO .AND. ALPHI.EQ.ZERO ) THEN
93: *
94: * H = I
95: *
96: TAU = ZERO
97: ELSE
98: *
99: * general case
100: *
101: BETA = -SIGN( DLAPY3( ALPHR, ALPHI, XNORM ), ALPHR )
102: SAFMIN = DLAMCH( 'S' ) / DLAMCH( 'E' )
103: RSAFMN = ONE / SAFMIN
104: *
105: KNT = 0
106: IF( ABS( BETA ).LT.SAFMIN ) THEN
107: *
108: * XNORM, BETA may be inaccurate; scale X and recompute them
109: *
110: 10 CONTINUE
111: KNT = KNT + 1
112: CALL ZDSCAL( N-1, RSAFMN, X, INCX )
113: BETA = BETA*RSAFMN
114: ALPHI = ALPHI*RSAFMN
115: ALPHR = ALPHR*RSAFMN
116: IF( ABS( BETA ).LT.SAFMIN )
117: $ GO TO 10
118: *
119: * New BETA is at most 1, at least SAFMIN
120: *
121: XNORM = DZNRM2( N-1, X, INCX )
122: ALPHA = DCMPLX( ALPHR, ALPHI )
123: BETA = -SIGN( DLAPY3( ALPHR, ALPHI, XNORM ), ALPHR )
124: END IF
125: TAU = DCMPLX( ( BETA-ALPHR ) / BETA, -ALPHI / BETA )
126: ALPHA = ZLADIV( DCMPLX( ONE ), ALPHA-BETA )
127: CALL ZSCAL( N-1, ALPHA, X, INCX )
128: *
129: * If ALPHA is subnormal, it may lose relative accuracy
130: *
131: DO 20 J = 1, KNT
132: BETA = BETA*SAFMIN
133: 20 CONTINUE
134: ALPHA = BETA
135: END IF
136: *
137: RETURN
138: *
139: * End of ZLARFG
140: *
141: END
CVSweb interface <joel.bertrand@systella.fr>