PLASMA  2.4.5
PLASMA - Parallel Linear Algebra for Scalable Multi-core Architectures
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups
ddrvge.f
Go to the documentation of this file.
1  SUBROUTINE ddrvge( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
2  $ a, afac, asav, b, bsav, x, xact, s, work,
3  $ rwork, iwork, nout )
4 *
5  include 'plasmaf.h'
6 *
7 * -- LAPACK test routine (version 3.1) --
8 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
9 * November 2006
10 *
11 * .. Scalar Arguments ..
12  LOGICAL tsterr
13  INTEGER nmax, nn, nout, nrhs
14  DOUBLE PRECISION thresh
15 * ..
16 * .. Array Arguments ..
17  LOGICAL dotype( * )
18  INTEGER iwork( * ), nval( * )
19  DOUBLE PRECISION a( * ), afac( * ), asav( * ), b( * ),
20  $ bsav( * ), rwork( * ), s( * ), work( * ),
21  $ x( * ), xact( * )
22 * ..
23 *
24 * Purpose
25 * =======
26 *
27 * DDRVGE tests the driver routines DGESV and -SVX.
28 *
29 * Arguments
30 * =========
31 *
32 * DOTYPE (input) LOGICAL array, dimension (NTYPES)
33 * The matrix types to be used for testing. Matrices of type j
34 * (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
35 * .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
36 *
37 * NN (input) INTEGER
38 * The number of values of N contained in the vector NVAL.
39 *
40 * NVAL (input) INTEGER array, dimension (NN)
41 * The values of the matrix column dimension N.
42 *
43 * NRHS (input) INTEGER
44 * The number of right hand side vectors to be generated for
45 * each linear system.
46 *
47 * THRESH (input) DOUBLE PRECISION
48 * The threshold value for the test ratios. A result is
49 * included in the output file if RESULT >= THRESH. To have
50 * every test ratio printed, use THRESH = 0.
51 *
52 * TSTERR (input) LOGICAL
53 * Flag that indicates whether error exits are to be tested.
54 *
55 * NMAX (input) INTEGER
56 * The maximum value permitted for N, used in dimensioning the
57 * work arrays.
58 *
59 * A (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
60 *
61 * AFAC (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
62 *
63 * ASAV (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
64 *
65 * B (workspace) DOUBLE PRECISION array, dimension (NMAX*NRHS)
66 *
67 * BSAV (workspace) DOUBLE PRECISION array, dimension (NMAX*NRHS)
68 *
69 * X (workspace) DOUBLE PRECISION array, dimension (NMAX*NRHS)
70 *
71 * XACT (workspace) DOUBLE PRECISION array, dimension (NMAX*NRHS)
72 *
73 * S (workspace) DOUBLE PRECISION array, dimension (2*NMAX)
74 *
75 * WORK (workspace) DOUBLE PRECISION array, dimension
76 * (NMAX*max(3,NRHS))
77 *
78 * RWORK (workspace) DOUBLE PRECISION array, dimension (2*NRHS+NMAX)
79 *
80 * IWORK (workspace) INTEGER array, dimension (2*NMAX)
81 *
82 * NOUT (input) INTEGER
83 * The unit number for output.
84 *
85 * =====================================================================
86 *
87 * .. Parameters ..
88  DOUBLE PRECISION one, zero
89  parameter( one = 1.0d+0, zero = 0.0d+0 )
90  INTEGER ntypes
91  parameter( ntypes = 11 )
92  INTEGER ntests
93  parameter( ntests = 7 )
94 * ONLY NOTRANS SUPPORTED !!!
95  INTEGER ntran
96  parameter( ntran = 1 )
97 * ..
98 * .. Local Scalars ..
99  LOGICAL equil, nofact, prefac, trfcon, zerot
100  CHARACTER dist, equed, fact, trans, type, xtype
101  CHARACTER*3 path
102  INTEGER i, iequed, ifact, imat, in, info, ioff, itran,
103  $ izero, k, k1, kl, ku, lda, lwork, mode, n, nb,
104  $ nbmin, nerrs, nfact, nfail, nimat, nrun, nt, ib
105  INTEGER hl( 2 ), hpiv( 2 )
106  INTEGER plasma_trans
107  DOUBLE PRECISION ainvnm, amax, anorm, anormi, anormo, cndnum,
108  $ colcnd, rcond, rcondc, rcondi, rcondo, roldc,
109  $ roldi, roldo, rowcnd, rpvgrw
110 * ..
111 * .. Local Arrays ..
112  CHARACTER equeds( 4 ), facts( 3 ), transs( ntran )
113  INTEGER iseed( 4 ), iseedy( 4 ), plasma_transs( ntran )
114  DOUBLE PRECISION result( ntests )
115 * ..
116 * .. External Functions ..
117  LOGICAL lsame
118  DOUBLE PRECISION dget06, dlamch, dlange, dlantr
119  EXTERNAL lsame, dget06, dlamch, dlange, dlantr
120 * ..
121 * .. External Subroutines ..
122  EXTERNAL aladhd, alaerh, alasvm, derrvx, dgeequ, dgesv,
123  $ dgesvx, dget02, dget04, dgetrf,
124  $ dgetri, dlacpy, dlaqge, dlarhs, dlaset, dlatb4,
125  $ dlatms, xlaenv
126 * ..
127 * .. Intrinsic Functions ..
128  INTRINSIC abs, max
129 * ..
130 * .. Scalars in Common ..
131  LOGICAL lerr, ok
132  CHARACTER*32 srnamt
133  INTEGER infot, nunit
134 * ..
135 * .. Common blocks ..
136  common / infoc / infot, nunit, ok, lerr
137  common / srnamc / srnamt
138 * ..
139 * .. Data statements ..
140  DATA iseedy / 1988, 1989, 1990, 1991 /
141 * DATA TRANSS / 'N', 'T', 'C' /
142  DATA transs / 'N' /
143  DATA plasma_transs / plasmanotrans /
144  DATA facts / 'F', 'N', 'E' /
145  DATA equeds / 'N', 'R', 'C', 'B' /
146 * ..
147 * .. Executable Statements ..
148 *
149 * Initialize constants and the random number seed.
150 *
151  path( 1: 1 ) = 'Double precision'
152  path( 2: 3 ) = 'GE'
153  rcondo = zero
154  rcondi = zero
155  nrun = 0
156  nfail = 0
157  nerrs = 0
158  DO 10 i = 1, 4
159  iseed( i ) = iseedy( i )
160  10 continue
161 *
162 * Test the error exits
163 *
164  IF( tsterr )
165  $ CALL derrvx( path, nout )
166  infot = 0
167 *
168 * Set the block size and minimum block size for testing.
169 *
170  nb = 128
171  ib = 32
172  nbmin = 32
173  CALL xlaenv( 1, nb )
174  CALL xlaenv( 2, nbmin )
175  CALL plasma_set( plasma_tile_size, nb, info )
176  CALL plasma_set( plasma_inner_block_size, ib, info )
177 *
178 * Do for each value of N in NVAL
179 *
180  DO 90 in = 1, nn
181  n = nval( in )
182  lda = max( n, 1 )
183  xtype = 'N'
184  nimat = ntypes
185  IF( n.LE.0 )
186  $ nimat = 1
187 *
188 * ALLOCATE L and IPIV
189 *
190 c$$$ CALL PLASMA_ALLOC_WORKSPACE_DGETRF_INCPIV(
191 c$$$ $ N, N, HL, HPIV, INFO )
192 *
193  DO 80 imat = 1, nimat
194 *
195 * Do the tests only if DOTYPE( IMAT ) is true.
196 *
197  IF( .NOT.dotype( imat ) )
198  $ go to 80
199 *
200 * Skip types 5, 6, or 7 if the matrix size is too small.
201 *
202  zerot = imat.GE.5 .AND. imat.LE.7
203  IF( zerot .AND. n.LT.imat-4 )
204  $ go to 80
205 *
206 * Set up parameters with DLATB4 and generate a test matrix
207 * with DLATMS.
208 *
209  CALL dlatb4( path, imat, n, n, type, kl, ku, anorm, mode,
210  $ cndnum, dist )
211  rcondc = one / cndnum
212 *
213  srnamt = 'DLATMS'
214  CALL dlatms( n, n, dist, iseed, type, rwork, mode, cndnum,
215  $ anorm, kl, ku, 'No packing', a, lda, work,
216  $ info )
217 *
218 * Check error code from DLATMS.
219 *
220  IF( info.NE.0 ) THEN
221  CALL alaerh( path, 'DLATMS', info, 0, ' ', n, n, -1, -1,
222  $ -1, imat, nfail, nerrs, nout )
223  go to 80
224  END IF
225 *
226 * For types 5-7, zero one or more columns of the matrix to
227 * test that INFO is returned correctly.
228 *
229  IF( zerot ) THEN
230  IF( imat.EQ.5 ) THEN
231  izero = 1
232  ELSE IF( imat.EQ.6 ) THEN
233  izero = n
234  ELSE
235  izero = n / 2 + 1
236  END IF
237  ioff = ( izero-1 )*lda
238  IF( imat.LT.7 ) THEN
239  DO 20 i = 1, n
240  a( ioff+i ) = zero
241  20 continue
242  ELSE
243  CALL dlaset( 'Full', n, n-izero+1, zero, zero,
244  $ a( ioff+1 ), lda )
245  END IF
246  ELSE
247  izero = 0
248  END IF
249 *
250 * Save a copy of the matrix A in ASAV.
251 *
252  CALL dlacpy( 'Full', n, n, a, lda, asav, lda )
253 *
254  DO 70 iequed = 1, 4
255  equed = equeds( iequed )
256  IF( iequed.EQ.1 ) THEN
257  nfact = 3
258  ELSE
259  nfact = 1
260  END IF
261 *
262  DO 60 ifact = 1, nfact
263  fact = facts( ifact )
264  prefac = lsame( fact, 'F' )
265  nofact = lsame( fact, 'N' )
266  equil = lsame( fact, 'E' )
267 *
268  IF( zerot ) THEN
269  IF( prefac )
270  $ go to 60
271  rcondo = zero
272  rcondi = zero
273 *
274  ELSE IF( .NOT.nofact ) THEN
275 *
276 * Compute the condition number for comparison with
277 * the value returned by DGESVX (FACT = 'N' reuses
278 * the condition number from the previous iteration
279 * with FACT = 'F').
280 *
281  CALL dlacpy( 'Full', n, n, asav, lda, afac, lda )
282  IF( equil .OR. iequed.GT.1 ) THEN
283 *
284 * Compute row and column scale factors to
285 * equilibrate the matrix A.
286 *
287  CALL dgeequ( n, n, afac, lda, s, s( n+1 ),
288  $ rowcnd, colcnd, amax, info )
289  IF( info.EQ.0 .AND. n.GT.0 ) THEN
290  IF( lsame( equed, 'R' ) ) THEN
291  rowcnd = zero
292  colcnd = one
293  ELSE IF( lsame( equed, 'C' ) ) THEN
294  rowcnd = one
295  colcnd = zero
296  ELSE IF( lsame( equed, 'B' ) ) THEN
297  rowcnd = zero
298  colcnd = zero
299  END IF
300 *
301 * Equilibrate the matrix.
302 *
303  CALL dlaqge( n, n, afac, lda, s, s( n+1 ),
304  $ rowcnd, colcnd, amax, equed )
305  END IF
306  END IF
307 *
308 * Save the condition number of the non-equilibrated
309 * system for use in DGET04.
310 *
311  IF( equil ) THEN
312  roldo = rcondo
313  roldi = rcondi
314  END IF
315 *
316 * Compute the 1-norm and infinity-norm of A.
317 *
318  anormo = dlange( '1', n, n, afac, lda, rwork )
319  anormi = dlange( 'I', n, n, afac, lda, rwork )
320 *
321 * Factor the matrix A.
322 *
323 c$$$ CALL PLASMA_DGETRF_INCPIV( N, N, AFAC, LDA,
324 c$$$ $ HL, HPIV, INFO )
325  CALL plasma_dgetrf( n, n, afac, lda,
326  $ iwork, info )
327  END IF
328 *
329  DO 50 itran = 1, ntran
330 *
331 * Do for each value of TRANS.
332 *
333  trans = transs( itran )
334  plasma_trans = plasma_transs( itran )
335  IF( itran.EQ.1 ) THEN
336  rcondc = rcondo
337  ELSE
338  rcondc = rcondi
339  END IF
340 *
341 * Restore the matrix A.
342 *
343  CALL dlacpy( 'Full', n, n, asav, lda, a, lda )
344 *
345 * Form an exact solution and set the right hand side.
346 *
347  srnamt = 'DLARHS'
348  CALL dlarhs( path, xtype, 'Full', trans, n, n, kl,
349  $ ku, nrhs, a, lda, xact, lda, b, lda,
350  $ iseed, info )
351  xtype = 'C'
352  CALL dlacpy( 'Full', n, nrhs, b, lda, bsav, lda )
353 *
354  IF( nofact .AND. itran.EQ.1 ) THEN
355 *
356 * --- Test DGESV ---
357 *
358 * Compute the LU factorization of the matrix and
359 * solve the system.
360 *
361  CALL dlacpy( 'Full', n, n, a, lda, afac, lda )
362  CALL dlacpy( 'Full', n, nrhs, b, lda, x, lda )
363 *
364  srnamt = 'DGESV '
365 c$$$ CALL PLASMA_DGESV_INCPIV( N, NRHS, AFAC, LDA,
366 c$$$ $ HL, HPIV, X, LDA, INFO )
367  CALL plasma_dgesv( n, nrhs, afac, lda,
368  $ iwork, x, lda, info )
369 *
370 * Check error code from DGESV .
371 *
372  IF( info.NE.izero )
373  $ CALL alaerh( path, 'DGESV ', info, izero,
374  $ ' ', n, n, -1, -1, nrhs, imat,
375  $ nfail, nerrs, nout )
376  IF( izero.EQ.0 ) THEN
377 *
378 * Compute residual of the computed solution.
379 *
380  CALL dlacpy( 'Full', n, nrhs, b, lda, work,
381  $ lda )
382  CALL dget02( 'No transpose', n, n, nrhs, a,
383  $ lda, x, lda, work, lda, rwork,
384  $ result( 1 ) )
385 *
386 * Check solution from generated exact solution.
387 *
388  CALL dget04( n, nrhs, x, lda, xact, lda,
389  $ rcondc, result( 2 ) )
390  nt = 2
391  END IF
392 *
393 * Print information about the tests that did not
394 * pass the threshold.
395 *
396  DO 30 k = 1, nt
397  IF( result( k ).GE.thresh ) THEN
398  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
399  $ CALL aladhd( nout, path )
400  WRITE( nout, fmt = 9999 )'DGESV ', n, ib,
401  $ imat, k, result( k )
402  nfail = nfail + 1
403  END IF
404  30 continue
405  nrun = nrun + nt
406  END IF
407  50 continue
408  60 continue
409  70 continue
410  80 continue
411 *
412 * DEALLOCATE HL and HPIV
413 *
414 c$$$ CALL PLASMA_DEALLOC_HANDLE( HL, INFO )
415 c$$$ CALL PLASMA_DEALLOC_HANDLE( HPIV, INFO )
416  90 continue
417 *
418 * Print a summary of the results.
419 *
420  CALL alasvm( path, nout, nfail, nrun, nerrs )
421 *
422  9999 format( 1x, a, ', N =', i5,', NB=', i5, ', type ', i2,
423  $ ', test(', i2, ') =', g12.5 )
424  9998 format( 1x, a, ', FACT=''', a1, ''', TRANS=''', a1, ''', N=', i5,
425  $ ', type ', i2, ', test(', i1, ')=', g12.5 )
426  9997 format( 1x, a, ', FACT=''', a1, ''', TRANS=''', a1, ''', N=', i5,
427  $ ', EQUED=''', a1, ''', type ', i2, ', test(', i1, ')=',
428  $ g12.5 )
429  return
430 *
431 * End of DDRVGE
432 *
433  END