MAGMA  1.2.0
MatrixAlgebraonGPUandMulticoreArchitectures
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups
dchkeq.f
Go to the documentation of this file.
1  SUBROUTINE dchkeq( THRESH, NOUT )
2 *
3 * -- LAPACK test routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8  INTEGER nout
9  DOUBLE PRECISION thresh
10 * ..
11 *
12 * Purpose
13 * =======
14 *
15 * DCHKEQ tests DGEEQU, DGBEQU, DPOEQU, DPPEQU and DPBEQU
16 *
17 * Arguments
18 * =========
19 *
20 * THRESH (input) DOUBLE PRECISION
21 * Threshold for testing routines. Should be between 2 and 10.
22 *
23 * NOUT (input) INTEGER
24 * The unit number for output.
25 *
26 * =====================================================================
27 *
28 * .. Parameters ..
29  DOUBLE PRECISION zero, one, ten
30  parameter( zero = 0.0d0, one = 1.0d+0, ten = 1.0d1 )
31  INTEGER nsz, nszb
32  parameter( nsz = 5, nszb = 3*nsz-2 )
33  INTEGER nszp, npow
34  parameter( nszp = ( nsz*( nsz+1 ) ) / 2,
35  $ npow = 2*nsz+1 )
36 * ..
37 * .. Local Scalars ..
38  LOGICAL ok
39  CHARACTER*3 path
40  INTEGER i, info, j, kl, ku, m, n
41  DOUBLE PRECISION ccond, eps, norm, ratio, rcmax, rcmin, rcond
42 * ..
43 * .. Local Arrays ..
44  DOUBLE PRECISION a( nsz, nsz ), ab( nszb, nsz ), ap( nszp ),
45  $ c( nsz ), pow( npow ), r( nsz ), reslts( 5 ),
46  $ rpow( npow )
47 * ..
48 * .. External Functions ..
49  DOUBLE PRECISION dlamch
50  EXTERNAL dlamch
51 * ..
52 * .. External Subroutines ..
53  EXTERNAL dgbequ, dgeequ, dpbequ, dpoequ, dppequ
54 * ..
55 * .. Intrinsic Functions ..
56  INTRINSIC abs, max, min
57 * ..
58 * .. Executable Statements ..
59 *
60  path( 1: 1 ) = 'Double precision'
61  path( 2: 3 ) = 'EQ'
62 *
63  eps = dlamch( 'P' )
64  DO 10 i = 1, 5
65  reslts( i ) = zero
66  10 continue
67  DO 20 i = 1, npow
68  pow( i ) = ten**( i-1 )
69  rpow( i ) = one / pow( i )
70  20 continue
71 *
72 * Test DGEEQU
73 *
74  DO 80 n = 0, nsz
75  DO 70 m = 0, nsz
76 *
77  DO 40 j = 1, nsz
78  DO 30 i = 1, nsz
79  IF( i.LE.m .AND. j.LE.n ) THEN
80  a( i, j ) = pow( i+j+1 )*( -1 )**( i+j )
81  ELSE
82  a( i, j ) = zero
83  END IF
84  30 continue
85  40 continue
86 *
87  CALL dgeequ( m, n, a, nsz, r, c, rcond, ccond, norm, info )
88 *
89  IF( info.NE.0 ) THEN
90  reslts( 1 ) = one
91  ELSE
92  IF( n.NE.0 .AND. m.NE.0 ) THEN
93  reslts( 1 ) = max( reslts( 1 ),
94  $ abs( ( rcond-rpow( m ) ) / rpow( m ) ) )
95  reslts( 1 ) = max( reslts( 1 ),
96  $ abs( ( ccond-rpow( n ) ) / rpow( n ) ) )
97  reslts( 1 ) = max( reslts( 1 ),
98  $ abs( ( norm-pow( n+m+1 ) ) / pow( n+m+
99  $ 1 ) ) )
100  DO 50 i = 1, m
101  reslts( 1 ) = max( reslts( 1 ),
102  $ abs( ( r( i )-rpow( i+n+1 ) ) /
103  $ rpow( i+n+1 ) ) )
104  50 continue
105  DO 60 j = 1, n
106  reslts( 1 ) = max( reslts( 1 ),
107  $ abs( ( c( j )-pow( n-j+1 ) ) /
108  $ pow( n-j+1 ) ) )
109  60 continue
110  END IF
111  END IF
112 *
113  70 continue
114  80 continue
115 *
116 * Test with zero rows and columns
117 *
118  DO 90 j = 1, nsz
119  a( max( nsz-1, 1 ), j ) = zero
120  90 continue
121  CALL dgeequ( nsz, nsz, a, nsz, r, c, rcond, ccond, norm, info )
122  IF( info.NE.max( nsz-1, 1 ) )
123  $ reslts( 1 ) = one
124 *
125  DO 100 j = 1, nsz
126  a( max( nsz-1, 1 ), j ) = one
127  100 continue
128  DO 110 i = 1, nsz
129  a( i, max( nsz-1, 1 ) ) = zero
130  110 continue
131  CALL dgeequ( nsz, nsz, a, nsz, r, c, rcond, ccond, norm, info )
132  IF( info.NE.nsz+max( nsz-1, 1 ) )
133  $ reslts( 1 ) = one
134  reslts( 1 ) = reslts( 1 ) / eps
135 *
136 * Test DGBEQU
137 *
138  DO 250 n = 0, nsz
139  DO 240 m = 0, nsz
140  DO 230 kl = 0, max( m-1, 0 )
141  DO 220 ku = 0, max( n-1, 0 )
142 *
143  DO 130 j = 1, nsz
144  DO 120 i = 1, nszb
145  ab( i, j ) = zero
146  120 continue
147  130 continue
148  DO 150 j = 1, n
149  DO 140 i = 1, m
150  IF( i.LE.min( m, j+kl ) .AND. i.GE.
151  $ max( 1, j-ku ) .AND. j.LE.n ) THEN
152  ab( ku+1+i-j, j ) = pow( i+j+1 )*
153  $ ( -1 )**( i+j )
154  END IF
155  140 continue
156  150 continue
157 *
158  CALL dgbequ( m, n, kl, ku, ab, nszb, r, c, rcond,
159  $ ccond, norm, info )
160 *
161  IF( info.NE.0 ) THEN
162  IF( .NOT.( ( n+kl.LT.m .AND. info.EQ.n+kl+1 ) .OR.
163  $ ( m+ku.LT.n .AND. info.EQ.2*m+ku+1 ) ) ) THEN
164  reslts( 2 ) = one
165  END IF
166  ELSE
167  IF( n.NE.0 .AND. m.NE.0 ) THEN
168 *
169  rcmin = r( 1 )
170  rcmax = r( 1 )
171  DO 160 i = 1, m
172  rcmin = min( rcmin, r( i ) )
173  rcmax = max( rcmax, r( i ) )
174  160 continue
175  ratio = rcmin / rcmax
176  reslts( 2 ) = max( reslts( 2 ),
177  $ abs( ( rcond-ratio ) / ratio ) )
178 *
179  rcmin = c( 1 )
180  rcmax = c( 1 )
181  DO 170 j = 1, n
182  rcmin = min( rcmin, c( j ) )
183  rcmax = max( rcmax, c( j ) )
184  170 continue
185  ratio = rcmin / rcmax
186  reslts( 2 ) = max( reslts( 2 ),
187  $ abs( ( ccond-ratio ) / ratio ) )
188 *
189  reslts( 2 ) = max( reslts( 2 ),
190  $ abs( ( norm-pow( n+m+1 ) ) /
191  $ pow( n+m+1 ) ) )
192  DO 190 i = 1, m
193  rcmax = zero
194  DO 180 j = 1, n
195  IF( i.LE.j+kl .AND. i.GE.j-ku ) THEN
196  ratio = abs( r( i )*pow( i+j+1 )*
197  $ c( j ) )
198  rcmax = max( rcmax, ratio )
199  END IF
200  180 continue
201  reslts( 2 ) = max( reslts( 2 ),
202  $ abs( one-rcmax ) )
203  190 continue
204 *
205  DO 210 j = 1, n
206  rcmax = zero
207  DO 200 i = 1, m
208  IF( i.LE.j+kl .AND. i.GE.j-ku ) THEN
209  ratio = abs( r( i )*pow( i+j+1 )*
210  $ c( j ) )
211  rcmax = max( rcmax, ratio )
212  END IF
213  200 continue
214  reslts( 2 ) = max( reslts( 2 ),
215  $ abs( one-rcmax ) )
216  210 continue
217  END IF
218  END IF
219 *
220  220 continue
221  230 continue
222  240 continue
223  250 continue
224  reslts( 2 ) = reslts( 2 ) / eps
225 *
226 * Test DPOEQU
227 *
228  DO 290 n = 0, nsz
229 *
230  DO 270 i = 1, nsz
231  DO 260 j = 1, nsz
232  IF( i.LE.n .AND. j.EQ.i ) THEN
233  a( i, j ) = pow( i+j+1 )*( -1 )**( i+j )
234  ELSE
235  a( i, j ) = zero
236  END IF
237  260 continue
238  270 continue
239 *
240  CALL dpoequ( n, a, nsz, r, rcond, norm, info )
241 *
242  IF( info.NE.0 ) THEN
243  reslts( 3 ) = one
244  ELSE
245  IF( n.NE.0 ) THEN
246  reslts( 3 ) = max( reslts( 3 ),
247  $ abs( ( rcond-rpow( n ) ) / rpow( n ) ) )
248  reslts( 3 ) = max( reslts( 3 ),
249  $ abs( ( norm-pow( 2*n+1 ) ) / pow( 2*n+
250  $ 1 ) ) )
251  DO 280 i = 1, n
252  reslts( 3 ) = max( reslts( 3 ),
253  $ abs( ( r( i )-rpow( i+1 ) ) / rpow( i+
254  $ 1 ) ) )
255  280 continue
256  END IF
257  END IF
258  290 continue
259  a( max( nsz-1, 1 ), max( nsz-1, 1 ) ) = -one
260  CALL dpoequ( nsz, a, nsz, r, rcond, norm, info )
261  IF( info.NE.max( nsz-1, 1 ) )
262  $ reslts( 3 ) = one
263  reslts( 3 ) = reslts( 3 ) / eps
264 *
265 * Test DPPEQU
266 *
267  DO 360 n = 0, nsz
268 *
269 * Upper triangular packed storage
270 *
271  DO 300 i = 1, ( n*( n+1 ) ) / 2
272  ap( i ) = zero
273  300 continue
274  DO 310 i = 1, n
275  ap( ( i*( i+1 ) ) / 2 ) = pow( 2*i+1 )
276  310 continue
277 *
278  CALL dppequ( 'U', n, ap, r, rcond, norm, info )
279 *
280  IF( info.NE.0 ) THEN
281  reslts( 4 ) = one
282  ELSE
283  IF( n.NE.0 ) THEN
284  reslts( 4 ) = max( reslts( 4 ),
285  $ abs( ( rcond-rpow( n ) ) / rpow( n ) ) )
286  reslts( 4 ) = max( reslts( 4 ),
287  $ abs( ( norm-pow( 2*n+1 ) ) / pow( 2*n+
288  $ 1 ) ) )
289  DO 320 i = 1, n
290  reslts( 4 ) = max( reslts( 4 ),
291  $ abs( ( r( i )-rpow( i+1 ) ) / rpow( i+
292  $ 1 ) ) )
293  320 continue
294  END IF
295  END IF
296 *
297 * Lower triangular packed storage
298 *
299  DO 330 i = 1, ( n*( n+1 ) ) / 2
300  ap( i ) = zero
301  330 continue
302  j = 1
303  DO 340 i = 1, n
304  ap( j ) = pow( 2*i+1 )
305  j = j + ( n-i+1 )
306  340 continue
307 *
308  CALL dppequ( 'L', n, ap, r, rcond, norm, info )
309 *
310  IF( info.NE.0 ) THEN
311  reslts( 4 ) = one
312  ELSE
313  IF( n.NE.0 ) THEN
314  reslts( 4 ) = max( reslts( 4 ),
315  $ abs( ( rcond-rpow( n ) ) / rpow( n ) ) )
316  reslts( 4 ) = max( reslts( 4 ),
317  $ abs( ( norm-pow( 2*n+1 ) ) / pow( 2*n+
318  $ 1 ) ) )
319  DO 350 i = 1, n
320  reslts( 4 ) = max( reslts( 4 ),
321  $ abs( ( r( i )-rpow( i+1 ) ) / rpow( i+
322  $ 1 ) ) )
323  350 continue
324  END IF
325  END IF
326 *
327  360 continue
328  i = ( nsz*( nsz+1 ) ) / 2 - 2
329  ap( i ) = -one
330  CALL dppequ( 'L', nsz, ap, r, rcond, norm, info )
331  IF( info.NE.max( nsz-1, 1 ) )
332  $ reslts( 4 ) = one
333  reslts( 4 ) = reslts( 4 ) / eps
334 *
335 * Test DPBEQU
336 *
337  DO 460 n = 0, nsz
338  DO 450 kl = 0, max( n-1, 0 )
339 *
340 * Test upper triangular storage
341 *
342  DO 380 j = 1, nsz
343  DO 370 i = 1, nszb
344  ab( i, j ) = zero
345  370 continue
346  380 continue
347  DO 390 j = 1, n
348  ab( kl+1, j ) = pow( 2*j+1 )
349  390 continue
350 *
351  CALL dpbequ( 'U', n, kl, ab, nszb, r, rcond, norm, info )
352 *
353  IF( info.NE.0 ) THEN
354  reslts( 5 ) = one
355  ELSE
356  IF( n.NE.0 ) THEN
357  reslts( 5 ) = max( reslts( 5 ),
358  $ abs( ( rcond-rpow( n ) ) / rpow( n ) ) )
359  reslts( 5 ) = max( reslts( 5 ),
360  $ abs( ( norm-pow( 2*n+1 ) ) / pow( 2*n+
361  $ 1 ) ) )
362  DO 400 i = 1, n
363  reslts( 5 ) = max( reslts( 5 ),
364  $ abs( ( r( i )-rpow( i+1 ) ) /
365  $ rpow( i+1 ) ) )
366  400 continue
367  END IF
368  END IF
369  IF( n.NE.0 ) THEN
370  ab( kl+1, max( n-1, 1 ) ) = -one
371  CALL dpbequ( 'U', n, kl, ab, nszb, r, rcond, norm, info )
372  IF( info.NE.max( n-1, 1 ) )
373  $ reslts( 5 ) = one
374  END IF
375 *
376 * Test lower triangular storage
377 *
378  DO 420 j = 1, nsz
379  DO 410 i = 1, nszb
380  ab( i, j ) = zero
381  410 continue
382  420 continue
383  DO 430 j = 1, n
384  ab( 1, j ) = pow( 2*j+1 )
385  430 continue
386 *
387  CALL dpbequ( 'L', n, kl, ab, nszb, r, rcond, norm, info )
388 *
389  IF( info.NE.0 ) THEN
390  reslts( 5 ) = one
391  ELSE
392  IF( n.NE.0 ) THEN
393  reslts( 5 ) = max( reslts( 5 ),
394  $ abs( ( rcond-rpow( n ) ) / rpow( n ) ) )
395  reslts( 5 ) = max( reslts( 5 ),
396  $ abs( ( norm-pow( 2*n+1 ) ) / pow( 2*n+
397  $ 1 ) ) )
398  DO 440 i = 1, n
399  reslts( 5 ) = max( reslts( 5 ),
400  $ abs( ( r( i )-rpow( i+1 ) ) /
401  $ rpow( i+1 ) ) )
402  440 continue
403  END IF
404  END IF
405  IF( n.NE.0 ) THEN
406  ab( 1, max( n-1, 1 ) ) = -one
407  CALL dpbequ( 'L', n, kl, ab, nszb, r, rcond, norm, info )
408  IF( info.NE.max( n-1, 1 ) )
409  $ reslts( 5 ) = one
410  END IF
411  450 continue
412  460 continue
413  reslts( 5 ) = reslts( 5 ) / eps
414  ok = ( reslts( 1 ).LE.thresh ) .AND.
415  $ ( reslts( 2 ).LE.thresh ) .AND.
416  $ ( reslts( 3 ).LE.thresh ) .AND.
417  $ ( reslts( 4 ).LE.thresh ) .AND. ( reslts( 5 ).LE.thresh )
418  WRITE( nout, fmt = * )
419  IF( ok ) THEN
420  WRITE( nout, fmt = 9999 )path
421  ELSE
422  IF( reslts( 1 ).GT.thresh )
423  $ WRITE( nout, fmt = 9998 )reslts( 1 ), thresh
424  IF( reslts( 2 ).GT.thresh )
425  $ WRITE( nout, fmt = 9997 )reslts( 2 ), thresh
426  IF( reslts( 3 ).GT.thresh )
427  $ WRITE( nout, fmt = 9996 )reslts( 3 ), thresh
428  IF( reslts( 4 ).GT.thresh )
429  $ WRITE( nout, fmt = 9995 )reslts( 4 ), thresh
430  IF( reslts( 5 ).GT.thresh )
431  $ WRITE( nout, fmt = 9994 )reslts( 5 ), thresh
432  END IF
433  9999 format( 1x, 'All tests for ', a3,
434  $ ' routines passed the threshold' )
435  9998 format( ' DGEEQU failed test with value ', d10.3, ' exceeding',
436  $ ' threshold ', d10.3 )
437  9997 format( ' DGBEQU failed test with value ', d10.3, ' exceeding',
438  $ ' threshold ', d10.3 )
439  9996 format( ' DPOEQU failed test with value ', d10.3, ' exceeding',
440  $ ' threshold ', d10.3 )
441  9995 format( ' DPPEQU failed test with value ', d10.3, ' exceeding',
442  $ ' threshold ', d10.3 )
443  9994 format( ' DPBEQU failed test with value ', d10.3, ' exceeding',
444  $ ' threshold ', d10.3 )
445  return
446 *
447 * End of DCHKEQ
448 *
449  END