LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine cdrvrfp ( integer  NOUT,
integer  NN,
integer, dimension( nn )  NVAL,
integer  NNS,
integer, dimension( nns )  NSVAL,
integer  NNT,
integer, dimension( nnt )  NTVAL,
real  THRESH,
complex, dimension( * )  A,
complex, dimension( * )  ASAV,
complex, dimension( * )  AFAC,
complex, dimension( * )  AINV,
complex, dimension( * )  B,
complex, dimension( * )  BSAV,
complex, dimension( * )  XACT,
complex, dimension( * )  X,
complex, dimension( * )  ARF,
complex, dimension( * )  ARFINV,
complex, dimension( * )  C_WORK_CLATMS,
complex, dimension( * )  C_WORK_CPOT02,
complex, dimension( * )  C_WORK_CPOT03,
real, dimension( * )  S_WORK_CLATMS,
real, dimension( * )  S_WORK_CLANHE,
real, dimension( * )  S_WORK_CPOT01,
real, dimension( * )  S_WORK_CPOT02,
real, dimension( * )  S_WORK_CPOT03 
)

CDRVRFP

Purpose: 
 CDRVRFP tests the LAPACK RFP routines:
     CPFTRF, CPFTRS, and CPFTRI.

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

 The tests are performed in Full Format, convertion back and forth from
 full format to RFP format are performed using the routines CTRTTF and
 CTFTTR.

 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 CPFTRF, 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 CPFTRF 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 REAL
                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 array, dimension (NMAX*NMAX)
[out]ASAV
          ASAV is COMPLEX array, dimension (NMAX*NMAX)
[out]AFAC
          AFAC is COMPLEX array, dimension (NMAX*NMAX)
[out]AINV
          AINV is COMPLEX array, dimension (NMAX*NMAX)
[out]B
          B is COMPLEX array, dimension (NMAX*MAXRHS)
[out]BSAV
          BSAV is COMPLEX array, dimension (NMAX*MAXRHS)
[out]XACT
          XACT is COMPLEX array, dimension (NMAX*MAXRHS)
[out]X
          X is COMPLEX array, dimension (NMAX*MAXRHS)
[out]ARF
          ARF is COMPLEX array, dimension ((NMAX*(NMAX+1))/2)
[out]ARFINV
          ARFINV is COMPLEX array, dimension ((NMAX*(NMAX+1))/2)
[out]C_WORK_CLATMS
          C_WORK_CLATMS is COMPLEX array, dimension ( 3*NMAX )
[out]C_WORK_CPOT02
          C_WORK_CPOT02 is COMPLEX array, dimension ( NMAX*MAXRHS )
[out]C_WORK_CPOT03
          C_WORK_CPOT03 is COMPLEX array, dimension ( NMAX*NMAX )
[out]S_WORK_CLATMS
          S_WORK_CLATMS is REAL array, dimension ( NMAX )
[out]S_WORK_CLANHE
          S_WORK_CLANHE is REAL array, dimension ( NMAX )
[out]S_WORK_CPOT01
          S_WORK_CPOT01 is REAL array, dimension ( NMAX )
[out]S_WORK_CPOT02
          S_WORK_CPOT02 is REAL array, dimension ( NMAX )
[out]S_WORK_CPOT03
          S_WORK_CPOT03 is REAL array, dimension ( NMAX )
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2013

Definition at line 246 of file cdrvrfp.f.

246 *
247 * -- LAPACK test routine (version 3.5.0) --
248 * -- LAPACK is a software package provided by Univ. of Tennessee, --
249 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
250 * November 2013
251 *
252 * .. Scalar Arguments ..
253  INTEGER nn, nns, nnt, nout
254  REAL thresh
255 * ..
256 * .. Array Arguments ..
257  INTEGER nval( nn ), nsval( nns ), ntval( nnt )
258  COMPLEX a( * )
259  COMPLEX ainv( * )
260  COMPLEX asav( * )
261  COMPLEX b( * )
262  COMPLEX bsav( * )
263  COMPLEX afac( * )
264  COMPLEX arf( * )
265  COMPLEX arfinv( * )
266  COMPLEX xact( * )
267  COMPLEX x( * )
268  COMPLEX c_work_clatms( * )
269  COMPLEX c_work_cpot02( * )
270  COMPLEX c_work_cpot03( * )
271  REAL s_work_clatms( * )
272  REAL s_work_clanhe( * )
273  REAL s_work_cpot01( * )
274  REAL s_work_cpot02( * )
275  REAL s_work_cpot03( * )
276 * ..
277 *
278 * =====================================================================
279 *
280 * .. Parameters ..
281  REAL one, zero
282  parameter ( one = 1.0e+0, zero = 0.0e+0 )
283  INTEGER ntests
284  parameter ( ntests = 4 )
285 * ..
286 * .. Local Scalars ..
287  LOGICAL zerot
288  INTEGER i, info, iuplo, lda, ldb, imat, nerrs, nfail,
289  + nrhs, nrun, izero, ioff, k, nt, n, iform, iin,
290  + iit, iis
291  CHARACTER dist, ctype, uplo, cform
292  INTEGER kl, ku, mode
293  REAL anorm, ainvnm, cndnum, rcondc
294 * ..
295 * .. Local Arrays ..
296  CHARACTER uplos( 2 ), forms( 2 )
297  INTEGER iseed( 4 ), iseedy( 4 )
298  REAL result( ntests )
299 * ..
300 * .. External Functions ..
301  REAL clanhe
302  EXTERNAL clanhe
303 * ..
304 * .. External Subroutines ..
305  EXTERNAL aladhd, alaerh, alasvm, cget04, ctfttr, clacpy,
306  + claipd, clarhs, clatb4, clatms, cpftri, cpftrf,
307  + cpftrs, cpot01, cpot02, cpot03, cpotri, cpotrf,
308  + ctrttf
309 * ..
310 * .. Scalars in Common ..
311  CHARACTER*32 srnamt
312 * ..
313 * .. Common blocks ..
314  COMMON / srnamc / srnamt
315 * ..
316 * .. Data statements ..
317  DATA iseedy / 1988, 1989, 1990, 1991 /
318  DATA uplos / 'U', 'L' /
319  DATA forms / 'N', 'C' /
320 * ..
321 * .. Executable Statements ..
322 *
323 * Initialize constants and the random number seed.
324 *
325  nrun = 0
326  nfail = 0
327  nerrs = 0
328  DO 10 i = 1, 4
329  iseed( i ) = iseedy( i )
330  10 CONTINUE
331 *
332  DO 130 iin = 1, nn
333 *
334  n = nval( iin )
335  lda = max( n, 1 )
336  ldb = max( n, 1 )
337 *
338  DO 980 iis = 1, nns
339 *
340  nrhs = nsval( iis )
341 *
342  DO 120 iit = 1, nnt
343 *
344  imat = ntval( iit )
345 *
346 * If N.EQ.0, only consider the first type
347 *
348  IF( n.EQ.0 .AND. iit.GE.1 ) GO TO 120
349 *
350 * Skip types 3, 4, or 5 if the matrix size is too small.
351 *
352  IF( imat.EQ.4 .AND. n.LE.1 ) GO TO 120
353  IF( imat.EQ.5 .AND. n.LE.2 ) GO TO 120
354 *
355 * Do first for UPLO = 'U', then for UPLO = 'L'
356 *
357  DO 110 iuplo = 1, 2
358  uplo = uplos( iuplo )
359 *
360 * Do first for CFORM = 'N', then for CFORM = 'C'
361 *
362  DO 100 iform = 1, 2
363  cform = forms( iform )
364 *
365 * Set up parameters with CLATB4 and generate a test
366 * matrix with CLATMS.
367 *
368  CALL clatb4( 'CPO', imat, n, n, ctype, kl, ku,
369  + anorm, mode, cndnum, dist )
370 *
371  srnamt = 'CLATMS'
372  CALL clatms( n, n, dist, iseed, ctype,
373  + s_work_clatms,
374  + mode, cndnum, anorm, kl, ku, uplo, a,
375  + lda, c_work_clatms, info )
376 *
377 * Check error code from CLATMS.
378 *
379  IF( info.NE.0 ) THEN
380  CALL alaerh( 'CPF', 'CLATMS', info, 0, uplo, n,
381  + n, -1, -1, -1, iit, nfail, nerrs,
382  + nout )
383  GO TO 100
384  END IF
385 *
386 * For types 3-5, zero one row and column of the matrix to
387 * test that INFO is returned correctly.
388 *
389  zerot = imat.GE.3 .AND. imat.LE.5
390  IF( zerot ) THEN
391  IF( iit.EQ.3 ) THEN
392  izero = 1
393  ELSE IF( iit.EQ.4 ) THEN
394  izero = n
395  ELSE
396  izero = n / 2 + 1
397  END IF
398  ioff = ( izero-1 )*lda
399 *
400 * Set row and column IZERO of A to 0.
401 *
402  IF( iuplo.EQ.1 ) THEN
403  DO 20 i = 1, izero - 1
404  a( ioff+i ) = zero
405  20 CONTINUE
406  ioff = ioff + izero
407  DO 30 i = izero, n
408  a( ioff ) = zero
409  ioff = ioff + lda
410  30 CONTINUE
411  ELSE
412  ioff = izero
413  DO 40 i = 1, izero - 1
414  a( ioff ) = zero
415  ioff = ioff + lda
416  40 CONTINUE
417  ioff = ioff - izero
418  DO 50 i = izero, n
419  a( ioff+i ) = zero
420  50 CONTINUE
421  END IF
422  ELSE
423  izero = 0
424  END IF
425 *
426 * Set the imaginary part of the diagonals.
427 *
428  CALL claipd( n, a, lda+1, 0 )
429 *
430 * Save a copy of the matrix A in ASAV.
431 *
432  CALL clacpy( uplo, n, n, a, lda, asav, lda )
433 *
434 * Compute the condition number of A (RCONDC).
435 *
436  IF( zerot ) THEN
437  rcondc = zero
438  ELSE
439 *
440 * Compute the 1-norm of A.
441 *
442  anorm = clanhe( '1', uplo, n, a, lda,
443  + s_work_clanhe )
444 *
445 * Factor the matrix A.
446 *
447  CALL cpotrf( uplo, n, a, lda, info )
448 *
449 * Form the inverse of A.
450 *
451  CALL cpotri( uplo, n, a, lda, info )
452 *
453 * Compute the 1-norm condition number of A.
454 *
455  IF ( n .NE. 0 ) THEN
456  ainvnm = clanhe( '1', uplo, n, a, lda,
457  + s_work_clanhe )
458  rcondc = ( one / anorm ) / ainvnm
459 *
460 * Restore the matrix A.
461 *
462  CALL clacpy( uplo, n, n, asav, lda, a, lda )
463  END IF
464 
465 *
466  END IF
467 *
468 * Form an exact solution and set the right hand side.
469 *
470  srnamt = 'CLARHS'
471  CALL clarhs( 'CPO', 'N', uplo, ' ', n, n, kl, ku,
472  + nrhs, a, lda, xact, lda, b, lda,
473  + iseed, info )
474  CALL clacpy( 'Full', n, nrhs, b, lda, bsav, lda )
475 *
476 * Compute the L*L' or U'*U factorization of the
477 * matrix and solve the system.
478 *
479  CALL clacpy( uplo, n, n, a, lda, afac, lda )
480  CALL clacpy( 'Full', n, nrhs, b, ldb, x, ldb )
481 *
482  srnamt = 'CTRTTF'
483  CALL ctrttf( cform, uplo, n, afac, lda, arf, info )
484  srnamt = 'CPFTRF'
485  CALL cpftrf( cform, uplo, n, arf, info )
486 *
487 * Check error code from CPFTRF.
488 *
489  IF( info.NE.izero ) THEN
490 *
491 * LANGOU: there is a small hick here: IZERO should
492 * always be INFO however if INFO is ZERO, ALAERH does not
493 * complain.
494 *
495  CALL alaerh( 'CPF', 'CPFSV ', info, izero,
496  + uplo, n, n, -1, -1, nrhs, iit,
497  + nfail, nerrs, nout )
498  GO TO 100
499  END IF
500 *
501 * Skip the tests if INFO is not 0.
502 *
503  IF( info.NE.0 ) THEN
504  GO TO 100
505  END IF
506 *
507  srnamt = 'CPFTRS'
508  CALL cpftrs( cform, uplo, n, nrhs, arf, x, ldb,
509  + info )
510 *
511  srnamt = 'CTFTTR'
512  CALL ctfttr( cform, uplo, n, arf, afac, lda, info )
513 *
514 * Reconstruct matrix from factors and compute
515 * residual.
516 *
517  CALL clacpy( uplo, n, n, afac, lda, asav, lda )
518  CALL cpot01( uplo, n, a, lda, afac, lda,
519  + s_work_cpot01, result( 1 ) )
520  CALL clacpy( uplo, n, n, asav, lda, afac, lda )
521 *
522 * Form the inverse and compute the residual.
523 *
524  IF(mod(n,2).EQ.0)THEN
525  CALL clacpy( 'A', n+1, n/2, arf, n+1, arfinv,
526  + n+1 )
527  ELSE
528  CALL clacpy( 'A', n, (n+1)/2, arf, n, arfinv,
529  + n )
530  END IF
531 *
532  srnamt = 'CPFTRI'
533  CALL cpftri( cform, uplo, n, arfinv , info )
534 *
535  srnamt = 'CTFTTR'
536  CALL ctfttr( cform, uplo, n, arfinv, ainv, lda,
537  + info )
538 *
539 * Check error code from CPFTRI.
540 *
541  IF( info.NE.0 )
542  + CALL alaerh( 'CPO', 'CPFTRI', info, 0, uplo, n,
543  + n, -1, -1, -1, imat, nfail, nerrs,
544  + nout )
545 *
546  CALL cpot03( uplo, n, a, lda, ainv, lda,
547  + c_work_cpot03, lda, s_work_cpot03,
548  + rcondc, result( 2 ) )
549 *
550 * Compute residual of the computed solution.
551 *
552  CALL clacpy( 'Full', n, nrhs, b, lda,
553  + c_work_cpot02, lda )
554  CALL cpot02( uplo, n, nrhs, a, lda, x, lda,
555  + c_work_cpot02, lda, s_work_cpot02,
556  + result( 3 ) )
557 *
558 * Check solution from generated exact solution.
559 *
560  CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
561  + result( 4 ) )
562  nt = 4
563 *
564 * Print information about the tests that did not
565 * pass the threshold.
566 *
567  DO 60 k = 1, nt
568  IF( result( k ).GE.thresh ) THEN
569  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
570  + CALL aladhd( nout, 'CPF' )
571  WRITE( nout, fmt = 9999 )'CPFSV ', uplo,
572  + n, iit, k, result( k )
573  nfail = nfail + 1
574  END IF
575  60 CONTINUE
576  nrun = nrun + nt
577  100 CONTINUE
578  110 CONTINUE
579  120 CONTINUE
580  980 CONTINUE
581  130 CONTINUE
582 *
583 * Print a summary of the results.
584 *
585  CALL alasvm( 'CPF', nout, nfail, nrun, nerrs )
586 *
587  9999 FORMAT( 1x, a6, ', UPLO=''', a1, ''', N =', i5, ', type ', i1,
588  + ', test(', i1, ')=', g12.5 )
589 *
590  RETURN
591 *
592 * End of CDRVRFP
593 *
alasvm
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:75
alaerh
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
Definition: alaerh.f:149
clarhs
subroutine clarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
CLARHS
Definition: clarhs.f:211
cpotri
subroutine cpotri(UPLO, N, A, LDA, INFO)
CPOTRI
Definition: cpotri.f:97
claipd
subroutine claipd(N, A, INDA, VINDA)
CLAIPD
Definition: claipd.f:85
cpot01
subroutine cpot01(UPLO, N, A, LDA, AFAC, LDAFAC, RWORK, RESID)
CPOT01
Definition: cpot01.f:108
clanhe
real function clanhe(NORM, UPLO, N, A, LDA, WORK)
CLANHE returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex Hermitian matrix.
Definition: clanhe.f:126
cpftrs
subroutine cpftrs(TRANSR, UPLO, N, NRHS, A, B, LDB, INFO)
CPFTRS
Definition: cpftrs.f:222
cpot02
subroutine cpot02(UPLO, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
CPOT02
Definition: cpot02.f:129
ctfttr
subroutine ctfttr(TRANSR, UPLO, N, ARF, A, LDA, INFO)
CTFTTR copies a triangular matrix from the rectangular full packed format (TF) to the standard full f...
Definition: ctfttr.f:218
clatms
subroutine clatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
CLATMS
Definition: clatms.f:334
aladhd
subroutine aladhd(IOUNIT, PATH)
ALADHD
Definition: aladhd.f:80
cpotrf
subroutine cpotrf(UPLO, N, A, LDA, INFO)
CPOTRF
Definition: cpotrf.f:109
clacpy
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
cpot03
subroutine cpot03(UPLO, N, A, LDA, AINV, LDAINV, WORK, LDWORK, RWORK, RCOND, RESID)
CPOT03
Definition: cpot03.f:128
cpftri
subroutine cpftri(TRANSR, UPLO, N, A, INFO)
CPFTRI
Definition: cpftri.f:214
clatb4
subroutine clatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
CLATB4
Definition: clatb4.f:123
cpftrf
subroutine cpftrf(TRANSR, UPLO, N, A, INFO)
CPFTRF
Definition: cpftrf.f:213
cget04
subroutine cget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
CGET04
Definition: cget04.f:104
ctrttf
subroutine ctrttf(TRANSR, UPLO, N, A, LDA, ARF, INFO)
CTRTTF copies a triangular matrix from the standard full format (TR) to the rectangular full packed f...
Definition: ctrttf.f:218

Here is the call graph for this function:

Here is the caller graph for this function: