LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine dget39 ( double precision  RMAX,
integer  LMAX,
integer  NINFO,
integer  KNT 
)

DGET39

Purpose: 
 DGET39 tests DLAQTR, a routine for solving the real or
 special complex quasi upper triangular system

      op(T)*p = scale*c,
 or
      op(T + iB)*(p+iq) = scale*(c+id),

 in real arithmetic. T is upper quasi-triangular.
 If it is complex, then the first diagonal block of T must be
 1 by 1, B has the special structure

                B = [ b(1) b(2) ... b(n) ]
                    [       w            ]
                    [           w        ]
                    [              .     ]
                    [                 w  ]

 op(A) = A or A', where A' denotes the conjugate transpose of
 the matrix A.

 On input, X = [ c ].  On output, X = [ p ].
               [ d ]                  [ q ]

 Scale is an output less than or equal to 1, chosen to avoid
 overflow in X.
 This subroutine is specially designed for the condition number
 estimation in the eigenproblem routine DTRSNA.

 The test code verifies that the following residual is order 1:

      ||(T+i*B)*(x1+i*x2) - scale*(d1+i*d2)||
    -----------------------------------------
        max(ulp*(||T||+||B||)*(||x1||+||x2||),
            (||T||+||B||)*smlnum/ulp,
            smlnum)

 (The (||T||+||B||)*smlnum/ulp term accounts for possible
  (gradual or nongradual) underflow in x1 and x2.)
Parameters
[out]RMAX
          RMAX is DOUBLE PRECISION
          Value of the largest test ratio.
[out]LMAX
          LMAX is INTEGER
          Example number where largest test ratio achieved.
[out]NINFO
          NINFO is INTEGER
          Number of examples where INFO is nonzero.
[out]KNT
          KNT is INTEGER
          Total number of examples tested.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 105 of file dget39.f.

105 *
106 * -- LAPACK test routine (version 3.4.0) --
107 * -- LAPACK is a software package provided by Univ. of Tennessee, --
108 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
109 * November 2011
110 *
111 * .. Scalar Arguments ..
112  INTEGER knt, lmax, ninfo
113  DOUBLE PRECISION rmax
114 * ..
115 *
116 * =====================================================================
117 *
118 * .. Parameters ..
119  INTEGER ldt, ldt2
120  parameter ( ldt = 10, ldt2 = 2*ldt )
121  DOUBLE PRECISION zero, one
122  parameter ( zero = 0.0d0, one = 1.0d0 )
123 * ..
124 * .. Local Scalars ..
125  INTEGER i, info, ivm1, ivm2, ivm3, ivm4, ivm5, j, k, n,
126  $ ndim
127  DOUBLE PRECISION bignum, domin, dumm, eps, norm, normtb, resid,
128  $ scale, smlnum, w, xnorm
129 * ..
130 * .. External Functions ..
131  INTEGER idamax
132  DOUBLE PRECISION dasum, ddot, dlamch, dlange
133  EXTERNAL idamax, dasum, ddot, dlamch, dlange
134 * ..
135 * .. External Subroutines ..
136  EXTERNAL dcopy, dgemv, dlabad, dlaqtr
137 * ..
138 * .. Intrinsic Functions ..
139  INTRINSIC abs, cos, dble, max, sin, sqrt
140 * ..
141 * .. Local Arrays ..
142  INTEGER idim( 6 ), ival( 5, 5, 6 )
143  DOUBLE PRECISION b( ldt ), d( ldt2 ), dum( 1 ), t( ldt, ldt ),
144  $ vm1( 5 ), vm2( 5 ), vm3( 5 ), vm4( 5 ),
145  $ vm5( 3 ), work( ldt ), x( ldt2 ), y( ldt2 )
146 * ..
147 * .. Data statements ..
148  DATA idim / 4, 5*5 /
149  DATA ival / 3, 4*0, 1, 1, -1, 0, 0, 3, 2, 1, 0, 0,
150  $ 4, 3, 2, 2, 0, 5*0, 1, 4*0, 2, 2, 3*0, 3, 3, 4,
151  $ 0, 0, 4, 2, 2, 3, 0, 4*1, 5, 1, 4*0, 2, 4, -2,
152  $ 0, 0, 3, 3, 4, 0, 0, 4, 2, 2, 3, 0, 5*1, 1,
153  $ 4*0, 2, 1, -1, 0, 0, 9, 8, 1, 0, 0, 4, 9, 1, 2,
154  $ -1, 5*2, 9, 4*0, 6, 4, 0, 0, 0, 3, 2, 1, 1, 0,
155  $ 5, 1, -1, 1, 0, 5*2, 4, 4*0, 2, 2, 0, 0, 0, 1,
156  $ 4, 4, 0, 0, 2, 4, 2, 2, -1, 5*2 /
157 * ..
158 * .. Executable Statements ..
159 *
160 * Get machine parameters
161 *
162  eps = dlamch( 'P' )
163  smlnum = dlamch( 'S' )
164  bignum = one / smlnum
165  CALL dlabad( smlnum, bignum )
166 *
167 * Set up test case parameters
168 *
169  vm1( 1 ) = one
170  vm1( 2 ) = sqrt( smlnum )
171  vm1( 3 ) = sqrt( vm1( 2 ) )
172  vm1( 4 ) = sqrt( bignum )
173  vm1( 5 ) = sqrt( vm1( 4 ) )
174 *
175  vm2( 1 ) = one
176  vm2( 2 ) = sqrt( smlnum )
177  vm2( 3 ) = sqrt( vm2( 2 ) )
178  vm2( 4 ) = sqrt( bignum )
179  vm2( 5 ) = sqrt( vm2( 4 ) )
180 *
181  vm3( 1 ) = one
182  vm3( 2 ) = sqrt( smlnum )
183  vm3( 3 ) = sqrt( vm3( 2 ) )
184  vm3( 4 ) = sqrt( bignum )
185  vm3( 5 ) = sqrt( vm3( 4 ) )
186 *
187  vm4( 1 ) = one
188  vm4( 2 ) = sqrt( smlnum )
189  vm4( 3 ) = sqrt( vm4( 2 ) )
190  vm4( 4 ) = sqrt( bignum )
191  vm4( 5 ) = sqrt( vm4( 4 ) )
192 *
193  vm5( 1 ) = one
194  vm5( 2 ) = eps
195  vm5( 3 ) = sqrt( smlnum )
196 *
197 * Initalization
198 *
199  knt = 0
200  rmax = zero
201  ninfo = 0
202  smlnum = smlnum / eps
203 *
204 * Begin test loop
205 *
206  DO 140 ivm5 = 1, 3
207  DO 130 ivm4 = 1, 5
208  DO 120 ivm3 = 1, 5
209  DO 110 ivm2 = 1, 5
210  DO 100 ivm1 = 1, 5
211  DO 90 ndim = 1, 6
212 *
213  n = idim( ndim )
214  DO 20 i = 1, n
215  DO 10 j = 1, n
216  t( i, j ) = dble( ival( i, j, ndim ) )*
217  $ vm1( ivm1 )
218  IF( i.GE.j )
219  $ t( i, j ) = t( i, j )*vm5( ivm5 )
220  10 CONTINUE
221  20 CONTINUE
222 *
223  w = one*vm2( ivm2 )
224 *
225  DO 30 i = 1, n
226  b( i ) = cos( dble( i ) )*vm3( ivm3 )
227  30 CONTINUE
228 *
229  DO 40 i = 1, 2*n
230  d( i ) = sin( dble( i ) )*vm4( ivm4 )
231  40 CONTINUE
232 *
233  norm = dlange( '1', n, n, t, ldt, work )
234  k = idamax( n, b, 1 )
235  normtb = norm + abs( b( k ) ) + abs( w )
236 *
237  CALL dcopy( n, d, 1, x, 1 )
238  knt = knt + 1
239  CALL dlaqtr( .false., .true., n, t, ldt, dum,
240  $ dumm, scale, x, work, info )
241  IF( info.NE.0 )
242  $ ninfo = ninfo + 1
243 *
244 * || T*x - scale*d || /
245 * max(ulp*||T||*||x||,smlnum/ulp*||T||,smlnum)
246 *
247  CALL dcopy( n, d, 1, y, 1 )
248  CALL dgemv( 'No transpose', n, n, one, t, ldt,
249  $ x, 1, -scale, y, 1 )
250  xnorm = dasum( n, x, 1 )
251  resid = dasum( n, y, 1 )
252  domin = max( smlnum, ( smlnum / eps )*norm,
253  $ ( norm*eps )*xnorm )
254  resid = resid / domin
255  IF( resid.GT.rmax ) THEN
256  rmax = resid
257  lmax = knt
258  END IF
259 *
260  CALL dcopy( n, d, 1, x, 1 )
261  knt = knt + 1
262  CALL dlaqtr( .true., .true., n, t, ldt, dum,
263  $ dumm, scale, x, work, info )
264  IF( info.NE.0 )
265  $ ninfo = ninfo + 1
266 *
267 * || T*x - scale*d || /
268 * max(ulp*||T||*||x||,smlnum/ulp*||T||,smlnum)
269 *
270  CALL dcopy( n, d, 1, y, 1 )
271  CALL dgemv( 'Transpose', n, n, one, t, ldt, x,
272  $ 1, -scale, y, 1 )
273  xnorm = dasum( n, x, 1 )
274  resid = dasum( n, y, 1 )
275  domin = max( smlnum, ( smlnum / eps )*norm,
276  $ ( norm*eps )*xnorm )
277  resid = resid / domin
278  IF( resid.GT.rmax ) THEN
279  rmax = resid
280  lmax = knt
281  END IF
282 *
283  CALL dcopy( 2*n, d, 1, x, 1 )
284  knt = knt + 1
285  CALL dlaqtr( .false., .false., n, t, ldt, b, w,
286  $ scale, x, work, info )
287  IF( info.NE.0 )
288  $ ninfo = ninfo + 1
289 *
290 * ||(T+i*B)*(x1+i*x2) - scale*(d1+i*d2)|| /
291 * max(ulp*(||T||+||B||)*(||x1||+||x2||),
292 * smlnum/ulp * (||T||+||B||), smlnum )
293 *
294 *
295  CALL dcopy( 2*n, d, 1, y, 1 )
296  y( 1 ) = ddot( n, b, 1, x( 1+n ), 1 ) +
297  $ scale*y( 1 )
298  DO 50 i = 2, n
299  y( i ) = w*x( i+n ) + scale*y( i )
300  50 CONTINUE
301  CALL dgemv( 'No transpose', n, n, one, t, ldt,
302  $ x, 1, -one, y, 1 )
303 *
304  y( 1+n ) = ddot( n, b, 1, x, 1 ) -
305  $ scale*y( 1+n )
306  DO 60 i = 2, n
307  y( i+n ) = w*x( i ) - scale*y( i+n )
308  60 CONTINUE
309  CALL dgemv( 'No transpose', n, n, one, t, ldt,
310  $ x( 1+n ), 1, one, y( 1+n ), 1 )
311 *
312  resid = dasum( 2*n, y, 1 )
313  domin = max( smlnum, ( smlnum / eps )*normtb,
314  $ eps*( normtb*dasum( 2*n, x, 1 ) ) )
315  resid = resid / domin
316  IF( resid.GT.rmax ) THEN
317  rmax = resid
318  lmax = knt
319  END IF
320 *
321  CALL dcopy( 2*n, d, 1, x, 1 )
322  knt = knt + 1
323  CALL dlaqtr( .true., .false., n, t, ldt, b, w,
324  $ scale, x, work, info )
325  IF( info.NE.0 )
326  $ ninfo = ninfo + 1
327 *
328 * ||(T+i*B)*(x1+i*x2) - scale*(d1+i*d2)|| /
329 * max(ulp*(||T||+||B||)*(||x1||+||x2||),
330 * smlnum/ulp * (||T||+||B||), smlnum )
331 *
332  CALL dcopy( 2*n, d, 1, y, 1 )
333  y( 1 ) = b( 1 )*x( 1+n ) - scale*y( 1 )
334  DO 70 i = 2, n
335  y( i ) = b( i )*x( 1+n ) + w*x( i+n ) -
336  $ scale*y( i )
337  70 CONTINUE
338  CALL dgemv( 'Transpose', n, n, one, t, ldt, x,
339  $ 1, one, y, 1 )
340 *
341  y( 1+n ) = b( 1 )*x( 1 ) + scale*y( 1+n )
342  DO 80 i = 2, n
343  y( i+n ) = b( i )*x( 1 ) + w*x( i ) +
344  $ scale*y( i+n )
345  80 CONTINUE
346  CALL dgemv( 'Transpose', n, n, one, t, ldt,
347  $ x( 1+n ), 1, -one, y( 1+n ), 1 )
348 *
349  resid = dasum( 2*n, y, 1 )
350  domin = max( smlnum, ( smlnum / eps )*normtb,
351  $ eps*( normtb*dasum( 2*n, x, 1 ) ) )
352  resid = resid / domin
353  IF( resid.GT.rmax ) THEN
354  rmax = resid
355  lmax = knt
356  END IF
357 *
358  90 CONTINUE
359  100 CONTINUE
360  110 CONTINUE
361  120 CONTINUE
362  130 CONTINUE
363  140 CONTINUE
364 *
365  RETURN
366 *
367 * End of DGET39
368 *
dcopy
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:53
idamax
integer function idamax(N, DX, INCX)
IDAMAX
Definition: idamax.f:53
dlamch
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
ddot
double precision function ddot(N, DX, INCX, DY, INCY)
DDOT
Definition: ddot.f:53
dgemv
subroutine dgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
DGEMV
Definition: dgemv.f:158
dlabad
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:76
dlange
double precision function dlange(NORM, M, N, A, LDA, WORK)
DLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: dlange.f:116
dlaqtr
subroutine dlaqtr(LTRAN, LREAL, N, T, LDT, B, W, SCALE, X, WORK, INFO)
DLAQTR solves a real quasi-triangular system of equations, or a complex quasi-triangular system of sp...
Definition: dlaqtr.f:167
dasum
double precision function dasum(N, DX, INCX)
DASUM
Definition: dasum.f:53

Here is the call graph for this function:

Here is the caller graph for this function: