LAPACK  3.6.1 LAPACK: Linear Algebra PACKage
 subroutine slarfgp ( integer N, real ALPHA, real, dimension( * ) X, integer INCX, real TAU )

SLARFGP generates an elementary reflector (Householder matrix) with non-negative beta.

Purpose:
``` SLARFGP generates a real elementary reflector H of order n, such
that

H * ( alpha ) = ( beta ),   H**T * H = I.
(   x   )   (   0  )

where alpha and beta are scalars, beta is non-negative, and x is
an (n-1)-element real vector.  H is represented in the form

H = I - tau * ( 1 ) * ( 1 v**T ) ,
( v )

where tau is a real scalar and v is a real (n-1)-element
vector.

If the elements of x are all zero, then tau = 0 and H is taken to be
the unit matrix.```
Parameters
 [in] N ``` N is INTEGER The order of the elementary reflector.``` [in,out] ALPHA ``` ALPHA is REAL On entry, the value alpha. On exit, it is overwritten with the value beta.``` [in,out] X ``` X is REAL array, dimension (1+(N-2)*abs(INCX)) On entry, the vector x. On exit, it is overwritten with the vector v.``` [in] INCX ``` INCX is INTEGER The increment between elements of X. INCX > 0.``` [out] TAU ``` TAU is REAL The value tau.```
Date
November 2015

Definition at line 106 of file slarfgp.f.

106 *
107 * -- LAPACK auxiliary routine (version 3.6.0) --
108 * -- LAPACK is a software package provided by Univ. of Tennessee, --
109 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
110 * November 2015
111 *
112 * .. Scalar Arguments ..
113  INTEGER incx, n
114  REAL alpha, tau
115 * ..
116 * .. Array Arguments ..
117  REAL x( * )
118 * ..
119 *
120 * =====================================================================
121 *
122 * .. Parameters ..
123  REAL two, one, zero
124  parameter ( two = 2.0e+0, one = 1.0e+0, zero = 0.0e+0 )
125 * ..
126 * .. Local Scalars ..
127  INTEGER j, knt
128  REAL beta, bignum, savealpha, smlnum, xnorm
129 * ..
130 * .. External Functions ..
131  REAL slamch, slapy2, snrm2
132  EXTERNAL slamch, slapy2, snrm2
133 * ..
134 * .. Intrinsic Functions ..
135  INTRINSIC abs, sign
136 * ..
137 * .. External Subroutines ..
138  EXTERNAL sscal
139 * ..
140 * .. Executable Statements ..
141 *
142  IF( n.LE.0 ) THEN
143  tau = zero
144  RETURN
145  END IF
146 *
147  xnorm = snrm2( n-1, x, incx )
148 *
149  IF( xnorm.EQ.zero ) THEN
150 *
151 * H = [+/-1, 0; I], sign chosen so ALPHA >= 0.
152 *
153  IF( alpha.GE.zero ) THEN
154 * When TAU.eq.ZERO, the vector is special-cased to be
155 * all zeros in the application routines. We do not need
156 * to clear it.
157  tau = zero
158  ELSE
159 * However, the application routines rely on explicit
160 * zero checks when TAU.ne.ZERO, and we must clear X.
161  tau = two
162  DO j = 1, n-1
163  x( 1 + (j-1)*incx ) = 0
164  END DO
165  alpha = -alpha
166  END IF
167  ELSE
168 *
169 * general case
170 *
171  beta = sign( slapy2( alpha, xnorm ), alpha )
172  smlnum = slamch( 'S' ) / slamch( 'E' )
173  knt = 0
174  IF( abs( beta ).LT.smlnum ) THEN
175 *
176 * XNORM, BETA may be inaccurate; scale X and recompute them
177 *
178  bignum = one / smlnum
179  10 CONTINUE
180  knt = knt + 1
181  CALL sscal( n-1, bignum, x, incx )
182  beta = beta*bignum
183  alpha = alpha*bignum
184  IF( abs( beta ).LT.smlnum )
185  \$ GO TO 10
186 *
187 * New BETA is at most 1, at least SMLNUM
188 *
189  xnorm = snrm2( n-1, x, incx )
190  beta = sign( slapy2( alpha, xnorm ), alpha )
191  END IF
192  savealpha = alpha
193  alpha = alpha + beta
194  IF( beta.LT.zero ) THEN
195  beta = -beta
196  tau = -alpha / beta
197  ELSE
198  alpha = xnorm * (xnorm/alpha)
199  tau = alpha / beta
200  alpha = -alpha
201  END IF
202 *
203  IF ( abs(tau).LE.smlnum ) THEN
204 *
205 * In the case where the computed TAU ends up being a denormalized number,
206 * it loses relative accuracy. This is a BIG problem. Solution: flush TAU
207 * to ZERO. This explains the next IF statement.
208 *
209 * (Bug report provided by Pat Quillen from MathWorks on Jul 29, 2009.)
210 * (Thanks Pat. Thanks MathWorks.)
211 *
212  IF( savealpha.GE.zero ) THEN
213  tau = zero
214  ELSE
215  tau = two
216  DO j = 1, n-1
217  x( 1 + (j-1)*incx ) = 0
218  END DO
219  beta = -savealpha
220  END IF
221 *
222  ELSE
223 *
224 * This is the general case.
225 *
226  CALL sscal( n-1, one / alpha, x, incx )
227 *
228  END IF
229 *
230 * If BETA is subnormal, it may lose relative accuracy
231 *
232  DO 20 j = 1, knt
233  beta = beta*smlnum
234  20 CONTINUE
235  alpha = beta
236  END IF
237 *
238  RETURN
239 *
240 * End of SLARFGP
241 *
real function slapy2(X, Y)
SLAPY2 returns sqrt(x2+y2).
Definition: slapy2.f:65
real function snrm2(N, X, INCX)
SNRM2
Definition: snrm2.f:56
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:55
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69

Here is the call graph for this function:

Here is the caller graph for this function: