LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ cdrvrfp()

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, conversion 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.

Definition at line 238 of file cdrvrfp.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  REAL THRESH
252 * ..
253 * .. Array Arguments ..
254  INTEGER NVAL( NN ), NSVAL( NNS ), NTVAL( NNT )
255  COMPLEX A( * )
256  COMPLEX AINV( * )
257  COMPLEX ASAV( * )
258  COMPLEX B( * )
259  COMPLEX BSAV( * )
260  COMPLEX AFAC( * )
261  COMPLEX ARF( * )
262  COMPLEX ARFINV( * )
263  COMPLEX XACT( * )
264  COMPLEX X( * )
265  COMPLEX C_WORK_CLATMS( * )
266  COMPLEX C_WORK_CPOT02( * )
267  COMPLEX C_WORK_CPOT03( * )
268  REAL S_WORK_CLATMS( * )
269  REAL S_WORK_CLANHE( * )
270  REAL S_WORK_CPOT01( * )
271  REAL S_WORK_CPOT02( * )
272  REAL S_WORK_CPOT03( * )
273 * ..
274 *
275 * =====================================================================
276 *
277 * .. Parameters ..
278  REAL ONE, ZERO
279  parameter( one = 1.0e+0, zero = 0.0e+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  REAL ANORM, AINVNM, CNDNUM, RCONDC
291 * ..
292 * .. Local Arrays ..
293  CHARACTER UPLOS( 2 ), FORMS( 2 )
294  INTEGER ISEED( 4 ), ISEEDY( 4 )
295  REAL RESULT( NTESTS )
296 * ..
297 * .. External Functions ..
298  REAL CLANHE
299  EXTERNAL clanhe
300 * ..
301 * .. External Subroutines ..
302  EXTERNAL aladhd, alaerh, alasvm, cget04, ctfttr, clacpy,
305  + ctrttf
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 CLATB4 and generate a test
363 * matrix with CLATMS.
364 *
365  CALL clatb4( 'CPO', imat, n, n, ctype, kl, ku,
366  + anorm, mode, cndnum, dist )
367 *
368  srnamt = 'CLATMS'
369  CALL clatms( n, n, dist, iseed, ctype,
370  + s_work_clatms,
371  + mode, cndnum, anorm, kl, ku, uplo, a,
372  + lda, c_work_clatms, info )
373 *
374 * Check error code from CLATMS.
375 *
376  IF( info.NE.0 ) THEN
377  CALL alaerh( 'CPF', 'CLATMS', 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 claipd( n, a, lda+1, 0 )
426 *
427 * Save a copy of the matrix A in ASAV.
428 *
429  CALL clacpy( 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 = clanhe( '1', uplo, n, a, lda,
440  + s_work_clanhe )
441 *
442 * Factor the matrix A.
443 *
444  CALL cpotrf( uplo, n, a, lda, info )
445 *
446 * Form the inverse of A.
447 *
448  CALL cpotri( 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 = clanhe( '1', uplo, n, a, lda,
455  + s_work_clanhe )
456  rcondc = ( one / anorm ) / ainvnm
457 *
458 * Restore the matrix A.
459 *
460  CALL clacpy( 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 = 'CLARHS'
468  CALL clarhs( 'CPO', 'N', uplo, ' ', n, n, kl, ku,
469  + nrhs, a, lda, xact, lda, b, lda,
470  + iseed, info )
471  CALL clacpy( '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 clacpy( uplo, n, n, a, lda, afac, lda )
477  CALL clacpy( 'Full', n, nrhs, b, ldb, x, ldb )
478 *
479  srnamt = 'CTRTTF'
480  CALL ctrttf( cform, uplo, n, afac, lda, arf, info )
481  srnamt = 'CPFTRF'
482  CALL cpftrf( cform, uplo, n, arf, info )
483 *
484 * Check error code from CPFTRF.
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( 'CPF', 'CPFSV ', 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 = 'CPFTRS'
505  CALL cpftrs( cform, uplo, n, nrhs, arf, x, ldb,
506  + info )
507 *
508  srnamt = 'CTFTTR'
509  CALL ctfttr( cform, uplo, n, arf, afac, lda, info )
510 *
511 * Reconstruct matrix from factors and compute
512 * residual.
513 *
514  CALL clacpy( uplo, n, n, afac, lda, asav, lda )
515  CALL cpot01( uplo, n, a, lda, afac, lda,
516  + s_work_cpot01, result( 1 ) )
517  CALL clacpy( 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 clacpy( 'A', n+1, n/2, arf, n+1, arfinv,
523  + n+1 )
524  ELSE
525  CALL clacpy( 'A', n, (n+1)/2, arf, n, arfinv,
526  + n )
527  END IF
528 *
529  srnamt = 'CPFTRI'
530  CALL cpftri( cform, uplo, n, arfinv , info )
531 *
532  srnamt = 'CTFTTR'
533  CALL ctfttr( cform, uplo, n, arfinv, ainv, lda,
534  + info )
535 *
536 * Check error code from CPFTRI.
537 *
538  IF( info.NE.0 )
539  + CALL alaerh( 'CPO', 'CPFTRI', info, 0, uplo, n,
540  + n, -1, -1, -1, imat, nfail, nerrs,
541  + nout )
542 *
543  CALL cpot03( uplo, n, a, lda, ainv, lda,
544  + c_work_cpot03, lda, s_work_cpot03,
545  + rcondc, result( 2 ) )
546 *
547 * Compute residual of the computed solution.
548 *
549  CALL clacpy( 'Full', n, nrhs, b, lda,
550  + c_work_cpot02, lda )
551  CALL cpot02( uplo, n, nrhs, a, lda, x, lda,
552  + c_work_cpot02, lda, s_work_cpot02,
553  + result( 3 ) )
554 *
555 * Check solution from generated exact solution.
556 *
557  CALL cget04( 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, 'CPF' )
568  WRITE( nout, fmt = 9999 )'CPFSV ', 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( 'CPF', 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 CDRVRFP
590 *
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:73
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 clarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
CLARHS
Definition: clarhs.f:208
subroutine clatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
CLATB4
Definition: clatb4.f:121
subroutine cget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
CGET04
Definition: cget04.f:102
subroutine cpot03(UPLO, N, A, LDA, AINV, LDAINV, WORK, LDWORK, RWORK, RCOND, RESID)
CPOT03
Definition: cpot03.f:126
subroutine claipd(N, A, INDA, VINDA)
CLAIPD
Definition: claipd.f:83
subroutine cpot01(UPLO, N, A, LDA, AFAC, LDAFAC, RWORK, RESID)
CPOT01
Definition: cpot01.f:106
subroutine cpot02(UPLO, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
CPOT02
Definition: cpot02.f:127
subroutine clatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
CLATMS
Definition: clatms.f:332
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,...
Definition: clanhe.f:124
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 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:216
subroutine cpftrs(TRANSR, UPLO, N, NRHS, A, B, LDB, INFO)
CPFTRS
Definition: cpftrs.f:220
subroutine cpftri(TRANSR, UPLO, N, A, INFO)
CPFTRI
Definition: cpftri.f:212
subroutine cpftrf(TRANSR, UPLO, N, A, INFO)
CPFTRF
Definition: cpftrf.f:211
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:216
subroutine cpotrf(UPLO, N, A, LDA, INFO)
CPOTRF
Definition: cpotrf.f:107
subroutine cpotri(UPLO, N, A, LDA, INFO)
CPOTRI
Definition: cpotri.f:95
Here is the call graph for this function:
Here is the caller graph for this function: