LAPACK  3.9.0
LAPACK: Linear Algebra PACKage
sdrvsy.f
Go to the documentation of this file.
1 *> \brief \b SDRVSY
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 SDRVSY( 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 A( * ), AFAC( * ), AINV( * ), B( * ),
24 * $ RWORK( * ), WORK( * ), X( * ), XACT( * )
25 * ..
26 *
27 *
28 *> \par Purpose:
29 * =============
30 *>
31 *> \verbatim
32 *>
33 *> SDRVSY tests the driver routines SSYSV and -SVX.
34 *> \endverbatim
35 *
36 * Arguments:
37 * ==========
38 *
39 *> \param[in] DOTYPE
40 *> \verbatim
41 *> DOTYPE is LOGICAL array, dimension (NTYPES)
42 *> The matrix types to be used for testing. Matrices of type j
43 *> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
44 *> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
45 *> \endverbatim
46 *>
47 *> \param[in] NN
48 *> \verbatim
49 *> NN is INTEGER
50 *> The number of values of N contained in the vector NVAL.
51 *> \endverbatim
52 *>
53 *> \param[in] NVAL
54 *> \verbatim
55 *> NVAL is INTEGER array, dimension (NN)
56 *> The values of the matrix dimension N.
57 *> \endverbatim
58 *>
59 *> \param[in] NRHS
60 *> \verbatim
61 *> NRHS is INTEGER
62 *> The number of right hand side vectors to be generated for
63 *> each linear system.
64 *> \endverbatim
65 *>
66 *> \param[in] THRESH
67 *> \verbatim
68 *> THRESH is REAL
69 *> The threshold value for the test ratios. A result is
70 *> included in the output file if RESULT >= THRESH. To have
71 *> every test ratio printed, use THRESH = 0.
72 *> \endverbatim
73 *>
74 *> \param[in] TSTERR
75 *> \verbatim
76 *> TSTERR is LOGICAL
77 *> Flag that indicates whether error exits are to be tested.
78 *> \endverbatim
79 *>
80 *> \param[in] NMAX
81 *> \verbatim
82 *> NMAX is INTEGER
83 *> The maximum value permitted for N, used in dimensioning the
84 *> work arrays.
85 *> \endverbatim
86 *>
87 *> \param[out] A
88 *> \verbatim
89 *> A is REAL array, dimension (NMAX*NMAX)
90 *> \endverbatim
91 *>
92 *> \param[out] AFAC
93 *> \verbatim
94 *> AFAC is REAL array, dimension (NMAX*NMAX)
95 *> \endverbatim
96 *>
97 *> \param[out] AINV
98 *> \verbatim
99 *> AINV is REAL array, dimension (NMAX*NMAX)
100 *> \endverbatim
101 *>
102 *> \param[out] B
103 *> \verbatim
104 *> B is REAL array, dimension (NMAX*NRHS)
105 *> \endverbatim
106 *>
107 *> \param[out] X
108 *> \verbatim
109 *> X is REAL array, dimension (NMAX*NRHS)
110 *> \endverbatim
111 *>
112 *> \param[out] XACT
113 *> \verbatim
114 *> XACT is REAL array, dimension (NMAX*NRHS)
115 *> \endverbatim
116 *>
117 *> \param[out] WORK
118 *> \verbatim
119 *> WORK is REAL array, dimension (NMAX*max(2,NRHS))
120 *> \endverbatim
121 *>
122 *> \param[out] RWORK
123 *> \verbatim
124 *> RWORK is REAL array, dimension (NMAX+2*NRHS)
125 *> \endverbatim
126 *>
127 *> \param[out] IWORK
128 *> \verbatim
129 *> IWORK is INTEGER array, dimension (2*NMAX)
130 *> \endverbatim
131 *>
132 *> \param[in] NOUT
133 *> \verbatim
134 *> NOUT is INTEGER
135 *> The unit number for output.
136 *> \endverbatim
137 *
138 * Authors:
139 * ========
140 *
141 *> \author Univ. of Tennessee
142 *> \author Univ. of California Berkeley
143 *> \author Univ. of Colorado Denver
144 *> \author NAG Ltd.
145 *
146 *> \date November 2013
147 *
148 *> \ingroup single_lin
149 *
150 * =====================================================================
151  SUBROUTINE sdrvsy( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
152  $ A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK,
153  $ NOUT )
154 *
155 * -- LAPACK test routine (version 3.5.0) --
156 * -- LAPACK is a software package provided by Univ. of Tennessee, --
157 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
158 * November 2013
159 *
160 * .. Scalar Arguments ..
161  LOGICAL tsterr
162  INTEGER nmax, nn, nout, nrhs
163  REAL thresh
164 * ..
165 * .. Array Arguments ..
166  LOGICAL dotype( * )
167  INTEGER IWORK( * ), NVAL( * )
168  REAL A( * ), AFAC( * ), AINV( * ), B( * ),
169  $ rwork( * ), work( * ), x( * ), xact( * )
170 * ..
171 *
172 * =====================================================================
173 *
174 * .. Parameters ..
175  REAL ONE, ZERO
176  PARAMETER ( ONE = 1.0e+0, zero = 0.0e+0 )
177  INTEGER NTYPES, NTESTS
178  parameter( ntypes = 10, ntests = 6 )
179  INTEGER NFACT
180  parameter( nfact = 2 )
181 * ..
182 * .. Local Scalars ..
183  LOGICAL ZEROT
184  CHARACTER DIST, FACT, TYPE, UPLO, XTYPE
185  CHARACTER*3 PATH
186  INTEGER I, I1, I2, IFACT, IMAT, IN, INFO, IOFF, IUPLO,
187  $ izero, j, k, k1, kl, ku, lda, lwork, mode, n,
188  $ nb, nbmin, nerrs, nfail, nimat, nrun, nt
189  REAL AINVNM, ANORM, CNDNUM, RCOND, RCONDC
190 * ..
191 * .. Local Arrays ..
192  CHARACTER FACTS( NFACT ), UPLOS( 2 )
193  INTEGER ISEED( 4 ), ISEEDY( 4 )
194  REAL RESULT( NTESTS )
195 * ..
196 * .. External Functions ..
197  REAL SGET06, SLANSY
198  EXTERNAL SGET06, SLANSY
199 * ..
200 * .. External Subroutines ..
201  EXTERNAL aladhd, alaerh, alasvm, serrvx, sget04, slacpy,
202  $ slarhs, slaset, slatb4, slatms, spot02, spot05,
203  $ ssysv, ssysvx, ssyt01, ssytrf, ssytri2, 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 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 ) = 'Single precision'
226  path( 2: 3 ) = 'SY'
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 serrvx( 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 SLATB4 and generate a test matrix
277 * with SLATMS.
278 *
279  CALL slatb4( path, imat, n, n, TYPE, kl, ku, anorm, mode,
280  $ cndnum, dist )
281 *
282  srnamt = 'SLATMS'
283  CALL slatms( n, n, dist, iseed, TYPE, rwork, mode,
284  $ cndnum, anorm, kl, ku, uplo, a, lda, work,
285  $ info )
286 *
287 * Check error code from SLATMS.
288 *
289  IF( info.NE.0 ) THEN
290  CALL alaerh( path, 'SLATMS', 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  DO 150 ifact = 1, nfact
363 *
364 * Do first for FACT = 'F', then for other values.
365 *
366  fact = facts( ifact )
367 *
368 * Compute the condition number for comparison with
369 * the value returned by SSYSVX.
370 *
371  IF( zerot ) THEN
372  IF( ifact.EQ.1 )
373  $ GO TO 150
374  rcondc = zero
375 *
376  ELSE IF( ifact.EQ.1 ) THEN
377 *
378 * Compute the 1-norm of A.
379 *
380  anorm = slansy( '1', uplo, n, a, lda, rwork )
381 *
382 * Factor the matrix A.
383 *
384  CALL slacpy( uplo, n, n, a, lda, afac, lda )
385  CALL ssytrf( uplo, n, afac, lda, iwork, work,
386  $ lwork, info )
387 *
388 * Compute inv(A) and take its norm.
389 *
390  CALL slacpy( uplo, n, n, afac, lda, ainv, lda )
391  lwork = (n+nb+1)*(nb+3)
392  CALL ssytri2( uplo, n, ainv, lda, iwork, work,
393  $ lwork, info )
394  ainvnm = slansy( '1', uplo, n, ainv, lda, rwork )
395 *
396 * Compute the 1-norm condition number of A.
397 *
398  IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
399  rcondc = one
400  ELSE
401  rcondc = ( one / anorm ) / ainvnm
402  END IF
403  END IF
404 *
405 * Form an exact solution and set the right hand side.
406 *
407  srnamt = 'SLARHS'
408  CALL slarhs( path, xtype, uplo, ' ', n, n, kl, ku,
409  $ nrhs, a, lda, xact, lda, b, lda, iseed,
410  $ info )
411  xtype = 'C'
412 *
413 * --- Test SSYSV ---
414 *
415  IF( ifact.EQ.2 ) THEN
416  CALL slacpy( uplo, n, n, a, lda, afac, lda )
417  CALL slacpy( 'Full', n, nrhs, b, lda, x, lda )
418 *
419 * Factor the matrix and solve the system using SSYSV.
420 *
421  srnamt = 'SSYSV '
422  CALL ssysv( uplo, n, nrhs, afac, lda, iwork, x,
423  $ lda, work, lwork, info )
424 *
425 * Adjust the expected value of INFO to account for
426 * pivoting.
427 *
428  k = izero
429  IF( k.GT.0 ) THEN
430  100 CONTINUE
431  IF( iwork( k ).LT.0 ) THEN
432  IF( iwork( k ).NE.-k ) THEN
433  k = -iwork( k )
434  GO TO 100
435  END IF
436  ELSE IF( iwork( k ).NE.k ) THEN
437  k = iwork( k )
438  GO TO 100
439  END IF
440  END IF
441 *
442 * Check error code from SSYSV .
443 *
444  IF( info.NE.k ) THEN
445  CALL alaerh( path, 'SSYSV ', info, k, uplo, n,
446  $ n, -1, -1, nrhs, imat, nfail,
447  $ nerrs, nout )
448  GO TO 120
449  ELSE IF( info.NE.0 ) THEN
450  GO TO 120
451  END IF
452 *
453 * Reconstruct matrix from factors and compute
454 * residual.
455 *
456  CALL ssyt01( uplo, n, a, lda, afac, lda, iwork,
457  $ ainv, lda, rwork, result( 1 ) )
458 *
459 * Compute residual of the computed solution.
460 *
461  CALL slacpy( 'Full', n, nrhs, b, lda, work, lda )
462  CALL spot02( uplo, n, nrhs, a, lda, x, lda, work,
463  $ lda, rwork, result( 2 ) )
464 *
465 * Check solution from generated exact solution.
466 *
467  CALL sget04( n, nrhs, x, lda, xact, lda, rcondc,
468  $ result( 3 ) )
469  nt = 3
470 *
471 * Print information about the tests that did not pass
472 * the threshold.
473 *
474  DO 110 k = 1, nt
475  IF( result( k ).GE.thresh ) THEN
476  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
477  $ CALL aladhd( nout, path )
478  WRITE( nout, fmt = 9999 )'SSYSV ', uplo, n,
479  $ imat, k, result( k )
480  nfail = nfail + 1
481  END IF
482  110 CONTINUE
483  nrun = nrun + nt
484  120 CONTINUE
485  END IF
486 *
487 * --- Test SSYSVX ---
488 *
489  IF( ifact.EQ.2 )
490  $ CALL slaset( uplo, n, n, zero, zero, afac, lda )
491  CALL slaset( 'Full', n, nrhs, zero, zero, x, lda )
492 *
493 * Solve the system and compute the condition number and
494 * error bounds using SSYSVX.
495 *
496  srnamt = 'SSYSVX'
497  CALL ssysvx( fact, uplo, n, nrhs, a, lda, afac, lda,
498  $ iwork, b, lda, x, lda, rcond, rwork,
499  $ rwork( nrhs+1 ), work, lwork,
500  $ iwork( n+1 ), info )
501 *
502 * Adjust the expected value of INFO to account for
503 * pivoting.
504 *
505  k = izero
506  IF( k.GT.0 ) THEN
507  130 CONTINUE
508  IF( iwork( k ).LT.0 ) THEN
509  IF( iwork( k ).NE.-k ) THEN
510  k = -iwork( k )
511  GO TO 130
512  END IF
513  ELSE IF( iwork( k ).NE.k ) THEN
514  k = iwork( k )
515  GO TO 130
516  END IF
517  END IF
518 *
519 * Check the error code from SSYSVX.
520 *
521  IF( info.NE.k ) THEN
522  CALL alaerh( path, 'SSYSVX', info, k, fact // uplo,
523  $ n, n, -1, -1, nrhs, imat, nfail,
524  $ nerrs, nout )
525  GO TO 150
526  END IF
527 *
528  IF( info.EQ.0 ) THEN
529  IF( ifact.GE.2 ) THEN
530 *
531 * Reconstruct matrix from factors and compute
532 * residual.
533 *
534  CALL ssyt01( uplo, n, a, lda, afac, lda, iwork,
535  $ ainv, lda, rwork( 2*nrhs+1 ),
536  $ result( 1 ) )
537  k1 = 1
538  ELSE
539  k1 = 2
540  END IF
541 *
542 * Compute residual of the computed solution.
543 *
544  CALL slacpy( 'Full', n, nrhs, b, lda, work, lda )
545  CALL spot02( uplo, n, nrhs, a, lda, x, lda, work,
546  $ lda, rwork( 2*nrhs+1 ), result( 2 ) )
547 *
548 * Check solution from generated exact solution.
549 *
550  CALL sget04( n, nrhs, x, lda, xact, lda, rcondc,
551  $ result( 3 ) )
552 *
553 * Check the error bounds from iterative refinement.
554 *
555  CALL spot05( uplo, n, nrhs, a, lda, b, lda, x, lda,
556  $ xact, lda, rwork, rwork( nrhs+1 ),
557  $ result( 4 ) )
558  ELSE
559  k1 = 6
560  END IF
561 *
562 * Compare RCOND from SSYSVX with the computed value
563 * in RCONDC.
564 *
565  result( 6 ) = sget06( rcond, rcondc )
566 *
567 * Print information about the tests that did not pass
568 * the threshold.
569 *
570  DO 140 k = k1, 6
571  IF( result( k ).GE.thresh ) THEN
572  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
573  $ CALL aladhd( nout, path )
574  WRITE( nout, fmt = 9998 )'SSYSVX', fact, uplo,
575  $ n, imat, k, result( k )
576  nfail = nfail + 1
577  END IF
578  140 CONTINUE
579  nrun = nrun + 7 - k1
580 *
581  150 CONTINUE
582 *
583  160 CONTINUE
584  170 CONTINUE
585  180 CONTINUE
586 *
587 * Print a summary of the results.
588 *
589  CALL alasvm( path, nout, nfail, nrun, nerrs )
590 *
591  9999 FORMAT( 1x, a, ', UPLO=''', a1, ''', N =', i5, ', type ', i2,
592  $ ', test ', i2, ', ratio =', g12.5 )
593  9998 FORMAT( 1x, a, ', FACT=''', a1, ''', UPLO=''', a1, ''', N =', i5,
594  $ ', type ', i2, ', test ', i2, ', ratio =', g12.5 )
595  RETURN
596 *
597 * End of SDRVSY
598 *
599  END
slarhs
subroutine slarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
SLARHS
Definition: slarhs.f:206
ssysvx
subroutine ssysvx(FACT, UPLO, N, NRHS, A, LDA, AF, LDAF, IPIV, B, LDB, X, LDX, RCOND, FERR, BERR, WORK, LWORK, IWORK, INFO)
SSYSVX computes the solution to system of linear equations A * X = B for SY matrices
Definition: ssysvx.f:286
ssytri2
subroutine ssytri2(UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO)
SSYTRI2
Definition: ssytri2.f:129
alasvm
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:75
spot02
subroutine spot02(UPLO, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
SPOT02
Definition: spot02.f:129
ssytrf
subroutine ssytrf(UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO)
SSYTRF
Definition: ssytrf.f:184
ssysv
subroutine ssysv(UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK, LWORK, INFO)
SSYSV computes the solution to system of linear equations A * X = B for SY matrices
Definition: ssysv.f:173
slacpy
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:105
ssyt01
subroutine ssyt01(UPLO, N, A, LDA, AFAC, LDAFAC, IPIV, C, LDC, RWORK, RESID)
SSYT01
Definition: ssyt01.f:126
aladhd
subroutine aladhd(IOUNIT, PATH)
ALADHD
Definition: aladhd.f:92
slaset
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: slaset.f:112
sget04
subroutine sget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
SGET04
Definition: sget04.f:104
sdrvsy
subroutine sdrvsy(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
SDRVSY
Definition: sdrvsy.f:154
alaerh
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
Definition: alaerh.f:149
slatms
subroutine slatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
SLATMS
Definition: slatms.f:323
slatb4
subroutine slatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
SLATB4
Definition: slatb4.f:122
spot05
subroutine spot05(UPLO, N, NRHS, A, LDA, B, LDB, X, LDX, XACT, LDXACT, FERR, BERR, RESLTS)
SPOT05
Definition: spot05.f:166
xlaenv
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
Definition: xlaenv.f:83
serrvx
subroutine serrvx(PATH, NUNIT)
SERRVX
Definition: serrvx.f:57