LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
sget31.f
Go to the documentation of this file.
1 *> \brief \b SGET31
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 SGET31( RMAX, LMAX, NINFO, KNT )
12 *
13 * .. Scalar Arguments ..
14 * INTEGER KNT, LMAX
15 * REAL RMAX
16 * ..
17 * .. Array Arguments ..
18 * INTEGER NINFO( 2 )
19 * ..
20 *
21 *
22 *> \par Purpose:
23 * =============
24 *>
25 *> \verbatim
26 *>
27 *> SGET31 tests SLALN2, a routine for solving
28 *>
29 *> (ca A - w D)X = sB
30 *>
31 *> where A is an NA by NA matrix (NA=1 or 2 only), w is a real (NW=1) or
32 *> complex (NW=2) constant, ca is a real constant, D is an NA by NA real
33 *> diagonal matrix, and B is an NA by NW matrix (when NW=2 the second
34 *> column of B contains the imaginary part of the solution). The code
35 *> returns X and s, where s is a scale factor, less than or equal to 1,
36 *> which is chosen to avoid overflow in X.
37 *>
38 *> If any singular values of ca A-w D are less than another input
39 *> parameter SMIN, they are perturbed up to SMIN.
40 *>
41 *> The test condition is that the scaled residual
42 *>
43 *> norm( (ca A-w D)*X - s*B ) /
44 *> ( max( ulp*norm(ca A-w D), SMIN )*norm(X) )
45 *>
46 *> should be on the order of 1. Here, ulp is the machine precision.
47 *> Also, it is verified that SCALE is less than or equal to 1, and that
48 *> XNORM = infinity-norm(X).
49 *> \endverbatim
50 *
51 * Arguments:
52 * ==========
53 *
54 *> \param[out] RMAX
55 *> \verbatim
56 *> RMAX is REAL
57 *> Value of the largest test ratio.
58 *> \endverbatim
59 *>
60 *> \param[out] LMAX
61 *> \verbatim
62 *> LMAX is INTEGER
63 *> Example number where largest test ratio achieved.
64 *> \endverbatim
65 *>
66 *> \param[out] NINFO
67 *> \verbatim
68 *> NINFO is INTEGER array, dimension (3)
69 *> NINFO(1) = number of examples with INFO less than 0
70 *> NINFO(2) = number of examples with INFO greater than 0
71 *> \endverbatim
72 *>
73 *> \param[out] KNT
74 *> \verbatim
75 *> KNT is INTEGER
76 *> Total number of examples tested.
77 *> \endverbatim
78 *
79 * Authors:
80 * ========
81 *
82 *> \author Univ. of Tennessee
83 *> \author Univ. of California Berkeley
84 *> \author Univ. of Colorado Denver
85 *> \author NAG Ltd.
86 *
87 *> \ingroup single_eig
88 *
89 * =====================================================================
90  SUBROUTINE sget31( RMAX, LMAX, NINFO, KNT )
91 *
92 * -- LAPACK test routine --
93 * -- LAPACK is a software package provided by Univ. of Tennessee, --
94 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
95 *
96 * .. Scalar Arguments ..
97  INTEGER KNT, LMAX
98  REAL RMAX
99 * ..
100 * .. Array Arguments ..
101  INTEGER NINFO( 2 )
102 * ..
103 *
104 * =====================================================================
105 *
106 * .. Parameters ..
107  REAL ZERO, HALF, ONE
108  parameter( zero = 0.0e0, half = 0.5e0, one = 1.0e0 )
109  REAL TWO, THREE, FOUR
110  parameter( two = 2.0e0, three = 3.0e0, four = 4.0e0 )
111  REAL SEVEN, TEN
112  parameter( seven = 7.0e0, ten = 10.0e0 )
113  REAL TWNONE
114  parameter( twnone = 21.0e0 )
115 * ..
116 * .. Local Scalars ..
117  INTEGER IA, IB, ICA, ID1, ID2, INFO, ISMIN, ITRANS,
118  $ IWI, IWR, NA, NW
119  REAL BIGNUM, CA, D1, D2, DEN, EPS, RES, SCALE, SMIN,
120  $ SMLNUM, TMP, UNFL, WI, WR, XNORM
121 * ..
122 * .. Local Arrays ..
123  LOGICAL LTRANS( 0: 1 )
124  REAL A( 2, 2 ), B( 2, 2 ), VAB( 3 ), VCA( 5 ),
125  $ VDD( 4 ), VSMIN( 4 ), VWI( 4 ), VWR( 4 ),
126  $ X( 2, 2 )
127 * ..
128 * .. External Functions ..
129  REAL SLAMCH
130  EXTERNAL slamch
131 * ..
132 * .. External Subroutines ..
133  EXTERNAL slabad, slaln2
134 * ..
135 * .. Intrinsic Functions ..
136  INTRINSIC abs, max, sqrt
137 * ..
138 * .. Data statements ..
139  DATA ltrans / .false., .true. /
140 * ..
141 * .. Executable Statements ..
142 *
143 * Get machine parameters
144 *
145  eps = slamch( 'P' )
146  unfl = slamch( 'U' )
147  smlnum = slamch( 'S' ) / eps
148  bignum = one / smlnum
149  CALL slabad( smlnum, bignum )
150 *
151 * Set up test case parameters
152 *
153  vsmin( 1 ) = smlnum
154  vsmin( 2 ) = eps
155  vsmin( 3 ) = one / ( ten*ten )
156  vsmin( 4 ) = one / eps
157  vab( 1 ) = sqrt( smlnum )
158  vab( 2 ) = one
159  vab( 3 ) = sqrt( bignum )
160  vwr( 1 ) = zero
161  vwr( 2 ) = half
162  vwr( 3 ) = two
163  vwr( 4 ) = one
164  vwi( 1 ) = smlnum
165  vwi( 2 ) = eps
166  vwi( 3 ) = one
167  vwi( 4 ) = two
168  vdd( 1 ) = sqrt( smlnum )
169  vdd( 2 ) = one
170  vdd( 3 ) = two
171  vdd( 4 ) = sqrt( bignum )
172  vca( 1 ) = zero
173  vca( 2 ) = sqrt( smlnum )
174  vca( 3 ) = eps
175  vca( 4 ) = half
176  vca( 5 ) = one
177 *
178  knt = 0
179  ninfo( 1 ) = 0
180  ninfo( 2 ) = 0
181  lmax = 0
182  rmax = zero
183 *
184 * Begin test loop
185 *
186  DO 190 id1 = 1, 4
187  d1 = vdd( id1 )
188  DO 180 id2 = 1, 4
189  d2 = vdd( id2 )
190  DO 170 ica = 1, 5
191  ca = vca( ica )
192  DO 160 itrans = 0, 1
193  DO 150 ismin = 1, 4
194  smin = vsmin( ismin )
195 *
196  na = 1
197  nw = 1
198  DO 30 ia = 1, 3
199  a( 1, 1 ) = vab( ia )
200  DO 20 ib = 1, 3
201  b( 1, 1 ) = vab( ib )
202  DO 10 iwr = 1, 4
203  IF( d1.EQ.one .AND. d2.EQ.one .AND. ca.EQ.
204  $ one ) THEN
205  wr = vwr( iwr )*a( 1, 1 )
206  ELSE
207  wr = vwr( iwr )
208  END IF
209  wi = zero
210  CALL slaln2( ltrans( itrans ), na, nw,
211  $ smin, ca, a, 2, d1, d2, b, 2,
212  $ wr, wi, x, 2, scale, xnorm,
213  $ info )
214  IF( info.LT.0 )
215  $ ninfo( 1 ) = ninfo( 1 ) + 1
216  IF( info.GT.0 )
217  $ ninfo( 2 ) = ninfo( 2 ) + 1
218  res = abs( ( ca*a( 1, 1 )-wr*d1 )*
219  $ x( 1, 1 )-scale*b( 1, 1 ) )
220  IF( info.EQ.0 ) THEN
221  den = max( eps*( abs( ( ca*a( 1,
222  $ 1 )-wr*d1 )*x( 1, 1 ) ) ),
223  $ smlnum )
224  ELSE
225  den = max( smin*abs( x( 1, 1 ) ),
226  $ smlnum )
227  END IF
228  res = res / den
229  IF( abs( x( 1, 1 ) ).LT.unfl .AND.
230  $ abs( b( 1, 1 ) ).LE.smlnum*
231  $ abs( ca*a( 1, 1 )-wr*d1 ) )res = zero
232  IF( scale.GT.one )
233  $ res = res + one / eps
234  res = res + abs( xnorm-abs( x( 1, 1 ) ) )
235  $ / max( smlnum, xnorm ) / eps
236  IF( info.NE.0 .AND. info.NE.1 )
237  $ res = res + one / eps
238  knt = knt + 1
239  IF( res.GT.rmax ) THEN
240  lmax = knt
241  rmax = res
242  END IF
243  10 CONTINUE
244  20 CONTINUE
245  30 CONTINUE
246 *
247  na = 1
248  nw = 2
249  DO 70 ia = 1, 3
250  a( 1, 1 ) = vab( ia )
251  DO 60 ib = 1, 3
252  b( 1, 1 ) = vab( ib )
253  b( 1, 2 ) = -half*vab( ib )
254  DO 50 iwr = 1, 4
255  IF( d1.EQ.one .AND. d2.EQ.one .AND. ca.EQ.
256  $ one ) THEN
257  wr = vwr( iwr )*a( 1, 1 )
258  ELSE
259  wr = vwr( iwr )
260  END IF
261  DO 40 iwi = 1, 4
262  IF( d1.EQ.one .AND. d2.EQ.one .AND.
263  $ ca.EQ.one ) THEN
264  wi = vwi( iwi )*a( 1, 1 )
265  ELSE
266  wi = vwi( iwi )
267  END IF
268  CALL slaln2( ltrans( itrans ), na, nw,
269  $ smin, ca, a, 2, d1, d2, b,
270  $ 2, wr, wi, x, 2, scale,
271  $ xnorm, info )
272  IF( info.LT.0 )
273  $ ninfo( 1 ) = ninfo( 1 ) + 1
274  IF( info.GT.0 )
275  $ ninfo( 2 ) = ninfo( 2 ) + 1
276  res = abs( ( ca*a( 1, 1 )-wr*d1 )*
277  $ x( 1, 1 )+( wi*d1 )*x( 1, 2 )-
278  $ scale*b( 1, 1 ) )
279  res = res + abs( ( -wi*d1 )*x( 1, 1 )+
280  $ ( ca*a( 1, 1 )-wr*d1 )*x( 1, 2 )-
281  $ scale*b( 1, 2 ) )
282  IF( info.EQ.0 ) THEN
283  den = max( eps*( max( abs( ca*a( 1,
284  $ 1 )-wr*d1 ), abs( d1*wi ) )*
285  $ ( abs( x( 1, 1 ) )+abs( x( 1,
286  $ 2 ) ) ) ), smlnum )
287  ELSE
288  den = max( smin*( abs( x( 1,
289  $ 1 ) )+abs( x( 1, 2 ) ) ),
290  $ smlnum )
291  END IF
292  res = res / den
293  IF( abs( x( 1, 1 ) ).LT.unfl .AND.
294  $ abs( x( 1, 2 ) ).LT.unfl .AND.
295  $ abs( b( 1, 1 ) ).LE.smlnum*
296  $ abs( ca*a( 1, 1 )-wr*d1 ) )
297  $ res = zero
298  IF( scale.GT.one )
299  $ res = res + one / eps
300  res = res + abs( xnorm-
301  $ abs( x( 1, 1 ) )-
302  $ abs( x( 1, 2 ) ) ) /
303  $ max( smlnum, xnorm ) / eps
304  IF( info.NE.0 .AND. info.NE.1 )
305  $ res = res + one / eps
306  knt = knt + 1
307  IF( res.GT.rmax ) THEN
308  lmax = knt
309  rmax = res
310  END IF
311  40 CONTINUE
312  50 CONTINUE
313  60 CONTINUE
314  70 CONTINUE
315 *
316  na = 2
317  nw = 1
318  DO 100 ia = 1, 3
319  a( 1, 1 ) = vab( ia )
320  a( 1, 2 ) = -three*vab( ia )
321  a( 2, 1 ) = -seven*vab( ia )
322  a( 2, 2 ) = twnone*vab( ia )
323  DO 90 ib = 1, 3
324  b( 1, 1 ) = vab( ib )
325  b( 2, 1 ) = -two*vab( ib )
326  DO 80 iwr = 1, 4
327  IF( d1.EQ.one .AND. d2.EQ.one .AND. ca.EQ.
328  $ one ) THEN
329  wr = vwr( iwr )*a( 1, 1 )
330  ELSE
331  wr = vwr( iwr )
332  END IF
333  wi = zero
334  CALL slaln2( ltrans( itrans ), na, nw,
335  $ smin, ca, a, 2, d1, d2, b, 2,
336  $ wr, wi, x, 2, scale, xnorm,
337  $ info )
338  IF( info.LT.0 )
339  $ ninfo( 1 ) = ninfo( 1 ) + 1
340  IF( info.GT.0 )
341  $ ninfo( 2 ) = ninfo( 2 ) + 1
342  IF( itrans.EQ.1 ) THEN
343  tmp = a( 1, 2 )
344  a( 1, 2 ) = a( 2, 1 )
345  a( 2, 1 ) = tmp
346  END IF
347  res = abs( ( ca*a( 1, 1 )-wr*d1 )*
348  $ x( 1, 1 )+( ca*a( 1, 2 ) )*
349  $ x( 2, 1 )-scale*b( 1, 1 ) )
350  res = res + abs( ( ca*a( 2, 1 ) )*
351  $ x( 1, 1 )+( ca*a( 2, 2 )-wr*d2 )*
352  $ x( 2, 1 )-scale*b( 2, 1 ) )
353  IF( info.EQ.0 ) THEN
354  den = max( eps*( max( abs( ca*a( 1,
355  $ 1 )-wr*d1 )+abs( ca*a( 1, 2 ) ),
356  $ abs( ca*a( 2, 1 ) )+abs( ca*a( 2,
357  $ 2 )-wr*d2 ) )*max( abs( x( 1,
358  $ 1 ) ), abs( x( 2, 1 ) ) ) ),
359  $ smlnum )
360  ELSE
361  den = max( eps*( max( smin / eps,
362  $ max( abs( ca*a( 1,
363  $ 1 )-wr*d1 )+abs( ca*a( 1, 2 ) ),
364  $ abs( ca*a( 2, 1 ) )+abs( ca*a( 2,
365  $ 2 )-wr*d2 ) ) )*max( abs( x( 1,
366  $ 1 ) ), abs( x( 2, 1 ) ) ) ),
367  $ smlnum )
368  END IF
369  res = res / den
370  IF( abs( x( 1, 1 ) ).LT.unfl .AND.
371  $ abs( x( 2, 1 ) ).LT.unfl .AND.
372  $ abs( b( 1, 1 ) )+abs( b( 2, 1 ) ).LE.
373  $ smlnum*( abs( ca*a( 1,
374  $ 1 )-wr*d1 )+abs( ca*a( 1,
375  $ 2 ) )+abs( ca*a( 2,
376  $ 1 ) )+abs( ca*a( 2, 2 )-wr*d2 ) ) )
377  $ res = zero
378  IF( scale.GT.one )
379  $ res = res + one / eps
380  res = res + abs( xnorm-
381  $ max( abs( x( 1, 1 ) ), abs( x( 2,
382  $ 1 ) ) ) ) / max( smlnum, xnorm ) /
383  $ eps
384  IF( info.NE.0 .AND. info.NE.1 )
385  $ res = res + one / eps
386  knt = knt + 1
387  IF( res.GT.rmax ) THEN
388  lmax = knt
389  rmax = res
390  END IF
391  80 CONTINUE
392  90 CONTINUE
393  100 CONTINUE
394 *
395  na = 2
396  nw = 2
397  DO 140 ia = 1, 3
398  a( 1, 1 ) = vab( ia )*two
399  a( 1, 2 ) = -three*vab( ia )
400  a( 2, 1 ) = -seven*vab( ia )
401  a( 2, 2 ) = twnone*vab( ia )
402  DO 130 ib = 1, 3
403  b( 1, 1 ) = vab( ib )
404  b( 2, 1 ) = -two*vab( ib )
405  b( 1, 2 ) = four*vab( ib )
406  b( 2, 2 ) = -seven*vab( ib )
407  DO 120 iwr = 1, 4
408  IF( d1.EQ.one .AND. d2.EQ.one .AND. ca.EQ.
409  $ one ) THEN
410  wr = vwr( iwr )*a( 1, 1 )
411  ELSE
412  wr = vwr( iwr )
413  END IF
414  DO 110 iwi = 1, 4
415  IF( d1.EQ.one .AND. d2.EQ.one .AND.
416  $ ca.EQ.one ) THEN
417  wi = vwi( iwi )*a( 1, 1 )
418  ELSE
419  wi = vwi( iwi )
420  END IF
421  CALL slaln2( ltrans( itrans ), na, nw,
422  $ smin, ca, a, 2, d1, d2, b,
423  $ 2, wr, wi, x, 2, scale,
424  $ xnorm, info )
425  IF( info.LT.0 )
426  $ ninfo( 1 ) = ninfo( 1 ) + 1
427  IF( info.GT.0 )
428  $ ninfo( 2 ) = ninfo( 2 ) + 1
429  IF( itrans.EQ.1 ) THEN
430  tmp = a( 1, 2 )
431  a( 1, 2 ) = a( 2, 1 )
432  a( 2, 1 ) = tmp
433  END IF
434  res = abs( ( ca*a( 1, 1 )-wr*d1 )*
435  $ x( 1, 1 )+( ca*a( 1, 2 ) )*
436  $ x( 2, 1 )+( wi*d1 )*x( 1, 2 )-
437  $ scale*b( 1, 1 ) )
438  res = res + abs( ( ca*a( 1,
439  $ 1 )-wr*d1 )*x( 1, 2 )+
440  $ ( ca*a( 1, 2 ) )*x( 2, 2 )-
441  $ ( wi*d1 )*x( 1, 1 )-scale*
442  $ b( 1, 2 ) )
443  res = res + abs( ( ca*a( 2, 1 ) )*
444  $ x( 1, 1 )+( ca*a( 2, 2 )-wr*d2 )*
445  $ x( 2, 1 )+( wi*d2 )*x( 2, 2 )-
446  $ scale*b( 2, 1 ) )
447  res = res + abs( ( ca*a( 2, 1 ) )*
448  $ x( 1, 2 )+( ca*a( 2, 2 )-wr*d2 )*
449  $ x( 2, 2 )-( wi*d2 )*x( 2, 1 )-
450  $ scale*b( 2, 2 ) )
451  IF( info.EQ.0 ) THEN
452  den = max( eps*( max( abs( ca*a( 1,
453  $ 1 )-wr*d1 )+abs( ca*a( 1,
454  $ 2 ) )+abs( wi*d1 ),
455  $ abs( ca*a( 2,
456  $ 1 ) )+abs( ca*a( 2,
457  $ 2 )-wr*d2 )+abs( wi*d2 ) )*
458  $ max( abs( x( 1,
459  $ 1 ) )+abs( x( 2, 1 ) ),
460  $ abs( x( 1, 2 ) )+abs( x( 2,
461  $ 2 ) ) ) ), smlnum )
462  ELSE
463  den = max( eps*( max( smin / eps,
464  $ max( abs( ca*a( 1,
465  $ 1 )-wr*d1 )+abs( ca*a( 1,
466  $ 2 ) )+abs( wi*d1 ),
467  $ abs( ca*a( 2,
468  $ 1 ) )+abs( ca*a( 2,
469  $ 2 )-wr*d2 )+abs( wi*d2 ) ) )*
470  $ max( abs( x( 1,
471  $ 1 ) )+abs( x( 2, 1 ) ),
472  $ abs( x( 1, 2 ) )+abs( x( 2,
473  $ 2 ) ) ) ), smlnum )
474  END IF
475  res = res / den
476  IF( abs( x( 1, 1 ) ).LT.unfl .AND.
477  $ abs( x( 2, 1 ) ).LT.unfl .AND.
478  $ abs( x( 1, 2 ) ).LT.unfl .AND.
479  $ abs( x( 2, 2 ) ).LT.unfl .AND.
480  $ abs( b( 1, 1 ) )+
481  $ abs( b( 2, 1 ) ).LE.smlnum*
482  $ ( abs( ca*a( 1, 1 )-wr*d1 )+
483  $ abs( ca*a( 1, 2 ) )+abs( ca*a( 2,
484  $ 1 ) )+abs( ca*a( 2,
485  $ 2 )-wr*d2 )+abs( wi*d2 )+abs( wi*
486  $ d1 ) ) )res = zero
487  IF( scale.GT.one )
488  $ res = res + one / eps
489  res = res + abs( xnorm-
490  $ max( abs( x( 1, 1 ) )+abs( x( 1,
491  $ 2 ) ), abs( x( 2,
492  $ 1 ) )+abs( x( 2, 2 ) ) ) ) /
493  $ max( smlnum, xnorm ) / eps
494  IF( info.NE.0 .AND. info.NE.1 )
495  $ res = res + one / eps
496  knt = knt + 1
497  IF( res.GT.rmax ) THEN
498  lmax = knt
499  rmax = res
500  END IF
501  110 CONTINUE
502  120 CONTINUE
503  130 CONTINUE
504  140 CONTINUE
505  150 CONTINUE
506  160 CONTINUE
507  170 CONTINUE
508  180 CONTINUE
509  190 CONTINUE
510 *
511  RETURN
512 *
513 * End of SGET31
514 *
515  END
subroutine slabad(SMALL, LARGE)
SLABAD
Definition: slabad.f:74
subroutine slaln2(LTRANS, NA, NW, SMIN, CA, A, LDA, D1, D2, B, LDB, WR, WI, X, LDX, SCALE, XNORM, INFO)
SLALN2 solves a 1-by-1 or 2-by-2 linear system of equations of the specified form.
Definition: slaln2.f:218
subroutine sget31(RMAX, LMAX, NINFO, KNT)
SGET31
Definition: sget31.f:91