LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dchkgt.f
Go to the documentation of this file.
1*> \brief \b DCHKGT
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8* Definition:
9* ===========
10*
11* SUBROUTINE DCHKGT( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
12* A, AF, B, X, XACT, WORK, RWORK, IWORK, NOUT )
13*
14* .. Scalar Arguments ..
15* LOGICAL TSTERR
16* INTEGER NN, NNS, NOUT
17* DOUBLE PRECISION THRESH
18* ..
19* .. Array Arguments ..
20* LOGICAL DOTYPE( * )
21* INTEGER IWORK( * ), NSVAL( * ), NVAL( * )
22* DOUBLE PRECISION A( * ), AF( * ), B( * ), RWORK( * ), WORK( * ),
23* $ X( * ), XACT( * )
24* ..
25*
26*
27*> \par Purpose:
28* =============
29*>
30*> \verbatim
31*>
32*> DCHKGT tests DGTTRF, -TRS, -RFS, and -CON
33*> \endverbatim
34*
35* Arguments:
36* ==========
37*
38*> \param[in] DOTYPE
39*> \verbatim
40*> DOTYPE is LOGICAL array, dimension (NTYPES)
41*> The matrix types to be used for testing. Matrices of type j
42*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
43*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
44*> \endverbatim
45*>
46*> \param[in] NN
47*> \verbatim
48*> NN is INTEGER
49*> The number of values of N contained in the vector NVAL.
50*> \endverbatim
51*>
52*> \param[in] NVAL
53*> \verbatim
54*> NVAL is INTEGER array, dimension (NN)
55*> The values of the matrix dimension N.
56*> \endverbatim
57*>
58*> \param[in] NNS
59*> \verbatim
60*> NNS is INTEGER
61*> The number of values of NRHS contained in the vector NSVAL.
62*> \endverbatim
63*>
64*> \param[in] NSVAL
65*> \verbatim
66*> NSVAL is INTEGER array, dimension (NNS)
67*> The values of the number of right hand sides NRHS.
68*> \endverbatim
69*>
70*> \param[in] THRESH
71*> \verbatim
72*> THRESH is DOUBLE PRECISION
73*> The threshold value for the test ratios. A result is
74*> included in the output file if RESULT >= THRESH. To have
75*> every test ratio printed, use THRESH = 0.
76*> \endverbatim
77*>
78*> \param[in] TSTERR
79*> \verbatim
80*> TSTERR is LOGICAL
81*> Flag that indicates whether error exits are to be tested.
82*> \endverbatim
83*>
84*> \param[out] A
85*> \verbatim
86*> A is DOUBLE PRECISION array, dimension (NMAX*4)
87*> \endverbatim
88*>
89*> \param[out] AF
90*> \verbatim
91*> AF is DOUBLE PRECISION array, dimension (NMAX*4)
92*> \endverbatim
93*>
94*> \param[out] B
95*> \verbatim
96*> B is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
97*> where NSMAX is the largest entry in NSVAL.
98*> \endverbatim
99*>
100*> \param[out] X
101*> \verbatim
102*> X is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
103*> \endverbatim
104*>
105*> \param[out] XACT
106*> \verbatim
107*> XACT is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
108*> \endverbatim
109*>
110*> \param[out] WORK
111*> \verbatim
112*> WORK is DOUBLE PRECISION array, dimension
113*> (NMAX*max(3,NSMAX))
114*> \endverbatim
115*>
116*> \param[out] RWORK
117*> \verbatim
118*> RWORK is DOUBLE PRECISION array, dimension
119*> (max(NMAX,2*NSMAX))
120*> \endverbatim
121*>
122*> \param[out] IWORK
123*> \verbatim
124*> IWORK is INTEGER array, dimension (2*NMAX)
125*> \endverbatim
126*>
127*> \param[in] NOUT
128*> \verbatim
129*> NOUT is INTEGER
130*> The unit number for output.
131*> \endverbatim
132*
133* Authors:
134* ========
135*
136*> \author Univ. of Tennessee
137*> \author Univ. of California Berkeley
138*> \author Univ. of Colorado Denver
139*> \author NAG Ltd.
140*
141*> \ingroup double_lin
142*
143* =====================================================================
144 SUBROUTINE dchkgt( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
145 $ A, AF, B, X, XACT, WORK, RWORK, IWORK, NOUT )
146*
147* -- LAPACK test routine --
148* -- LAPACK is a software package provided by Univ. of Tennessee, --
149* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
150*
151* .. Scalar Arguments ..
152 LOGICAL TSTERR
153 INTEGER NN, NNS, NOUT
154 DOUBLE PRECISION THRESH
155* ..
156* .. Array Arguments ..
157 LOGICAL DOTYPE( * )
158 INTEGER IWORK( * ), NSVAL( * ), NVAL( * )
159 DOUBLE PRECISION A( * ), AF( * ), B( * ), RWORK( * ), WORK( * ),
160 $ x( * ), xact( * )
161* ..
162*
163* =====================================================================
164*
165* .. Parameters ..
166 DOUBLE PRECISION ONE, ZERO
167 parameter( one = 1.0d+0, zero = 0.0d+0 )
168 INTEGER NTYPES
169 parameter( ntypes = 12 )
170 INTEGER NTESTS
171 parameter( ntests = 7 )
172* ..
173* .. Local Scalars ..
174 LOGICAL TRFCON, ZEROT
175 CHARACTER DIST, NORM, TRANS, TYPE
176 CHARACTER*3 PATH
177 INTEGER I, IMAT, IN, INFO, IRHS, ITRAN, IX, IZERO, J,
178 $ k, kl, koff, ku, lda, m, mode, n, nerrs, nfail,
179 $ nimat, nrhs, nrun
180 DOUBLE PRECISION AINVNM, ANORM, COND, RCOND, RCONDC, RCONDI,
181 $ rcondo
182* ..
183* .. Local Arrays ..
184 CHARACTER TRANSS( 3 )
185 INTEGER ISEED( 4 ), ISEEDY( 4 )
186 DOUBLE PRECISION RESULT( NTESTS ), Z( 3 )
187* ..
188* .. External Functions ..
189 DOUBLE PRECISION DASUM, DGET06, DLANGT
190 EXTERNAL dasum, dget06, dlangt
191* ..
192* .. External Subroutines ..
193 EXTERNAL alaerh, alahd, alasum, dcopy, derrge, dget04,
196 $ dscal
197* ..
198* .. Intrinsic Functions ..
199 INTRINSIC max
200* ..
201* .. Scalars in Common ..
202 LOGICAL LERR, OK
203 CHARACTER*32 SRNAMT
204 INTEGER INFOT, NUNIT
205* ..
206* .. Common blocks ..
207 COMMON / infoc / infot, nunit, ok, lerr
208 COMMON / srnamc / srnamt
209* ..
210* .. Data statements ..
211 DATA iseedy / 0, 0, 0, 1 / , transs / 'N', 'T',
212 $ 'C' /
213* ..
214* .. Executable Statements ..
215*
216 path( 1: 1 ) = 'Double precision'
217 path( 2: 3 ) = 'GT'
218 nrun = 0
219 nfail = 0
220 nerrs = 0
221 DO 10 i = 1, 4
222 iseed( i ) = iseedy( i )
223 10 CONTINUE
224*
225* Test the error exits
226*
227 IF( tsterr )
228 $ CALL derrge( path, nout )
229 infot = 0
230*
231 DO 110 in = 1, nn
232*
233* Do for each value of N in NVAL.
234*
235 n = nval( in )
236 m = max( n-1, 0 )
237 lda = max( 1, n )
238 nimat = ntypes
239 IF( n.LE.0 )
240 $ nimat = 1
241*
242 DO 100 imat = 1, nimat
243*
244* Do the tests only if DOTYPE( IMAT ) is true.
245*
246 IF( .NOT.dotype( imat ) )
247 $ GO TO 100
248*
249* Set up parameters with DLATB4.
250*
251 CALL dlatb4( path, imat, n, n, TYPE, kl, ku, anorm, mode,
252 $ cond, dist )
253*
254 zerot = imat.GE.8 .AND. imat.LE.10
255 IF( imat.LE.6 ) THEN
256*
257* Types 1-6: generate matrices of known condition number.
258*
259 koff = max( 2-ku, 3-max( 1, n ) )
260 srnamt = 'DLATMS'
261 CALL dlatms( n, n, dist, iseed, TYPE, rwork, mode, cond,
262 $ anorm, kl, ku, 'Z', af( koff ), 3, work,
263 $ info )
264*
265* Check the error code from DLATMS.
266*
267 IF( info.NE.0 ) THEN
268 CALL alaerh( path, 'DLATMS', info, 0, ' ', n, n, kl,
269 $ ku, -1, imat, nfail, nerrs, nout )
270 GO TO 100
271 END IF
272 izero = 0
273*
274 IF( n.GT.1 ) THEN
275 CALL dcopy( n-1, af( 4 ), 3, a, 1 )
276 CALL dcopy( n-1, af( 3 ), 3, a( n+m+1 ), 1 )
277 END IF
278 CALL dcopy( n, af( 2 ), 3, a( m+1 ), 1 )
279 ELSE
280*
281* Types 7-12: generate tridiagonal matrices with
282* unknown condition numbers.
283*
284 IF( .NOT.zerot .OR. .NOT.dotype( 7 ) ) THEN
285*
286* Generate a matrix with elements from [-1,1].
287*
288 CALL dlarnv( 2, iseed, n+2*m, a )
289 IF( anorm.NE.one )
290 $ CALL dscal( n+2*m, anorm, a, 1 )
291 ELSE IF( izero.GT.0 ) THEN
292*
293* Reuse the last matrix by copying back the zeroed out
294* elements.
295*
296 IF( izero.EQ.1 ) THEN
297 a( n ) = z( 2 )
298 IF( n.GT.1 )
299 $ a( 1 ) = z( 3 )
300 ELSE IF( izero.EQ.n ) THEN
301 a( 3*n-2 ) = z( 1 )
302 a( 2*n-1 ) = z( 2 )
303 ELSE
304 a( 2*n-2+izero ) = z( 1 )
305 a( n-1+izero ) = z( 2 )
306 a( izero ) = z( 3 )
307 END IF
308 END IF
309*
310* If IMAT > 7, set one column of the matrix to 0.
311*
312 IF( .NOT.zerot ) THEN
313 izero = 0
314 ELSE IF( imat.EQ.8 ) THEN
315 izero = 1
316 z( 2 ) = a( n )
317 a( n ) = zero
318 IF( n.GT.1 ) THEN
319 z( 3 ) = a( 1 )
320 a( 1 ) = zero
321 END IF
322 ELSE IF( imat.EQ.9 ) THEN
323 izero = n
324 z( 1 ) = a( 3*n-2 )
325 z( 2 ) = a( 2*n-1 )
326 a( 3*n-2 ) = zero
327 a( 2*n-1 ) = zero
328 ELSE
329 izero = ( n+1 ) / 2
330 DO 20 i = izero, n - 1
331 a( 2*n-2+i ) = zero
332 a( n-1+i ) = zero
333 a( i ) = zero
334 20 CONTINUE
335 a( 3*n-2 ) = zero
336 a( 2*n-1 ) = zero
337 END IF
338 END IF
339*
340*+ TEST 1
341* Factor A as L*U and compute the ratio
342* norm(L*U - A) / (n * norm(A) * EPS )
343*
344 CALL dcopy( n+2*m, a, 1, af, 1 )
345 srnamt = 'DGTTRF'
346 CALL dgttrf( n, af, af( m+1 ), af( n+m+1 ), af( n+2*m+1 ),
347 $ iwork, info )
348*
349* Check error code from DGTTRF.
350*
351 IF( info.NE.izero )
352 $ CALL alaerh( path, 'DGTTRF', info, izero, ' ', n, n, 1,
353 $ 1, -1, imat, nfail, nerrs, nout )
354 trfcon = info.NE.0
355*
356 CALL dgtt01( n, a, a( m+1 ), a( n+m+1 ), af, af( m+1 ),
357 $ af( n+m+1 ), af( n+2*m+1 ), iwork, work, lda,
358 $ rwork, result( 1 ) )
359*
360* Print the test ratio if it is .GE. THRESH.
361*
362 IF( result( 1 ).GE.thresh ) THEN
363 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
364 $ CALL alahd( nout, path )
365 WRITE( nout, fmt = 9999 )n, imat, 1, result( 1 )
366 nfail = nfail + 1
367 END IF
368 nrun = nrun + 1
369*
370 DO 50 itran = 1, 2
371 trans = transs( itran )
372 IF( itran.EQ.1 ) THEN
373 norm = 'O'
374 ELSE
375 norm = 'I'
376 END IF
377 anorm = dlangt( norm, n, a, a( m+1 ), a( n+m+1 ) )
378*
379 IF( .NOT.trfcon ) THEN
380*
381* Use DGTTRS to solve for one column at a time of inv(A)
382* or inv(A^T), computing the maximum column sum as we
383* go.
384*
385 ainvnm = zero
386 DO 40 i = 1, n
387 DO 30 j = 1, n
388 x( j ) = zero
389 30 CONTINUE
390 x( i ) = one
391 CALL dgttrs( trans, n, 1, af, af( m+1 ),
392 $ af( n+m+1 ), af( n+2*m+1 ), iwork, x,
393 $ lda, info )
394 ainvnm = max( ainvnm, dasum( n, x, 1 ) )
395 40 CONTINUE
396*
397* Compute RCONDC = 1 / (norm(A) * norm(inv(A))
398*
399 IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
400 rcondc = one
401 ELSE
402 rcondc = ( one / anorm ) / ainvnm
403 END IF
404 IF( itran.EQ.1 ) THEN
405 rcondo = rcondc
406 ELSE
407 rcondi = rcondc
408 END IF
409 ELSE
410 rcondc = zero
411 END IF
412*
413*+ TEST 7
414* Estimate the reciprocal of the condition number of the
415* matrix.
416*
417 srnamt = 'DGTCON'
418 CALL dgtcon( norm, n, af, af( m+1 ), af( n+m+1 ),
419 $ af( n+2*m+1 ), iwork, anorm, rcond, work,
420 $ iwork( n+1 ), info )
421*
422* Check error code from DGTCON.
423*
424 IF( info.NE.0 )
425 $ CALL alaerh( path, 'DGTCON', info, 0, norm, n, n, -1,
426 $ -1, -1, imat, nfail, nerrs, nout )
427*
428 result( 7 ) = dget06( rcond, rcondc )
429*
430* Print the test ratio if it is .GE. THRESH.
431*
432 IF( result( 7 ).GE.thresh ) THEN
433 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
434 $ CALL alahd( nout, path )
435 WRITE( nout, fmt = 9997 )norm, n, imat, 7,
436 $ result( 7 )
437 nfail = nfail + 1
438 END IF
439 nrun = nrun + 1
440 50 CONTINUE
441*
442* Skip the remaining tests if the matrix is singular.
443*
444 IF( trfcon )
445 $ GO TO 100
446*
447 DO 90 irhs = 1, nns
448 nrhs = nsval( irhs )
449*
450* Generate NRHS random solution vectors.
451*
452 ix = 1
453 DO 60 j = 1, nrhs
454 CALL dlarnv( 2, iseed, n, xact( ix ) )
455 ix = ix + lda
456 60 CONTINUE
457*
458 DO 80 itran = 1, 3
459 trans = transs( itran )
460 IF( itran.EQ.1 ) THEN
461 rcondc = rcondo
462 ELSE
463 rcondc = rcondi
464 END IF
465*
466* Set the right hand side.
467*
468 CALL dlagtm( trans, n, nrhs, one, a, a( m+1 ),
469 $ a( n+m+1 ), xact, lda, zero, b, lda )
470*
471*+ TEST 2
472* Solve op(A) * X = B and compute the residual.
473*
474 CALL dlacpy( 'Full', n, nrhs, b, lda, x, lda )
475 srnamt = 'DGTTRS'
476 CALL dgttrs( trans, n, nrhs, af, af( m+1 ),
477 $ af( n+m+1 ), af( n+2*m+1 ), iwork, x,
478 $ lda, info )
479*
480* Check error code from DGTTRS.
481*
482 IF( info.NE.0 )
483 $ CALL alaerh( path, 'DGTTRS', info, 0, trans, n, n,
484 $ -1, -1, nrhs, imat, nfail, nerrs,
485 $ nout )
486*
487 CALL dlacpy( 'Full', n, nrhs, b, lda, work, lda )
488 CALL dgtt02( trans, n, nrhs, a, a( m+1 ), a( n+m+1 ),
489 $ x, lda, work, lda, result( 2 ) )
490*
491*+ TEST 3
492* Check solution from generated exact solution.
493*
494 CALL dget04( n, nrhs, x, lda, xact, lda, rcondc,
495 $ result( 3 ) )
496*
497*+ TESTS 4, 5, and 6
498* Use iterative refinement to improve the solution.
499*
500 srnamt = 'DGTRFS'
501 CALL dgtrfs( trans, n, nrhs, a, a( m+1 ), a( n+m+1 ),
502 $ af, af( m+1 ), af( n+m+1 ),
503 $ af( n+2*m+1 ), iwork, b, lda, x, lda,
504 $ rwork, rwork( nrhs+1 ), work,
505 $ iwork( n+1 ), info )
506*
507* Check error code from DGTRFS.
508*
509 IF( info.NE.0 )
510 $ CALL alaerh( path, 'DGTRFS', info, 0, trans, n, n,
511 $ -1, -1, nrhs, imat, nfail, nerrs,
512 $ nout )
513*
514 CALL dget04( n, nrhs, x, lda, xact, lda, rcondc,
515 $ result( 4 ) )
516 CALL dgtt05( trans, n, nrhs, a, a( m+1 ), a( n+m+1 ),
517 $ b, lda, x, lda, xact, lda, rwork,
518 $ rwork( nrhs+1 ), result( 5 ) )
519*
520* Print information about the tests that did not pass
521* the threshold.
522*
523 DO 70 k = 2, 6
524 IF( result( k ).GE.thresh ) THEN
525 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
526 $ CALL alahd( nout, path )
527 WRITE( nout, fmt = 9998 )trans, n, nrhs, imat,
528 $ k, result( k )
529 nfail = nfail + 1
530 END IF
531 70 CONTINUE
532 nrun = nrun + 5
533 80 CONTINUE
534 90 CONTINUE
535*
536 100 CONTINUE
537 110 CONTINUE
538*
539* Print a summary of the results.
540*
541 CALL alasum( path, nout, nfail, nrun, nerrs )
542*
543 9999 FORMAT( 12x, 'N =', i5, ',', 10x, ' type ', i2, ', test(', i2,
544 $ ') = ', g12.5 )
545 9998 FORMAT( ' TRANS=''', a1, ''', N =', i5, ', NRHS=', i3, ', type ',
546 $ i2, ', test(', i2, ') = ', g12.5 )
547 9997 FORMAT( ' NORM =''', a1, ''', N =', i5, ',', 10x, ' type ', i2,
548 $ ', test(', i2, ') = ', g12.5 )
549 RETURN
550*
551* End of DCHKGT
552*
553 END
subroutine alasum(type, nout, nfail, nrun, nerrs)
ALASUM
Definition alasum.f:73
subroutine alaerh(path, subnam, info, infoe, opts, m, n, kl, ku, n5, imat, nfail, nerrs, nout)
ALAERH
Definition alaerh.f:147
subroutine alahd(iounit, path)
ALAHD
Definition alahd.f:107
subroutine dchkgt(dotype, nn, nval, nns, nsval, thresh, tsterr, a, af, b, x, xact, work, rwork, iwork, nout)
DCHKGT
Definition dchkgt.f:146
subroutine derrge(path, nunit)
DERRGE
Definition derrge.f:55
subroutine dget04(n, nrhs, x, ldx, xact, ldxact, rcond, resid)
DGET04
Definition dget04.f:102
subroutine dgtt01(n, dl, d, du, dlf, df, duf, du2, ipiv, work, ldwork, rwork, resid)
DGTT01
Definition dgtt01.f:134
subroutine dgtt02(trans, n, nrhs, dl, d, du, x, ldx, b, ldb, resid)
DGTT02
Definition dgtt02.f:125
subroutine dgtt05(trans, n, nrhs, dl, d, du, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
DGTT05
Definition dgtt05.f:165
subroutine dlatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
DLATB4
Definition dlatb4.f:120
subroutine dlatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
DLATMS
Definition dlatms.f:321
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
Definition dcopy.f:82
subroutine dgtcon(norm, n, dl, d, du, du2, ipiv, anorm, rcond, work, iwork, info)
DGTCON
Definition dgtcon.f:146
subroutine dgtrfs(trans, n, nrhs, dl, d, du, dlf, df, duf, du2, ipiv, b, ldb, x, ldx, ferr, berr, work, iwork, info)
DGTRFS
Definition dgtrfs.f:209
subroutine dgttrf(n, dl, d, du, du2, ipiv, info)
DGTTRF
Definition dgttrf.f:124
subroutine dgttrs(trans, n, nrhs, dl, d, du, du2, ipiv, b, ldb, info)
DGTTRS
Definition dgttrs.f:138
subroutine dlacpy(uplo, m, n, a, lda, b, ldb)
DLACPY copies all or part of one two-dimensional array to another.
Definition dlacpy.f:103
subroutine dlagtm(trans, n, nrhs, alpha, dl, d, du, x, ldx, beta, b, ldb)
DLAGTM performs a matrix-matrix product of the form C = αAB+βC, where A is a tridiagonal matrix,...
Definition dlagtm.f:145
subroutine dlarnv(idist, iseed, n, x)
DLARNV returns a vector of random numbers from a uniform or normal distribution.
Definition dlarnv.f:97
subroutine dscal(n, da, dx, incx)
DSCAL
Definition dscal.f:79