LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ zdrvrfp()

subroutine zdrvrfp ( integer  nout,
integer  nn,
integer, dimension( nn )  nval,
integer  nns,
integer, dimension( nns )  nsval,
integer  nnt,
integer, dimension( nnt )  ntval,
double precision  thresh,
complex*16, dimension( * )  a,
complex*16, dimension( * )  asav,
complex*16, dimension( * )  afac,
complex*16, dimension( * )  ainv,
complex*16, dimension( * )  b,
complex*16, dimension( * )  bsav,
complex*16, dimension( * )  xact,
complex*16, dimension( * )  x,
complex*16, dimension( * )  arf,
complex*16, dimension( * )  arfinv,
complex*16, dimension( * )  z_work_zlatms,
complex*16, dimension( * )  z_work_zpot02,
complex*16, dimension( * )  z_work_zpot03,
double precision, dimension( * )  d_work_zlatms,
double precision, dimension( * )  d_work_zlanhe,
double precision, dimension( * )  d_work_zpot01,
double precision, dimension( * )  d_work_zpot02,
double precision, dimension( * )  d_work_zpot03 
)

ZDRVRFP

Purpose:
 ZDRVRFP tests the LAPACK RFP routines:
     ZPFTRF, ZPFTRS, and ZPFTRI.

 This testing routine follow the same tests as ZDRVPO (test for the full
 format Symmetric Positive Definite solver).

 The tests are performed in Full Format, conversion back and forth from
 full format to RFP format are performed using the routines ZTRTTF and
 ZTFTTR.

 First, a specific matrix A of size N is created. There is nine types of
 different matrixes possible.
  1. Diagonal                        6. Random, CNDNUM = sqrt(0.1/EPS)
  2. Random, CNDNUM = 2              7. Random, CNDNUM = 0.1/EPS
 *3. First row and column zero       8. Scaled near underflow
 *4. Last row and column zero        9. Scaled near overflow
 *5. Middle row and column zero
 (* - tests error exits from ZPFTRF, no test ratios are computed)
 A solution XACT of size N-by-NRHS is created and the associated right
 hand side B as well. Then ZPFTRF is called to compute L (or U), the
 Cholesky factor of A. Then L (or U) is used to solve the linear system
 of equations AX = B. This gives X. Then L (or U) is used to compute the
 inverse of A, AINV. The following four tests are then performed:
 (1) norm( L*L' - A ) / ( N * norm(A) * EPS ) or
     norm( U'*U - A ) / ( N * norm(A) * EPS ),
 (2) norm(B - A*X) / ( norm(A) * norm(X) * EPS ),
 (3) norm( I - A*AINV ) / ( N * norm(A) * norm(AINV) * EPS ),
 (4) ( norm(X-XACT) * RCOND ) / ( norm(XACT) * EPS ),
 where EPS is the machine precision, RCOND the condition number of A, and
 norm( . ) the 1-norm for (1,2,3) and the inf-norm for (4).
 Errors occur when INFO parameter is not as expected. Failures occur when
 a test ratios is greater than THRES.
Parameters
[in]NOUT
          NOUT is INTEGER
                The unit number for output.
[in]NN
          NN is INTEGER
                The number of values of N contained in the vector NVAL.
[in]NVAL
          NVAL is INTEGER array, dimension (NN)
                The values of the matrix dimension N.
[in]NNS
          NNS is INTEGER
                The number of values of NRHS contained in the vector NSVAL.
[in]NSVAL
          NSVAL is INTEGER array, dimension (NNS)
                The values of the number of right-hand sides NRHS.
[in]NNT
          NNT is INTEGER
                The number of values of MATRIX TYPE contained in the vector NTVAL.
[in]NTVAL
          NTVAL is INTEGER array, dimension (NNT)
                The values of matrix type (between 0 and 9 for PO/PP/PF matrices).
[in]THRESH
          THRESH is DOUBLE PRECISION
                The threshold value for the test ratios.  A result is
                included in the output file if RESULT >= THRESH.  To have
                every test ratio printed, use THRESH = 0.
[out]A
          A is COMPLEX*16 array, dimension (NMAX*NMAX)
[out]ASAV
          ASAV is COMPLEX*16 array, dimension (NMAX*NMAX)
[out]AFAC
          AFAC is COMPLEX*16 array, dimension (NMAX*NMAX)
[out]AINV
          AINV is COMPLEX*16 array, dimension (NMAX*NMAX)
[out]B
          B is COMPLEX*16 array, dimension (NMAX*MAXRHS)
[out]BSAV
          BSAV is COMPLEX*16 array, dimension (NMAX*MAXRHS)
[out]XACT
          XACT is COMPLEX*16 array, dimension (NMAX*MAXRHS)
[out]X
          X is COMPLEX*16 array, dimension (NMAX*MAXRHS)
[out]ARF
          ARF is COMPLEX*16 array, dimension ((NMAX*(NMAX+1))/2)
[out]ARFINV
          ARFINV is COMPLEX*16 array, dimension ((NMAX*(NMAX+1))/2)
[out]Z_WORK_ZLATMS
          Z_WORK_ZLATMS is COMPLEX*16 array, dimension ( 3*NMAX )
[out]Z_WORK_ZPOT02
          Z_WORK_ZPOT02 is COMPLEX*16 array, dimension ( NMAX*MAXRHS )
[out]Z_WORK_ZPOT03
          Z_WORK_ZPOT03 is COMPLEX*16 array, dimension ( NMAX*NMAX )
[out]D_WORK_ZLATMS
          D_WORK_ZLATMS is DOUBLE PRECISION array, dimension ( NMAX )
[out]D_WORK_ZLANHE
          D_WORK_ZLANHE is DOUBLE PRECISION array, dimension ( NMAX )
[out]D_WORK_ZPOT01
          D_WORK_ZPOT01 is DOUBLE PRECISION array, dimension ( NMAX )
[out]D_WORK_ZPOT02
          D_WORK_ZPOT02 is DOUBLE PRECISION array, dimension ( NMAX )
[out]D_WORK_ZPOT03
          D_WORK_ZPOT03 is DOUBLE PRECISION array, dimension ( NMAX )
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 238 of file zdrvrfp.f.

244*
245* -- LAPACK test routine --
246* -- LAPACK is a software package provided by Univ. of Tennessee, --
247* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
248*
249* .. Scalar Arguments ..
250 INTEGER NN, NNS, NNT, NOUT
251 DOUBLE PRECISION THRESH
252* ..
253* .. Array Arguments ..
254 INTEGER NVAL( NN ), NSVAL( NNS ), NTVAL( NNT )
255 COMPLEX*16 A( * )
256 COMPLEX*16 AINV( * )
257 COMPLEX*16 ASAV( * )
258 COMPLEX*16 B( * )
259 COMPLEX*16 BSAV( * )
260 COMPLEX*16 AFAC( * )
261 COMPLEX*16 ARF( * )
262 COMPLEX*16 ARFINV( * )
263 COMPLEX*16 XACT( * )
264 COMPLEX*16 X( * )
265 COMPLEX*16 Z_WORK_ZLATMS( * )
266 COMPLEX*16 Z_WORK_ZPOT02( * )
267 COMPLEX*16 Z_WORK_ZPOT03( * )
268 DOUBLE PRECISION D_WORK_ZLATMS( * )
269 DOUBLE PRECISION D_WORK_ZLANHE( * )
270 DOUBLE PRECISION D_WORK_ZPOT01( * )
271 DOUBLE PRECISION D_WORK_ZPOT02( * )
272 DOUBLE PRECISION D_WORK_ZPOT03( * )
273* ..
274*
275* =====================================================================
276*
277* .. Parameters ..
278 DOUBLE PRECISION ONE, ZERO
279 parameter( one = 1.0d+0, zero = 0.0d+0 )
280 INTEGER NTESTS
281 parameter( ntests = 4 )
282* ..
283* .. Local Scalars ..
284 LOGICAL ZEROT
285 INTEGER I, INFO, IUPLO, LDA, LDB, IMAT, NERRS, NFAIL,
286 + NRHS, NRUN, IZERO, IOFF, K, NT, N, IFORM, IIN,
287 + IIT, IIS
288 CHARACTER DIST, CTYPE, UPLO, CFORM
289 INTEGER KL, KU, MODE
290 DOUBLE PRECISION ANORM, AINVNM, CNDNUM, RCONDC
291* ..
292* .. Local Arrays ..
293 CHARACTER UPLOS( 2 ), FORMS( 2 )
294 INTEGER ISEED( 4 ), ISEEDY( 4 )
295 DOUBLE PRECISION RESULT( NTESTS )
296* ..
297* .. External Functions ..
298 DOUBLE PRECISION ZLANHE
299 EXTERNAL zlanhe
300* ..
301* .. External Subroutines ..
302 EXTERNAL aladhd, alaerh, alasvm, zget04, ztfttr, zlacpy,
303 + zlaipd, zlarhs, zlatb4, zlatms, zpftri, zpftrf,
304 + zpftrs, zpot01, zpot02, zpot03, zpotri, zpotrf,
305 + ztrttf
306* ..
307* .. Scalars in Common ..
308 CHARACTER*32 SRNAMT
309* ..
310* .. Common blocks ..
311 COMMON / srnamc / srnamt
312* ..
313* .. Data statements ..
314 DATA iseedy / 1988, 1989, 1990, 1991 /
315 DATA uplos / 'U', 'L' /
316 DATA forms / 'N', 'C' /
317* ..
318* .. Executable Statements ..
319*
320* Initialize constants and the random number seed.
321*
322 nrun = 0
323 nfail = 0
324 nerrs = 0
325 DO 10 i = 1, 4
326 iseed( i ) = iseedy( i )
327 10 CONTINUE
328*
329 DO 130 iin = 1, nn
330*
331 n = nval( iin )
332 lda = max( n, 1 )
333 ldb = max( n, 1 )
334*
335 DO 980 iis = 1, nns
336*
337 nrhs = nsval( iis )
338*
339 DO 120 iit = 1, nnt
340*
341 imat = ntval( iit )
342*
343* If N.EQ.0, only consider the first type
344*
345 IF( n.EQ.0 .AND. iit.GE.1 ) GO TO 120
346*
347* Skip types 3, 4, or 5 if the matrix size is too small.
348*
349 IF( imat.EQ.4 .AND. n.LE.1 ) GO TO 120
350 IF( imat.EQ.5 .AND. n.LE.2 ) GO TO 120
351*
352* Do first for UPLO = 'U', then for UPLO = 'L'
353*
354 DO 110 iuplo = 1, 2
355 uplo = uplos( iuplo )
356*
357* Do first for CFORM = 'N', then for CFORM = 'C'
358*
359 DO 100 iform = 1, 2
360 cform = forms( iform )
361*
362* Set up parameters with ZLATB4 and generate a test
363* matrix with ZLATMS.
364*
365 CALL zlatb4( 'ZPO', imat, n, n, ctype, kl, ku,
366 + anorm, mode, cndnum, dist )
367*
368 srnamt = 'ZLATMS'
369 CALL zlatms( n, n, dist, iseed, ctype,
370 + d_work_zlatms,
371 + mode, cndnum, anorm, kl, ku, uplo, a,
372 + lda, z_work_zlatms, info )
373*
374* Check error code from ZLATMS.
375*
376 IF( info.NE.0 ) THEN
377 CALL alaerh( 'ZPF', 'ZLATMS', info, 0, uplo, n,
378 + n, -1, -1, -1, iit, nfail, nerrs,
379 + nout )
380 GO TO 100
381 END IF
382*
383* For types 3-5, zero one row and column of the matrix to
384* test that INFO is returned correctly.
385*
386 zerot = imat.GE.3 .AND. imat.LE.5
387 IF( zerot ) THEN
388 IF( iit.EQ.3 ) THEN
389 izero = 1
390 ELSE IF( iit.EQ.4 ) THEN
391 izero = n
392 ELSE
393 izero = n / 2 + 1
394 END IF
395 ioff = ( izero-1 )*lda
396*
397* Set row and column IZERO of A to 0.
398*
399 IF( iuplo.EQ.1 ) THEN
400 DO 20 i = 1, izero - 1
401 a( ioff+i ) = zero
402 20 CONTINUE
403 ioff = ioff + izero
404 DO 30 i = izero, n
405 a( ioff ) = zero
406 ioff = ioff + lda
407 30 CONTINUE
408 ELSE
409 ioff = izero
410 DO 40 i = 1, izero - 1
411 a( ioff ) = zero
412 ioff = ioff + lda
413 40 CONTINUE
414 ioff = ioff - izero
415 DO 50 i = izero, n
416 a( ioff+i ) = zero
417 50 CONTINUE
418 END IF
419 ELSE
420 izero = 0
421 END IF
422*
423* Set the imaginary part of the diagonals.
424*
425 CALL zlaipd( n, a, lda+1, 0 )
426*
427* Save a copy of the matrix A in ASAV.
428*
429 CALL zlacpy( uplo, n, n, a, lda, asav, lda )
430*
431* Compute the condition number of A (RCONDC).
432*
433 IF( zerot ) THEN
434 rcondc = zero
435 ELSE
436*
437* Compute the 1-norm of A.
438*
439 anorm = zlanhe( '1', uplo, n, a, lda,
440 + d_work_zlanhe )
441*
442* Factor the matrix A.
443*
444 CALL zpotrf( uplo, n, a, lda, info )
445*
446* Form the inverse of A.
447*
448 CALL zpotri( uplo, n, a, lda, info )
449
450 IF ( n .NE. 0 ) THEN
451*
452* Compute the 1-norm condition number of A.
453*
454 ainvnm = zlanhe( '1', uplo, n, a, lda,
455 + d_work_zlanhe )
456 rcondc = ( one / anorm ) / ainvnm
457*
458* Restore the matrix A.
459*
460 CALL zlacpy( uplo, n, n, asav, lda, a, lda )
461 END IF
462*
463 END IF
464*
465* Form an exact solution and set the right hand side.
466*
467 srnamt = 'ZLARHS'
468 CALL zlarhs( 'ZPO', 'N', uplo, ' ', n, n, kl, ku,
469 + nrhs, a, lda, xact, lda, b, lda,
470 + iseed, info )
471 CALL zlacpy( 'Full', n, nrhs, b, lda, bsav, lda )
472*
473* Compute the L*L' or U'*U factorization of the
474* matrix and solve the system.
475*
476 CALL zlacpy( uplo, n, n, a, lda, afac, lda )
477 CALL zlacpy( 'Full', n, nrhs, b, ldb, x, ldb )
478*
479 srnamt = 'ZTRTTF'
480 CALL ztrttf( cform, uplo, n, afac, lda, arf, info )
481 srnamt = 'ZPFTRF'
482 CALL zpftrf( cform, uplo, n, arf, info )
483*
484* Check error code from ZPFTRF.
485*
486 IF( info.NE.izero ) THEN
487*
488* LANGOU: there is a small hick here: IZERO should
489* always be INFO however if INFO is ZERO, ALAERH does not
490* complain.
491*
492 CALL alaerh( 'ZPF', 'ZPFSV ', info, izero,
493 + uplo, n, n, -1, -1, nrhs, iit,
494 + nfail, nerrs, nout )
495 GO TO 100
496 END IF
497*
498* Skip the tests if INFO is not 0.
499*
500 IF( info.NE.0 ) THEN
501 GO TO 100
502 END IF
503*
504 srnamt = 'ZPFTRS'
505 CALL zpftrs( cform, uplo, n, nrhs, arf, x, ldb,
506 + info )
507*
508 srnamt = 'ZTFTTR'
509 CALL ztfttr( cform, uplo, n, arf, afac, lda, info )
510*
511* Reconstruct matrix from factors and compute
512* residual.
513*
514 CALL zlacpy( uplo, n, n, afac, lda, asav, lda )
515 CALL zpot01( uplo, n, a, lda, afac, lda,
516 + d_work_zpot01, result( 1 ) )
517 CALL zlacpy( uplo, n, n, asav, lda, afac, lda )
518*
519* Form the inverse and compute the residual.
520*
521 IF(mod(n,2).EQ.0)THEN
522 CALL zlacpy( 'A', n+1, n/2, arf, n+1, arfinv,
523 + n+1 )
524 ELSE
525 CALL zlacpy( 'A', n, (n+1)/2, arf, n, arfinv,
526 + n )
527 END IF
528*
529 srnamt = 'ZPFTRI'
530 CALL zpftri( cform, uplo, n, arfinv , info )
531*
532 srnamt = 'ZTFTTR'
533 CALL ztfttr( cform, uplo, n, arfinv, ainv, lda,
534 + info )
535*
536* Check error code from ZPFTRI.
537*
538 IF( info.NE.0 )
539 + CALL alaerh( 'ZPO', 'ZPFTRI', info, 0, uplo, n,
540 + n, -1, -1, -1, imat, nfail, nerrs,
541 + nout )
542*
543 CALL zpot03( uplo, n, a, lda, ainv, lda,
544 + z_work_zpot03, lda, d_work_zpot03,
545 + rcondc, result( 2 ) )
546*
547* Compute residual of the computed solution.
548*
549 CALL zlacpy( 'Full', n, nrhs, b, lda,
550 + z_work_zpot02, lda )
551 CALL zpot02( uplo, n, nrhs, a, lda, x, lda,
552 + z_work_zpot02, lda, d_work_zpot02,
553 + result( 3 ) )
554*
555* Check solution from generated exact solution.
556*
557 CALL zget04( n, nrhs, x, lda, xact, lda, rcondc,
558 + result( 4 ) )
559 nt = 4
560*
561* Print information about the tests that did not
562* pass the threshold.
563*
564 DO 60 k = 1, nt
565 IF( result( k ).GE.thresh ) THEN
566 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
567 + CALL aladhd( nout, 'ZPF' )
568 WRITE( nout, fmt = 9999 )'ZPFSV ', uplo,
569 + n, iit, k, result( k )
570 nfail = nfail + 1
571 END IF
572 60 CONTINUE
573 nrun = nrun + nt
574 100 CONTINUE
575 110 CONTINUE
576 120 CONTINUE
577 980 CONTINUE
578 130 CONTINUE
579*
580* Print a summary of the results.
581*
582 CALL alasvm( 'ZPF', nout, nfail, nrun, nerrs )
583*
584 9999 FORMAT( 1x, a6, ', UPLO=''', a1, ''', N =', i5, ', type ', i1,
585 + ', test(', i1, ')=', g12.5 )
586*
587 RETURN
588*
589* End of ZDRVRFP
590*
alasvm
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
Definition alasvm.f:73
zlarhs
subroutine zlarhs(path, xtype, uplo, trans, m, n, kl, ku, nrhs, a, lda, x, ldx, b, ldb, iseed, info)
ZLARHS
Definition zlarhs.f:208
aladhd
subroutine aladhd(iounit, path)
ALADHD
Definition aladhd.f:90
alaerh
subroutine alaerh(path, subnam, info, infoe, opts, m, n, kl, ku, n5, imat, nfail, nerrs, nout)
ALAERH
Definition alaerh.f:147
zlacpy
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:103
zlanhe
double precision function zlanhe(norm, uplo, n, a, lda, work)
ZLANHE returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition zlanhe.f:124
zpftrf
subroutine zpftrf(transr, uplo, n, a, info)
ZPFTRF
Definition zpftrf.f:211
zpftri
subroutine zpftri(transr, uplo, n, a, info)
ZPFTRI
Definition zpftri.f:212
zpftrs
subroutine zpftrs(transr, uplo, n, nrhs, a, b, ldb, info)
ZPFTRS
Definition zpftrs.f:220
zpotrf
subroutine zpotrf(uplo, n, a, lda, info)
ZPOTRF
Definition zpotrf.f:107
zpotri
subroutine zpotri(uplo, n, a, lda, info)
ZPOTRI
Definition zpotri.f:95
ztfttr
subroutine ztfttr(transr, uplo, n, arf, a, lda, info)
ZTFTTR copies a triangular matrix from the rectangular full packed format (TF) to the standard full f...
Definition ztfttr.f:216
ztrttf
subroutine ztrttf(transr, uplo, n, a, lda, arf, info)
ZTRTTF copies a triangular matrix from the standard full format (TR) to the rectangular full packed f...
Definition ztrttf.f:216
zget04
subroutine zget04(n, nrhs, x, ldx, xact, ldxact, rcond, resid)
ZGET04
Definition zget04.f:102
zlaipd
subroutine zlaipd(n, a, inda, vinda)
ZLAIPD
Definition zlaipd.f:83
zlatb4
subroutine zlatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
ZLATB4
Definition zlatb4.f:121
zlatms
subroutine zlatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
ZLATMS
Definition zlatms.f:332
zpot01
subroutine zpot01(uplo, n, a, lda, afac, ldafac, rwork, resid)
ZPOT01
Definition zpot01.f:106
zpot02
subroutine zpot02(uplo, n, nrhs, a, lda, x, ldx, b, ldb, rwork, resid)
ZPOT02
Definition zpot02.f:127
zpot03
subroutine zpot03(uplo, n, a, lda, ainv, ldainv, work, ldwork, rwork, rcond, resid)
ZPOT03
Definition zpot03.f:126
Here is the call graph for this function:
Here is the caller graph for this function: