124 parameter( two = 2.0e+0, one = 1.0e+0, zero = 0.0e+0 )
128 REAL alphi, alphr, beta, bignum, smlnum, xnorm
137 INTRINSIC abs, aimag, cmplx,
REAL, sign
150 alphr =
REAL( alpha )
151 alphi = aimag( alpha )
153 IF( xnorm.EQ.zero )
THEN
157 IF( alphi.EQ.zero )
THEN
158 IF( alphr.GE.zero )
THEN
168 x( 1 + (
j-1)*incx ) = zero
174 xnorm =
slapy2( alphr, alphi )
175 tau = cmplx( one - alphr / xnorm, -alphi / xnorm )
177 x( 1 + (
j-1)*incx ) = zero
185 beta = sign(
slapy3( alphr, alphi, xnorm ), alphr )
187 bignum = one / smlnum
190 IF( abs( beta ).LT.smlnum )
THEN
196 CALL
csscal( n-1, bignum,
x, incx )
200 IF( abs( beta ).LT.smlnum )
206 alpha = cmplx( alphr, alphi )
207 beta = sign(
slapy3( alphr, alphi, xnorm ), alphr )
211 IF( beta.LT.zero )
THEN
215 alphr = alphi * (alphi/
REAL( alpha ))
216 alphr = alphr + xnorm * (xnorm/
REAL( alpha ))
217 tau = cmplx( alphr/beta, -alphi/beta )
218 alpha = cmplx( -alphr, alphi )
220 alpha =
cladiv( cmplx( one ), alpha )
222 IF ( abs(tau).LE.smlnum )
THEN
231 alphr =
REAL( savealpha )
232 alphi = aimag( savealpha )
233 IF( alphi.EQ.zero )
THEN
234 IF( alphr.GE.zero )
THEN
239 x( 1 + (
j-1)*incx ) = zero
244 xnorm =
slapy2( alphr, alphi )
245 tau = cmplx( one - alphr / xnorm, -alphi / xnorm )
247 x( 1 + (
j-1)*incx ) = zero
256 CALL
cscal( n-1, alpha,
x, incx )
subroutine cscal(N, CA, CX, INCX)
CSCAL
COMPLEX function cladiv(X, Y)
CLADIV performs complex division in real arithmetic, avoiding unnecessary overflow.
REAL function slamch(CMACH)
SLAMCH
REAL function slapy2(X, Y)
SLAPY2 returns sqrt(x2+y2).
REAL function slapy3(X, Y, Z)
SLAPY3 returns sqrt(x2+y2+z2).
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine clarfgp(N, ALPHA, X, INCX, TAU)
CLARFGP generates an elementary reflector (Householder matrix) with non-negatibe beta.
REAL function scnrm2(N, X, INCX)
SCNRM2
subroutine csscal(N, SA, CX, INCX)
CSSCAL