LAPACK  3.6.1 LAPACK: Linear Algebra PACKage
dlarfg.f
Go to the documentation of this file.
1 *> \brief \b DLARFG generates an elementary reflector (Householder matrix).
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlarfg.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlarfg.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlarfg.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE DLARFG( N, ALPHA, X, INCX, TAU )
22 *
23 * .. Scalar Arguments ..
24 * INTEGER INCX, N
25 * DOUBLE PRECISION ALPHA, TAU
26 * ..
27 * .. Array Arguments ..
28 * DOUBLE PRECISION X( * )
29 * ..
30 *
31 *
32 *> \par Purpose:
33 * =============
34 *>
35 *> \verbatim
36 *>
37 *> DLARFG generates a real elementary reflector H of order n, such
38 *> that
39 *>
40 *> H * ( alpha ) = ( beta ), H**T * H = I.
41 *> ( x ) ( 0 )
42 *>
43 *> where alpha and beta are scalars, and x is an (n-1)-element real
44 *> vector. H is represented in the form
45 *>
46 *> H = I - tau * ( 1 ) * ( 1 v**T ) ,
47 *> ( v )
48 *>
49 *> where tau is a real scalar and v is a real (n-1)-element
50 *> vector.
51 *>
52 *> If the elements of x are all zero, then tau = 0 and H is taken to be
53 *> the unit matrix.
54 *>
55 *> Otherwise 1 <= tau <= 2.
56 *> \endverbatim
57 *
58 * Arguments:
59 * ==========
60 *
61 *> \param[in] N
62 *> \verbatim
63 *> N is INTEGER
64 *> The order of the elementary reflector.
65 *> \endverbatim
66 *>
67 *> \param[in,out] ALPHA
68 *> \verbatim
69 *> ALPHA is DOUBLE PRECISION
70 *> On entry, the value alpha.
71 *> On exit, it is overwritten with the value beta.
72 *> \endverbatim
73 *>
74 *> \param[in,out] X
75 *> \verbatim
76 *> X is DOUBLE PRECISION array, dimension
77 *> (1+(N-2)*abs(INCX))
78 *> On entry, the vector x.
79 *> On exit, it is overwritten with the vector v.
80 *> \endverbatim
81 *>
82 *> \param[in] INCX
83 *> \verbatim
84 *> INCX is INTEGER
85 *> The increment between elements of X. INCX > 0.
86 *> \endverbatim
87 *>
88 *> \param[out] TAU
89 *> \verbatim
90 *> TAU is DOUBLE PRECISION
91 *> The value tau.
92 *> \endverbatim
93 *
94 * Authors:
95 * ========
96 *
97 *> \author Univ. of Tennessee
98 *> \author Univ. of California Berkeley
99 *> \author Univ. of Colorado Denver
100 *> \author NAG Ltd.
101 *
102 *> \date September 2012
103 *
104 *> \ingroup doubleOTHERauxiliary
105 *
106 * =====================================================================
107  SUBROUTINE dlarfg( N, ALPHA, X, INCX, TAU )
108 *
109 * -- LAPACK auxiliary routine (version 3.4.2) --
110 * -- LAPACK is a software package provided by Univ. of Tennessee, --
111 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
112 * September 2012
113 *
114 * .. Scalar Arguments ..
115  INTEGER INCX, N
116  DOUBLE PRECISION ALPHA, TAU
117 * ..
118 * .. Array Arguments ..
119  DOUBLE PRECISION X( * )
120 * ..
121 *
122 * =====================================================================
123 *
124 * .. Parameters ..
125  DOUBLE PRECISION ONE, ZERO
126  parameter ( one = 1.0d+0, zero = 0.0d+0 )
127 * ..
128 * .. Local Scalars ..
129  INTEGER J, KNT
130  DOUBLE PRECISION BETA, RSAFMN, SAFMIN, XNORM
131 * ..
132 * .. External Functions ..
133  DOUBLE PRECISION DLAMCH, DLAPY2, DNRM2
134  EXTERNAL dlamch, dlapy2, dnrm2
135 * ..
136 * .. Intrinsic Functions ..
137  INTRINSIC abs, sign
138 * ..
139 * .. External Subroutines ..
140  EXTERNAL dscal
141 * ..
142 * .. Executable Statements ..
143 *
144  IF( n.LE.1 ) THEN
145  tau = zero
146  RETURN
147  END IF
148 *
149  xnorm = dnrm2( n-1, x, incx )
150 *
151  IF( xnorm.EQ.zero ) THEN
152 *
153 * H = I
154 *
155  tau = zero
156  ELSE
157 *
158 * general case
159 *
160  beta = -sign( dlapy2( alpha, xnorm ), alpha )
161  safmin = dlamch( 'S' ) / dlamch( 'E' )
162  knt = 0
163  IF( abs( beta ).LT.safmin ) THEN
164 *
165 * XNORM, BETA may be inaccurate; scale X and recompute them
166 *
167  rsafmn = one / safmin
168  10 CONTINUE
169  knt = knt + 1
170  CALL dscal( n-1, rsafmn, x, incx )
171  beta = beta*rsafmn
172  alpha = alpha*rsafmn
173  IF( abs( beta ).LT.safmin )
174  \$ GO TO 10
175 *
176 * New BETA is at most 1, at least SAFMIN
177 *
178  xnorm = dnrm2( n-1, x, incx )
179  beta = -sign( dlapy2( alpha, xnorm ), alpha )
180  END IF
181  tau = ( beta-alpha ) / beta
182  CALL dscal( n-1, one / ( alpha-beta ), x, incx )
183 *
184 * If ALPHA is subnormal, it may lose relative accuracy
185 *
186  DO 20 j = 1, knt
187  beta = beta*safmin
188  20 CONTINUE
189  alpha = beta
190  END IF
191 *
192  RETURN
193 *
194 * End of DLARFG
195 *
196  END
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:55
subroutine dlarfg(N, ALPHA, X, INCX, TAU)
DLARFG generates an elementary reflector (Householder matrix).
Definition: dlarfg.f:108