LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
alahd.f
Go to the documentation of this file.
1 *> \brief \b ALAHD
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 ALAHD( IOUNIT, PATH )
12 *
13 * .. Scalar Arguments ..
14 * CHARACTER*3 PATH
15 * INTEGER IOUNIT
16 * ..
17 *
18 *
19 *> \par Purpose:
20 * =============
21 *>
22 *> \verbatim
23 *>
24 *> ALAHD prints header information for the different test paths.
25 *> \endverbatim
26 *
27 * Arguments:
28 * ==========
29 *
30 *> \param[in] IOUNIT
31 *> \verbatim
32 *> IOUNIT is INTEGER
33 *> The unit number to which the header information should be
34 *> printed.
35 *> \endverbatim
36 *>
37 *> \param[in] PATH
38 *> \verbatim
39 *> PATH is CHARACTER*3
40 *> The name of the path for which the header information is to
41 *> be printed. Current paths are
42 *> _GE: General matrices
43 *> _GB: General band
44 *> _GT: General Tridiagonal
45 *> _PO: Symmetric or Hermitian positive definite
46 *> _PS: Symmetric or Hermitian positive semi-definite
47 *> _PP: Symmetric or Hermitian positive definite packed
48 *> _PB: Symmetric or Hermitian positive definite band
49 *> _PT: Symmetric or Hermitian positive definite tridiagonal
50 *> _SY: Symmetric indefinite
51 *> _SP: Symmetric indefinite packed
52 *> _HE: (complex) Hermitian indefinite
53 *> _HP: (complex) Hermitian indefinite packed
54 *> _TR: Triangular
55 *> _TP: Triangular packed
56 *> _TB: Triangular band
57 *> _QR: QR (general matrices)
58 *> _LQ: LQ (general matrices)
59 *> _QL: QL (general matrices)
60 *> _RQ: RQ (general matrices)
61 *> _QP: QR with column pivoting
62 *> _TZ: Trapezoidal
63 *> _LS: Least Squares driver routines
64 *> _LU: LU variants
65 *> _CH: Cholesky variants
66 *> _QS: QR variants
67 *> _QT: QRT (general matrices)
68 *> _QX: QRT (triangular-pentagonal matrices)
69 *> The first character must be one of S, D, C, or Z (C or Z only
70 *> if complex).
71 *> \endverbatim
72 *
73 * Authors:
74 * ========
75 *
76 *> \author Univ. of Tennessee
77 *> \author Univ. of California Berkeley
78 *> \author Univ. of Colorado Denver
79 *> \author NAG Ltd.
80 *
81 *> \date April 2012
82 *
83 *> \ingroup aux_lin
84 *
85 * =====================================================================
86  SUBROUTINE alahd( IOUNIT, PATH )
87 *
88 * -- LAPACK test routine (version 3.4.1) --
89 * -- LAPACK is a software package provided by Univ. of Tennessee, --
90 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
91 * April 2012
92 *
93 * .. Scalar Arguments ..
94  CHARACTER*3 path
95  INTEGER iounit
96 * ..
97 *
98 * =====================================================================
99 *
100 * .. Local Scalars ..
101  LOGICAL corz, sord
102  CHARACTER c1, c3
103  CHARACTER*2 p2
104  CHARACTER*4 eigcnm
105  CHARACTER*32 subnam
106  CHARACTER*9 sym
107 * ..
108 * .. External Functions ..
109  LOGICAL lsame, lsamen
110  EXTERNAL lsame, lsamen
111 * ..
112 * .. Intrinsic Functions ..
113  INTRINSIC len_trim
114 * ..
115 * .. Executable Statements ..
116 *
117  IF( iounit.LE.0 )
118  $ return
119  c1 = path( 1: 1 )
120  c3 = path( 3: 3 )
121  p2 = path( 2: 3 )
122  sord = lsame( c1, 'S' ) .OR. lsame( c1, 'D' )
123  corz = lsame( c1, 'C' ) .OR. lsame( c1, 'Z' )
124  IF( .NOT.( sord .OR. corz ) )
125  $ return
126 *
127  IF( lsamen( 2, p2, 'GE' ) ) THEN
128 *
129 * GE: General dense
130 *
131  WRITE( iounit, fmt = 9999 )path
132  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
133  WRITE( iounit, fmt = 9979 )
134  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
135  WRITE( iounit, fmt = 9962 )1
136  WRITE( iounit, fmt = 9961 )2
137  WRITE( iounit, fmt = 9960 )3
138  WRITE( iounit, fmt = 9959 )4
139  WRITE( iounit, fmt = 9958 )5
140  WRITE( iounit, fmt = 9957 )6
141  WRITE( iounit, fmt = 9956 )7
142  WRITE( iounit, fmt = 9955 )8
143  WRITE( iounit, fmt = '( '' Messages:'' )' )
144 *
145  ELSE IF( lsamen( 2, p2, 'GB' ) ) THEN
146 *
147 * GB: General band
148 *
149  WRITE( iounit, fmt = 9998 )path
150  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
151  WRITE( iounit, fmt = 9978 )
152  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
153  WRITE( iounit, fmt = 9962 )1
154  WRITE( iounit, fmt = 9960 )2
155  WRITE( iounit, fmt = 9959 )3
156  WRITE( iounit, fmt = 9958 )4
157  WRITE( iounit, fmt = 9957 )5
158  WRITE( iounit, fmt = 9956 )6
159  WRITE( iounit, fmt = 9955 )7
160  WRITE( iounit, fmt = '( '' Messages:'' )' )
161 *
162  ELSE IF( lsamen( 2, p2, 'GT' ) ) THEN
163 *
164 * GT: General tridiagonal
165 *
166  WRITE( iounit, fmt = 9997 )path
167  WRITE( iounit, fmt = 9977 )
168  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
169  WRITE( iounit, fmt = 9962 )1
170  WRITE( iounit, fmt = 9960 )2
171  WRITE( iounit, fmt = 9959 )3
172  WRITE( iounit, fmt = 9958 )4
173  WRITE( iounit, fmt = 9957 )5
174  WRITE( iounit, fmt = 9956 )6
175  WRITE( iounit, fmt = 9955 )7
176  WRITE( iounit, fmt = '( '' Messages:'' )' )
177 *
178  ELSE IF( lsamen( 2, p2, 'PO' ) .OR. lsamen( 2, p2, 'PP' ) ) THEN
179 *
180 * PO: Positive definite full
181 * PP: Positive definite packed
182 *
183  IF( sord ) THEN
184  sym = 'Symmetric'
185  ELSE
186  sym = 'Hermitian'
187  END IF
188  IF( lsame( c3, 'O' ) ) THEN
189  WRITE( iounit, fmt = 9996 )path, sym
190  ELSE
191  WRITE( iounit, fmt = 9995 )path, sym
192  END IF
193  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
194  WRITE( iounit, fmt = 9975 )path
195  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
196  WRITE( iounit, fmt = 9954 )1
197  WRITE( iounit, fmt = 9961 )2
198  WRITE( iounit, fmt = 9960 )3
199  WRITE( iounit, fmt = 9959 )4
200  WRITE( iounit, fmt = 9958 )5
201  WRITE( iounit, fmt = 9957 )6
202  WRITE( iounit, fmt = 9956 )7
203  WRITE( iounit, fmt = 9955 )8
204  WRITE( iounit, fmt = '( '' Messages:'' )' )
205 *
206  ELSE IF( lsamen( 2, p2, 'PS' ) ) THEN
207 *
208 * PS: Positive semi-definite full
209 *
210  IF( sord ) THEN
211  sym = 'Symmetric'
212  ELSE
213  sym = 'Hermitian'
214  END IF
215  IF( lsame( c1, 'S' ) .OR. lsame( c1, 'C' ) ) THEN
216  eigcnm = '1E04'
217  ELSE
218  eigcnm = '1D12'
219  END IF
220  WRITE( iounit, fmt = 9995 )path, sym
221  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
222  WRITE( iounit, fmt = 8973 )eigcnm, eigcnm, eigcnm
223  WRITE( iounit, fmt = '( '' Difference:'' )' )
224  WRITE( iounit, fmt = 8972 )c1
225  WRITE( iounit, fmt = '( '' Test ratio:'' )' )
226  WRITE( iounit, fmt = 8950 )
227  WRITE( iounit, fmt = '( '' Messages:'' )' )
228  ELSE IF( lsamen( 2, p2, 'PB' ) ) THEN
229 *
230 * PB: Positive definite band
231 *
232  IF( sord ) THEN
233  WRITE( iounit, fmt = 9994 )path, 'Symmetric'
234  ELSE
235  WRITE( iounit, fmt = 9994 )path, 'Hermitian'
236  END IF
237  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
238  WRITE( iounit, fmt = 9973 )path
239  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
240  WRITE( iounit, fmt = 9954 )1
241  WRITE( iounit, fmt = 9960 )2
242  WRITE( iounit, fmt = 9959 )3
243  WRITE( iounit, fmt = 9958 )4
244  WRITE( iounit, fmt = 9957 )5
245  WRITE( iounit, fmt = 9956 )6
246  WRITE( iounit, fmt = 9955 )7
247  WRITE( iounit, fmt = '( '' Messages:'' )' )
248 *
249  ELSE IF( lsamen( 2, p2, 'PT' ) ) THEN
250 *
251 * PT: Positive definite tridiagonal
252 *
253  IF( sord ) THEN
254  WRITE( iounit, fmt = 9993 )path, 'Symmetric'
255  ELSE
256  WRITE( iounit, fmt = 9993 )path, 'Hermitian'
257  END IF
258  WRITE( iounit, fmt = 9976 )
259  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
260  WRITE( iounit, fmt = 9952 )1
261  WRITE( iounit, fmt = 9960 )2
262  WRITE( iounit, fmt = 9959 )3
263  WRITE( iounit, fmt = 9958 )4
264  WRITE( iounit, fmt = 9957 )5
265  WRITE( iounit, fmt = 9956 )6
266  WRITE( iounit, fmt = 9955 )7
267  WRITE( iounit, fmt = '( '' Messages:'' )' )
268 *
269  ELSE IF( lsamen( 2, p2, 'SY' ) ) THEN
270 *
271 * SY: Symmetric indefinite full
272 *
273  IF( lsame( c3, 'Y' ) ) THEN
274  WRITE( iounit, fmt = 9992 )path, 'Symmetric'
275  ELSE
276  WRITE( iounit, fmt = 9991 )path, 'Symmetric'
277  END IF
278  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
279  IF( sord ) THEN
280  WRITE( iounit, fmt = 9972 )
281  ELSE
282  WRITE( iounit, fmt = 9971 )
283  END IF
284  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
285  WRITE( iounit, fmt = 9953 )1
286  WRITE( iounit, fmt = 9961 )2
287  WRITE( iounit, fmt = 9960 )3
288  WRITE( iounit, fmt = 9960 )4
289  WRITE( iounit, fmt = 9959 )5
290  WRITE( iounit, fmt = 9958 )6
291  WRITE( iounit, fmt = 9956 )7
292  WRITE( iounit, fmt = 9957 )8
293  WRITE( iounit, fmt = 9955 )9
294  WRITE( iounit, fmt = '( '' Messages:'' )' )
295 *
296  ELSE IF( lsamen( 2, p2, 'SP' ) ) THEN
297 *
298 * SP: Symmetric indefinite packed
299 *
300  IF( lsame( c3, 'Y' ) ) THEN
301  WRITE( iounit, fmt = 9992 )path, 'Symmetric'
302  ELSE
303  WRITE( iounit, fmt = 9991 )path, 'Symmetric'
304  END IF
305  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
306  IF( sord ) THEN
307  WRITE( iounit, fmt = 9972 )
308  ELSE
309  WRITE( iounit, fmt = 9971 )
310  END IF
311  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
312  WRITE( iounit, fmt = 9953 )1
313  WRITE( iounit, fmt = 9961 )2
314  WRITE( iounit, fmt = 9960 )3
315  WRITE( iounit, fmt = 9959 )4
316  WRITE( iounit, fmt = 9958 )5
317  WRITE( iounit, fmt = 9956 )6
318  WRITE( iounit, fmt = 9957 )7
319  WRITE( iounit, fmt = 9955 )8
320  WRITE( iounit, fmt = '( '' Messages:'' )' )
321 *
322  ELSE IF( lsamen( 2, p2, 'HE' ) ) THEN
323 *
324 * HE: Hermitian indefinite full
325 *
326  IF( lsame( c3, 'E' ) ) THEN
327  WRITE( iounit, fmt = 9992 )path, 'Hermitian'
328  ELSE
329  WRITE( iounit, fmt = 9991 )path, 'Hermitian'
330  END IF
331  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
332  IF( sord ) THEN
333  WRITE( iounit, fmt = 9972 )
334  ELSE
335  WRITE( iounit, fmt = 9971 )
336  END IF
337  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
338  WRITE( iounit, fmt = 9953 )1
339  WRITE( iounit, fmt = 9961 )2
340  WRITE( iounit, fmt = 9960 )3
341  WRITE( iounit, fmt = 9960 )4
342  WRITE( iounit, fmt = 9959 )5
343  WRITE( iounit, fmt = 9958 )6
344  WRITE( iounit, fmt = 9956 )7
345  WRITE( iounit, fmt = 9957 )8
346  WRITE( iounit, fmt = 9955 )9
347  WRITE( iounit, fmt = '( '' Messages:'' )' )
348 *
349  ELSE IF( lsamen( 2, p2, 'HP' ) ) THEN
350 *
351 * HP: Hermitian indefinite packed
352 *
353  IF( lsame( c3, 'E' ) ) THEN
354  WRITE( iounit, fmt = 9992 )path, 'Hermitian'
355  ELSE
356  WRITE( iounit, fmt = 9991 )path, 'Hermitian'
357  END IF
358  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
359  WRITE( iounit, fmt = 9972 )
360  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
361  WRITE( iounit, fmt = 9953 )1
362  WRITE( iounit, fmt = 9961 )2
363  WRITE( iounit, fmt = 9960 )3
364  WRITE( iounit, fmt = 9959 )4
365  WRITE( iounit, fmt = 9958 )5
366  WRITE( iounit, fmt = 9956 )6
367  WRITE( iounit, fmt = 9957 )7
368  WRITE( iounit, fmt = 9955 )8
369  WRITE( iounit, fmt = '( '' Messages:'' )' )
370 *
371  ELSE IF( lsamen( 2, p2, 'TR' ) .OR. lsamen( 2, p2, 'TP' ) ) THEN
372 *
373 * TR: Triangular full
374 * TP: Triangular packed
375 *
376  IF( lsame( c3, 'R' ) ) THEN
377  WRITE( iounit, fmt = 9990 )path
378  subnam = path( 1: 1 ) // 'LATRS'
379  ELSE
380  WRITE( iounit, fmt = 9989 )path
381  subnam = path( 1: 1 ) // 'LATPS'
382  END IF
383  WRITE( iounit, fmt = 9966 )path
384  WRITE( iounit, fmt = 9965 )subnam(1:len_trim( subnam ))
385  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
386  WRITE( iounit, fmt = 9961 )1
387  WRITE( iounit, fmt = 9960 )2
388  WRITE( iounit, fmt = 9959 )3
389  WRITE( iounit, fmt = 9958 )4
390  WRITE( iounit, fmt = 9957 )5
391  WRITE( iounit, fmt = 9956 )6
392  WRITE( iounit, fmt = 9955 )7
393  WRITE( iounit, fmt = 9951 )subnam(1:len_trim( subnam )), 8
394  WRITE( iounit, fmt = '( '' Messages:'' )' )
395 *
396  ELSE IF( lsamen( 2, p2, 'TB' ) ) THEN
397 *
398 * TB: Triangular band
399 *
400  WRITE( iounit, fmt = 9988 )path
401  subnam = path( 1: 1 ) // 'LATBS'
402  WRITE( iounit, fmt = 9964 )path
403  WRITE( iounit, fmt = 9963 )subnam(1:len_trim( subnam ))
404  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
405  WRITE( iounit, fmt = 9960 )1
406  WRITE( iounit, fmt = 9959 )2
407  WRITE( iounit, fmt = 9958 )3
408  WRITE( iounit, fmt = 9957 )4
409  WRITE( iounit, fmt = 9956 )5
410  WRITE( iounit, fmt = 9955 )6
411  WRITE( iounit, fmt = 9951 )subnam(1:len_trim( subnam )), 7
412  WRITE( iounit, fmt = '( '' Messages:'' )' )
413 *
414  ELSE IF( lsamen( 2, p2, 'QR' ) ) THEN
415 *
416 * QR decomposition of rectangular matrices
417 *
418  WRITE( iounit, fmt = 9987 )path, 'QR'
419  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
420  WRITE( iounit, fmt = 9970 )
421  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
422  WRITE( iounit, fmt = 9950 )1
423  WRITE( iounit, fmt = 6950 )8
424  WRITE( iounit, fmt = 9946 )2
425  WRITE( iounit, fmt = 9944 )3, 'M'
426  WRITE( iounit, fmt = 9943 )4, 'M'
427  WRITE( iounit, fmt = 9942 )5, 'M'
428  WRITE( iounit, fmt = 9941 )6, 'M'
429  WRITE( iounit, fmt = 9960 )7
430  WRITE( iounit, fmt = 6660 )9
431  WRITE( iounit, fmt = '( '' Messages:'' )' )
432 *
433  ELSE IF( lsamen( 2, p2, 'LQ' ) ) THEN
434 *
435 * LQ decomposition of rectangular matrices
436 *
437  WRITE( iounit, fmt = 9987 )path, 'LQ'
438  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
439  WRITE( iounit, fmt = 9970 )
440  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
441  WRITE( iounit, fmt = 9949 )1
442  WRITE( iounit, fmt = 9945 )2
443  WRITE( iounit, fmt = 9944 )3, 'N'
444  WRITE( iounit, fmt = 9943 )4, 'N'
445  WRITE( iounit, fmt = 9942 )5, 'N'
446  WRITE( iounit, fmt = 9941 )6, 'N'
447  WRITE( iounit, fmt = 9960 )7
448  WRITE( iounit, fmt = '( '' Messages:'' )' )
449 *
450  ELSE IF( lsamen( 2, p2, 'QL' ) ) THEN
451 *
452 * QL decomposition of rectangular matrices
453 *
454  WRITE( iounit, fmt = 9987 )path, 'QL'
455  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
456  WRITE( iounit, fmt = 9970 )
457  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
458  WRITE( iounit, fmt = 9948 )1
459  WRITE( iounit, fmt = 9946 )2
460  WRITE( iounit, fmt = 9944 )3, 'M'
461  WRITE( iounit, fmt = 9943 )4, 'M'
462  WRITE( iounit, fmt = 9942 )5, 'M'
463  WRITE( iounit, fmt = 9941 )6, 'M'
464  WRITE( iounit, fmt = 9960 )7
465  WRITE( iounit, fmt = '( '' Messages:'' )' )
466 *
467  ELSE IF( lsamen( 2, p2, 'RQ' ) ) THEN
468 *
469 * RQ decomposition of rectangular matrices
470 *
471  WRITE( iounit, fmt = 9987 )path, 'RQ'
472  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
473  WRITE( iounit, fmt = 9970 )
474  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
475  WRITE( iounit, fmt = 9947 )1
476  WRITE( iounit, fmt = 9945 )2
477  WRITE( iounit, fmt = 9944 )3, 'N'
478  WRITE( iounit, fmt = 9943 )4, 'N'
479  WRITE( iounit, fmt = 9942 )5, 'N'
480  WRITE( iounit, fmt = 9941 )6, 'N'
481  WRITE( iounit, fmt = 9960 )7
482  WRITE( iounit, fmt = '( '' Messages:'' )' )
483 *
484  ELSE IF( lsamen( 2, p2, 'QP' ) ) THEN
485 *
486 * QR decomposition with column pivoting
487 *
488  WRITE( iounit, fmt = 9986 )path
489  WRITE( iounit, fmt = 9969 )
490  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
491  WRITE( iounit, fmt = 9940 )1
492  WRITE( iounit, fmt = 9939 )2
493  WRITE( iounit, fmt = 9938 )3
494  WRITE( iounit, fmt = '( '' Messages:'' )' )
495 *
496  ELSE IF( lsamen( 2, p2, 'TZ' ) ) THEN
497 *
498 * TZ: Trapezoidal
499 *
500  WRITE( iounit, fmt = 9985 )path
501  WRITE( iounit, fmt = 9968 )
502  WRITE( iounit, fmt = 9929 )c1, c1
503  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
504  WRITE( iounit, fmt = 9940 )1
505  WRITE( iounit, fmt = 9937 )2
506  WRITE( iounit, fmt = 9938 )3
507  WRITE( iounit, fmt = 9940 )4
508  WRITE( iounit, fmt = 9937 )5
509  WRITE( iounit, fmt = 9938 )6
510  WRITE( iounit, fmt = '( '' Messages:'' )' )
511 *
512  ELSE IF( lsamen( 2, p2, 'LS' ) ) THEN
513 *
514 * LS: Least Squares driver routines for
515 * LS, LSD, LSS, LSX and LSY.
516 *
517  WRITE( iounit, fmt = 9984 )path
518  WRITE( iounit, fmt = 9967 )
519  WRITE( iounit, fmt = 9921 )c1, c1, c1, c1, c1
520  WRITE( iounit, fmt = 9935 )1
521  WRITE( iounit, fmt = 9931 )2
522  WRITE( iounit, fmt = 9933 )3
523  WRITE( iounit, fmt = 9935 )4
524  WRITE( iounit, fmt = 9934 )5
525  WRITE( iounit, fmt = 9932 )6
526  WRITE( iounit, fmt = 9920 )
527  WRITE( iounit, fmt = '( '' Messages:'' )' )
528 *
529  ELSE IF( lsamen( 2, p2, 'LU' ) ) THEN
530 *
531 * LU factorization variants
532 *
533  WRITE( iounit, fmt = 9983 )path
534  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
535  WRITE( iounit, fmt = 9979 )
536  WRITE( iounit, fmt = '( '' Test ratio:'' )' )
537  WRITE( iounit, fmt = 9962 )1
538  WRITE( iounit, fmt = '( '' Messages:'' )' )
539 *
540  ELSE IF( lsamen( 2, p2, 'CH' ) ) THEN
541 *
542 * Cholesky factorization variants
543 *
544  WRITE( iounit, fmt = 9982 )path
545  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
546  WRITE( iounit, fmt = 9974 )
547  WRITE( iounit, fmt = '( '' Test ratio:'' )' )
548  WRITE( iounit, fmt = 9954 )1
549  WRITE( iounit, fmt = '( '' Messages:'' )' )
550 *
551  ELSE IF( lsamen( 2, p2, 'QS' ) ) THEN
552 *
553 * QR factorization variants
554 *
555  WRITE( iounit, fmt = 9981 )path
556  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
557  WRITE( iounit, fmt = 9970 )
558  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
559 *
560  ELSE IF( lsamen( 2, p2, 'QT' ) ) THEN
561 *
562 * QRT (general matrices)
563 *
564  WRITE( iounit, fmt = 8000 ) path
565  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
566  WRITE( iounit, fmt = 8011 ) 1
567  WRITE( iounit, fmt = 8012 ) 2
568  WRITE( iounit, fmt = 8013 ) 3
569  WRITE( iounit, fmt = 8014 ) 4
570  WRITE( iounit, fmt = 8015 ) 5
571  WRITE( iounit, fmt = 8016 ) 6
572 *
573  ELSE IF( lsamen( 2, p2, 'QX' ) ) THEN
574 *
575 * QRT (triangular-pentagonal)
576 *
577  WRITE( iounit, fmt = 8001 ) path
578  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
579  WRITE( iounit, fmt = 8017 ) 1
580  WRITE( iounit, fmt = 8018 ) 2
581  WRITE( iounit, fmt = 8019 ) 3
582  WRITE( iounit, fmt = 8020 ) 4
583  WRITE( iounit, fmt = 8021 ) 5
584  WRITE( iounit, fmt = 8022 ) 6
585 *
586  ELSE
587 *
588 * Print error message if no header is available.
589 *
590  WRITE( iounit, fmt = 9980 )path
591  END IF
592 *
593 * First line of header
594 *
595  9999 format( / 1x, a3, ': General dense matrices' )
596  9998 format( / 1x, a3, ': General band matrices' )
597  9997 format( / 1x, a3, ': General tridiagonal' )
598  9996 format( / 1x, a3, ': ', a9, ' positive definite matrices' )
599  9995 format( / 1x, a3, ': ', a9, ' positive definite packed matrices'
600  $ )
601  9994 format( / 1x, a3, ': ', a9, ' positive definite band matrices' )
602  9993 format( / 1x, a3, ': ', a9, ' positive definite tridiagonal' )
603  9992 format( / 1x, a3, ': ', a9, ' indefinite matrices' )
604  9991 format( / 1x, a3, ': ', a9, ' indefinite packed matrices' )
605  9990 format( / 1x, a3, ': Triangular matrices' )
606  9989 format( / 1x, a3, ': Triangular packed matrices' )
607  9988 format( / 1x, a3, ': Triangular band matrices' )
608  9987 format( / 1x, a3, ': ', a2, ' factorization of general matrices'
609  $ )
610  9986 format( / 1x, a3, ': QR factorization with column pivoting' )
611  9985 format( / 1x, a3, ': RQ factorization of trapezoidal matrix' )
612  9984 format( / 1x, a3, ': Least squares driver routines' )
613  9983 format( / 1x, a3, ': LU factorization variants' )
614  9982 format( / 1x, a3, ': Cholesky factorization variants' )
615  9981 format( / 1x, a3, ': QR factorization variants' )
616  9980 format( / 1x, a3, ': No header available' )
617  8000 format( / 1x, a3, ': QRT factorization for general matrices' )
618  8001 format( / 1x, a3, ': QRT factorization for ',
619  $ 'triangular-pentagonal matrices' )
620 
621 *
622 * GE matrix types
623 *
624  9979 format( 4x, '1. Diagonal', 24x, '7. Last n/2 columns zero', / 4x,
625  $ '2. Upper triangular', 16x,
626  $ '8. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
627  $ '3. Lower triangular', 16x, '9. Random, CNDNUM = 0.1/EPS',
628  $ / 4x, '4. Random, CNDNUM = 2', 13x,
629  $ '10. Scaled near underflow', / 4x, '5. First column zero',
630  $ 14x, '11. Scaled near overflow', / 4x,
631  $ '6. Last column zero' )
632 *
633 * GB matrix types
634 *
635  9978 format( 4x, '1. Random, CNDNUM = 2', 14x,
636  $ '5. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
637  $ '2. First column zero', 15x, '6. Random, CNDNUM = .01/EPS',
638  $ / 4x, '3. Last column zero', 16x,
639  $ '7. Scaled near underflow', / 4x,
640  $ '4. Last n/2 columns zero', 11x, '8. Scaled near overflow' )
641 *
642 * GT matrix types
643 *
644  9977 format( ' Matrix types (1-6 have specified condition numbers):',
645  $ / 4x, '1. Diagonal', 24x, '7. Random, unspecified CNDNUM',
646  $ / 4x, '2. Random, CNDNUM = 2', 14x, '8. First column zero',
647  $ / 4x, '3. Random, CNDNUM = sqrt(0.1/EPS)', 2x,
648  $ '9. Last column zero', / 4x, '4. Random, CNDNUM = 0.1/EPS',
649  $ 7x, '10. Last n/2 columns zero', / 4x,
650  $ '5. Scaled near underflow', 10x,
651  $ '11. Scaled near underflow', / 4x,
652  $ '6. Scaled near overflow', 11x, '12. Scaled near overflow' )
653 *
654 * PT matrix types
655 *
656  9976 format( ' Matrix types (1-6 have specified condition numbers):',
657  $ / 4x, '1. Diagonal', 24x, '7. Random, unspecified CNDNUM',
658  $ / 4x, '2. Random, CNDNUM = 2', 14x,
659  $ '8. First row and column zero', / 4x,
660  $ '3. Random, CNDNUM = sqrt(0.1/EPS)', 2x,
661  $ '9. Last row and column zero', / 4x,
662  $ '4. Random, CNDNUM = 0.1/EPS', 7x,
663  $ '10. Middle row and column zero', / 4x,
664  $ '5. Scaled near underflow', 10x,
665  $ '11. Scaled near underflow', / 4x,
666  $ '6. Scaled near overflow', 11x, '12. Scaled near overflow' )
667 *
668 * PO, PP matrix types
669 *
670  9975 format( 4x, '1. Diagonal', 24x,
671  $ '6. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
672  $ '2. Random, CNDNUM = 2', 14x, '7. Random, CNDNUM = 0.1/EPS',
673  $ / 3x, '*3. First row and column zero', 7x,
674  $ '8. Scaled near underflow', / 3x,
675  $ '*4. Last row and column zero', 8x,
676  $ '9. Scaled near overflow', / 3x,
677  $ '*5. Middle row and column zero', / 3x,
678  $ '(* - tests error exits from ', a3,
679  $ 'TRF, no test ratios are computed)' )
680 *
681 * CH matrix types
682 *
683  9974 format( 4x, '1. Diagonal', 24x,
684  $ '6. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
685  $ '2. Random, CNDNUM = 2', 14x, '7. Random, CNDNUM = 0.1/EPS',
686  $ / 3x, '*3. First row and column zero', 7x,
687  $ '8. Scaled near underflow', / 3x,
688  $ '*4. Last row and column zero', 8x,
689  $ '9. Scaled near overflow', / 3x,
690  $ '*5. Middle row and column zero', / 3x,
691  $ '(* - tests error exits, no test ratios are computed)' )
692 *
693 * PS matrix types
694 *
695  8973 format( 4x, '1. Diagonal', / 4x, '2. Random, CNDNUM = 2', 14x,
696  $ / 3x, '*3. Nonzero eigenvalues of: D(1:RANK-1)=1 and ',
697  $ 'D(RANK) = 1.0/', a4, / 3x,
698  $ '*4. Nonzero eigenvalues of: D(1)=1 and ',
699  $ ' D(2:RANK) = 1.0/', a4, / 3x,
700  $ '*5. Nonzero eigenvalues of: D(I) = ', a4,
701  $ '**(-(I-1)/(RANK-1)) ', ' I=1:RANK', / 4x,
702  $ '6. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
703  $ '7. Random, CNDNUM = 0.1/EPS', / 4x,
704  $ '8. Scaled near underflow', / 4x, '9. Scaled near overflow',
705  $ / 3x, '(* - Semi-definite tests )' )
706  8972 format( 3x, 'RANK minus computed rank, returned by ', a, 'PSTRF' )
707 *
708 * PB matrix types
709 *
710  9973 format( 4x, '1. Random, CNDNUM = 2', 14x,
711  $ '5. Random, CNDNUM = sqrt(0.1/EPS)', / 3x,
712  $ '*2. First row and column zero', 7x,
713  $ '6. Random, CNDNUM = 0.1/EPS', / 3x,
714  $ '*3. Last row and column zero', 8x,
715  $ '7. Scaled near underflow', / 3x,
716  $ '*4. Middle row and column zero', 6x,
717  $ '8. Scaled near overflow', / 3x,
718  $ '(* - tests error exits from ', a3,
719  $ 'TRF, no test ratios are computed)' )
720 *
721 * SSY, SSP, CHE, CHP matrix types
722 *
723  9972 format( 4x, '1. Diagonal', 24x,
724  $ '6. Last n/2 rows and columns zero', / 4x,
725  $ '2. Random, CNDNUM = 2', 14x,
726  $ '7. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
727  $ '3. First row and column zero', 7x,
728  $ '8. Random, CNDNUM = 0.1/EPS', / 4x,
729  $ '4. Last row and column zero', 8x,
730  $ '9. Scaled near underflow', / 4x,
731  $ '5. Middle row and column zero', 5x,
732  $ '10. Scaled near overflow' )
733 *
734 * CSY, CSP matrix types
735 *
736  9971 format( 4x, '1. Diagonal', 24x,
737  $ '7. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
738  $ '2. Random, CNDNUM = 2', 14x, '8. Random, CNDNUM = 0.1/EPS',
739  $ / 4x, '3. First row and column zero', 7x,
740  $ '9. Scaled near underflow', / 4x,
741  $ '4. Last row and column zero', 7x,
742  $ '10. Scaled near overflow', / 4x,
743  $ '5. Middle row and column zero', 5x,
744  $ '11. Block diagonal matrix', / 4x,
745  $ '6. Last n/2 rows and columns zero' )
746 *
747 * QR matrix types
748 *
749  9970 format( 4x, '1. Diagonal', 24x,
750  $ '5. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
751  $ '2. Upper triangular', 16x, '6. Random, CNDNUM = 0.1/EPS',
752  $ / 4x, '3. Lower triangular', 16x,
753  $ '7. Scaled near underflow', / 4x, '4. Random, CNDNUM = 2',
754  $ 14x, '8. Scaled near overflow' )
755 *
756 * QP matrix types
757 *
758  9969 format( ' Matrix types (2-6 have condition 1/EPS):', / 4x,
759  $ '1. Zero matrix', 21x, '4. First n/2 columns fixed', / 4x,
760  $ '2. One small eigenvalue', 12x, '5. Last n/2 columns fixed',
761  $ / 4x, '3. Geometric distribution', 10x,
762  $ '6. Every second column fixed' )
763 *
764 * TZ matrix types
765 *
766  9968 format( ' Matrix types (2-3 have condition 1/EPS):', / 4x,
767  $ '1. Zero matrix', / 4x, '2. One small eigenvalue', / 4x,
768  $ '3. Geometric distribution' )
769 *
770 * LS matrix types
771 *
772  9967 format( ' Matrix types (1-3: full rank, 4-6: rank deficient):',
773  $ / 4x, '1 and 4. Normal scaling', / 4x,
774  $ '2 and 5. Scaled near overflow', / 4x,
775  $ '3 and 6. Scaled near underflow' )
776 *
777 * TR, TP matrix types
778 *
779  9966 format( ' Matrix types for ', a3, ' routines:', / 4x,
780  $ '1. Diagonal', 24x, '6. Scaled near overflow', / 4x,
781  $ '2. Random, CNDNUM = 2', 14x, '7. Identity', / 4x,
782  $ '3. Random, CNDNUM = sqrt(0.1/EPS) ',
783  $ '8. Unit triangular, CNDNUM = 2', / 4x,
784  $ '4. Random, CNDNUM = 0.1/EPS', 8x,
785  $ '9. Unit, CNDNUM = sqrt(0.1/EPS)', / 4x,
786  $ '5. Scaled near underflow', 10x,
787  $ '10. Unit, CNDNUM = 0.1/EPS' )
788  9965 format( ' Special types for testing ', a, ':', / 3x,
789  $ '11. Matrix elements are O(1), large right hand side', / 3x,
790  $ '12. First diagonal causes overflow,',
791  $ ' offdiagonal column norms < 1', / 3x,
792  $ '13. First diagonal causes overflow,',
793  $ ' offdiagonal column norms > 1', / 3x,
794  $ '14. Growth factor underflows, solution does not overflow',
795  $ / 3x, '15. Small diagonal causes gradual overflow', / 3x,
796  $ '16. One zero diagonal element', / 3x,
797  $ '17. Large offdiagonals cause overflow when adding a column'
798  $ , / 3x, '18. Unit triangular with large right hand side' )
799 *
800 * TB matrix types
801 *
802  9964 format( ' Matrix types for ', a3, ' routines:', / 4x,
803  $ '1. Random, CNDNUM = 2', 14x, '6. Identity', / 4x,
804  $ '2. Random, CNDNUM = sqrt(0.1/EPS) ',
805  $ '7. Unit triangular, CNDNUM = 2', / 4x,
806  $ '3. Random, CNDNUM = 0.1/EPS', 8x,
807  $ '8. Unit, CNDNUM = sqrt(0.1/EPS)', / 4x,
808  $ '4. Scaled near underflow', 11x,
809  $ '9. Unit, CNDNUM = 0.1/EPS', / 4x,
810  $ '5. Scaled near overflow' )
811  9963 format( ' Special types for testing ', a, ':', / 3x,
812  $ '10. Matrix elements are O(1), large right hand side', / 3x,
813  $ '11. First diagonal causes overflow,',
814  $ ' offdiagonal column norms < 1', / 3x,
815  $ '12. First diagonal causes overflow,',
816  $ ' offdiagonal column norms > 1', / 3x,
817  $ '13. Growth factor underflows, solution does not overflow',
818  $ / 3x, '14. Small diagonal causes gradual overflow', / 3x,
819  $ '15. One zero diagonal element', / 3x,
820  $ '16. Large offdiagonals cause overflow when adding a column'
821  $ , / 3x, '17. Unit triangular with large right hand side' )
822 *
823 * Test ratios
824 *
825  9962 format( 3x, i2, ': norm( L * U - A ) / ( N * norm(A) * EPS )' )
826  9961 format( 3x, i2, ': norm( I - A*AINV ) / ',
827  $ '( N * norm(A) * norm(AINV) * EPS )' )
828  9960 format( 3x, i2, ': norm( B - A * X ) / ',
829  $ '( norm(A) * norm(X) * EPS )' )
830  6660 format( 3x, i2, ': diagonal is not non-negative')
831  9959 format( 3x, i2, ': norm( X - XACT ) / ',
832  $ '( norm(XACT) * CNDNUM * EPS )' )
833  9958 format( 3x, i2, ': norm( X - XACT ) / ',
834  $ '( norm(XACT) * CNDNUM * EPS ), refined' )
835  9957 format( 3x, i2, ': norm( X - XACT ) / ',
836  $ '( norm(XACT) * (error bound) )' )
837  9956 format( 3x, i2, ': (backward error) / EPS' )
838  9955 format( 3x, i2, ': RCOND * CNDNUM - 1.0' )
839  9954 format( 3x, i2, ': norm( U'' * U - A ) / ( N * norm(A) * EPS )',
840  $ ', or', / 7x, 'norm( L * L'' - A ) / ( N * norm(A) * EPS )'
841  $ )
842  8950 format( 3x,
843  $ 'norm( P * U'' * U * P'' - A ) / ( N * norm(A) * EPS )',
844  $ ', or', / 3x,
845  $ 'norm( P * L * L'' * P'' - A ) / ( N * norm(A) * EPS )' )
846  9953 format( 3x, i2, ': norm( U*D*U'' - A ) / ( N * norm(A) * EPS )',
847  $ ', or', / 7x, 'norm( L*D*L'' - A ) / ( N * norm(A) * EPS )'
848  $ )
849  9952 format( 3x, i2, ': norm( U''*D*U - A ) / ( N * norm(A) * EPS )',
850  $ ', or', / 7x, 'norm( L*D*L'' - A ) / ( N * norm(A) * EPS )'
851  $ )
852  9951 format( ' Test ratio for ', a, ':', / 3x, i2,
853  $ ': norm( s*b - A*x ) / ( norm(A) * norm(x) * EPS )' )
854  9950 format( 3x, i2, ': norm( R - Q'' * A ) / ( M * norm(A) * EPS )' )
855  6950 format( 3x, i2, ': norm( R - Q'
856 ' * A ) / ( M * norm(A) * EPS ) $ [RFPG]' )
857  9949 format( 3x, i2, ': norm( L - A * Q'' ) / ( N * norm(A) * EPS )' )
858  9948 format( 3x, i2, ': norm( L - Q'' * A ) / ( M * norm(A) * EPS )' )
859  9947 format( 3x, i2, ': norm( R - A * Q'' ) / ( N * norm(A) * EPS )' )
860  9946 format( 3x, i2, ': norm( I - Q''*Q ) / ( M * EPS )' )
861  9945 format( 3x, i2, ': norm( I - Q*Q'' ) / ( N * EPS )' )
862  9944 format( 3x, i2, ': norm( Q*C - Q*C ) / ', '( ', a1,
863  $ ' * norm(C) * EPS )' )
864  9943 format( 3x, i2, ': norm( C*Q - C*Q ) / ', '( ', a1,
865  $ ' * norm(C) * EPS )' )
866  9942 format( 3x, i2, ': norm( Q''*C - Q''*C )/ ', '( ', a1,
867  $ ' * norm(C) * EPS )' )
868  9941 format( 3x, i2, ': norm( C*Q'' - C*Q'' )/ ', '( ', a1,
869  $ ' * norm(C) * EPS )' )
870  9940 format( 3x, i2, ': norm(svd(A) - svd(R)) / ',
871  $ '( M * norm(svd(R)) * EPS )' )
872  9939 format( 3x, i2, ': norm( A*P - Q*R ) / ( M * norm(A) * EPS )'
873  $ )
874  9938 format( 3x, i2, ': norm( I - Q''*Q ) / ( M * EPS )' )
875  9937 format( 3x, i2, ': norm( A - R*Q ) / ( M * norm(A) * EPS )'
876  $ )
877  9935 format( 3x, i2, ': norm( B - A * X ) / ',
878  $ '( max(M,N) * norm(A) * norm(X) * EPS )' )
879  9934 format( 3x, i2, ': norm( (A*X-B)'' *A ) / ',
880  $ '( max(M,N,NRHS) * norm(A) * norm(B) * EPS )' )
881  9933 format( 3x, i2, ': norm(svd(A)-svd(R)) / ',
882  $ '( min(M,N) * norm(svd(R)) * EPS )' )
883  9932 format( 3x, i2, ': Check if X is in the row space of A or A''' )
884  9931 format( 3x, i2, ': norm( (A*X-B)'' *A ) / ',
885  $ '( max(M,N,NRHS) * norm(A) * norm(B) * EPS )', / 7x,
886  $ 'if TRANS=''N'.GE.' and MN or TRANS=''T'.LT.' and MN, ',
887  $ 'otherwise', / 7x,
888  $ 'check if X is in the row space of A or A'' ',
889  $ '(overdetermined case)' )
890  9929 format( ' Test ratios (1-3: ', a1, 'TZRQF, 4-6: ', a1,
891  $ 'TZRZF):' )
892  9920 format( 3x, ' 7-10: same as 3-6', 3x, ' 11-14: same as 3-6',
893  $ 3x, ' 15-18: same as 3-6' )
894  9921 format( ' Test ratios:', / ' (1-2: ', a1, 'GELS, 3-6: ', a1,
895  $ 'GELSX, 7-10: ', a1, 'GELSY, 11-14: ', a1, 'GELSS, 15-18: ',
896  $ a1, 'GELSD)' )
897  8011 format(3x,i2,': norm( R - Q''*A ) / ( M * norm(A) * EPS )' )
898  8012 format(3x,i2,': norm( I - Q''*Q ) / ( M * EPS )' )
899  8013 format(3x,i2,': norm( Q*C - Q*C ) / ( M * norm(C) * EPS )' )
900  8014 format(3x,i2,': norm( Q''*C - Q''*C ) / ( M * norm(C) * EPS )')
901  8015 format(3x,i2,': norm( C*Q - C*Q ) / ( M * norm(C) * EPS )' )
902  8016 format(3x,i2,': norm( C*Q'' - C*Q'' ) / ( M * norm(C) * EPS )')
903  8017 format(3x,i2,': norm( R - Q''*A ) / ( (M+N) * norm(A) * EPS )' )
904  8018 format(3x,i2,': norm( I - Q''*Q ) / ( (M+N) * EPS )' )
905  8019 format(3x,i2,': norm( Q*C - Q*C ) / ( (M+N) * norm(C) * EPS )' )
906  8020 format(3x,i2,
907  $ ': norm( Q''*C - Q''*C ) / ( (M+N) * norm(C) * EPS )')
908  8021 format(3x,i2,': norm( C*Q - C*Q ) / ( (M+N) * norm(C) * EPS )' )
909  8022 format(3x,i2,
910  $ ': norm( C*Q'' - C*Q'' ) / ( (M+N) * norm(C) * EPS )')
911 *
912  return
913 *
914 * End of ALAHD
915 *
916  END