LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
cdrvhe.f
Go to the documentation of this file.
1*> \brief \b CDRVHE
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 CDRVHE( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
12* A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK,
13* NOUT )
14*
15* .. Scalar Arguments ..
16* LOGICAL TSTERR
17* INTEGER NMAX, NN, NOUT, NRHS
18* REAL THRESH
19* ..
20* .. Array Arguments ..
21* LOGICAL DOTYPE( * )
22* INTEGER IWORK( * ), NVAL( * )
23* REAL RWORK( * )
24* COMPLEX A( * ), AFAC( * ), AINV( * ), B( * ),
25* $ WORK( * ), X( * ), XACT( * )
26* ..
27*
28*
29*> \par Purpose:
30* =============
31*>
32*> \verbatim
33*>
34*> CDRVHE tests the driver routines CHESV and -SVX.
35*> \endverbatim
36*
37* Arguments:
38* ==========
39*
40*> \param[in] DOTYPE
41*> \verbatim
42*> DOTYPE is LOGICAL array, dimension (NTYPES)
43*> The matrix types to be used for testing. Matrices of type j
44*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
45*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
46*> \endverbatim
47*>
48*> \param[in] NN
49*> \verbatim
50*> NN is INTEGER
51*> The number of values of N contained in the vector NVAL.
52*> \endverbatim
53*>
54*> \param[in] NVAL
55*> \verbatim
56*> NVAL is INTEGER array, dimension (NN)
57*> The values of the matrix dimension N.
58*> \endverbatim
59*>
60*> \param[in] NRHS
61*> \verbatim
62*> NRHS is INTEGER
63*> The number of right hand side vectors to be generated for
64*> each linear system.
65*> \endverbatim
66*>
67*> \param[in] THRESH
68*> \verbatim
69*> THRESH is REAL
70*> The threshold value for the test ratios. A result is
71*> included in the output file if RESULT >= THRESH. To have
72*> every test ratio printed, use THRESH = 0.
73*> \endverbatim
74*>
75*> \param[in] TSTERR
76*> \verbatim
77*> TSTERR is LOGICAL
78*> Flag that indicates whether error exits are to be tested.
79*> \endverbatim
80*>
81*> \param[in] NMAX
82*> \verbatim
83*> NMAX is INTEGER
84*> The maximum value permitted for N, used in dimensioning the
85*> work arrays.
86*> \endverbatim
87*>
88*> \param[out] A
89*> \verbatim
90*> A is COMPLEX array, dimension (NMAX*NMAX)
91*> \endverbatim
92*>
93*> \param[out] AFAC
94*> \verbatim
95*> AFAC is COMPLEX array, dimension (NMAX*NMAX)
96*> \endverbatim
97*>
98*> \param[out] AINV
99*> \verbatim
100*> AINV is COMPLEX array, dimension (NMAX*NMAX)
101*> \endverbatim
102*>
103*> \param[out] B
104*> \verbatim
105*> B is COMPLEX array, dimension (NMAX*NRHS)
106*> \endverbatim
107*>
108*> \param[out] X
109*> \verbatim
110*> X is COMPLEX array, dimension (NMAX*NRHS)
111*> \endverbatim
112*>
113*> \param[out] XACT
114*> \verbatim
115*> XACT is COMPLEX array, dimension (NMAX*NRHS)
116*> \endverbatim
117*>
118*> \param[out] WORK
119*> \verbatim
120*> WORK is COMPLEX array, dimension (NMAX*max(2,NRHS))
121*> \endverbatim
122*>
123*> \param[out] RWORK
124*> \verbatim
125*> RWORK is REAL array, dimension (NMAX+2*NRHS)
126*> \endverbatim
127*>
128*> \param[out] IWORK
129*> \verbatim
130*> IWORK is INTEGER array, dimension (NMAX)
131*> \endverbatim
132*>
133*> \param[in] NOUT
134*> \verbatim
135*> NOUT is INTEGER
136*> The unit number for output.
137*> \endverbatim
138*
139* Authors:
140* ========
141*
142*> \author Univ. of Tennessee
143*> \author Univ. of California Berkeley
144*> \author Univ. of Colorado Denver
145*> \author NAG Ltd.
146*
147*> \ingroup complex_lin
148*
149* =====================================================================
150 SUBROUTINE cdrvhe( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
151 $ A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK,
152 $ NOUT )
153*
154* -- LAPACK test routine --
155* -- LAPACK is a software package provided by Univ. of Tennessee, --
156* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
157*
158* .. Scalar Arguments ..
159 LOGICAL TSTERR
160 INTEGER NMAX, NN, NOUT, NRHS
161 REAL THRESH
162* ..
163* .. Array Arguments ..
164 LOGICAL DOTYPE( * )
165 INTEGER IWORK( * ), NVAL( * )
166 REAL RWORK( * )
167 COMPLEX A( * ), AFAC( * ), AINV( * ), B( * ),
168 $ work( * ), x( * ), xact( * )
169* ..
170*
171* =====================================================================
172*
173* .. Parameters ..
174 REAL ONE, ZERO
175 PARAMETER ( ONE = 1.0e+0, zero = 0.0e+0 )
176 INTEGER NTYPES, NTESTS
177 parameter( ntypes = 10, ntests = 6 )
178 INTEGER NFACT
179 parameter( nfact = 2 )
180* ..
181* .. Local Scalars ..
182 LOGICAL ZEROT
183 CHARACTER DIST, FACT, TYPE, UPLO, XTYPE
184 CHARACTER*3 PATH
185 INTEGER I, I1, I2, IFACT, IMAT, IN, INFO, IOFF, IUPLO,
186 $ izero, j, k, k1, kl, ku, lda, lwork, mode, n,
187 $ nb, nbmin, nerrs, nfail, nimat, nrun, nt
188 REAL AINVNM, ANORM, CNDNUM, RCOND, RCONDC
189* ..
190* .. Local Arrays ..
191 CHARACTER FACTS( NFACT ), UPLOS( 2 )
192 INTEGER ISEED( 4 ), ISEEDY( 4 )
193 REAL RESULT( NTESTS )
194* ..
195* .. External Functions ..
196 REAL CLANHE, SGET06
197 EXTERNAL CLANHE, SGET06
198* ..
199* .. External Subroutines ..
200 EXTERNAL aladhd, alaerh, alasvm, cerrvx, cget04, chesv,
203 $ cpot05, xlaenv
204* ..
205* .. Scalars in Common ..
206 LOGICAL LERR, OK
207 CHARACTER*32 SRNAMT
208 INTEGER INFOT, NUNIT
209* ..
210* .. Common blocks ..
211 COMMON / infoc / infot, nunit, ok, lerr
212 COMMON / srnamc / srnamt
213* ..
214* .. Intrinsic Functions ..
215 INTRINSIC cmplx, max, min
216* ..
217* .. Data statements ..
218 DATA iseedy / 1988, 1989, 1990, 1991 /
219 DATA uplos / 'U', 'L' / , facts / 'F', 'N' /
220* ..
221* .. Executable Statements ..
222*
223* Initialize constants and the random number seed.
224*
225 path( 1: 1 ) = 'Complex precision'
226 path( 2: 3 ) = 'HE'
227 nrun = 0
228 nfail = 0
229 nerrs = 0
230 DO 10 i = 1, 4
231 iseed( i ) = iseedy( i )
232 10 CONTINUE
233 lwork = max( 2*nmax, nmax*nrhs )
234*
235* Test the error exits
236*
237 IF( tsterr )
238 $ CALL cerrvx( path, nout )
239 infot = 0
240*
241* Set the block size and minimum block size for testing.
242*
243 nb = 1
244 nbmin = 2
245 CALL xlaenv( 1, nb )
246 CALL xlaenv( 2, nbmin )
247*
248* Do for each value of N in NVAL
249*
250 DO 180 in = 1, nn
251 n = nval( in )
252 lda = max( n, 1 )
253 xtype = 'N'
254 nimat = ntypes
255 IF( n.LE.0 )
256 $ nimat = 1
257*
258 DO 170 imat = 1, nimat
259*
260* Do the tests only if DOTYPE( IMAT ) is true.
261*
262 IF( .NOT.dotype( imat ) )
263 $ GO TO 170
264*
265* Skip types 3, 4, 5, or 6 if the matrix size is too small.
266*
267 zerot = imat.GE.3 .AND. imat.LE.6
268 IF( zerot .AND. n.LT.imat-2 )
269 $ GO TO 170
270*
271* Do first for UPLO = 'U', then for UPLO = 'L'
272*
273 DO 160 iuplo = 1, 2
274 uplo = uplos( iuplo )
275*
276* Set up parameters with CLATB4 and generate a test matrix
277* with CLATMS.
278*
279 CALL clatb4( path, imat, n, n, TYPE, kl, ku, anorm, mode,
280 $ cndnum, dist )
281*
282 srnamt = 'CLATMS'
283 CALL clatms( n, n, dist, iseed, TYPE, rwork, mode,
284 $ cndnum, anorm, kl, ku, uplo, a, lda, work,
285 $ info )
286*
287* Check error code from CLATMS.
288*
289 IF( info.NE.0 ) THEN
290 CALL alaerh( path, 'CLATMS', info, 0, uplo, n, n, -1,
291 $ -1, -1, imat, nfail, nerrs, nout )
292 GO TO 160
293 END IF
294*
295* For types 3-6, zero one or more rows and columns of the
296* matrix to test that INFO is returned correctly.
297*
298 IF( zerot ) THEN
299 IF( imat.EQ.3 ) THEN
300 izero = 1
301 ELSE IF( imat.EQ.4 ) THEN
302 izero = n
303 ELSE
304 izero = n / 2 + 1
305 END IF
306*
307 IF( imat.LT.6 ) THEN
308*
309* Set row and column IZERO to zero.
310*
311 IF( iuplo.EQ.1 ) THEN
312 ioff = ( izero-1 )*lda
313 DO 20 i = 1, izero - 1
314 a( ioff+i ) = zero
315 20 CONTINUE
316 ioff = ioff + izero
317 DO 30 i = izero, n
318 a( ioff ) = zero
319 ioff = ioff + lda
320 30 CONTINUE
321 ELSE
322 ioff = izero
323 DO 40 i = 1, izero - 1
324 a( ioff ) = zero
325 ioff = ioff + lda
326 40 CONTINUE
327 ioff = ioff - izero
328 DO 50 i = izero, n
329 a( ioff+i ) = zero
330 50 CONTINUE
331 END IF
332 ELSE
333 ioff = 0
334 IF( iuplo.EQ.1 ) THEN
335*
336* Set the first IZERO rows and columns to zero.
337*
338 DO 70 j = 1, n
339 i2 = min( j, izero )
340 DO 60 i = 1, i2
341 a( ioff+i ) = zero
342 60 CONTINUE
343 ioff = ioff + lda
344 70 CONTINUE
345 ELSE
346*
347* Set the last IZERO rows and columns to zero.
348*
349 DO 90 j = 1, n
350 i1 = max( j, izero )
351 DO 80 i = i1, n
352 a( ioff+i ) = zero
353 80 CONTINUE
354 ioff = ioff + lda
355 90 CONTINUE
356 END IF
357 END IF
358 ELSE
359 izero = 0
360 END IF
361*
362* Set the imaginary part of the diagonals.
363*
364 CALL claipd( n, a, lda+1, 0 )
365*
366 DO 150 ifact = 1, nfact
367*
368* Do first for FACT = 'F', then for other values.
369*
370 fact = facts( ifact )
371*
372* Compute the condition number for comparison with
373* the value returned by CHESVX.
374*
375 IF( zerot ) THEN
376 IF( ifact.EQ.1 )
377 $ GO TO 150
378 rcondc = zero
379*
380 ELSE IF( ifact.EQ.1 ) THEN
381*
382* Compute the 1-norm of A.
383*
384 anorm = clanhe( '1', uplo, n, a, lda, rwork )
385*
386* Factor the matrix A.
387*
388 CALL clacpy( uplo, n, n, a, lda, afac, lda )
389 CALL chetrf( uplo, n, afac, lda, iwork, work,
390 $ lwork, info )
391*
392* Compute inv(A) and take its norm.
393*
394 CALL clacpy( uplo, n, n, afac, lda, ainv, lda )
395 lwork = (n+nb+1)*(nb+3)
396 CALL chetri2( uplo, n, ainv, lda, iwork, work,
397 $ lwork, info )
398 ainvnm = clanhe( '1', uplo, n, ainv, lda, rwork )
399*
400* Compute the 1-norm condition number of A.
401*
402 IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
403 rcondc = one
404 ELSE
405 rcondc = ( one / anorm ) / ainvnm
406 END IF
407 END IF
408*
409* Form an exact solution and set the right hand side.
410*
411 srnamt = 'CLARHS'
412 CALL clarhs( path, xtype, uplo, ' ', n, n, kl, ku,
413 $ nrhs, a, lda, xact, lda, b, lda, iseed,
414 $ info )
415 xtype = 'C'
416*
417* --- Test CHESV ---
418*
419 IF( ifact.EQ.2 ) THEN
420 CALL clacpy( uplo, n, n, a, lda, afac, lda )
421 CALL clacpy( 'Full', n, nrhs, b, lda, x, lda )
422*
423* Factor the matrix and solve the system using CHESV.
424*
425 srnamt = 'CHESV '
426 CALL chesv( uplo, n, nrhs, afac, lda, iwork, x,
427 $ lda, work, lwork, info )
428*
429* Adjust the expected value of INFO to account for
430* pivoting.
431*
432 k = izero
433 IF( k.GT.0 ) THEN
434 100 CONTINUE
435 IF( iwork( k ).LT.0 ) THEN
436 IF( iwork( k ).NE.-k ) THEN
437 k = -iwork( k )
438 GO TO 100
439 END IF
440 ELSE IF( iwork( k ).NE.k ) THEN
441 k = iwork( k )
442 GO TO 100
443 END IF
444 END IF
445*
446* Check error code from CHESV .
447*
448 IF( info.NE.k ) THEN
449 CALL alaerh( path, 'CHESV ', info, k, uplo, n,
450 $ n, -1, -1, nrhs, imat, nfail,
451 $ nerrs, nout )
452 GO TO 120
453 ELSE IF( info.NE.0 ) THEN
454 GO TO 120
455 END IF
456*
457* Reconstruct matrix from factors and compute
458* residual.
459*
460 CALL chet01( uplo, n, a, lda, afac, lda, iwork,
461 $ ainv, lda, rwork, result( 1 ) )
462*
463* Compute residual of the computed solution.
464*
465 CALL clacpy( 'Full', n, nrhs, b, lda, work, lda )
466 CALL cpot02( uplo, n, nrhs, a, lda, x, lda, work,
467 $ lda, rwork, result( 2 ) )
468*
469* Check solution from generated exact solution.
470*
471 CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
472 $ result( 3 ) )
473 nt = 3
474*
475* Print information about the tests that did not pass
476* the threshold.
477*
478 DO 110 k = 1, nt
479 IF( result( k ).GE.thresh ) THEN
480 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
481 $ CALL aladhd( nout, path )
482 WRITE( nout, fmt = 9999 )'CHESV ', uplo, n,
483 $ imat, k, result( k )
484 nfail = nfail + 1
485 END IF
486 110 CONTINUE
487 nrun = nrun + nt
488 120 CONTINUE
489 END IF
490*
491* --- Test CHESVX ---
492*
493 IF( ifact.EQ.2 )
494 $ CALL claset( uplo, n, n, cmplx( zero ),
495 $ cmplx( zero ), afac, lda )
496 CALL claset( 'Full', n, nrhs, cmplx( zero ),
497 $ cmplx( zero ), x, lda )
498*
499* Solve the system and compute the condition number and
500* error bounds using CHESVX.
501*
502 srnamt = 'CHESVX'
503 CALL chesvx( fact, uplo, n, nrhs, a, lda, afac, lda,
504 $ iwork, b, lda, x, lda, rcond, rwork,
505 $ rwork( nrhs+1 ), work, lwork,
506 $ rwork( 2*nrhs+1 ), info )
507*
508* Adjust the expected value of INFO to account for
509* pivoting.
510*
511 k = izero
512 IF( k.GT.0 ) THEN
513 130 CONTINUE
514 IF( iwork( k ).LT.0 ) THEN
515 IF( iwork( k ).NE.-k ) THEN
516 k = -iwork( k )
517 GO TO 130
518 END IF
519 ELSE IF( iwork( k ).NE.k ) THEN
520 k = iwork( k )
521 GO TO 130
522 END IF
523 END IF
524*
525* Check the error code from CHESVX.
526*
527 IF( info.NE.k ) THEN
528 CALL alaerh( path, 'CHESVX', info, k, fact // uplo,
529 $ n, n, -1, -1, nrhs, imat, nfail,
530 $ nerrs, nout )
531 GO TO 150
532 END IF
533*
534 IF( info.EQ.0 ) THEN
535 IF( ifact.GE.2 ) THEN
536*
537* Reconstruct matrix from factors and compute
538* residual.
539*
540 CALL chet01( uplo, n, a, lda, afac, lda, iwork,
541 $ ainv, lda, rwork( 2*nrhs+1 ),
542 $ result( 1 ) )
543 k1 = 1
544 ELSE
545 k1 = 2
546 END IF
547*
548* Compute residual of the computed solution.
549*
550 CALL clacpy( 'Full', n, nrhs, b, lda, work, lda )
551 CALL cpot02( uplo, n, nrhs, a, lda, x, lda, work,
552 $ lda, rwork( 2*nrhs+1 ), result( 2 ) )
553*
554* Check solution from generated exact solution.
555*
556 CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
557 $ result( 3 ) )
558*
559* Check the error bounds from iterative refinement.
560*
561 CALL cpot05( uplo, n, nrhs, a, lda, b, lda, x, lda,
562 $ xact, lda, rwork, rwork( nrhs+1 ),
563 $ result( 4 ) )
564 ELSE
565 k1 = 6
566 END IF
567*
568* Compare RCOND from CHESVX with the computed value
569* in RCONDC.
570*
571 result( 6 ) = sget06( rcond, rcondc )
572*
573* Print information about the tests that did not pass
574* the threshold.
575*
576 DO 140 k = k1, 6
577 IF( result( k ).GE.thresh ) THEN
578 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
579 $ CALL aladhd( nout, path )
580 WRITE( nout, fmt = 9998 )'CHESVX', fact, uplo,
581 $ n, imat, k, result( k )
582 nfail = nfail + 1
583 END IF
584 140 CONTINUE
585 nrun = nrun + 7 - k1
586*
587 150 CONTINUE
588*
589 160 CONTINUE
590 170 CONTINUE
591 180 CONTINUE
592*
593* Print a summary of the results.
594*
595 CALL alasvm( path, nout, nfail, nrun, nerrs )
596*
597 9999 FORMAT( 1x, a, ', UPLO=''', a1, ''', N =', i5, ', type ', i2,
598 $ ', test ', i2, ', ratio =', g12.5 )
599 9998 FORMAT( 1x, a, ', FACT=''', a1, ''', UPLO=''', a1, ''', N =', i5,
600 $ ', type ', i2, ', test ', i2, ', ratio =', g12.5 )
601 RETURN
602*
603* End of CDRVHE
604*
605 END
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
Definition alasvm.f:73
subroutine clarhs(path, xtype, uplo, trans, m, n, kl, ku, nrhs, a, lda, x, ldx, b, ldb, iseed, info)
CLARHS
Definition clarhs.f:208
subroutine xlaenv(ispec, nvalue)
XLAENV
Definition xlaenv.f:81
subroutine aladhd(iounit, path)
ALADHD
Definition aladhd.f:90
subroutine alaerh(path, subnam, info, infoe, opts, m, n, kl, ku, n5, imat, nfail, nerrs, nout)
ALAERH
Definition alaerh.f:147
subroutine cdrvhe(dotype, nn, nval, nrhs, thresh, tsterr, nmax, a, afac, ainv, b, x, xact, work, rwork, iwork, nout)
CDRVHE
Definition cdrvhe.f:153
subroutine cerrvx(path, nunit)
CERRVX
Definition cerrvx.f:55
subroutine cget04(n, nrhs, x, ldx, xact, ldxact, rcond, resid)
CGET04
Definition cget04.f:102
subroutine chet01(uplo, n, a, lda, afac, ldafac, ipiv, c, ldc, rwork, resid)
CHET01
Definition chet01.f:126
subroutine claipd(n, a, inda, vinda)
CLAIPD
Definition claipd.f:83
subroutine clatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
CLATB4
Definition clatb4.f:121
subroutine clatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
CLATMS
Definition clatms.f:332
subroutine cpot02(uplo, n, nrhs, a, lda, x, ldx, b, ldb, rwork, resid)
CPOT02
Definition cpot02.f:127
subroutine cpot05(uplo, n, nrhs, a, lda, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
CPOT05
Definition cpot05.f:165
subroutine chesv(uplo, n, nrhs, a, lda, ipiv, b, ldb, work, lwork, info)
CHESV computes the solution to system of linear equations A * X = B for HE matrices
Definition chesv.f:171
subroutine chesvx(fact, uplo, n, nrhs, a, lda, af, ldaf, ipiv, b, ldb, x, ldx, rcond, ferr, berr, work, lwork, rwork, info)
CHESVX computes the solution to system of linear equations A * X = B for HE matrices
Definition chesvx.f:285
subroutine chetrf(uplo, n, a, lda, ipiv, work, lwork, info)
CHETRF
Definition chetrf.f:177
subroutine chetri2(uplo, n, a, lda, ipiv, work, lwork, info)
CHETRI2
Definition chetri2.f:127
subroutine clacpy(uplo, m, n, a, lda, b, ldb)
CLACPY copies all or part of one two-dimensional array to another.
Definition clacpy.f:103
subroutine claset(uplo, m, n, alpha, beta, a, lda)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition claset.f:106