LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ sdrvrfp()

subroutine sdrvrfp ( integer  NOUT,
integer  NN,
integer, dimension( nn )  NVAL,
integer  NNS,
integer, dimension( nns )  NSVAL,
integer  NNT,
integer, dimension( nnt )  NTVAL,
real  THRESH,
real, dimension( * )  A,
real, dimension( * )  ASAV,
real, dimension( * )  AFAC,
real, dimension( * )  AINV,
real, dimension( * )  B,
real, dimension( * )  BSAV,
real, dimension( * )  XACT,
real, dimension( * )  X,
real, dimension( * )  ARF,
real, dimension( * )  ARFINV,
real, dimension( * )  S_WORK_SLATMS,
real, dimension( * )  S_WORK_SPOT01,
real, dimension( * )  S_TEMP_SPOT02,
real, dimension( * )  S_TEMP_SPOT03,
real, dimension( * )  S_WORK_SLANSY,
real, dimension( * )  S_WORK_SPOT02,
real, dimension( * )  S_WORK_SPOT03 
)

SDRVRFP

Purpose:
 SDRVRFP tests the LAPACK RFP routines:
     SPFTRF, SPFTRS, and SPFTRI.

 This testing routine follow the same tests as DDRVPO (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 STRTTF and
 STFTTR.

 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 SPFTRF, 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 SPFTRF 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 REAL array, dimension (NMAX*NMAX)
[out]ASAV
          ASAV is REAL array, dimension (NMAX*NMAX)
[out]AFAC
          AFAC is REAL array, dimension (NMAX*NMAX)
[out]AINV
          AINV is REAL array, dimension (NMAX*NMAX)
[out]B
          B is REAL array, dimension (NMAX*MAXRHS)
[out]BSAV
          BSAV is REAL array, dimension (NMAX*MAXRHS)
[out]XACT
          XACT is REAL array, dimension (NMAX*MAXRHS)
[out]X
          X is REAL array, dimension (NMAX*MAXRHS)
[out]ARF
          ARF is REAL array, dimension ((NMAX*(NMAX+1))/2)
[out]ARFINV
          ARFINV is REAL array, dimension ((NMAX*(NMAX+1))/2)
[out]S_WORK_SLATMS
          S_WORK_SLATMS is REAL array, dimension ( 3*NMAX )
[out]S_WORK_SPOT01
          S_WORK_SPOT01 is REAL array, dimension ( NMAX )
[out]S_TEMP_SPOT02
          S_TEMP_SPOT02 is REAL array, dimension ( NMAX*MAXRHS )
[out]S_TEMP_SPOT03
          S_TEMP_SPOT03 is REAL array, dimension ( NMAX*NMAX )
[out]S_WORK_SLANSY
          S_WORK_SLANSY is REAL array, dimension ( NMAX )
[out]S_WORK_SPOT02
          S_WORK_SPOT02 is REAL array, dimension ( NMAX )
[out]S_WORK_SPOT03
          S_WORK_SPOT03 is REAL array, dimension ( NMAX )
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 232 of file sdrvrfp.f.

238 *
239 * -- LAPACK test routine --
240 * -- LAPACK is a software package provided by Univ. of Tennessee, --
241 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
242 *
243 * .. Scalar Arguments ..
244  INTEGER NN, NNS, NNT, NOUT
245  REAL THRESH
246 * ..
247 * .. Array Arguments ..
248  INTEGER NVAL( NN ), NSVAL( NNS ), NTVAL( NNT )
249  REAL A( * )
250  REAL AINV( * )
251  REAL ASAV( * )
252  REAL B( * )
253  REAL BSAV( * )
254  REAL AFAC( * )
255  REAL ARF( * )
256  REAL ARFINV( * )
257  REAL XACT( * )
258  REAL X( * )
259  REAL S_WORK_SLATMS( * )
260  REAL S_WORK_SPOT01( * )
261  REAL S_TEMP_SPOT02( * )
262  REAL S_TEMP_SPOT03( * )
263  REAL S_WORK_SLANSY( * )
264  REAL S_WORK_SPOT02( * )
265  REAL S_WORK_SPOT03( * )
266 * ..
267 *
268 * =====================================================================
269 *
270 * .. Parameters ..
271  REAL ONE, ZERO
272  parameter( one = 1.0e+0, zero = 0.0e+0 )
273  INTEGER NTESTS
274  parameter( ntests = 4 )
275 * ..
276 * .. Local Scalars ..
277  LOGICAL ZEROT
278  INTEGER I, INFO, IUPLO, LDA, LDB, IMAT, NERRS, NFAIL,
279  + NRHS, NRUN, IZERO, IOFF, K, NT, N, IFORM, IIN,
280  + IIT, IIS
281  CHARACTER DIST, CTYPE, UPLO, CFORM
282  INTEGER KL, KU, MODE
283  REAL ANORM, AINVNM, CNDNUM, RCONDC
284 * ..
285 * .. Local Arrays ..
286  CHARACTER UPLOS( 2 ), FORMS( 2 )
287  INTEGER ISEED( 4 ), ISEEDY( 4 )
288  REAL RESULT( NTESTS )
289 * ..
290 * .. External Functions ..
291  REAL SLANSY
292  EXTERNAL slansy
293 * ..
294 * .. External Subroutines ..
295  EXTERNAL aladhd, alaerh, alasvm, sget04, stfttr, slacpy,
298 * ..
299 * .. Scalars in Common ..
300  CHARACTER*32 SRNAMT
301 * ..
302 * .. Common blocks ..
303  COMMON / srnamc / srnamt
304 * ..
305 * .. Data statements ..
306  DATA iseedy / 1988, 1989, 1990, 1991 /
307  DATA uplos / 'U', 'L' /
308  DATA forms / 'N', 'T' /
309 * ..
310 * .. Executable Statements ..
311 *
312 * Initialize constants and the random number seed.
313 *
314  nrun = 0
315  nfail = 0
316  nerrs = 0
317  DO 10 i = 1, 4
318  iseed( i ) = iseedy( i )
319  10 CONTINUE
320 *
321  DO 130 iin = 1, nn
322 *
323  n = nval( iin )
324  lda = max( n, 1 )
325  ldb = max( n, 1 )
326 *
327  DO 980 iis = 1, nns
328 *
329  nrhs = nsval( iis )
330 *
331  DO 120 iit = 1, nnt
332 *
333  imat = ntval( iit )
334 *
335 * If N.EQ.0, only consider the first type
336 *
337  IF( n.EQ.0 .AND. iit.GE.1 ) GO TO 120
338 *
339 * Skip types 3, 4, or 5 if the matrix size is too small.
340 *
341  IF( imat.EQ.4 .AND. n.LE.1 ) GO TO 120
342  IF( imat.EQ.5 .AND. n.LE.2 ) GO TO 120
343 *
344 * Do first for UPLO = 'U', then for UPLO = 'L'
345 *
346  DO 110 iuplo = 1, 2
347  uplo = uplos( iuplo )
348 *
349 * Do first for CFORM = 'N', then for CFORM = 'C'
350 *
351  DO 100 iform = 1, 2
352  cform = forms( iform )
353 *
354 * Set up parameters with SLATB4 and generate a test
355 * matrix with SLATMS.
356 *
357  CALL slatb4( 'SPO', imat, n, n, ctype, kl, ku,
358  + anorm, mode, cndnum, dist )
359 *
360  srnamt = 'SLATMS'
361  CALL slatms( n, n, dist, iseed, ctype,
362  + s_work_slatms,
363  + mode, cndnum, anorm, kl, ku, uplo, a,
364  + lda, s_work_slatms, info )
365 *
366 * Check error code from SLATMS.
367 *
368  IF( info.NE.0 ) THEN
369  CALL alaerh( 'SPF', 'SLATMS', info, 0, uplo, n,
370  + n, -1, -1, -1, iit, nfail, nerrs,
371  + nout )
372  GO TO 100
373  END IF
374 *
375 * For types 3-5, zero one row and column of the matrix to
376 * test that INFO is returned correctly.
377 *
378  zerot = imat.GE.3 .AND. imat.LE.5
379  IF( zerot ) THEN
380  IF( iit.EQ.3 ) THEN
381  izero = 1
382  ELSE IF( iit.EQ.4 ) THEN
383  izero = n
384  ELSE
385  izero = n / 2 + 1
386  END IF
387  ioff = ( izero-1 )*lda
388 *
389 * Set row and column IZERO of A to 0.
390 *
391  IF( iuplo.EQ.1 ) THEN
392  DO 20 i = 1, izero - 1
393  a( ioff+i ) = zero
394  20 CONTINUE
395  ioff = ioff + izero
396  DO 30 i = izero, n
397  a( ioff ) = zero
398  ioff = ioff + lda
399  30 CONTINUE
400  ELSE
401  ioff = izero
402  DO 40 i = 1, izero - 1
403  a( ioff ) = zero
404  ioff = ioff + lda
405  40 CONTINUE
406  ioff = ioff - izero
407  DO 50 i = izero, n
408  a( ioff+i ) = zero
409  50 CONTINUE
410  END IF
411  ELSE
412  izero = 0
413  END IF
414 *
415 * Save a copy of the matrix A in ASAV.
416 *
417  CALL slacpy( uplo, n, n, a, lda, asav, lda )
418 *
419 * Compute the condition number of A (RCONDC).
420 *
421  IF( zerot ) THEN
422  rcondc = zero
423  ELSE
424 *
425 * Compute the 1-norm of A.
426 *
427  anorm = slansy( '1', uplo, n, a, lda,
428  + s_work_slansy )
429 *
430 * Factor the matrix A.
431 *
432  CALL spotrf( uplo, n, a, lda, info )
433 *
434 * Form the inverse of A.
435 *
436  CALL spotri( uplo, n, a, lda, info )
437 
438  IF ( n .NE. 0 ) THEN
439 *
440 * Compute the 1-norm condition number of A.
441 *
442  ainvnm = slansy( '1', uplo, n, a, lda,
443  + s_work_slansy )
444  rcondc = ( one / anorm ) / ainvnm
445 *
446 * Restore the matrix A.
447 *
448  CALL slacpy( uplo, n, n, asav, lda, a, lda )
449  END IF
450 *
451  END IF
452 *
453 * Form an exact solution and set the right hand side.
454 *
455  srnamt = 'SLARHS'
456  CALL slarhs( 'SPO', 'N', uplo, ' ', n, n, kl, ku,
457  + nrhs, a, lda, xact, lda, b, lda,
458  + iseed, info )
459  CALL slacpy( 'Full', n, nrhs, b, lda, bsav, lda )
460 *
461 * Compute the L*L' or U'*U factorization of the
462 * matrix and solve the system.
463 *
464  CALL slacpy( uplo, n, n, a, lda, afac, lda )
465  CALL slacpy( 'Full', n, nrhs, b, ldb, x, ldb )
466 *
467  srnamt = 'STRTTF'
468  CALL strttf( cform, uplo, n, afac, lda, arf, info )
469  srnamt = 'SPFTRF'
470  CALL spftrf( cform, uplo, n, arf, info )
471 *
472 * Check error code from SPFTRF.
473 *
474  IF( info.NE.izero ) THEN
475 *
476 * LANGOU: there is a small hick here: IZERO should
477 * always be INFO however if INFO is ZERO, ALAERH does not
478 * complain.
479 *
480  CALL alaerh( 'SPF', 'SPFSV ', info, izero,
481  + uplo, n, n, -1, -1, nrhs, iit,
482  + nfail, nerrs, nout )
483  GO TO 100
484  END IF
485 *
486 * Skip the tests if INFO is not 0.
487 *
488  IF( info.NE.0 ) THEN
489  GO TO 100
490  END IF
491 *
492  srnamt = 'SPFTRS'
493  CALL spftrs( cform, uplo, n, nrhs, arf, x, ldb,
494  + info )
495 *
496  srnamt = 'STFTTR'
497  CALL stfttr( cform, uplo, n, arf, afac, lda, info )
498 *
499 * Reconstruct matrix from factors and compute
500 * residual.
501 *
502  CALL slacpy( uplo, n, n, afac, lda, asav, lda )
503  CALL spot01( uplo, n, a, lda, afac, lda,
504  + s_work_spot01, result( 1 ) )
505  CALL slacpy( uplo, n, n, asav, lda, afac, lda )
506 *
507 * Form the inverse and compute the residual.
508 *
509  IF(mod(n,2).EQ.0)THEN
510  CALL slacpy( 'A', n+1, n/2, arf, n+1, arfinv,
511  + n+1 )
512  ELSE
513  CALL slacpy( 'A', n, (n+1)/2, arf, n, arfinv,
514  + n )
515  END IF
516 *
517  srnamt = 'SPFTRI'
518  CALL spftri( cform, uplo, n, arfinv , info )
519 *
520  srnamt = 'STFTTR'
521  CALL stfttr( cform, uplo, n, arfinv, ainv, lda,
522  + info )
523 *
524 * Check error code from SPFTRI.
525 *
526  IF( info.NE.0 )
527  + CALL alaerh( 'SPO', 'SPFTRI', info, 0, uplo, n,
528  + n, -1, -1, -1, imat, nfail, nerrs,
529  + nout )
530 *
531  CALL spot03( uplo, n, a, lda, ainv, lda,
532  + s_temp_spot03, lda, s_work_spot03,
533  + rcondc, result( 2 ) )
534 *
535 * Compute residual of the computed solution.
536 *
537  CALL slacpy( 'Full', n, nrhs, b, lda,
538  + s_temp_spot02, lda )
539  CALL spot02( uplo, n, nrhs, a, lda, x, lda,
540  + s_temp_spot02, lda, s_work_spot02,
541  + result( 3 ) )
542 *
543 * Check solution from generated exact solution.
544 
545  CALL sget04( n, nrhs, x, lda, xact, lda, rcondc,
546  + result( 4 ) )
547  nt = 4
548 *
549 * Print information about the tests that did not
550 * pass the threshold.
551 *
552  DO 60 k = 1, nt
553  IF( result( k ).GE.thresh ) THEN
554  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
555  + CALL aladhd( nout, 'SPF' )
556  WRITE( nout, fmt = 9999 )'SPFSV ', uplo,
557  + n, iit, k, result( k )
558  nfail = nfail + 1
559  END IF
560  60 CONTINUE
561  nrun = nrun + nt
562  100 CONTINUE
563  110 CONTINUE
564  120 CONTINUE
565  980 CONTINUE
566  130 CONTINUE
567 *
568 * Print a summary of the results.
569 *
570  CALL alasvm( 'SPF', nout, nfail, nrun, nerrs )
571 *
572  9999 FORMAT( 1x, a6, ', UPLO=''', a1, ''', N =', i5, ', type ', i1,
573  + ', test(', i1, ')=', g12.5 )
574 *
575  RETURN
576 *
577 * End of SDRVRFP
578 *
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:103
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 slatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
SLATMS
Definition: slatms.f:321
subroutine spftrf(TRANSR, UPLO, N, A, INFO)
SPFTRF
Definition: spftrf.f:198
subroutine spftri(TRANSR, UPLO, N, A, INFO)
SPFTRI
Definition: spftri.f:191
subroutine strttf(TRANSR, UPLO, N, A, LDA, ARF, INFO)
STRTTF copies a triangular matrix from the standard full format (TR) to the rectangular full packed f...
Definition: strttf.f:194
subroutine stfttr(TRANSR, UPLO, N, ARF, A, LDA, INFO)
STFTTR copies a triangular matrix from the rectangular full packed format (TF) to the standard full f...
Definition: stfttr.f:196
subroutine spftrs(TRANSR, UPLO, N, NRHS, A, B, LDB, INFO)
SPFTRS
Definition: spftrs.f:199
subroutine spotri(UPLO, N, A, LDA, INFO)
SPOTRI
Definition: spotri.f:95
subroutine spotrf(UPLO, N, A, LDA, INFO)
SPOTRF
Definition: spotrf.f:107
real function slansy(NORM, UPLO, N, A, LDA, WORK)
SLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition: slansy.f:122
subroutine slarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
SLARHS
Definition: slarhs.f:205
subroutine slatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
SLATB4
Definition: slatb4.f:120
subroutine spot01(UPLO, N, A, LDA, AFAC, LDAFAC, RWORK, RESID)
SPOT01
Definition: spot01.f:104
subroutine spot03(UPLO, N, A, LDA, AINV, LDAINV, WORK, LDWORK, RWORK, RCOND, RESID)
SPOT03
Definition: spot03.f:125
subroutine sget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
SGET04
Definition: sget04.f:102
subroutine spot02(UPLO, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
SPOT02
Definition: spot02.f:127
Here is the call graph for this function:
Here is the caller graph for this function: