LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dgtt05.f
Go to the documentation of this file.
1*> \brief \b DGTT05
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 DGTT05( TRANS, N, NRHS, DL, D, DU, B, LDB, X, LDX,
12* XACT, LDXACT, FERR, BERR, RESLTS )
13*
14* .. Scalar Arguments ..
15* CHARACTER TRANS
16* INTEGER LDB, LDX, LDXACT, N, NRHS
17* ..
18* .. Array Arguments ..
19* DOUBLE PRECISION B( LDB, * ), BERR( * ), D( * ), DL( * ),
20* $ DU( * ), FERR( * ), RESLTS( * ), X( LDX, * ),
21* $ XACT( LDXACT, * )
22* ..
23*
24*
25*> \par Purpose:
26* =============
27*>
28*> \verbatim
29*>
30*> DGTT05 tests the error bounds from iterative refinement for the
31*> computed solution to a system of equations A*X = B, where A is a
32*> general tridiagonal matrix of order n and op(A) = A or A**T,
33*> depending on TRANS.
34*>
35*> RESLTS(1) = test of the error bound
36*> = norm(X - XACT) / ( norm(X) * FERR )
37*>
38*> A large value is returned if this ratio is not less than one.
39*>
40*> RESLTS(2) = residual from the iterative refinement routine
41*> = the maximum of BERR / ( NZ*EPS + (*) ), where
42*> (*) = NZ*UNFL / (min_i (abs(op(A))*abs(X) +abs(b))_i )
43*> and NZ = max. number of nonzeros in any row of A, plus 1
44*> \endverbatim
45*
46* Arguments:
47* ==========
48*
49*> \param[in] TRANS
50*> \verbatim
51*> TRANS is CHARACTER*1
52*> Specifies the form of the system of equations.
53*> = 'N': A * X = B (No transpose)
54*> = 'T': A**T * X = B (Transpose)
55*> = 'C': A**H * X = B (Conjugate transpose = Transpose)
56*> \endverbatim
57*>
58*> \param[in] N
59*> \verbatim
60*> N is INTEGER
61*> The number of rows of the matrices X and XACT. N >= 0.
62*> \endverbatim
63*>
64*> \param[in] NRHS
65*> \verbatim
66*> NRHS is INTEGER
67*> The number of columns of the matrices X and XACT. NRHS >= 0.
68*> \endverbatim
69*>
70*> \param[in] DL
71*> \verbatim
72*> DL is DOUBLE PRECISION array, dimension (N-1)
73*> The (n-1) sub-diagonal elements of A.
74*> \endverbatim
75*>
76*> \param[in] D
77*> \verbatim
78*> D is DOUBLE PRECISION array, dimension (N)
79*> The diagonal elements of A.
80*> \endverbatim
81*>
82*> \param[in] DU
83*> \verbatim
84*> DU is DOUBLE PRECISION array, dimension (N-1)
85*> The (n-1) super-diagonal elements of A.
86*> \endverbatim
87*>
88*> \param[in] B
89*> \verbatim
90*> B is DOUBLE PRECISION array, dimension (LDB,NRHS)
91*> The right hand side vectors for the system of linear
92*> equations.
93*> \endverbatim
94*>
95*> \param[in] LDB
96*> \verbatim
97*> LDB is INTEGER
98*> The leading dimension of the array B. LDB >= max(1,N).
99*> \endverbatim
100*>
101*> \param[in] X
102*> \verbatim
103*> X is DOUBLE PRECISION array, dimension (LDX,NRHS)
104*> The computed solution vectors. Each vector is stored as a
105*> column of the matrix X.
106*> \endverbatim
107*>
108*> \param[in] LDX
109*> \verbatim
110*> LDX is INTEGER
111*> The leading dimension of the array X. LDX >= max(1,N).
112*> \endverbatim
113*>
114*> \param[in] XACT
115*> \verbatim
116*> XACT is DOUBLE PRECISION array, dimension (LDX,NRHS)
117*> The exact solution vectors. Each vector is stored as a
118*> column of the matrix XACT.
119*> \endverbatim
120*>
121*> \param[in] LDXACT
122*> \verbatim
123*> LDXACT is INTEGER
124*> The leading dimension of the array XACT. LDXACT >= max(1,N).
125*> \endverbatim
126*>
127*> \param[in] FERR
128*> \verbatim
129*> FERR is DOUBLE PRECISION array, dimension (NRHS)
130*> The estimated forward error bounds for each solution vector
131*> X. If XTRUE is the true solution, FERR bounds the magnitude
132*> of the largest entry in (X - XTRUE) divided by the magnitude
133*> of the largest entry in X.
134*> \endverbatim
135*>
136*> \param[in] BERR
137*> \verbatim
138*> BERR is DOUBLE PRECISION array, dimension (NRHS)
139*> The componentwise relative backward error of each solution
140*> vector (i.e., the smallest relative change in any entry of A
141*> or B that makes X an exact solution).
142*> \endverbatim
143*>
144*> \param[out] RESLTS
145*> \verbatim
146*> RESLTS is DOUBLE PRECISION array, dimension (2)
147*> The maximum over the NRHS solution vectors of the ratios:
148*> RESLTS(1) = norm(X - XACT) / ( norm(X) * FERR )
149*> RESLTS(2) = BERR / ( NZ*EPS + (*) )
150*> \endverbatim
151*
152* Authors:
153* ========
154*
155*> \author Univ. of Tennessee
156*> \author Univ. of California Berkeley
157*> \author Univ. of Colorado Denver
158*> \author NAG Ltd.
159*
160*> \ingroup double_lin
161*
162* =====================================================================
163 SUBROUTINE dgtt05( TRANS, N, NRHS, DL, D, DU, B, LDB, X, LDX,
164 $ XACT, LDXACT, FERR, BERR, RESLTS )
165*
166* -- LAPACK test routine --
167* -- LAPACK is a software package provided by Univ. of Tennessee, --
168* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
169*
170* .. Scalar Arguments ..
171 CHARACTER TRANS
172 INTEGER LDB, LDX, LDXACT, N, NRHS
173* ..
174* .. Array Arguments ..
175 DOUBLE PRECISION B( LDB, * ), BERR( * ), D( * ), DL( * ),
176 $ du( * ), ferr( * ), reslts( * ), x( ldx, * ),
177 $ xact( ldxact, * )
178* ..
179*
180* =====================================================================
181*
182* .. Parameters ..
183 DOUBLE PRECISION ZERO, ONE
184 parameter( zero = 0.0d+0, one = 1.0d+0 )
185* ..
186* .. Local Scalars ..
187 LOGICAL NOTRAN
188 INTEGER I, IMAX, J, K, NZ
189 DOUBLE PRECISION AXBI, DIFF, EPS, ERRBND, OVFL, TMP, UNFL, XNORM
190* ..
191* .. External Functions ..
192 LOGICAL LSAME
193 INTEGER IDAMAX
194 DOUBLE PRECISION DLAMCH
195 EXTERNAL lsame, idamax, dlamch
196* ..
197* .. Intrinsic Functions ..
198 INTRINSIC abs, max, min
199* ..
200* .. Executable Statements ..
201*
202* Quick exit if N = 0 or NRHS = 0.
203*
204 IF( n.LE.0 .OR. nrhs.LE.0 ) THEN
205 reslts( 1 ) = zero
206 reslts( 2 ) = zero
207 RETURN
208 END IF
209*
210 eps = dlamch( 'Epsilon' )
211 unfl = dlamch( 'Safe minimum' )
212 ovfl = one / unfl
213 notran = lsame( trans, 'N' )
214 nz = 4
215*
216* Test 1: Compute the maximum of
217* norm(X - XACT) / ( norm(X) * FERR )
218* over all the vectors X and XACT using the infinity-norm.
219*
220 errbnd = zero
221 DO 30 j = 1, nrhs
222 imax = idamax( n, x( 1, j ), 1 )
223 xnorm = max( abs( x( imax, j ) ), unfl )
224 diff = zero
225 DO 10 i = 1, n
226 diff = max( diff, abs( x( i, j )-xact( i, j ) ) )
227 10 CONTINUE
228*
229 IF( xnorm.GT.one ) THEN
230 GO TO 20
231 ELSE IF( diff.LE.ovfl*xnorm ) THEN
232 GO TO 20
233 ELSE
234 errbnd = one / eps
235 GO TO 30
236 END IF
237*
238 20 CONTINUE
239 IF( diff / xnorm.LE.ferr( j ) ) THEN
240 errbnd = max( errbnd, ( diff / xnorm ) / ferr( j ) )
241 ELSE
242 errbnd = one / eps
243 END IF
244 30 CONTINUE
245 reslts( 1 ) = errbnd
246*
247* Test 2: Compute the maximum of BERR / ( NZ*EPS + (*) ), where
248* (*) = NZ*UNFL / (min_i (abs(op(A))*abs(X) +abs(b))_i )
249*
250 DO 60 k = 1, nrhs
251 IF( notran ) THEN
252 IF( n.EQ.1 ) THEN
253 axbi = abs( b( 1, k ) ) + abs( d( 1 )*x( 1, k ) )
254 ELSE
255 axbi = abs( b( 1, k ) ) + abs( d( 1 )*x( 1, k ) ) +
256 $ abs( du( 1 )*x( 2, k ) )
257 DO 40 i = 2, n - 1
258 tmp = abs( b( i, k ) ) + abs( dl( i-1 )*x( i-1, k ) )
259 $ + abs( d( i )*x( i, k ) ) +
260 $ abs( du( i )*x( i+1, k ) )
261 axbi = min( axbi, tmp )
262 40 CONTINUE
263 tmp = abs( b( n, k ) ) + abs( dl( n-1 )*x( n-1, k ) ) +
264 $ abs( d( n )*x( n, k ) )
265 axbi = min( axbi, tmp )
266 END IF
267 ELSE
268 IF( n.EQ.1 ) THEN
269 axbi = abs( b( 1, k ) ) + abs( d( 1 )*x( 1, k ) )
270 ELSE
271 axbi = abs( b( 1, k ) ) + abs( d( 1 )*x( 1, k ) ) +
272 $ abs( dl( 1 )*x( 2, k ) )
273 DO 50 i = 2, n - 1
274 tmp = abs( b( i, k ) ) + abs( du( i-1 )*x( i-1, k ) )
275 $ + abs( d( i )*x( i, k ) ) +
276 $ abs( dl( i )*x( i+1, k ) )
277 axbi = min( axbi, tmp )
278 50 CONTINUE
279 tmp = abs( b( n, k ) ) + abs( du( n-1 )*x( n-1, k ) ) +
280 $ abs( d( n )*x( n, k ) )
281 axbi = min( axbi, tmp )
282 END IF
283 END IF
284 tmp = berr( k ) / ( nz*eps+nz*unfl / max( axbi, nz*unfl ) )
285 IF( k.EQ.1 ) THEN
286 reslts( 2 ) = tmp
287 ELSE
288 reslts( 2 ) = max( reslts( 2 ), tmp )
289 END IF
290 60 CONTINUE
291*
292 RETURN
293*
294* End of DGTT05
295*
296 END
subroutine dgtt05(trans, n, nrhs, dl, d, du, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
DGTT05
Definition dgtt05.f:165