LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
cdrvgt.f
Go to the documentation of this file.
1 *> \brief \b CDRVGT
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 CDRVGT( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, A, AF,
12 * B, X, XACT, WORK, RWORK, IWORK, NOUT )
13 *
14 * .. Scalar Arguments ..
15 * LOGICAL TSTERR
16 * INTEGER NN, NOUT, NRHS
17 * REAL THRESH
18 * ..
19 * .. Array Arguments ..
20 * LOGICAL DOTYPE( * )
21 * INTEGER IWORK( * ), NVAL( * )
22 * REAL RWORK( * )
23 * COMPLEX A( * ), AF( * ), B( * ), WORK( * ), X( * ),
24 * $ XACT( * )
25 * ..
26 *
27 *
28 *> \par Purpose:
29 * =============
30 *>
31 *> \verbatim
32 *>
33 *> CDRVGT tests CGTSV 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 sides, NRHS >= 0.
63 *> \endverbatim
64 *>
65 *> \param[in] THRESH
66 *> \verbatim
67 *> THRESH is REAL
68 *> The threshold value for the test ratios. A result is
69 *> included in the output file if RESULT >= THRESH. To have
70 *> every test ratio printed, use THRESH = 0.
71 *> \endverbatim
72 *>
73 *> \param[in] TSTERR
74 *> \verbatim
75 *> TSTERR is LOGICAL
76 *> Flag that indicates whether error exits are to be tested.
77 *> \endverbatim
78 *>
79 *> \param[out] A
80 *> \verbatim
81 *> A is COMPLEX array, dimension (NMAX*4)
82 *> \endverbatim
83 *>
84 *> \param[out] AF
85 *> \verbatim
86 *> AF is COMPLEX array, dimension (NMAX*4)
87 *> \endverbatim
88 *>
89 *> \param[out] B
90 *> \verbatim
91 *> B is COMPLEX array, dimension (NMAX*NRHS)
92 *> \endverbatim
93 *>
94 *> \param[out] X
95 *> \verbatim
96 *> X is COMPLEX array, dimension (NMAX*NRHS)
97 *> \endverbatim
98 *>
99 *> \param[out] XACT
100 *> \verbatim
101 *> XACT is COMPLEX array, dimension (NMAX*NRHS)
102 *> \endverbatim
103 *>
104 *> \param[out] WORK
105 *> \verbatim
106 *> WORK is COMPLEX array, dimension
107 *> (NMAX*max(3,NRHS))
108 *> \endverbatim
109 *>
110 *> \param[out] RWORK
111 *> \verbatim
112 *> RWORK is REAL array, dimension (NMAX+2*NRHS)
113 *> \endverbatim
114 *>
115 *> \param[out] IWORK
116 *> \verbatim
117 *> IWORK is INTEGER array, dimension (2*NMAX)
118 *> \endverbatim
119 *>
120 *> \param[in] NOUT
121 *> \verbatim
122 *> NOUT is INTEGER
123 *> The unit number for output.
124 *> \endverbatim
125 *
126 * Authors:
127 * ========
128 *
129 *> \author Univ. of Tennessee
130 *> \author Univ. of California Berkeley
131 *> \author Univ. of Colorado Denver
132 *> \author NAG Ltd.
133 *
134 *> \date November 2011
135 *
136 *> \ingroup complex_lin
137 *
138 * =====================================================================
139  SUBROUTINE cdrvgt( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, A, AF,
140  $ b, x, xact, work, rwork, iwork, nout )
141 *
142 * -- LAPACK test routine (version 3.4.0) --
143 * -- LAPACK is a software package provided by Univ. of Tennessee, --
144 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
145 * November 2011
146 *
147 * .. Scalar Arguments ..
148  LOGICAL tsterr
149  INTEGER nn, nout, nrhs
150  REAL thresh
151 * ..
152 * .. Array Arguments ..
153  LOGICAL dotype( * )
154  INTEGER iwork( * ), nval( * )
155  REAL rwork( * )
156  COMPLEX a( * ), af( * ), b( * ), work( * ), x( * ),
157  $ xact( * )
158 * ..
159 *
160 * =====================================================================
161 *
162 * .. Parameters ..
163  REAL one, zero
164  parameter( one = 1.0e+0, zero = 0.0e+0 )
165  INTEGER ntypes
166  parameter( ntypes = 12 )
167  INTEGER ntests
168  parameter( ntests = 6 )
169 * ..
170 * .. Local Scalars ..
171  LOGICAL trfcon, zerot
172  CHARACTER dist, fact, trans, type
173  CHARACTER*3 path
174  INTEGER i, ifact, imat, in, info, itran, ix, izero, j,
175  $ k, k1, kl, koff, ku, lda, m, mode, n, nerrs,
176  $ nfail, nimat, nrun, nt
177  REAL ainvnm, anorm, anormi, anormo, cond, rcond,
178  $ rcondc, rcondi, rcondo
179 * ..
180 * .. Local Arrays ..
181  CHARACTER transs( 3 )
182  INTEGER iseed( 4 ), iseedy( 4 )
183  REAL result( ntests ), z( 3 )
184 * ..
185 * .. External Functions ..
186  REAL clangt, scasum, sget06
187  EXTERNAL clangt, scasum, sget06
188 * ..
189 * .. External Subroutines ..
190  EXTERNAL aladhd, alaerh, alasvm, ccopy, cerrvx, cget04,
191  $ cgtsv, cgtsvx, cgtt01, cgtt02, cgtt05, cgttrf,
192  $ cgttrs, clacpy, clagtm, clarnv, claset, clatb4,
193  $ clatms, csscal
194 * ..
195 * .. Intrinsic Functions ..
196  INTRINSIC cmplx, max
197 * ..
198 * .. Scalars in Common ..
199  LOGICAL lerr, ok
200  CHARACTER*32 srnamt
201  INTEGER infot, nunit
202 * ..
203 * .. Common blocks ..
204  common / infoc / infot, nunit, ok, lerr
205  common / srnamc / srnamt
206 * ..
207 * .. Data statements ..
208  DATA iseedy / 0, 0, 0, 1 / , transs / 'N', 'T',
209  $ 'C' /
210 * ..
211 * .. Executable Statements ..
212 *
213  path( 1: 1 ) = 'Complex precision'
214  path( 2: 3 ) = 'GT'
215  nrun = 0
216  nfail = 0
217  nerrs = 0
218  DO 10 i = 1, 4
219  iseed( i ) = iseedy( i )
220  10 continue
221 *
222 * Test the error exits
223 *
224  IF( tsterr )
225  $ CALL cerrvx( path, nout )
226  infot = 0
227 *
228  DO 140 in = 1, nn
229 *
230 * Do for each value of N in NVAL.
231 *
232  n = nval( in )
233  m = max( n-1, 0 )
234  lda = max( 1, n )
235  nimat = ntypes
236  IF( n.LE.0 )
237  $ nimat = 1
238 *
239  DO 130 imat = 1, nimat
240 *
241 * Do the tests only if DOTYPE( IMAT ) is true.
242 *
243  IF( .NOT.dotype( imat ) )
244  $ go to 130
245 *
246 * Set up parameters with CLATB4.
247 *
248  CALL clatb4( path, imat, n, n, type, kl, ku, anorm, mode,
249  $ cond, dist )
250 *
251  zerot = imat.GE.8 .AND. imat.LE.10
252  IF( imat.LE.6 ) THEN
253 *
254 * Types 1-6: generate matrices of known condition number.
255 *
256  koff = max( 2-ku, 3-max( 1, n ) )
257  srnamt = 'CLATMS'
258  CALL clatms( n, n, dist, iseed, type, rwork, mode, cond,
259  $ anorm, kl, ku, 'Z', af( koff ), 3, work,
260  $ info )
261 *
262 * Check the error code from CLATMS.
263 *
264  IF( info.NE.0 ) THEN
265  CALL alaerh( path, 'CLATMS', info, 0, ' ', n, n, kl,
266  $ ku, -1, imat, nfail, nerrs, nout )
267  go to 130
268  END IF
269  izero = 0
270 *
271  IF( n.GT.1 ) THEN
272  CALL ccopy( n-1, af( 4 ), 3, a, 1 )
273  CALL ccopy( n-1, af( 3 ), 3, a( n+m+1 ), 1 )
274  END IF
275  CALL ccopy( n, af( 2 ), 3, a( m+1 ), 1 )
276  ELSE
277 *
278 * Types 7-12: generate tridiagonal matrices with
279 * unknown condition numbers.
280 *
281  IF( .NOT.zerot .OR. .NOT.dotype( 7 ) ) THEN
282 *
283 * Generate a matrix with elements from [-1,1].
284 *
285  CALL clarnv( 2, iseed, n+2*m, a )
286  IF( anorm.NE.one )
287  $ CALL csscal( n+2*m, anorm, a, 1 )
288  ELSE IF( izero.GT.0 ) THEN
289 *
290 * Reuse the last matrix by copying back the zeroed out
291 * elements.
292 *
293  IF( izero.EQ.1 ) THEN
294  a( n ) = z( 2 )
295  IF( n.GT.1 )
296  $ a( 1 ) = z( 3 )
297  ELSE IF( izero.EQ.n ) THEN
298  a( 3*n-2 ) = z( 1 )
299  a( 2*n-1 ) = z( 2 )
300  ELSE
301  a( 2*n-2+izero ) = z( 1 )
302  a( n-1+izero ) = z( 2 )
303  a( izero ) = z( 3 )
304  END IF
305  END IF
306 *
307 * If IMAT > 7, set one column of the matrix to 0.
308 *
309  IF( .NOT.zerot ) THEN
310  izero = 0
311  ELSE IF( imat.EQ.8 ) THEN
312  izero = 1
313  z( 2 ) = a( n )
314  a( n ) = zero
315  IF( n.GT.1 ) THEN
316  z( 3 ) = a( 1 )
317  a( 1 ) = zero
318  END IF
319  ELSE IF( imat.EQ.9 ) THEN
320  izero = n
321  z( 1 ) = a( 3*n-2 )
322  z( 2 ) = a( 2*n-1 )
323  a( 3*n-2 ) = zero
324  a( 2*n-1 ) = zero
325  ELSE
326  izero = ( n+1 ) / 2
327  DO 20 i = izero, n - 1
328  a( 2*n-2+i ) = zero
329  a( n-1+i ) = zero
330  a( i ) = zero
331  20 continue
332  a( 3*n-2 ) = zero
333  a( 2*n-1 ) = zero
334  END IF
335  END IF
336 *
337  DO 120 ifact = 1, 2
338  IF( ifact.EQ.1 ) THEN
339  fact = 'F'
340  ELSE
341  fact = 'N'
342  END IF
343 *
344 * Compute the condition number for comparison with
345 * the value returned by CGTSVX.
346 *
347  IF( zerot ) THEN
348  IF( ifact.EQ.1 )
349  $ go to 120
350  rcondo = zero
351  rcondi = zero
352 *
353  ELSE IF( ifact.EQ.1 ) THEN
354  CALL ccopy( n+2*m, a, 1, af, 1 )
355 *
356 * Compute the 1-norm and infinity-norm of A.
357 *
358  anormo = clangt( '1', n, a, a( m+1 ), a( n+m+1 ) )
359  anormi = clangt( 'I', n, a, a( m+1 ), a( n+m+1 ) )
360 *
361 * Factor the matrix A.
362 *
363  CALL cgttrf( n, af, af( m+1 ), af( n+m+1 ),
364  $ af( n+2*m+1 ), iwork, info )
365 *
366 * Use CGTTRS to solve for one column at a time of
367 * inv(A), computing the maximum column sum as we go.
368 *
369  ainvnm = zero
370  DO 40 i = 1, n
371  DO 30 j = 1, n
372  x( j ) = zero
373  30 continue
374  x( i ) = one
375  CALL cgttrs( 'No transpose', n, 1, af, af( m+1 ),
376  $ af( n+m+1 ), af( n+2*m+1 ), iwork, x,
377  $ lda, info )
378  ainvnm = max( ainvnm, scasum( n, x, 1 ) )
379  40 continue
380 *
381 * Compute the 1-norm condition number of A.
382 *
383  IF( anormo.LE.zero .OR. ainvnm.LE.zero ) THEN
384  rcondo = one
385  ELSE
386  rcondo = ( one / anormo ) / ainvnm
387  END IF
388 *
389 * Use CGTTRS to solve for one column at a time of
390 * inv(A'), computing the maximum column sum as we go.
391 *
392  ainvnm = zero
393  DO 60 i = 1, n
394  DO 50 j = 1, n
395  x( j ) = zero
396  50 continue
397  x( i ) = one
398  CALL cgttrs( 'Conjugate transpose', n, 1, af,
399  $ af( m+1 ), af( n+m+1 ), af( n+2*m+1 ),
400  $ iwork, x, lda, info )
401  ainvnm = max( ainvnm, scasum( n, x, 1 ) )
402  60 continue
403 *
404 * Compute the infinity-norm condition number of A.
405 *
406  IF( anormi.LE.zero .OR. ainvnm.LE.zero ) THEN
407  rcondi = one
408  ELSE
409  rcondi = ( one / anormi ) / ainvnm
410  END IF
411  END IF
412 *
413  DO 110 itran = 1, 3
414  trans = transs( itran )
415  IF( itran.EQ.1 ) THEN
416  rcondc = rcondo
417  ELSE
418  rcondc = rcondi
419  END IF
420 *
421 * Generate NRHS random solution vectors.
422 *
423  ix = 1
424  DO 70 j = 1, nrhs
425  CALL clarnv( 2, iseed, n, xact( ix ) )
426  ix = ix + lda
427  70 continue
428 *
429 * Set the right hand side.
430 *
431  CALL clagtm( trans, n, nrhs, one, a, a( m+1 ),
432  $ a( n+m+1 ), xact, lda, zero, b, lda )
433 *
434  IF( ifact.EQ.2 .AND. itran.EQ.1 ) THEN
435 *
436 * --- Test CGTSV ---
437 *
438 * Solve the system using Gaussian elimination with
439 * partial pivoting.
440 *
441  CALL ccopy( n+2*m, a, 1, af, 1 )
442  CALL clacpy( 'Full', n, nrhs, b, lda, x, lda )
443 *
444  srnamt = 'CGTSV '
445  CALL cgtsv( n, nrhs, af, af( m+1 ), af( n+m+1 ), x,
446  $ lda, info )
447 *
448 * Check error code from CGTSV .
449 *
450  IF( info.NE.izero )
451  $ CALL alaerh( path, 'CGTSV ', info, izero, ' ',
452  $ n, n, 1, 1, nrhs, imat, nfail,
453  $ nerrs, nout )
454  nt = 1
455  IF( izero.EQ.0 ) THEN
456 *
457 * Check residual of computed solution.
458 *
459  CALL clacpy( 'Full', n, nrhs, b, lda, work,
460  $ lda )
461  CALL cgtt02( trans, n, nrhs, a, a( m+1 ),
462  $ a( n+m+1 ), x, lda, work, lda,
463  $ result( 2 ) )
464 *
465 * Check solution from generated exact solution.
466 *
467  CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
468  $ result( 3 ) )
469  nt = 3
470  END IF
471 *
472 * Print information about the tests that did not pass
473 * the threshold.
474 *
475  DO 80 k = 2, nt
476  IF( result( k ).GE.thresh ) THEN
477  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
478  $ CALL aladhd( nout, path )
479  WRITE( nout, fmt = 9999 )'CGTSV ', n, imat,
480  $ k, result( k )
481  nfail = nfail + 1
482  END IF
483  80 continue
484  nrun = nrun + nt - 1
485  END IF
486 *
487 * --- Test CGTSVX ---
488 *
489  IF( ifact.GT.1 ) THEN
490 *
491 * Initialize AF to zero.
492 *
493  DO 90 i = 1, 3*n - 2
494  af( i ) = zero
495  90 continue
496  END IF
497  CALL claset( 'Full', n, nrhs, cmplx( zero ),
498  $ cmplx( zero ), x, lda )
499 *
500 * Solve the system and compute the condition number and
501 * error bounds using CGTSVX.
502 *
503  srnamt = 'CGTSVX'
504  CALL cgtsvx( fact, trans, n, nrhs, a, a( m+1 ),
505  $ a( n+m+1 ), af, af( m+1 ), af( n+m+1 ),
506  $ af( n+2*m+1 ), iwork, b, lda, x, lda,
507  $ rcond, rwork, rwork( nrhs+1 ), work,
508  $ rwork( 2*nrhs+1 ), info )
509 *
510 * Check the error code from CGTSVX.
511 *
512  IF( info.NE.izero )
513  $ CALL alaerh( path, 'CGTSVX', info, izero,
514  $ fact // trans, n, n, 1, 1, nrhs, imat,
515  $ nfail, nerrs, nout )
516 *
517  IF( ifact.GE.2 ) THEN
518 *
519 * Reconstruct matrix from factors and compute
520 * residual.
521 *
522  CALL cgtt01( n, a, a( m+1 ), a( n+m+1 ), af,
523  $ af( m+1 ), af( n+m+1 ), af( n+2*m+1 ),
524  $ iwork, work, lda, rwork, result( 1 ) )
525  k1 = 1
526  ELSE
527  k1 = 2
528  END IF
529 *
530  IF( info.EQ.0 ) THEN
531  trfcon = .false.
532 *
533 * Check residual of computed solution.
534 *
535  CALL clacpy( 'Full', n, nrhs, b, lda, work, lda )
536  CALL cgtt02( trans, n, nrhs, a, a( m+1 ),
537  $ a( n+m+1 ), x, lda, work, lda,
538  $ result( 2 ) )
539 *
540 * Check solution from generated exact solution.
541 *
542  CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
543  $ result( 3 ) )
544 *
545 * Check the error bounds from iterative refinement.
546 *
547  CALL cgtt05( trans, n, nrhs, a, a( m+1 ),
548  $ a( n+m+1 ), b, lda, x, lda, xact, lda,
549  $ rwork, rwork( nrhs+1 ), result( 4 ) )
550  nt = 5
551  END IF
552 *
553 * Print information about the tests that did not pass
554 * the threshold.
555 *
556  DO 100 k = k1, nt
557  IF( result( k ).GE.thresh ) THEN
558  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
559  $ CALL aladhd( nout, path )
560  WRITE( nout, fmt = 9998 )'CGTSVX', fact, trans,
561  $ n, imat, k, result( k )
562  nfail = nfail + 1
563  END IF
564  100 continue
565 *
566 * Check the reciprocal of the condition number.
567 *
568  result( 6 ) = sget06( rcond, rcondc )
569  IF( result( 6 ).GE.thresh ) THEN
570  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
571  $ CALL aladhd( nout, path )
572  WRITE( nout, fmt = 9998 )'CGTSVX', fact, trans, n,
573  $ imat, k, result( k )
574  nfail = nfail + 1
575  END IF
576  nrun = nrun + nt - k1 + 2
577 *
578  110 continue
579  120 continue
580  130 continue
581  140 continue
582 *
583 * Print a summary of the results.
584 *
585  CALL alasvm( path, nout, nfail, nrun, nerrs )
586 *
587  9999 format( 1x, a, ', N =', i5, ', type ', i2, ', test ', i2,
588  $ ', ratio = ', g12.5 )
589  9998 format( 1x, a, ', FACT=''', a1, ''', TRANS=''', a1, ''', N =',
590  $ i5, ', type ', i2, ', test ', i2, ', ratio = ', g12.5 )
591  return
592 *
593 * End of CDRVGT
594 *
595  END