LAPACK  3.9.0
LAPACK: Linear Algebra PACKage
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 *> with partial (Bunch-Kaufman) pivoting
52 *> _SR: Symmetric indefinite,
53 *> with rook (bounded Bunch-Kaufman) pivoting
54 *> _SK: Symmetric indefinite,
55 *> with rook (bounded Bunch-Kaufman) pivoting
56 *> ( new storage format for factors:
57 *> L and diagonal of D is stored in A,
58 *> subdiagonal of D is stored in E )
59 *> _SP: Symmetric indefinite packed,
60 *> with partial (Bunch-Kaufman) pivoting
61 *> _HA: (complex) Hermitian ,
62 *> with Aasen Algorithm
63 *> _HE: (complex) Hermitian indefinite,
64 *> with partial (Bunch-Kaufman) pivoting
65 *> _HR: (complex) Hermitian indefinite,
66 *> with rook (bounded Bunch-Kaufman) pivoting
67 *> _HK: (complex) Hermitian indefinite,
68 *> with rook (bounded Bunch-Kaufman) pivoting
69 *> ( new storage format for factors:
70 *> L and diagonal of D is stored in A,
71 *> subdiagonal of D is stored in E )
72 *> _HP: (complex) Hermitian indefinite packed,
73 *> with partial (Bunch-Kaufman) pivoting
74 *> _TR: Triangular
75 *> _TP: Triangular packed
76 *> _TB: Triangular band
77 *> _QR: QR (general matrices)
78 *> _LQ: LQ (general matrices)
79 *> _QL: QL (general matrices)
80 *> _RQ: RQ (general matrices)
81 *> _QP: QR with column pivoting
82 *> _TZ: Trapezoidal
83 *> _LS: Least Squares driver routines
84 *> _LU: LU variants
85 *> _CH: Cholesky variants
86 *> _QS: QR variants
87 *> _QT: QRT (general matrices)
88 *> _QX: QRT (triangular-pentagonal matrices)
89 *> _TS: QR routines for tall-skinny and short-wide matrices
90 *> _HH: Householder reconstruction for tall-skinny matrices
91 *> The first character must be one of S, D, C, or Z (C or Z only
92 *> if complex).
93 *> \endverbatim
94 *
95 * Authors:
96 * ========
97 *
98 *> \author Univ. of Tennessee
99 *> \author Univ. of California Berkeley
100 *> \author Univ. of Colorado Denver
101 *> \author NAG Ltd.
102 *
103 *> \date June 2019
104 *
105 *> \ingroup aux_lin
106 *
107 * =====================================================================
108  SUBROUTINE alahd( IOUNIT, PATH )
109 *
110 * -- LAPACK test routine (version 3.9.0) --
111 * -- LAPACK is a software package provided by Univ. of Tennessee, --
112 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
113 * June 2019
114 *
115 * .. Scalar Arguments ..
116  CHARACTER*3 PATH
117  INTEGER IOUNIT
118 * ..
119 *
120 * =====================================================================
121 *
122 * .. Local Scalars ..
123  LOGICAL CORZ, SORD
124  CHARACTER C1, C3
125  CHARACTER*2 P2
126  CHARACTER*4 EIGCNM
127  CHARACTER*32 SUBNAM
128  CHARACTER*9 SYM
129 * ..
130 * .. External Functions ..
131  LOGICAL LSAME, LSAMEN
132  EXTERNAL lsame, lsamen
133 * ..
134 * .. Intrinsic Functions ..
135  INTRINSIC len_trim
136 * ..
137 * .. Executable Statements ..
138 *
139  IF( iounit.LE.0 )
140  $ RETURN
141  c1 = path( 1: 1 )
142  c3 = path( 3: 3 )
143  p2 = path( 2: 3 )
144  sord = lsame( c1, 'S' ) .OR. lsame( c1, 'D' )
145  corz = lsame( c1, 'C' ) .OR. lsame( c1, 'Z' )
146  IF( .NOT.( sord .OR. corz ) )
147  $ RETURN
148 *
149  IF( lsamen( 2, p2, 'GE' ) ) THEN
150 *
151 * GE: General dense
152 *
153  WRITE( iounit, fmt = 9999 )path
154  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
155  WRITE( iounit, fmt = 9979 )
156  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
157  WRITE( iounit, fmt = 9962 )1
158  WRITE( iounit, fmt = 9961 )2
159  WRITE( iounit, fmt = 9960 )3
160  WRITE( iounit, fmt = 9959 )4
161  WRITE( iounit, fmt = 9958 )5
162  WRITE( iounit, fmt = 9957 )6
163  WRITE( iounit, fmt = 9956 )7
164  WRITE( iounit, fmt = 9955 )8
165  WRITE( iounit, fmt = '( '' Messages:'' )' )
166 *
167  ELSE IF( lsamen( 2, p2, 'GB' ) ) THEN
168 *
169 * GB: General band
170 *
171  WRITE( iounit, fmt = 9998 )path
172  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
173  WRITE( iounit, fmt = 9978 )
174  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
175  WRITE( iounit, fmt = 9962 )1
176  WRITE( iounit, fmt = 9960 )2
177  WRITE( iounit, fmt = 9959 )3
178  WRITE( iounit, fmt = 9958 )4
179  WRITE( iounit, fmt = 9957 )5
180  WRITE( iounit, fmt = 9956 )6
181  WRITE( iounit, fmt = 9955 )7
182  WRITE( iounit, fmt = '( '' Messages:'' )' )
183 *
184  ELSE IF( lsamen( 2, p2, 'GT' ) ) THEN
185 *
186 * GT: General tridiagonal
187 *
188  WRITE( iounit, fmt = 9997 )path
189  WRITE( iounit, fmt = 9977 )
190  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
191  WRITE( iounit, fmt = 9962 )1
192  WRITE( iounit, fmt = 9960 )2
193  WRITE( iounit, fmt = 9959 )3
194  WRITE( iounit, fmt = 9958 )4
195  WRITE( iounit, fmt = 9957 )5
196  WRITE( iounit, fmt = 9956 )6
197  WRITE( iounit, fmt = 9955 )7
198  WRITE( iounit, fmt = '( '' Messages:'' )' )
199 *
200  ELSE IF( lsamen( 2, p2, 'PO' ) .OR. lsamen( 2, p2, 'PP' ) ) THEN
201 *
202 * PO: Positive definite full
203 * PP: Positive definite packed
204 *
205  IF( sord ) THEN
206  sym = 'Symmetric'
207  ELSE
208  sym = 'Hermitian'
209  END IF
210  IF( lsame( c3, 'O' ) ) THEN
211  WRITE( iounit, fmt = 9996 )path, sym
212  ELSE
213  WRITE( iounit, fmt = 9995 )path, sym
214  END IF
215  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
216  WRITE( iounit, fmt = 9975 )path
217  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
218  WRITE( iounit, fmt = 9954 )1
219  WRITE( iounit, fmt = 9961 )2
220  WRITE( iounit, fmt = 9960 )3
221  WRITE( iounit, fmt = 9959 )4
222  WRITE( iounit, fmt = 9958 )5
223  WRITE( iounit, fmt = 9957 )6
224  WRITE( iounit, fmt = 9956 )7
225  WRITE( iounit, fmt = 9955 )8
226  WRITE( iounit, fmt = '( '' Messages:'' )' )
227 *
228  ELSE IF( lsamen( 2, p2, 'PS' ) ) THEN
229 *
230 * PS: Positive semi-definite full
231 *
232  IF( sord ) THEN
233  sym = 'Symmetric'
234  ELSE
235  sym = 'Hermitian'
236  END IF
237  IF( lsame( c1, 'S' ) .OR. lsame( c1, 'C' ) ) THEN
238  eigcnm = '1E04'
239  ELSE
240  eigcnm = '1D12'
241  END IF
242  WRITE( iounit, fmt = 9995 )path, sym
243  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
244  WRITE( iounit, fmt = 8973 )eigcnm, eigcnm, eigcnm
245  WRITE( iounit, fmt = '( '' Difference:'' )' )
246  WRITE( iounit, fmt = 8972 )c1
247  WRITE( iounit, fmt = '( '' Test ratio:'' )' )
248  WRITE( iounit, fmt = 8950 )
249  WRITE( iounit, fmt = '( '' Messages:'' )' )
250  ELSE IF( lsamen( 2, p2, 'PB' ) ) THEN
251 *
252 * PB: Positive definite band
253 *
254  IF( sord ) THEN
255  WRITE( iounit, fmt = 9994 )path, 'Symmetric'
256  ELSE
257  WRITE( iounit, fmt = 9994 )path, 'Hermitian'
258  END IF
259  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
260  WRITE( iounit, fmt = 9973 )path
261  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
262  WRITE( iounit, fmt = 9954 )1
263  WRITE( iounit, fmt = 9960 )2
264  WRITE( iounit, fmt = 9959 )3
265  WRITE( iounit, fmt = 9958 )4
266  WRITE( iounit, fmt = 9957 )5
267  WRITE( iounit, fmt = 9956 )6
268  WRITE( iounit, fmt = 9955 )7
269  WRITE( iounit, fmt = '( '' Messages:'' )' )
270 *
271  ELSE IF( lsamen( 2, p2, 'PT' ) ) THEN
272 *
273 * PT: Positive definite tridiagonal
274 *
275  IF( sord ) THEN
276  WRITE( iounit, fmt = 9993 )path, 'Symmetric'
277  ELSE
278  WRITE( iounit, fmt = 9993 )path, 'Hermitian'
279  END IF
280  WRITE( iounit, fmt = 9976 )
281  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
282  WRITE( iounit, fmt = 9952 )1
283  WRITE( iounit, fmt = 9960 )2
284  WRITE( iounit, fmt = 9959 )3
285  WRITE( iounit, fmt = 9958 )4
286  WRITE( iounit, fmt = 9957 )5
287  WRITE( iounit, fmt = 9956 )6
288  WRITE( iounit, fmt = 9955 )7
289  WRITE( iounit, fmt = '( '' Messages:'' )' )
290 *
291  ELSE IF( lsamen( 2, p2, 'SY' ) ) THEN
292 *
293 * SY: Symmetric indefinite full,
294 * with partial (Bunch-Kaufman) pivoting algorithm
295 *
296  IF( lsame( c3, 'Y' ) ) THEN
297  WRITE( iounit, fmt = 9992 )path, 'Symmetric'
298  ELSE
299  WRITE( iounit, fmt = 9991 )path, 'Symmetric'
300  END IF
301  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
302  IF( sord ) THEN
303  WRITE( iounit, fmt = 9972 )
304  ELSE
305  WRITE( iounit, fmt = 9971 )
306  END IF
307  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
308  WRITE( iounit, fmt = 9953 )1
309  WRITE( iounit, fmt = 9961 )2
310  WRITE( iounit, fmt = 9960 )3
311  WRITE( iounit, fmt = 9960 )4
312  WRITE( iounit, fmt = 9959 )5
313  WRITE( iounit, fmt = 9958 )6
314  WRITE( iounit, fmt = 9956 )7
315  WRITE( iounit, fmt = 9957 )8
316  WRITE( iounit, fmt = 9955 )9
317  WRITE( iounit, fmt = '( '' Messages:'' )' )
318 *
319  ELSE IF( lsamen( 2, p2, 'SR' ) .OR. lsamen( 2, p2, 'SK') ) THEN
320 *
321 * SR: Symmetric indefinite full,
322 * with rook (bounded Bunch-Kaufman) pivoting algorithm
323 *
324 * SK: Symmetric indefinite full,
325 * with rook (bounded Bunch-Kaufman) pivoting algorithm,
326 * ( new storage format for factors:
327 * L and diagonal of D is stored in A,
328 * subdiagonal of D is stored in E )
329 *
330  WRITE( iounit, fmt = 9892 )path, 'Symmetric'
331 *
332  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
333  IF( sord ) THEN
334  WRITE( iounit, fmt = 9972 )
335  ELSE
336  WRITE( iounit, fmt = 9971 )
337  END IF
338 *
339  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
340  WRITE( iounit, fmt = 9953 )1
341  WRITE( iounit, fmt = 9961 )2
342  WRITE( iounit, fmt = 9927 )3
343  WRITE( iounit, fmt = 9928 )
344  WRITE( iounit, fmt = 9926 )4
345  WRITE( iounit, fmt = 9928 )
346  WRITE( iounit, fmt = 9960 )5
347  WRITE( iounit, fmt = 9959 )6
348  WRITE( iounit, fmt = 9955 )7
349  WRITE( iounit, fmt = '( '' Messages:'' )' )
350 *
351  ELSE IF( lsamen( 2, p2, 'SP' ) ) THEN
352 *
353 * SP: Symmetric indefinite packed,
354 * with partial (Bunch-Kaufman) pivoting algorithm
355 *
356  IF( lsame( c3, 'Y' ) ) THEN
357  WRITE( iounit, fmt = 9992 )path, 'Symmetric'
358  ELSE
359  WRITE( iounit, fmt = 9991 )path, 'Symmetric'
360  END IF
361  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
362  IF( sord ) THEN
363  WRITE( iounit, fmt = 9972 )
364  ELSE
365  WRITE( iounit, fmt = 9971 )
366  END IF
367  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
368  WRITE( iounit, fmt = 9953 )1
369  WRITE( iounit, fmt = 9961 )2
370  WRITE( iounit, fmt = 9960 )3
371  WRITE( iounit, fmt = 9959 )4
372  WRITE( iounit, fmt = 9958 )5
373  WRITE( iounit, fmt = 9956 )6
374  WRITE( iounit, fmt = 9957 )7
375  WRITE( iounit, fmt = 9955 )8
376  WRITE( iounit, fmt = '( '' Messages:'' )' )
377 *
378  ELSE IF( lsamen( 2, p2, 'HA' ) ) THEN
379 *
380 * HA: Hermitian,
381 * with Assen Algorithm
382 *
383  WRITE( iounit, fmt = 9992 )path, 'Hermitian'
384 *
385  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
386  WRITE( iounit, fmt = 9972 )
387 *
388  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
389  WRITE( iounit, fmt = 9953 )1
390  WRITE( iounit, fmt = 9961 )2
391  WRITE( iounit, fmt = 9960 )3
392  WRITE( iounit, fmt = 9960 )4
393  WRITE( iounit, fmt = 9959 )5
394  WRITE( iounit, fmt = 9958 )6
395  WRITE( iounit, fmt = 9956 )7
396  WRITE( iounit, fmt = 9957 )8
397  WRITE( iounit, fmt = 9955 )9
398  WRITE( iounit, fmt = '( '' Messages:'' )' )
399 *
400  ELSE IF( lsamen( 2, p2, 'HE' ) ) THEN
401 *
402 * HE: Hermitian indefinite full,
403 * with partial (Bunch-Kaufman) pivoting algorithm
404 *
405  WRITE( iounit, fmt = 9992 )path, 'Hermitian'
406 *
407  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
408  WRITE( iounit, fmt = 9972 )
409 *
410  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
411  WRITE( iounit, fmt = 9953 )1
412  WRITE( iounit, fmt = 9961 )2
413  WRITE( iounit, fmt = 9960 )3
414  WRITE( iounit, fmt = 9960 )4
415  WRITE( iounit, fmt = 9959 )5
416  WRITE( iounit, fmt = 9958 )6
417  WRITE( iounit, fmt = 9956 )7
418  WRITE( iounit, fmt = 9957 )8
419  WRITE( iounit, fmt = 9955 )9
420  WRITE( iounit, fmt = '( '' Messages:'' )' )
421 *
422  ELSE IF( lsamen( 2, p2, 'HR' ) .OR. lsamen( 2, p2, 'HR' ) ) THEN
423 *
424 * HR: Hermitian indefinite full,
425 * with rook (bounded Bunch-Kaufman) pivoting algorithm
426 *
427 * HK: Hermitian indefinite full,
428 * with rook (bounded Bunch-Kaufman) pivoting algorithm,
429 * ( new storage format for factors:
430 * L and diagonal of D is stored in A,
431 * subdiagonal of D is stored in E )
432 *
433  WRITE( iounit, fmt = 9892 )path, 'Hermitian'
434 *
435  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
436  WRITE( iounit, fmt = 9972 )
437 *
438  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
439  WRITE( iounit, fmt = 9953 )1
440  WRITE( iounit, fmt = 9961 )2
441  WRITE( iounit, fmt = 9927 )3
442  WRITE( iounit, fmt = 9928 )
443  WRITE( iounit, fmt = 9926 )4
444  WRITE( iounit, fmt = 9928 )
445  WRITE( iounit, fmt = 9960 )5
446  WRITE( iounit, fmt = 9959 )6
447  WRITE( iounit, fmt = 9955 )7
448  WRITE( iounit, fmt = '( '' Messages:'' )' )
449 *
450  ELSE IF( lsamen( 2, p2, 'HP' ) ) THEN
451 *
452 * HP: Hermitian indefinite packed,
453 * with partial (Bunch-Kaufman) pivoting algorithm
454 *
455  IF( lsame( c3, 'E' ) ) THEN
456  WRITE( iounit, fmt = 9992 )path, 'Hermitian'
457  ELSE
458  WRITE( iounit, fmt = 9991 )path, 'Hermitian'
459  END IF
460  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
461  WRITE( iounit, fmt = 9972 )
462  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
463  WRITE( iounit, fmt = 9953 )1
464  WRITE( iounit, fmt = 9961 )2
465  WRITE( iounit, fmt = 9960 )3
466  WRITE( iounit, fmt = 9959 )4
467  WRITE( iounit, fmt = 9958 )5
468  WRITE( iounit, fmt = 9956 )6
469  WRITE( iounit, fmt = 9957 )7
470  WRITE( iounit, fmt = 9955 )8
471  WRITE( iounit, fmt = '( '' Messages:'' )' )
472 *
473  ELSE IF( lsamen( 2, p2, 'TR' ) .OR. lsamen( 2, p2, 'TP' ) ) THEN
474 *
475 * TR: Triangular full
476 * TP: Triangular packed
477 *
478  IF( lsame( c3, 'R' ) ) THEN
479  WRITE( iounit, fmt = 9990 )path
480  subnam = path( 1: 1 ) // 'LATRS'
481  ELSE
482  WRITE( iounit, fmt = 9989 )path
483  subnam = path( 1: 1 ) // 'LATPS'
484  END IF
485  WRITE( iounit, fmt = 9966 )path
486  WRITE( iounit, fmt = 9965 )subnam(1:len_trim( subnam ))
487  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
488  WRITE( iounit, fmt = 9961 )1
489  WRITE( iounit, fmt = 9960 )2
490  WRITE( iounit, fmt = 9959 )3
491  WRITE( iounit, fmt = 9958 )4
492  WRITE( iounit, fmt = 9957 )5
493  WRITE( iounit, fmt = 9956 )6
494  WRITE( iounit, fmt = 9955 )7
495  WRITE( iounit, fmt = 9951 )subnam(1:len_trim( subnam )), 8
496  WRITE( iounit, fmt = '( '' Messages:'' )' )
497 *
498  ELSE IF( lsamen( 2, p2, 'TB' ) ) THEN
499 *
500 * TB: Triangular band
501 *
502  WRITE( iounit, fmt = 9988 )path
503  subnam = path( 1: 1 ) // 'LATBS'
504  WRITE( iounit, fmt = 9964 )path
505  WRITE( iounit, fmt = 9963 )subnam(1:len_trim( subnam ))
506  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
507  WRITE( iounit, fmt = 9960 )1
508  WRITE( iounit, fmt = 9959 )2
509  WRITE( iounit, fmt = 9958 )3
510  WRITE( iounit, fmt = 9957 )4
511  WRITE( iounit, fmt = 9956 )5
512  WRITE( iounit, fmt = 9955 )6
513  WRITE( iounit, fmt = 9951 )subnam(1:len_trim( subnam )), 7
514  WRITE( iounit, fmt = '( '' Messages:'' )' )
515 *
516  ELSE IF( lsamen( 2, p2, 'QR' ) ) THEN
517 *
518 * QR decomposition of rectangular matrices
519 *
520  WRITE( iounit, fmt = 9987 )path, 'QR'
521  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
522  WRITE( iounit, fmt = 9970 )
523  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
524  WRITE( iounit, fmt = 9950 )1
525  WRITE( iounit, fmt = 6950 )8
526  WRITE( iounit, fmt = 9946 )2
527  WRITE( iounit, fmt = 9944 )3, 'M'
528  WRITE( iounit, fmt = 9943 )4, 'M'
529  WRITE( iounit, fmt = 9942 )5, 'M'
530  WRITE( iounit, fmt = 9941 )6, 'M'
531  WRITE( iounit, fmt = 9960 )7
532  WRITE( iounit, fmt = 6660 )9
533  WRITE( iounit, fmt = '( '' Messages:'' )' )
534 *
535  ELSE IF( lsamen( 2, p2, 'LQ' ) ) THEN
536 *
537 * LQ decomposition of rectangular matrices
538 *
539  WRITE( iounit, fmt = 9987 )path, 'LQ'
540  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
541  WRITE( iounit, fmt = 9970 )
542  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
543  WRITE( iounit, fmt = 9949 )1
544  WRITE( iounit, fmt = 9945 )2
545  WRITE( iounit, fmt = 9944 )3, 'N'
546  WRITE( iounit, fmt = 9943 )4, 'N'
547  WRITE( iounit, fmt = 9942 )5, 'N'
548  WRITE( iounit, fmt = 9941 )6, 'N'
549  WRITE( iounit, fmt = 9960 )7
550  WRITE( iounit, fmt = '( '' Messages:'' )' )
551 *
552  ELSE IF( lsamen( 2, p2, 'QL' ) ) THEN
553 *
554 * QL decomposition of rectangular matrices
555 *
556  WRITE( iounit, fmt = 9987 )path, 'QL'
557  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
558  WRITE( iounit, fmt = 9970 )
559  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
560  WRITE( iounit, fmt = 9948 )1
561  WRITE( iounit, fmt = 9946 )2
562  WRITE( iounit, fmt = 9944 )3, 'M'
563  WRITE( iounit, fmt = 9943 )4, 'M'
564  WRITE( iounit, fmt = 9942 )5, 'M'
565  WRITE( iounit, fmt = 9941 )6, 'M'
566  WRITE( iounit, fmt = 9960 )7
567  WRITE( iounit, fmt = '( '' Messages:'' )' )
568 *
569  ELSE IF( lsamen( 2, p2, 'RQ' ) ) THEN
570 *
571 * RQ decomposition of rectangular matrices
572 *
573  WRITE( iounit, fmt = 9987 )path, 'RQ'
574  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
575  WRITE( iounit, fmt = 9970 )
576  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
577  WRITE( iounit, fmt = 9947 )1
578  WRITE( iounit, fmt = 9945 )2
579  WRITE( iounit, fmt = 9944 )3, 'N'
580  WRITE( iounit, fmt = 9943 )4, 'N'
581  WRITE( iounit, fmt = 9942 )5, 'N'
582  WRITE( iounit, fmt = 9941 )6, 'N'
583  WRITE( iounit, fmt = 9960 )7
584  WRITE( iounit, fmt = '( '' Messages:'' )' )
585 *
586  ELSE IF( lsamen( 2, p2, 'QP' ) ) THEN
587 *
588 * QR decomposition with column pivoting
589 *
590  WRITE( iounit, fmt = 9986 )path
591  WRITE( iounit, fmt = 9969 )
592  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
593  WRITE( iounit, fmt = 9940 )1
594  WRITE( iounit, fmt = 9939 )2
595  WRITE( iounit, fmt = 9938 )3
596  WRITE( iounit, fmt = '( '' Messages:'' )' )
597 *
598  ELSE IF( lsamen( 2, p2, 'TZ' ) ) THEN
599 *
600 * TZ: Trapezoidal
601 *
602  WRITE( iounit, fmt = 9985 )path
603  WRITE( iounit, fmt = 9968 )
604  WRITE( iounit, fmt = 9929 )c1
605  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
606  WRITE( iounit, fmt = 9940 )1
607  WRITE( iounit, fmt = 9937 )2
608  WRITE( iounit, fmt = 9938 )3
609  WRITE( iounit, fmt = '( '' Messages:'' )' )
610 *
611  ELSE IF( lsamen( 2, p2, 'LS' ) ) THEN
612 *
613 * LS: Least Squares driver routines for
614 * LS, LSD, LSS, LSX and LSY.
615 *
616  WRITE( iounit, fmt = 9984 )path
617  WRITE( iounit, fmt = 9967 )
618  WRITE( iounit, fmt = 9921 )c1, c1, c1, c1
619  WRITE( iounit, fmt = 9935 )1
620  WRITE( iounit, fmt = 9931 )2
621  WRITE( iounit, fmt = 9933 )3
622  WRITE( iounit, fmt = 9935 )4
623  WRITE( iounit, fmt = 9934 )5
624  WRITE( iounit, fmt = 9932 )6
625  WRITE( iounit, fmt = 9920 )
626  WRITE( iounit, fmt = '( '' Messages:'' )' )
627 *
628  ELSE IF( lsamen( 2, p2, 'LU' ) ) THEN
629 *
630 * LU factorization variants
631 *
632  WRITE( iounit, fmt = 9983 )path
633  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
634  WRITE( iounit, fmt = 9979 )
635  WRITE( iounit, fmt = '( '' Test ratio:'' )' )
636  WRITE( iounit, fmt = 9962 )1
637  WRITE( iounit, fmt = '( '' Messages:'' )' )
638 *
639  ELSE IF( lsamen( 2, p2, 'CH' ) ) THEN
640 *
641 * Cholesky factorization variants
642 *
643  WRITE( iounit, fmt = 9982 )path
644  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
645  WRITE( iounit, fmt = 9974 )
646  WRITE( iounit, fmt = '( '' Test ratio:'' )' )
647  WRITE( iounit, fmt = 9954 )1
648  WRITE( iounit, fmt = '( '' Messages:'' )' )
649 *
650  ELSE IF( lsamen( 2, p2, 'QS' ) ) THEN
651 *
652 * QR factorization variants
653 *
654  WRITE( iounit, fmt = 9981 )path
655  WRITE( iounit, fmt = '( '' Matrix types:'' )' )
656  WRITE( iounit, fmt = 9970 )
657  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
658 *
659  ELSE IF( lsamen( 2, p2, 'QT' ) ) THEN
660 *
661 * QRT (general matrices)
662 *
663  WRITE( iounit, fmt = 8000 ) path
664  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
665  WRITE( iounit, fmt = 8011 ) 1
666  WRITE( iounit, fmt = 8012 ) 2
667  WRITE( iounit, fmt = 8013 ) 3
668  WRITE( iounit, fmt = 8014 ) 4
669  WRITE( iounit, fmt = 8015 ) 5
670  WRITE( iounit, fmt = 8016 ) 6
671 *
672  ELSE IF( lsamen( 2, p2, 'QX' ) ) THEN
673 *
674 * QRT (triangular-pentagonal)
675 *
676  WRITE( iounit, fmt = 8001 ) path
677  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
678  WRITE( iounit, fmt = 8017 ) 1
679  WRITE( iounit, fmt = 8018 ) 2
680  WRITE( iounit, fmt = 8019 ) 3
681  WRITE( iounit, fmt = 8020 ) 4
682  WRITE( iounit, fmt = 8021 ) 5
683  WRITE( iounit, fmt = 8022 ) 6
684 *
685  ELSE IF( lsamen( 2, p2, 'TQ' ) ) THEN
686 *
687 * QRT (triangular-pentagonal)
688 *
689  WRITE( iounit, fmt = 8002 ) path
690  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
691  WRITE( iounit, fmt = 8023 ) 1
692  WRITE( iounit, fmt = 8024 ) 2
693  WRITE( iounit, fmt = 8025 ) 3
694  WRITE( iounit, fmt = 8026 ) 4
695  WRITE( iounit, fmt = 8027 ) 5
696  WRITE( iounit, fmt = 8028 ) 6
697 *
698  ELSE IF( lsamen( 2, p2, 'XQ' ) ) THEN
699 *
700 * QRT (triangular-pentagonal)
701 *
702  WRITE( iounit, fmt = 8003 ) path
703  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
704  WRITE( iounit, fmt = 8029 ) 1
705  WRITE( iounit, fmt = 8030 ) 2
706  WRITE( iounit, fmt = 8031 ) 3
707  WRITE( iounit, fmt = 8032 ) 4
708  WRITE( iounit, fmt = 8033 ) 5
709  WRITE( iounit, fmt = 8034 ) 6
710 *
711  ELSE IF( lsamen( 2, p2, 'TS' ) ) THEN
712 *
713 * TS: QR routines for tall-skinny and short-wide matrices
714 *
715  WRITE( iounit, fmt = 8004 ) path
716  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
717  WRITE( iounit, fmt = 8035 ) 1
718  WRITE( iounit, fmt = 8036 ) 2
719  WRITE( iounit, fmt = 8037 ) 3
720  WRITE( iounit, fmt = 8038 ) 4
721  WRITE( iounit, fmt = 8039 ) 5
722  WRITE( iounit, fmt = 8040 ) 6
723 *
724  ELSE IF( lsamen( 2, p2, 'HH' ) ) THEN
725 *
726 * HH: Householder reconstruction for tall-skinny matrices
727 *
728  WRITE( iounit, fmt = 8005 ) path
729  WRITE( iounit, fmt = '( '' Test ratios:'' )' )
730  WRITE( iounit, fmt = 8050 ) 1
731  WRITE( iounit, fmt = 8051 ) 2
732  WRITE( iounit, fmt = 8052 ) 3
733  WRITE( iounit, fmt = 8053 ) 4
734  WRITE( iounit, fmt = 8054 ) 5
735  WRITE( iounit, fmt = 8055 ) 6
736 *
737  ELSE
738 *
739 * Print error message if no header is available.
740 *
741  WRITE( iounit, fmt = 9980 )path
742  END IF
743 *
744 * First line of header
745 *
746  9999 FORMAT( / 1x, a3, ': General dense matrices' )
747  9998 FORMAT( / 1x, a3, ': General band matrices' )
748  9997 FORMAT( / 1x, a3, ': General tridiagonal' )
749  9996 FORMAT( / 1x, a3, ': ', a9, ' positive definite matrices' )
750  9995 FORMAT( / 1x, a3, ': ', a9, ' positive definite packed matrices'
751  $ )
752  9994 FORMAT( / 1x, a3, ': ', a9, ' positive definite band matrices' )
753  9993 FORMAT( / 1x, a3, ': ', a9, ' positive definite tridiagonal' )
754  9992 FORMAT( / 1x, a3, ': ', a9, ' indefinite matrices',
755  $ ', partial (Bunch-Kaufman) pivoting' )
756  9991 FORMAT( / 1x, a3, ': ', a9, ' indefinite packed matrices',
757  $ ', partial (Bunch-Kaufman) pivoting' )
758  9892 FORMAT( / 1x, a3, ': ', a9, ' indefinite matrices',
759  $ ', "rook" (bounded Bunch-Kaufman) pivoting' )
760  9891 FORMAT( / 1x, a3, ': ', a9, ' indefinite packed matrices',
761  $ ', "rook" (bounded Bunch-Kaufman) pivoting' )
762  9990 FORMAT( / 1x, a3, ': Triangular matrices' )
763  9989 FORMAT( / 1x, a3, ': Triangular packed matrices' )
764  9988 FORMAT( / 1x, a3, ': Triangular band matrices' )
765  9987 FORMAT( / 1x, a3, ': ', a2, ' factorization of general matrices'
766  $ )
767  9986 FORMAT( / 1x, a3, ': QR factorization with column pivoting' )
768  9985 FORMAT( / 1x, a3, ': RQ factorization of trapezoidal matrix' )
769  9984 FORMAT( / 1x, a3, ': Least squares driver routines' )
770  9983 FORMAT( / 1x, a3, ': LU factorization variants' )
771  9982 FORMAT( / 1x, a3, ': Cholesky factorization variants' )
772  9981 FORMAT( / 1x, a3, ': QR factorization variants' )
773  9980 FORMAT( / 1x, a3, ': No header available' )
774  8000 FORMAT( / 1x, a3, ': QRT factorization for general matrices' )
775  8001 FORMAT( / 1x, a3, ': QRT factorization for ',
776  $ 'triangular-pentagonal matrices' )
777  8002 FORMAT( / 1x, a3, ': LQT factorization for general matrices' )
778  8003 FORMAT( / 1x, a3, ': LQT factorization for ',
779  $ 'triangular-pentagonal matrices' )
780  8004 FORMAT( / 1x, a3, ': TS factorization for ',
781  $ 'tall-skinny or short-wide matrices' )
782  8005 FORMAT( / 1x, a3, ': Householder recostruction from TSQR',
783  $ ' factorization output ', /,' for tall-skinny matrices.' )
784 *
785 * GE matrix types
786 *
787  9979 FORMAT( 4x, '1. Diagonal', 24x, '7. Last n/2 columns zero', / 4x,
788  $ '2. Upper triangular', 16x,
789  $ '8. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
790  $ '3. Lower triangular', 16x, '9. Random, CNDNUM = 0.1/EPS',
791  $ / 4x, '4. Random, CNDNUM = 2', 13x,
792  $ '10. Scaled near underflow', / 4x, '5. First column zero',
793  $ 14x, '11. Scaled near overflow', / 4x,
794  $ '6. Last column zero' )
795 *
796 * GB matrix types
797 *
798  9978 FORMAT( 4x, '1. Random, CNDNUM = 2', 14x,
799  $ '5. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
800  $ '2. First column zero', 15x, '6. Random, CNDNUM = .01/EPS',
801  $ / 4x, '3. Last column zero', 16x,
802  $ '7. Scaled near underflow', / 4x,
803  $ '4. Last n/2 columns zero', 11x, '8. Scaled near overflow' )
804 *
805 * GT matrix types
806 *
807  9977 FORMAT( ' Matrix types (1-6 have specified condition numbers):',
808  $ / 4x, '1. Diagonal', 24x, '7. Random, unspecified CNDNUM',
809  $ / 4x, '2. Random, CNDNUM = 2', 14x, '8. First column zero',
810  $ / 4x, '3. Random, CNDNUM = sqrt(0.1/EPS)', 2x,
811  $ '9. Last column zero', / 4x, '4. Random, CNDNUM = 0.1/EPS',
812  $ 7x, '10. Last n/2 columns zero', / 4x,
813  $ '5. Scaled near underflow', 10x,
814  $ '11. Scaled near underflow', / 4x,
815  $ '6. Scaled near overflow', 11x, '12. Scaled near overflow' )
816 *
817 * PT matrix types
818 *
819  9976 FORMAT( ' Matrix types (1-6 have specified condition numbers):',
820  $ / 4x, '1. Diagonal', 24x, '7. Random, unspecified CNDNUM',
821  $ / 4x, '2. Random, CNDNUM = 2', 14x,
822  $ '8. First row and column zero', / 4x,
823  $ '3. Random, CNDNUM = sqrt(0.1/EPS)', 2x,
824  $ '9. Last row and column zero', / 4x,
825  $ '4. Random, CNDNUM = 0.1/EPS', 7x,
826  $ '10. Middle row and column zero', / 4x,
827  $ '5. Scaled near underflow', 10x,
828  $ '11. Scaled near underflow', / 4x,
829  $ '6. Scaled near overflow', 11x, '12. Scaled near overflow' )
830 *
831 * PO, PP matrix types
832 *
833  9975 FORMAT( 4x, '1. Diagonal', 24x,
834  $ '6. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
835  $ '2. Random, CNDNUM = 2', 14x, '7. Random, CNDNUM = 0.1/EPS',
836  $ / 3x, '*3. First row and column zero', 7x,
837  $ '8. Scaled near underflow', / 3x,
838  $ '*4. Last row and column zero', 8x,
839  $ '9. Scaled near overflow', / 3x,
840  $ '*5. Middle row and column zero', / 3x,
841  $ '(* - tests error exits from ', a3,
842  $ 'TRF, no test ratios are computed)' )
843 *
844 * CH matrix types
845 *
846  9974 FORMAT( 4x, '1. Diagonal', 24x,
847  $ '6. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
848  $ '2. Random, CNDNUM = 2', 14x, '7. Random, CNDNUM = 0.1/EPS',
849  $ / 3x, '*3. First row and column zero', 7x,
850  $ '8. Scaled near underflow', / 3x,
851  $ '*4. Last row and column zero', 8x,
852  $ '9. Scaled near overflow', / 3x,
853  $ '*5. Middle row and column zero', / 3x,
854  $ '(* - tests error exits, no test ratios are computed)' )
855 *
856 * PS matrix types
857 *
858  8973 FORMAT( 4x, '1. Diagonal', / 4x, '2. Random, CNDNUM = 2', 14x,
859  $ / 3x, '*3. Nonzero eigenvalues of: D(1:RANK-1)=1 and ',
860  $ 'D(RANK) = 1.0/', a4, / 3x,
861  $ '*4. Nonzero eigenvalues of: D(1)=1 and ',
862  $ ' D(2:RANK) = 1.0/', a4, / 3x,
863  $ '*5. Nonzero eigenvalues of: D(I) = ', a4,
864  $ '**(-(I-1)/(RANK-1)) ', ' I=1:RANK', / 4x,
865  $ '6. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
866  $ '7. Random, CNDNUM = 0.1/EPS', / 4x,
867  $ '8. Scaled near underflow', / 4x, '9. Scaled near overflow',
868  $ / 3x, '(* - Semi-definite tests )' )
869  8972 FORMAT( 3x, 'RANK minus computed rank, returned by ', a, 'PSTRF' )
870 *
871 * PB matrix types
872 *
873  9973 FORMAT( 4x, '1. Random, CNDNUM = 2', 14x,
874  $ '5. Random, CNDNUM = sqrt(0.1/EPS)', / 3x,
875  $ '*2. First row and column zero', 7x,
876  $ '6. Random, CNDNUM = 0.1/EPS', / 3x,
877  $ '*3. Last row and column zero', 8x,
878  $ '7. Scaled near underflow', / 3x,
879  $ '*4. Middle row and column zero', 6x,
880  $ '8. Scaled near overflow', / 3x,
881  $ '(* - tests error exits from ', a3,
882  $ 'TRF, no test ratios are computed)' )
883 *
884 * SSY, SSR, SSP, CHE, CHR, CHP matrix types
885 *
886  9972 FORMAT( 4x, '1. Diagonal', 24x,
887  $ '6. Last n/2 rows and columns zero', / 4x,
888  $ '2. Random, CNDNUM = 2', 14x,
889  $ '7. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
890  $ '3. First row and column zero', 7x,
891  $ '8. Random, CNDNUM = 0.1/EPS', / 4x,
892  $ '4. Last row and column zero', 8x,
893  $ '9. Scaled near underflow', / 4x,
894  $ '5. Middle row and column zero', 5x,
895  $ '10. Scaled near overflow' )
896 *
897 * CSY, CSR, CSP matrix types
898 *
899  9971 FORMAT( 4x, '1. Diagonal', 24x,
900  $ '7. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
901  $ '2. Random, CNDNUM = 2', 14x, '8. Random, CNDNUM = 0.1/EPS',
902  $ / 4x, '3. First row and column zero', 7x,
903  $ '9. Scaled near underflow', / 4x,
904  $ '4. Last row and column zero', 7x,
905  $ '10. Scaled near overflow', / 4x,
906  $ '5. Middle row and column zero', 5x,
907  $ '11. Block diagonal matrix', / 4x,
908  $ '6. Last n/2 rows and columns zero' )
909 *
910 * QR matrix types
911 *
912  9970 FORMAT( 4x, '1. Diagonal', 24x,
913  $ '5. Random, CNDNUM = sqrt(0.1/EPS)', / 4x,
914  $ '2. Upper triangular', 16x, '6. Random, CNDNUM = 0.1/EPS',
915  $ / 4x, '3. Lower triangular', 16x,
916  $ '7. Scaled near underflow', / 4x, '4. Random, CNDNUM = 2',
917  $ 14x, '8. Scaled near overflow' )
918 *
919 * QP matrix types
920 *
921  9969 FORMAT( ' Matrix types (2-6 have condition 1/EPS):', / 4x,
922  $ '1. Zero matrix', 21x, '4. First n/2 columns fixed', / 4x,
923  $ '2. One small eigenvalue', 12x, '5. Last n/2 columns fixed',
924  $ / 4x, '3. Geometric distribution', 10x,
925  $ '6. Every second column fixed' )
926 *
927 * TZ matrix types
928 *
929  9968 FORMAT( ' Matrix types (2-3 have condition 1/EPS):', / 4x,
930  $ '1. Zero matrix', / 4x, '2. One small eigenvalue', / 4x,
931  $ '3. Geometric distribution' )
932 *
933 * LS matrix types
934 *
935  9967 FORMAT( ' Matrix types (1-3: full rank, 4-6: rank deficient):',
936  $ / 4x, '1 and 4. Normal scaling', / 4x,
937  $ '2 and 5. Scaled near overflow', / 4x,
938  $ '3 and 6. Scaled near underflow' )
939 *
940 * TR, TP matrix types
941 *
942  9966 FORMAT( ' Matrix types for ', a3, ' routines:', / 4x,
943  $ '1. Diagonal', 24x, '6. Scaled near overflow', / 4x,
944  $ '2. Random, CNDNUM = 2', 14x, '7. Identity', / 4x,
945  $ '3. Random, CNDNUM = sqrt(0.1/EPS) ',
946  $ '8. Unit triangular, CNDNUM = 2', / 4x,
947  $ '4. Random, CNDNUM = 0.1/EPS', 8x,
948  $ '9. Unit, CNDNUM = sqrt(0.1/EPS)', / 4x,
949  $ '5. Scaled near underflow', 10x,
950  $ '10. Unit, CNDNUM = 0.1/EPS' )
951  9965 FORMAT( ' Special types for testing ', a, ':', / 3x,
952  $ '11. Matrix elements are O(1), large right hand side', / 3x,
953  $ '12. First diagonal causes overflow,',
954  $ ' offdiagonal column norms < 1', / 3x,
955  $ '13. First diagonal causes overflow,',
956  $ ' offdiagonal column norms > 1', / 3x,
957  $ '14. Growth factor underflows, solution does not overflow',
958  $ / 3x, '15. Small diagonal causes gradual overflow', / 3x,
959  $ '16. One zero diagonal element', / 3x,
960  $ '17. Large offdiagonals cause overflow when adding a column'
961  $ , / 3x, '18. Unit triangular with large right hand side' )
962 *
963 * TB matrix types
964 *
965  9964 FORMAT( ' Matrix types for ', a3, ' routines:', / 4x,
966  $ '1. Random, CNDNUM = 2', 14x, '6. Identity', / 4x,
967  $ '2. Random, CNDNUM = sqrt(0.1/EPS) ',
968  $ '7. Unit triangular, CNDNUM = 2', / 4x,
969  $ '3. Random, CNDNUM = 0.1/EPS', 8x,
970  $ '8. Unit, CNDNUM = sqrt(0.1/EPS)', / 4x,
971  $ '4. Scaled near underflow', 11x,
972  $ '9. Unit, CNDNUM = 0.1/EPS', / 4x,
973  $ '5. Scaled near overflow' )
974  9963 FORMAT( ' Special types for testing ', a, ':', / 3x,
975  $ '10. Matrix elements are O(1), large right hand side', / 3x,
976  $ '11. First diagonal causes overflow,',
977  $ ' offdiagonal column norms < 1', / 3x,
978  $ '12. First diagonal causes overflow,',
979  $ ' offdiagonal column norms > 1', / 3x,
980  $ '13. Growth factor underflows, solution does not overflow',
981  $ / 3x, '14. Small diagonal causes gradual overflow', / 3x,
982  $ '15. One zero diagonal element', / 3x,
983  $ '16. Large offdiagonals cause overflow when adding a column'
984  $ , / 3x, '17. Unit triangular with large right hand side' )
985 *
986 * Test ratios
987 *
988  9962 FORMAT( 3x, i2, ': norm( L * U - A ) / ( N * norm(A) * EPS )' )
989  9961 FORMAT( 3x, i2, ': norm( I - A*AINV ) / ',
990  $ '( N * norm(A) * norm(AINV) * EPS )' )
991  9960 FORMAT( 3x, i2, ': norm( B - A * X ) / ',
992  $ '( norm(A) * norm(X) * EPS )' )
993  6660 FORMAT( 3x, i2, ': diagonal is not non-negative')
994  9959 FORMAT( 3x, i2, ': norm( X - XACT ) / ',
995  $ '( norm(XACT) * CNDNUM * EPS )' )
996  9958 FORMAT( 3x, i2, ': norm( X - XACT ) / ',
997  $ '( norm(XACT) * CNDNUM * EPS ), refined' )
998  9957 FORMAT( 3x, i2, ': norm( X - XACT ) / ',
999  $ '( norm(XACT) * (error bound) )' )
1000  9956 FORMAT( 3x, i2, ': (backward error) / EPS' )
1001  9955 FORMAT( 3x, i2, ': RCOND * CNDNUM - 1.0' )
1002  9954 FORMAT( 3x, i2, ': norm( U'' * U - A ) / ( N * norm(A) * EPS )',
1003  $ ', or', / 7x, 'norm( L * L'' - A ) / ( N * norm(A) * EPS )'
1004  $ )
1005  8950 FORMAT( 3x,
1006  $ 'norm( P * U'' * U * P'' - A ) / ( N * norm(A) * EPS )',
1007  $ ', or', / 3x,
1008  $ 'norm( P * L * L'' * P'' - A ) / ( N * norm(A) * EPS )' )
1009  9953 FORMAT( 3x, i2, ': norm( U*D*U'' - A ) / ( N * norm(A) * EPS )',
1010  $ ', or', / 7x, 'norm( L*D*L'' - A ) / ( N * norm(A) * EPS )'
1011  $ )
1012  9952 FORMAT( 3x, i2, ': norm( U''*D*U - A ) / ( N * norm(A) * EPS )',
1013  $ ', or', / 7x, 'norm( L*D*L'' - A ) / ( N * norm(A) * EPS )'
1014  $ )
1015  9951 FORMAT( ' Test ratio for ', a, ':', / 3x, i2,
1016  $ ': norm( s*b - A*x ) / ( norm(A) * norm(x) * EPS )' )
1017  9950 FORMAT( 3x, i2, ': norm( R - Q'' * A ) / ( M * norm(A) * EPS )' )
1018  6950 FORMAT( 3x, i2, ': norm( R - Q'' * A ) / ( M * norm(A) * EPS )
1019  $ [RFPG]' )
1020  9949 FORMAT( 3x, i2, ': norm( L - A * Q'' ) / ( N * norm(A) * EPS )' )
1021  9948 FORMAT( 3x, i2, ': norm( L - Q'' * A ) / ( M * norm(A) * EPS )' )
1022  9947 FORMAT( 3x, i2, ': norm( R - A * Q'' ) / ( N * norm(A) * EPS )' )
1023  9946 FORMAT( 3x, i2, ': norm( I - Q''*Q ) / ( M * EPS )' )
1024  9945 FORMAT( 3x, i2, ': norm( I - Q*Q'' ) / ( N * EPS )' )
1025  9944 FORMAT( 3x, i2, ': norm( Q*C - Q*C ) / ', '( ', a1,
1026  $ ' * norm(C) * EPS )' )
1027  9943 FORMAT( 3x, i2, ': norm( C*Q - C*Q ) / ', '( ', a1,
1028  $ ' * norm(C) * EPS )' )
1029  9942 FORMAT( 3x, i2, ': norm( Q''*C - Q''*C )/ ', '( ', a1,
1030  $ ' * norm(C) * EPS )' )
1031  9941 FORMAT( 3x, i2, ': norm( C*Q'' - C*Q'' )/ ', '( ', a1,
1032  $ ' * norm(C) * EPS )' )
1033  9940 FORMAT( 3x, i2, ': norm(svd(A) - svd(R)) / ',
1034  $ '( M * norm(svd(R)) * EPS )' )
1035  9939 FORMAT( 3x, i2, ': norm( A*P - Q*R ) / ( M * norm(A) * EPS )'
1036  $ )
1037  9938 FORMAT( 3x, i2, ': norm( I - Q''*Q ) / ( M * EPS )' )
1038  9937 FORMAT( 3x, i2, ': norm( A - R*Q ) / ( M * norm(A) * EPS )'
1039  $ )
1040  9935 FORMAT( 3x, i2, ': norm( B - A * X ) / ',
1041  $ '( max(M,N) * norm(A) * norm(X) * EPS )' )
1042  9934 FORMAT( 3x, i2, ': norm( (A*X-B)'' *A ) / ',
1043  $ '( max(M,N,NRHS) * norm(A) * norm(B) * EPS )' )
1044  9933 FORMAT( 3x, i2, ': norm(svd(A)-svd(R)) / ',
1045  $ '( min(M,N) * norm(svd(R)) * EPS )' )
1046  9932 FORMAT( 3x, i2, ': Check if X is in the row space of A or A''' )
1047  9931 FORMAT( 3x, i2, ': norm( (A*X-B)'' *A ) / ',
1048  $ '( max(M,N,NRHS) * norm(A) * norm(B) * EPS )', / 7x,
1049  $ 'if TRANS=''N'.GE.' and MN or TRANS=''T'.LT.' and MN, ',
1050  $ 'otherwise', / 7x,
1051  $ 'check if X is in the row space of A or A'' ',
1052  $ '(overdetermined case)' )
1053  9929 FORMAT( ' Test ratios (1-3: ', a1, 'TZRZF):' )
1054  9920 FORMAT( 3x, ' 7-10: same as 3-6', 3x, ' 11-14: same as 3-6' )
1055  9921 FORMAT( ' Test ratios:', / ' (1-2: ', a1, 'GELS, 3-6: ', a1,
1056  $ 'GELSY, 7-10: ', a1, 'GELSS, 11-14: ', a1, 'GELSD, 15-16: ',
1057  $ a1, 'GETSLS)')
1058  9928 FORMAT( 7x, 'where ALPHA = ( 1 + SQRT( 17 ) ) / 8' )
1059  9927 FORMAT( 3x, i2, ': ABS( Largest element in L )', / 12x,
1060  $ ' - ( 1 / ( 1 - ALPHA ) ) + THRESH' )
1061  9926 FORMAT( 3x, i2, ': Largest 2-Norm of 2-by-2 pivots', / 12x,
1062  $ ' - ( ( 1 + ALPHA ) / ( 1 - ALPHA ) ) + THRESH' )
1063  8011 FORMAT(3x,i2,': norm( R - Q''*A ) / ( M * norm(A) * EPS )' )
1064  8012 FORMAT(3x,i2,': norm( I - Q''*Q ) / ( M * EPS )' )
1065  8013 FORMAT(3x,i2,': norm( Q*C - Q*C ) / ( M * norm(C) * EPS )' )
1066  8014 FORMAT(3x,i2,': norm( Q''*C - Q''*C ) / ( M * norm(C) * EPS )')
1067  8015 FORMAT(3x,i2,': norm( C*Q - C*Q ) / ( M * norm(C) * EPS )' )
1068  8016 FORMAT(3x,i2,': norm( C*Q'' - C*Q'' ) / ( M * norm(C) * EPS )')
1069  8017 FORMAT(3x,i2,': norm( R - Q''*A ) / ( (M+N) * norm(A) * EPS )' )
1070  8018 FORMAT(3x,i2,': norm( I - Q''*Q ) / ( (M+N) * EPS )' )
1071  8019 FORMAT(3x,i2,': norm( Q*C - Q*C ) / ( (M+N) * norm(C) * EPS )' )
1072  8020 FORMAT(3x,i2,
1073  $ ': norm( Q''*C - Q''*C ) / ( (M+N) * norm(C) * EPS )')
1074  8021 FORMAT(3x,i2,': norm( C*Q - C*Q ) / ( (M+N) * norm(C) * EPS )' )
1075  8022 FORMAT(3x,i2,
1076  $ ': norm( C*Q'' - C*Q'' ) / ( (M+N) * norm(C) * EPS )')
1077  8023 FORMAT(3x,i2,': norm( L - A*Q'' ) / ( (M+N) * norm(A) * EPS )' )
1078  8024 FORMAT(3x,i2,': norm( I - Q*Q'' ) / ( (M+N) * EPS )' )
1079  8025 FORMAT(3x,i2,': norm( Q*C - Q*C ) / ( (M+N) * norm(C) * EPS )' )
1080  8026 FORMAT(3x,i2,
1081  $ ': norm( Q''*C - Q''*C ) / ( (M+N) * norm(C) * EPS )')
1082  8027 FORMAT(3x,i2,': norm( C*Q - C*Q ) / ( (M+N) * norm(C) * EPS )' )
1083  8028 FORMAT(3x,i2,
1084  $ ': norm( C*Q'' - C*Q'' ) / ( (M+N) * norm(C) * EPS )')
1085  8029 FORMAT(3x,i2,': norm( L - A*Q'' ) / ( (M+N) * norm(A) * EPS )' )
1086  8030 FORMAT(3x,i2,': norm( I - Q*Q'' ) / ( (M+N) * EPS )' )
1087  8031 FORMAT(3x,i2,': norm( Q*C - Q*C ) / ( (M+N) * norm(C) * EPS )' )
1088  8032 FORMAT(3x,i2,
1089  $ ': norm( Q''*C - Q''*C ) / ( (M+N) * norm(C) * EPS )')
1090  8033 FORMAT(3x,i2,': norm( C*Q - C*Q ) / ( (M+N) * norm(C) * EPS )' )
1091  8034 FORMAT(3x,i2,
1092  $ ': norm( C*Q'' - C*Q'' ) / ( (M+N) * norm(C) * EPS )')
1093  8035 FORMAT(3x,i2,': norm( R - Q''*A ) / ( (M+N) * norm(A) * EPS )' )
1094  8036 FORMAT(3x,i2,': norm( I - Q''*Q ) / ( (M+N) * EPS )' )
1095  8037 FORMAT(3x,i2,': norm( Q*C - Q*C ) / ( (M+N) * norm(C) * EPS )' )
1096  8038 FORMAT(3x,i2,
1097  $ ': norm( Q''*C - Q''*C ) / ( (M+N) * norm(C) * EPS )')
1098  8039 FORMAT(3x,i2,': norm( C*Q - C*Q ) / ( (M+N) * norm(C) * EPS )' )
1099  8040 FORMAT(3x,i2,
1100  $ ': norm( C*Q'' - C*Q'' ) / ( (M+N) * norm(C) * EPS )')
1101 *
1102  8050 FORMAT(3x,i2,': norm( R - Q''*A ) / ( M * norm(A) * EPS )' )
1103  8051 FORMAT(3x,i2,': norm( I - Q''*Q ) / ( M * EPS )' )
1104  8052 FORMAT(3x,i2,': norm( Q*C - Q*C ) / ( M * norm(C) * EPS )' )
1105  8053 FORMAT(3x,i2,': norm( Q''*C - Q''*C ) / ( M * norm(C) * EPS )')
1106  8054 FORMAT(3x,i2,': norm( C*Q - C*Q ) / ( M * norm(C) * EPS )' )
1107  8055 FORMAT(3x,i2,': norm( C*Q'' - C*Q'' ) / ( M * norm(C) * EPS )')
1108 
1109 *
1110  RETURN
1111 *
1112 * End of ALAHD
1113 *
1114  END
alahd
subroutine alahd(IOUNIT, PATH)
ALAHD
Definition: alahd.f:109