104 SUBROUTINE sget39( RMAX, LMAX, NINFO, KNT )
112 INTEGER knt, lmax, ninfo
120 parameter( ldt = 10, ldt2 = 2*ldt )
122 parameter( zero = 0.0, one = 1.0 )
125 INTEGER i, info, ivm1, ivm2, ivm3, ivm4, ivm5,
j, k, n,
127 REAL bignum, domin, dumm, eps, norm, normtb, resid,
128 $ scale, smlnum, w, xnorm
139 INTRINSIC abs, cos, max,
REAL, sin, sqrt
142 INTEGER idim( 6 ), ival( 5, 5, 6 )
143 REAL b( ldt ), d( ldt2 ), dum( 1 ), t( ldt, ldt ),
144 $ vm1( 5 ), vm2( 5 ), vm3( 5 ), vm4( 5 ),
145 $ vm5( 3 ), work( ldt ),
x( ldt2 ), y( ldt2 )
149 DATA ival / 3, 4*0, 1, 1, -1, 0, 0, 3, 2, 1, 0, 0,
150 $ 4, 3, 2, 2, 0, 5*0, 1, 4*0, 2, 2, 3*0, 3, 3, 4,
151 $ 0, 0, 4, 2, 2, 3, 0, 4*1, 5, 1, 4*0, 2, 4, -2,
152 $ 0, 0, 3, 3, 4, 0, 0, 4, 2, 2, 3, 0, 5*1, 1,
153 $ 4*0, 2, 1, -1, 0, 0, 9, 8, 1, 0, 0, 4, 9, 1, 2,
154 $ -1, 5*2, 9, 4*0, 6, 4, 0, 0, 0, 3, 2, 1, 1, 0,
155 $ 5, 1, -1, 1, 0, 5*2, 4, 4*0, 2, 2, 0, 0, 0, 1,
156 $ 4, 4, 0, 0, 2, 4, 2, 2, -1, 5*2 /
164 bignum = one / smlnum
165 CALL
slabad( smlnum, bignum )
170 vm1( 2 ) = sqrt( smlnum )
171 vm1( 3 ) = sqrt( vm1( 2 ) )
172 vm1( 4 ) = sqrt( bignum )
173 vm1( 5 ) = sqrt( vm1( 4 ) )
176 vm2( 2 ) = sqrt( smlnum )
177 vm2( 3 ) = sqrt( vm2( 2 ) )
178 vm2( 4 ) = sqrt( bignum )
179 vm2( 5 ) = sqrt( vm2( 4 ) )
182 vm3( 2 ) = sqrt( smlnum )
183 vm3( 3 ) = sqrt( vm3( 2 ) )
184 vm3( 4 ) = sqrt( bignum )
185 vm3( 5 ) = sqrt( vm3( 4 ) )
188 vm4( 2 ) = sqrt( smlnum )
189 vm4( 3 ) = sqrt( vm4( 2 ) )
190 vm4( 4 ) = sqrt( bignum )
191 vm4( 5 ) = sqrt( vm4( 4 ) )
195 vm5( 3 ) = sqrt( smlnum )
202 smlnum = smlnum / eps
216 t( i,
j ) =
REAL( IVAL( I, J, NDIM ) )*
219 $ t( i,
j ) = t( i,
j )*vm5( ivm5 )
226 b( i ) = cos(
REAL( I ) )*vm3( ivm3 )
230 d( i ) = sin(
REAL( I ) )*vm4( ivm4 )
233 norm =
slange(
'1', n, n, t, ldt, work )
235 normtb = norm + abs(
b( k ) ) + abs( w )
237 CALL
scopy( n, d, 1,
x, 1 )
239 CALL
slaqtr( .false., .true., n, t, ldt, dum,
240 $ dumm, scale,
x, work, info )
247 CALL
scopy( n, d, 1, y, 1 )
248 CALL
sgemv(
'No transpose', n, n, one, t, ldt,
249 $
x, 1, -scale, y, 1 )
251 resid =
sasum( n, y, 1 )
252 domin = max( smlnum, ( smlnum / eps )*norm,
253 $ ( norm*eps )*xnorm )
254 resid = resid / domin
255 IF( resid.GT.rmax )
THEN
260 CALL
scopy( n, d, 1,
x, 1 )
262 CALL
slaqtr( .true., .true., n, t, ldt, dum,
263 $ dumm, scale,
x, work, info )
270 CALL
scopy( n, d, 1, y, 1 )
271 CALL
sgemv(
'Transpose', n, n, one, t, ldt,
x,
274 resid =
sasum( n, y, 1 )
275 domin = max( smlnum, ( smlnum / eps )*norm,
276 $ ( norm*eps )*xnorm )
277 resid = resid / domin
278 IF( resid.GT.rmax )
THEN
283 CALL
scopy( 2*n, d, 1,
x, 1 )
285 CALL
slaqtr( .false., .false., n, t, ldt,
b, w,
286 $ scale,
x, work, info )
295 CALL
scopy( 2*n, d, 1, y, 1 )
296 y( 1 ) =
sdot( n,
b, 1,
x( 1+n ), 1 ) +
299 y( i ) = w*
x( i+n ) + scale*y( i )
301 CALL
sgemv(
'No transpose', n, n, one, t, ldt,
304 y( 1+n ) =
sdot( n,
b, 1,
x, 1 ) -
307 y( i+n ) = w*
x( i ) - scale*y( i+n )
309 CALL
sgemv(
'No transpose', n, n, one, t, ldt,
310 $
x( 1+n ), 1, one, y( 1+n ), 1 )
312 resid =
sasum( 2*n, y, 1 )
313 domin = max( smlnum, ( smlnum / eps )*normtb,
314 $ eps*( normtb*
sasum( 2*n,
x, 1 ) ) )
315 resid = resid / domin
316 IF( resid.GT.rmax )
THEN
321 CALL
scopy( 2*n, d, 1,
x, 1 )
323 CALL
slaqtr( .true., .false., n, t, ldt,
b, w,
324 $ scale,
x, work, info )
332 CALL
scopy( 2*n, d, 1, y, 1 )
333 y( 1 ) =
b( 1 )*
x( 1+n ) - scale*y( 1 )
335 y( i ) =
b( i )*
x( 1+n ) + w*
x( i+n ) -
338 CALL
sgemv(
'Transpose', n, n, one, t, ldt,
x,
341 y( 1+n ) =
b( 1 )*
x( 1 ) + scale*y( 1+n )
343 y( i+n ) =
b( i )*
x( 1 ) + w*
x( i ) +
346 CALL
sgemv(
'Transpose', n, n, one, t, ldt,
347 $
x( 1+n ), 1, -one, y( 1+n ), 1 )
349 resid =
sasum( 2*n, y, 1 )
350 domin = max( smlnum, ( smlnum / eps )*normtb,
351 $ eps*( normtb*
sasum( 2*n,
x, 1 ) ) )
352 resid = resid / domin
353 IF( resid.GT.rmax )
THEN
INTEGER function isamax(N, SX, INCX)
ISAMAX
REAL function slamch(CMACH)
SLAMCH
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real b(3) integer i
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
subroutine sgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
SGEMV
REAL function slange(NORM, M, N, A, LDA, WORK)
SLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine slaqtr(LTRAN, LREAL, N, T, LDT, B, W, SCALE, X, WORK, INFO)
SLAQTR solves a real quasi-triangular system of equations, or a complex quasi-triangular system of sp...
subroutine slabad(SMALL, LARGE)
SLABAD
REAL function sasum(N, SX, INCX)
SASUM
subroutine sget39(RMAX, LMAX, NINFO, KNT)
SGET39
REAL function sdot(N, SX, INCX, SY, INCY)
SDOT