LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
ddrvpb.f
Go to the documentation of this file.
1*> \brief \b DDRVPB
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 DDRVPB( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
12* A, AFAC, ASAV, B, BSAV, X, XACT, S, WORK,
13* RWORK, IWORK, NOUT )
14*
15* .. Scalar Arguments ..
16* LOGICAL TSTERR
17* INTEGER NMAX, NN, NOUT, NRHS
18* DOUBLE PRECISION THRESH
19* ..
20* .. Array Arguments ..
21* LOGICAL DOTYPE( * )
22* INTEGER IWORK( * ), NVAL( * )
23* DOUBLE PRECISION A( * ), AFAC( * ), ASAV( * ), B( * ),
24* $ BSAV( * ), RWORK( * ), S( * ), WORK( * ),
25* $ X( * ), XACT( * )
26* ..
27*
28*
29*> \par Purpose:
30* =============
31*>
32*> \verbatim
33*>
34*> DDRVPB tests the driver routines DPBSV and -SVX.
35*> \endverbatim
36*
37* Arguments:
38* ==========
39*
40*> \param[in] DOTYPE
41*> \verbatim
42*> DOTYPE is LOGICAL array, dimension (NTYPES)
43*> The matrix types to be used for testing. Matrices of type j
44*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
45*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
46*> \endverbatim
47*>
48*> \param[in] NN
49*> \verbatim
50*> NN is INTEGER
51*> The number of values of N contained in the vector NVAL.
52*> \endverbatim
53*>
54*> \param[in] NVAL
55*> \verbatim
56*> NVAL is INTEGER array, dimension (NN)
57*> The values of the matrix dimension N.
58*> \endverbatim
59*>
60*> \param[in] NRHS
61*> \verbatim
62*> NRHS is INTEGER
63*> The number of right hand side vectors to be generated for
64*> each linear system.
65*> \endverbatim
66*>
67*> \param[in] THRESH
68*> \verbatim
69*> THRESH is DOUBLE PRECISION
70*> The threshold value for the test ratios. A result is
71*> included in the output file if RESULT >= THRESH. To have
72*> every test ratio printed, use THRESH = 0.
73*> \endverbatim
74*>
75*> \param[in] TSTERR
76*> \verbatim
77*> TSTERR is LOGICAL
78*> Flag that indicates whether error exits are to be tested.
79*> \endverbatim
80*>
81*> \param[in] NMAX
82*> \verbatim
83*> NMAX is INTEGER
84*> The maximum value permitted for N, used in dimensioning the
85*> work arrays.
86*> \endverbatim
87*>
88*> \param[out] A
89*> \verbatim
90*> A is DOUBLE PRECISION array, dimension (NMAX*NMAX)
91*> \endverbatim
92*>
93*> \param[out] AFAC
94*> \verbatim
95*> AFAC is DOUBLE PRECISION array, dimension (NMAX*NMAX)
96*> \endverbatim
97*>
98*> \param[out] ASAV
99*> \verbatim
100*> ASAV is DOUBLE PRECISION array, dimension (NMAX*NMAX)
101*> \endverbatim
102*>
103*> \param[out] B
104*> \verbatim
105*> B is DOUBLE PRECISION array, dimension (NMAX*NRHS)
106*> \endverbatim
107*>
108*> \param[out] BSAV
109*> \verbatim
110*> BSAV is DOUBLE PRECISION array, dimension (NMAX*NRHS)
111*> \endverbatim
112*>
113*> \param[out] X
114*> \verbatim
115*> X is DOUBLE PRECISION array, dimension (NMAX*NRHS)
116*> \endverbatim
117*>
118*> \param[out] XACT
119*> \verbatim
120*> XACT is DOUBLE PRECISION array, dimension (NMAX*NRHS)
121*> \endverbatim
122*>
123*> \param[out] S
124*> \verbatim
125*> S is DOUBLE PRECISION array, dimension (NMAX)
126*> \endverbatim
127*>
128*> \param[out] WORK
129*> \verbatim
130*> WORK is DOUBLE PRECISION array, dimension
131*> (NMAX*max(3,NRHS))
132*> \endverbatim
133*>
134*> \param[out] RWORK
135*> \verbatim
136*> RWORK is DOUBLE PRECISION array, dimension (NMAX+2*NRHS)
137*> \endverbatim
138*>
139*> \param[out] IWORK
140*> \verbatim
141*> IWORK is INTEGER array, dimension (NMAX)
142*> \endverbatim
143*>
144*> \param[in] NOUT
145*> \verbatim
146*> NOUT is INTEGER
147*> The unit number for output.
148*> \endverbatim
149*
150* Authors:
151* ========
152*
153*> \author Univ. of Tennessee
154*> \author Univ. of California Berkeley
155*> \author Univ. of Colorado Denver
156*> \author NAG Ltd.
157*
158*> \ingroup double_lin
159*
160* =====================================================================
161 SUBROUTINE ddrvpb( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
162 $ A, AFAC, ASAV, B, BSAV, X, XACT, S, WORK,
163 $ RWORK, IWORK, NOUT )
164*
165* -- LAPACK test routine --
166* -- LAPACK is a software package provided by Univ. of Tennessee, --
167* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
168*
169* .. Scalar Arguments ..
170 LOGICAL TSTERR
171 INTEGER NMAX, NN, NOUT, NRHS
172 DOUBLE PRECISION THRESH
173* ..
174* .. Array Arguments ..
175 LOGICAL DOTYPE( * )
176 INTEGER IWORK( * ), NVAL( * )
177 DOUBLE PRECISION A( * ), AFAC( * ), ASAV( * ), B( * ),
178 $ bsav( * ), rwork( * ), s( * ), work( * ),
179 $ x( * ), xact( * )
180* ..
181*
182* =====================================================================
183*
184* .. Parameters ..
185 DOUBLE PRECISION ONE, ZERO
186 PARAMETER ( ONE = 1.0d+0, zero = 0.0d+0 )
187 INTEGER NTYPES, NTESTS
188 parameter( ntypes = 8, ntests = 6 )
189 INTEGER NBW
190 parameter( nbw = 4 )
191* ..
192* .. Local Scalars ..
193 LOGICAL EQUIL, NOFACT, PREFAC, ZEROT
194 CHARACTER DIST, EQUED, FACT, PACKIT, TYPE, UPLO, XTYPE
195 CHARACTER*3 PATH
196 INTEGER I, I1, I2, IEQUED, IFACT, IKD, IMAT, IN, INFO,
197 $ ioff, iuplo, iw, izero, k, k1, kd, kl, koff,
198 $ ku, lda, ldab, mode, n, nb, nbmin, nerrs,
199 $ nfact, nfail, nimat, nkd, nrun, nt
200 DOUBLE PRECISION AINVNM, AMAX, ANORM, CNDNUM, RCOND, RCONDC,
201 $ ROLDC, SCOND
202* ..
203* .. Local Arrays ..
204 CHARACTER EQUEDS( 2 ), FACTS( 3 )
205 INTEGER ISEED( 4 ), ISEEDY( 4 ), KDVAL( NBW )
206 DOUBLE PRECISION RESULT( NTESTS )
207* ..
208* .. External Functions ..
209 LOGICAL LSAME
210 DOUBLE PRECISION DGET06, DLANGE, DLANSB
211 EXTERNAL lsame, dget06, dlange, dlansb
212* ..
213* .. External Subroutines ..
214 EXTERNAL aladhd, alaerh, alasvm, dcopy, derrvx, dget04,
218* ..
219* .. Intrinsic Functions ..
220 INTRINSIC max, min
221* ..
222* .. Scalars in Common ..
223 LOGICAL LERR, OK
224 CHARACTER*32 SRNAMT
225 INTEGER INFOT, NUNIT
226* ..
227* .. Common blocks ..
228 COMMON / infoc / infot, nunit, ok, lerr
229 COMMON / srnamc / srnamt
230* ..
231* .. Data statements ..
232 DATA iseedy / 1988, 1989, 1990, 1991 /
233 DATA facts / 'F', 'N', 'E' /
234 DATA equeds / 'N', 'Y' /
235* ..
236* .. Executable Statements ..
237*
238* Initialize constants and the random number seed.
239*
240 path( 1: 1 ) = 'Double precision'
241 path( 2: 3 ) = 'PB'
242 nrun = 0
243 nfail = 0
244 nerrs = 0
245 DO 10 i = 1, 4
246 iseed( i ) = iseedy( i )
247 10 CONTINUE
248*
249* Test the error exits
250*
251 IF( tsterr )
252 $ CALL derrvx( path, nout )
253 infot = 0
254 kdval( 1 ) = 0
255*
256* Set the block size and minimum block size for testing.
257*
258 nb = 1
259 nbmin = 2
260 CALL xlaenv( 1, nb )
261 CALL xlaenv( 2, nbmin )
262*
263* Do for each value of N in NVAL
264*
265 DO 110 in = 1, nn
266 n = nval( in )
267 lda = max( n, 1 )
268 xtype = 'N'
269*
270* Set limits on the number of loop iterations.
271*
272 nkd = max( 1, min( n, 4 ) )
273 nimat = ntypes
274 IF( n.EQ.0 )
275 $ nimat = 1
276*
277 kdval( 2 ) = n + ( n+1 ) / 4
278 kdval( 3 ) = ( 3*n-1 ) / 4
279 kdval( 4 ) = ( n+1 ) / 4
280*
281 DO 100 ikd = 1, nkd
282*
283* Do for KD = 0, (5*N+1)/4, (3N-1)/4, and (N+1)/4. This order
284* makes it easier to skip redundant values for small values
285* of N.
286*
287 kd = kdval( ikd )
288 ldab = kd + 1
289*
290* Do first for UPLO = 'U', then for UPLO = 'L'
291*
292 DO 90 iuplo = 1, 2
293 koff = 1
294 IF( iuplo.EQ.1 ) THEN
295 uplo = 'U'
296 packit = 'Q'
297 koff = max( 1, kd+2-n )
298 ELSE
299 uplo = 'L'
300 packit = 'B'
301 END IF
302*
303 DO 80 imat = 1, nimat
304*
305* Do the tests only if DOTYPE( IMAT ) is true.
306*
307 IF( .NOT.dotype( imat ) )
308 $ GO TO 80
309*
310* Skip types 2, 3, or 4 if the matrix size is too small.
311*
312 zerot = imat.GE.2 .AND. imat.LE.4
313 IF( zerot .AND. n.LT.imat-1 )
314 $ GO TO 80
315*
316 IF( .NOT.zerot .OR. .NOT.dotype( 1 ) ) THEN
317*
318* Set up parameters with DLATB4 and generate a test
319* matrix with DLATMS.
320*
321 CALL dlatb4( path, imat, n, n, TYPE, kl, ku, anorm,
322 $ mode, cndnum, dist )
323*
324 srnamt = 'DLATMS'
325 CALL dlatms( n, n, dist, iseed, TYPE, rwork, mode,
326 $ cndnum, anorm, kd, kd, packit,
327 $ a( koff ), ldab, work, info )
328*
329* Check error code from DLATMS.
330*
331 IF( info.NE.0 ) THEN
332 CALL alaerh( path, 'DLATMS', info, 0, uplo, n,
333 $ n, -1, -1, -1, imat, nfail, nerrs,
334 $ nout )
335 GO TO 80
336 END IF
337 ELSE IF( izero.GT.0 ) THEN
338*
339* Use the same matrix for types 3 and 4 as for type
340* 2 by copying back the zeroed out column,
341*
342 iw = 2*lda + 1
343 IF( iuplo.EQ.1 ) THEN
344 ioff = ( izero-1 )*ldab + kd + 1
345 CALL dcopy( izero-i1, work( iw ), 1,
346 $ a( ioff-izero+i1 ), 1 )
347 iw = iw + izero - i1
348 CALL dcopy( i2-izero+1, work( iw ), 1,
349 $ a( ioff ), max( ldab-1, 1 ) )
350 ELSE
351 ioff = ( i1-1 )*ldab + 1
352 CALL dcopy( izero-i1, work( iw ), 1,
353 $ a( ioff+izero-i1 ),
354 $ max( ldab-1, 1 ) )
355 ioff = ( izero-1 )*ldab + 1
356 iw = iw + izero - i1
357 CALL dcopy( i2-izero+1, work( iw ), 1,
358 $ a( ioff ), 1 )
359 END IF
360 END IF
361*
362* For types 2-4, zero one row and column of the matrix
363* to test that INFO is returned correctly.
364*
365 izero = 0
366 IF( zerot ) THEN
367 IF( imat.EQ.2 ) THEN
368 izero = 1
369 ELSE IF( imat.EQ.3 ) THEN
370 izero = n
371 ELSE
372 izero = n / 2 + 1
373 END IF
374*
375* Save the zeroed out row and column in WORK(*,3)
376*
377 iw = 2*lda
378 DO 20 i = 1, min( 2*kd+1, n )
379 work( iw+i ) = zero
380 20 CONTINUE
381 iw = iw + 1
382 i1 = max( izero-kd, 1 )
383 i2 = min( izero+kd, n )
384*
385 IF( iuplo.EQ.1 ) THEN
386 ioff = ( izero-1 )*ldab + kd + 1
387 CALL dswap( izero-i1, a( ioff-izero+i1 ), 1,
388 $ work( iw ), 1 )
389 iw = iw + izero - i1
390 CALL dswap( i2-izero+1, a( ioff ),
391 $ max( ldab-1, 1 ), work( iw ), 1 )
392 ELSE
393 ioff = ( i1-1 )*ldab + 1
394 CALL dswap( izero-i1, a( ioff+izero-i1 ),
395 $ max( ldab-1, 1 ), work( iw ), 1 )
396 ioff = ( izero-1 )*ldab + 1
397 iw = iw + izero - i1
398 CALL dswap( i2-izero+1, a( ioff ), 1,
399 $ work( iw ), 1 )
400 END IF
401 END IF
402*
403* Save a copy of the matrix A in ASAV.
404*
405 CALL dlacpy( 'Full', kd+1, n, a, ldab, asav, ldab )
406*
407 DO 70 iequed = 1, 2
408 equed = equeds( iequed )
409 IF( iequed.EQ.1 ) THEN
410 nfact = 3
411 ELSE
412 nfact = 1
413 END IF
414*
415 DO 60 ifact = 1, nfact
416 fact = facts( ifact )
417 prefac = lsame( fact, 'F' )
418 nofact = lsame( fact, 'N' )
419 equil = lsame( fact, 'E' )
420*
421 IF( zerot ) THEN
422 IF( prefac )
423 $ GO TO 60
424 rcondc = zero
425*
426 ELSE IF( .NOT.lsame( fact, 'N' ) ) THEN
427*
428* Compute the condition number for comparison
429* with the value returned by DPBSVX (FACT =
430* 'N' reuses the condition number from the
431* previous iteration with FACT = 'F').
432*
433 CALL dlacpy( 'Full', kd+1, n, asav, ldab,
434 $ afac, ldab )
435 IF( equil .OR. iequed.GT.1 ) THEN
436*
437* Compute row and column scale factors to
438* equilibrate the matrix A.
439*
440 CALL dpbequ( uplo, n, kd, afac, ldab, s,
441 $ scond, amax, info )
442 IF( info.EQ.0 .AND. n.GT.0 ) THEN
443 IF( iequed.GT.1 )
444 $ scond = zero
445*
446* Equilibrate the matrix.
447*
448 CALL dlaqsb( uplo, n, kd, afac, ldab,
449 $ s, scond, amax, equed )
450 END IF
451 END IF
452*
453* Save the condition number of the
454* non-equilibrated system for use in DGET04.
455*
456 IF( equil )
457 $ roldc = rcondc
458*
459* Compute the 1-norm of A.
460*
461 anorm = dlansb( '1', uplo, n, kd, afac, ldab,
462 $ rwork )
463*
464* Factor the matrix A.
465*
466 CALL dpbtrf( uplo, n, kd, afac, ldab, info )
467*
468* Form the inverse of A.
469*
470 CALL dlaset( 'Full', n, n, zero, one, a,
471 $ lda )
472 srnamt = 'DPBTRS'
473 CALL dpbtrs( uplo, n, kd, n, afac, ldab, a,
474 $ lda, info )
475*
476* Compute the 1-norm condition number of A.
477*
478 ainvnm = dlange( '1', n, n, a, lda, rwork )
479 IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
480 rcondc = one
481 ELSE
482 rcondc = ( one / anorm ) / ainvnm
483 END IF
484 END IF
485*
486* Restore the matrix A.
487*
488 CALL dlacpy( 'Full', kd+1, n, asav, ldab, a,
489 $ ldab )
490*
491* Form an exact solution and set the right hand
492* side.
493*
494 srnamt = 'DLARHS'
495 CALL dlarhs( path, xtype, uplo, ' ', n, n, kd,
496 $ kd, nrhs, a, ldab, xact, lda, b,
497 $ lda, iseed, info )
498 xtype = 'C'
499 CALL dlacpy( 'Full', n, nrhs, b, lda, bsav,
500 $ lda )
501*
502 IF( nofact ) THEN
503*
504* --- Test DPBSV ---
505*
506* Compute the L*L' or U'*U factorization of the
507* matrix and solve the system.
508*
509 CALL dlacpy( 'Full', kd+1, n, a, ldab, afac,
510 $ ldab )
511 CALL dlacpy( 'Full', n, nrhs, b, lda, x,
512 $ lda )
513*
514 srnamt = 'DPBSV '
515 CALL dpbsv( uplo, n, kd, nrhs, afac, ldab, x,
516 $ lda, info )
517*
518* Check error code from DPBSV .
519*
520 IF( info.NE.izero ) THEN
521 CALL alaerh( path, 'DPBSV ', info, izero,
522 $ uplo, n, n, kd, kd, nrhs,
523 $ imat, nfail, nerrs, nout )
524 GO TO 40
525 ELSE IF( info.NE.0 ) THEN
526 GO TO 40
527 END IF
528*
529* Reconstruct matrix from factors and compute
530* residual.
531*
532 CALL dpbt01( uplo, n, kd, a, ldab, afac,
533 $ ldab, rwork, result( 1 ) )
534*
535* Compute residual of the computed solution.
536*
537 CALL dlacpy( 'Full', n, nrhs, b, lda, work,
538 $ lda )
539 CALL dpbt02( uplo, n, kd, nrhs, a, ldab, x,
540 $ lda, work, lda, rwork,
541 $ result( 2 ) )
542*
543* Check solution from generated exact solution.
544*
545 CALL dget04( n, nrhs, x, lda, xact, lda,
546 $ rcondc, result( 3 ) )
547 nt = 3
548*
549* Print information about the tests that did
550* not pass the threshold.
551*
552 DO 30 k = 1, nt
553 IF( result( k ).GE.thresh ) THEN
554 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
555 $ CALL aladhd( nout, path )
556 WRITE( nout, fmt = 9999 )'DPBSV ',
557 $ uplo, n, kd, imat, k, result( k )
558 nfail = nfail + 1
559 END IF
560 30 CONTINUE
561 nrun = nrun + nt
562 40 CONTINUE
563 END IF
564*
565* --- Test DPBSVX ---
566*
567 IF( .NOT.prefac )
568 $ CALL dlaset( 'Full', kd+1, n, zero, zero,
569 $ afac, ldab )
570 CALL dlaset( 'Full', n, nrhs, zero, zero, x,
571 $ lda )
572 IF( iequed.GT.1 .AND. n.GT.0 ) THEN
573*
574* Equilibrate the matrix if FACT='F' and
575* EQUED='Y'
576*
577 CALL dlaqsb( uplo, n, kd, a, ldab, s, scond,
578 $ amax, equed )
579 END IF
580*
581* Solve the system and compute the condition
582* number and error bounds using DPBSVX.
583*
584 srnamt = 'DPBSVX'
585 CALL dpbsvx( fact, uplo, n, kd, nrhs, a, ldab,
586 $ afac, ldab, equed, s, b, lda, x,
587 $ lda, rcond, rwork, rwork( nrhs+1 ),
588 $ work, iwork, info )
589*
590* Check the error code from DPBSVX.
591*
592 IF( info.NE.izero ) THEN
593 CALL alaerh( path, 'DPBSVX', info, izero,
594 $ fact // uplo, n, n, kd, kd,
595 $ nrhs, imat, nfail, nerrs, nout )
596 GO TO 60
597 END IF
598*
599 IF( info.EQ.0 ) THEN
600 IF( .NOT.prefac ) THEN
601*
602* Reconstruct matrix from factors and
603* compute residual.
604*
605 CALL dpbt01( uplo, n, kd, a, ldab, afac,
606 $ ldab, rwork( 2*nrhs+1 ),
607 $ result( 1 ) )
608 k1 = 1
609 ELSE
610 k1 = 2
611 END IF
612*
613* Compute residual of the computed solution.
614*
615 CALL dlacpy( 'Full', n, nrhs, bsav, lda,
616 $ work, lda )
617 CALL dpbt02( uplo, n, kd, nrhs, asav, ldab,
618 $ x, lda, work, lda,
619 $ rwork( 2*nrhs+1 ), result( 2 ) )
620*
621* Check solution from generated exact solution.
622*
623 IF( nofact .OR. ( prefac .AND. lsame( equed,
624 $ 'N' ) ) ) THEN
625 CALL dget04( n, nrhs, x, lda, xact, lda,
626 $ rcondc, result( 3 ) )
627 ELSE
628 CALL dget04( n, nrhs, x, lda, xact, lda,
629 $ roldc, result( 3 ) )
630 END IF
631*
632* Check the error bounds from iterative
633* refinement.
634*
635 CALL dpbt05( uplo, n, kd, nrhs, asav, ldab,
636 $ b, lda, x, lda, xact, lda,
637 $ rwork, rwork( nrhs+1 ),
638 $ result( 4 ) )
639 ELSE
640 k1 = 6
641 END IF
642*
643* Compare RCOND from DPBSVX with the computed
644* value in RCONDC.
645*
646 result( 6 ) = dget06( rcond, rcondc )
647*
648* Print information about the tests that did not
649* pass the threshold.
650*
651 DO 50 k = k1, 6
652 IF( result( k ).GE.thresh ) THEN
653 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
654 $ CALL aladhd( nout, path )
655 IF( prefac ) THEN
656 WRITE( nout, fmt = 9997 )'DPBSVX',
657 $ fact, uplo, n, kd, equed, imat, k,
658 $ result( k )
659 ELSE
660 WRITE( nout, fmt = 9998 )'DPBSVX',
661 $ fact, uplo, n, kd, imat, k,
662 $ result( k )
663 END IF
664 nfail = nfail + 1
665 END IF
666 50 CONTINUE
667 nrun = nrun + 7 - k1
668 60 CONTINUE
669 70 CONTINUE
670 80 CONTINUE
671 90 CONTINUE
672 100 CONTINUE
673 110 CONTINUE
674*
675* Print a summary of the results.
676*
677 CALL alasvm( path, nout, nfail, nrun, nerrs )
678*
679 9999 FORMAT( 1x, a, ', UPLO=''', a1, ''', N =', i5, ', KD =', i5,
680 $ ', type ', i1, ', test(', i1, ')=', g12.5 )
681 9998 FORMAT( 1x, a, '( ''', a1, ''', ''', a1, ''', ', i5, ', ', i5,
682 $ ', ... ), type ', i1, ', test(', i1, ')=', g12.5 )
683 9997 FORMAT( 1x, a, '( ''', a1, ''', ''', a1, ''', ', i5, ', ', i5,
684 $ ', ... ), EQUED=''', a1, ''', type ', i1, ', test(', i1,
685 $ ')=', g12.5 )
686 RETURN
687*
688* End of DDRVPB
689*
690 END
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
Definition alasvm.f:73
subroutine dlarhs(path, xtype, uplo, trans, m, n, kl, ku, nrhs, a, lda, x, ldx, b, ldb, iseed, info)
DLARHS
Definition dlarhs.f:205
subroutine xlaenv(ispec, nvalue)
XLAENV
Definition xlaenv.f:81
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 ddrvpb(dotype, nn, nval, nrhs, thresh, tsterr, nmax, a, afac, asav, b, bsav, x, xact, s, work, rwork, iwork, nout)
DDRVPB
Definition ddrvpb.f:164
subroutine derrvx(path, nunit)
DERRVX
Definition derrvx.f:55
subroutine dget04(n, nrhs, x, ldx, xact, ldxact, rcond, resid)
DGET04
Definition dget04.f:102
subroutine dlatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
DLATB4
Definition dlatb4.f:120
subroutine dlatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
DLATMS
Definition dlatms.f:321
subroutine dpbt01(uplo, n, kd, a, lda, afac, ldafac, rwork, resid)
DPBT01
Definition dpbt01.f:119
subroutine dpbt02(uplo, n, kd, nrhs, a, lda, x, ldx, b, ldb, rwork, resid)
DPBT02
Definition dpbt02.f:136
subroutine dpbt05(uplo, n, kd, nrhs, ab, ldab, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
DPBT05
Definition dpbt05.f:171
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
Definition dcopy.f:82
subroutine dlacpy(uplo, m, n, a, lda, b, ldb)
DLACPY copies all or part of one two-dimensional array to another.
Definition dlacpy.f:103
subroutine dlaqsb(uplo, n, kd, ab, ldab, s, scond, amax, equed)
DLAQSB scales a symmetric/Hermitian band matrix, using scaling factors computed by spbequ.
Definition dlaqsb.f:140
subroutine dlaset(uplo, m, n, alpha, beta, a, lda)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition dlaset.f:110
subroutine dpbequ(uplo, n, kd, ab, ldab, s, scond, amax, info)
DPBEQU
Definition dpbequ.f:129
subroutine dpbsv(uplo, n, kd, nrhs, ab, ldab, b, ldb, info)
DPBSV computes the solution to system of linear equations A * X = B for OTHER matrices
Definition dpbsv.f:164
subroutine dpbsvx(fact, uplo, n, kd, nrhs, ab, ldab, afb, ldafb, equed, s, b, ldb, x, ldx, rcond, ferr, berr, work, iwork, info)
DPBSVX computes the solution to system of linear equations A * X = B for OTHER matrices
Definition dpbsvx.f:343
subroutine dpbtrf(uplo, n, kd, ab, ldab, info)
DPBTRF
Definition dpbtrf.f:142
subroutine dpbtrs(uplo, n, kd, nrhs, ab, ldab, b, ldb, info)
DPBTRS
Definition dpbtrs.f:121
subroutine dswap(n, dx, incx, dy, incy)
DSWAP
Definition dswap.f:82