LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches
dtrt02.f
Go to the documentation of this file.
1*> \brief \b DTRT02
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 DTRT02( UPLO, TRANS, DIAG, N, NRHS, A, LDA, X, LDX, B,
12* LDB, WORK, RESID )
13*
14* .. Scalar Arguments ..
15* CHARACTER DIAG, TRANS, UPLO
16* INTEGER LDA, LDB, LDX, N, NRHS
17* DOUBLE PRECISION RESID
18* ..
19* .. Array Arguments ..
20* DOUBLE PRECISION A( LDA, * ), B( LDB, * ), WORK( * ),
21* \$ X( LDX, * )
22* ..
23*
24*
25*> \par Purpose:
26* =============
27*>
28*> \verbatim
29*>
30*> DTRT02 computes the residual for the computed solution to a
31*> triangular system of linear equations op(A)*X = B, where A is a
32*> triangular matrix. The test ratio is the maximum over
33*> norm(b - op(A)*x) / ( ||op(A)||_1 * norm(x) * EPS ),
34*> where op(A) = A or A**T, b is the column of B, x is the solution
35*> vector, and EPS is the machine epsilon.
36*> The norm used is the 1-norm.
37*> \endverbatim
38*
39* Arguments:
40* ==========
41*
42*> \param[in] UPLO
43*> \verbatim
44*> UPLO is CHARACTER*1
45*> Specifies whether the matrix A is upper or lower triangular.
46*> = 'U': Upper triangular
47*> = 'L': Lower triangular
48*> \endverbatim
49*>
50*> \param[in] TRANS
51*> \verbatim
52*> TRANS is CHARACTER*1
53*> Specifies the operation applied to A.
54*> = 'N': A * X = B (No transpose)
55*> = 'T': A**T * X = B (Transpose)
56*> = 'C': A**H * X = B (Conjugate transpose = Transpose)
57*> \endverbatim
58*>
59*> \param[in] DIAG
60*> \verbatim
61*> DIAG is CHARACTER*1
62*> Specifies whether or not the matrix A is unit triangular.
63*> = 'N': Non-unit triangular
64*> = 'U': Unit triangular
65*> \endverbatim
66*>
67*> \param[in] N
68*> \verbatim
69*> N is INTEGER
70*> The order of the matrix A. N >= 0.
71*> \endverbatim
72*>
73*> \param[in] NRHS
74*> \verbatim
75*> NRHS is INTEGER
76*> The number of right hand sides, i.e., the number of columns
77*> of the matrices X and B. NRHS >= 0.
78*> \endverbatim
79*>
80*> \param[in] A
81*> \verbatim
82*> A is DOUBLE PRECISION array, dimension (LDA,N)
83*> The triangular matrix A. If UPLO = 'U', the leading n by n
84*> upper triangular part of the array A contains the upper
85*> triangular matrix, and the strictly lower triangular part of
86*> A is not referenced. If UPLO = 'L', the leading n by n lower
87*> triangular part of the array A contains the lower triangular
88*> matrix, and the strictly upper triangular part of A is not
89*> referenced. If DIAG = 'U', the diagonal elements of A are
90*> also not referenced and are assumed to be 1.
91*> \endverbatim
92*>
93*> \param[in] LDA
94*> \verbatim
95*> LDA is INTEGER
96*> The leading dimension of the array A. LDA >= max(1,N).
97*> \endverbatim
98*>
99*> \param[in] X
100*> \verbatim
101*> X is DOUBLE PRECISION array, dimension (LDX,NRHS)
102*> The computed solution vectors for the system of linear
103*> equations.
104*> \endverbatim
105*>
106*> \param[in] LDX
107*> \verbatim
108*> LDX is INTEGER
109*> The leading dimension of the array X. LDX >= max(1,N).
110*> \endverbatim
111*>
112*> \param[in] B
113*> \verbatim
114*> B is DOUBLE PRECISION array, dimension (LDB,NRHS)
115*> The right hand side vectors for the system of linear
116*> equations.
117*> \endverbatim
118*>
119*> \param[in] LDB
120*> \verbatim
121*> LDB is INTEGER
122*> The leading dimension of the array B. LDB >= max(1,N).
123*> \endverbatim
124*>
125*> \param[out] WORK
126*> \verbatim
127*> WORK is DOUBLE PRECISION array, dimension (N)
128*> \endverbatim
129*>
130*> \param[out] RESID
131*> \verbatim
132*> RESID is DOUBLE PRECISION
133*> The maximum over the number of right hand sides of
134*> norm(op(A)*X - B) / ( norm(op(A)) * norm(X) * EPS ).
135*> \endverbatim
136*
137* Authors:
138* ========
139*
140*> \author Univ. of Tennessee
141*> \author Univ. of California Berkeley
142*> \author Univ. of Colorado Denver
143*> \author NAG Ltd.
144*
145*> \ingroup double_lin
146*
147* =====================================================================
148 SUBROUTINE dtrt02( UPLO, TRANS, DIAG, N, NRHS, A, LDA, X, LDX, B,
149 \$ LDB, WORK, RESID )
150*
151* -- LAPACK test routine --
152* -- LAPACK is a software package provided by Univ. of Tennessee, --
153* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
154*
155* .. Scalar Arguments ..
156 CHARACTER DIAG, TRANS, UPLO
157 INTEGER LDA, LDB, LDX, N, NRHS
158 DOUBLE PRECISION RESID
159* ..
160* .. Array Arguments ..
161 DOUBLE PRECISION A( LDA, * ), B( LDB, * ), WORK( * ),
162 \$ x( ldx, * )
163* ..
164*
165* =====================================================================
166*
167* .. Parameters ..
168 DOUBLE PRECISION ZERO, ONE
169 parameter( zero = 0.0d+0, one = 1.0d+0 )
170* ..
171* .. Local Scalars ..
172 INTEGER J
173 DOUBLE PRECISION ANORM, BNORM, EPS, XNORM
174* ..
175* .. External Functions ..
176 LOGICAL LSAME
177 DOUBLE PRECISION DASUM, DLAMCH, DLANTR
178 EXTERNAL lsame, dasum, dlamch, dlantr
179* ..
180* .. External Subroutines ..
181 EXTERNAL daxpy, dcopy, dtrmv
182* ..
183* .. Intrinsic Functions ..
184 INTRINSIC max
185* ..
186* .. Executable Statements ..
187*
188* Quick exit if N = 0 or NRHS = 0
189*
190 IF( n.LE.0 .OR. nrhs.LE.0 ) THEN
191 resid = zero
192 RETURN
193 END IF
194*
195* Compute the 1-norm of op(A).
196*
197 IF( lsame( trans, 'N' ) ) THEN
198 anorm = dlantr( '1', uplo, diag, n, n, a, lda, work )
199 ELSE
200 anorm = dlantr( 'I', uplo, diag, n, n, a, lda, work )
201 END IF
202*
203* Exit with RESID = 1/EPS if ANORM = 0.
204*
205 eps = dlamch( 'Epsilon' )
206 IF( anorm.LE.zero ) THEN
207 resid = one / eps
208 RETURN
209 END IF
210*
211* Compute the maximum over the number of right hand sides of
212* norm(op(A)*X - B) / ( norm(op(A)) * norm(X) * EPS )
213*
214 resid = zero
215 DO 10 j = 1, nrhs
216 CALL dcopy( n, x( 1, j ), 1, work, 1 )
217 CALL dtrmv( uplo, trans, diag, n, a, lda, work, 1 )
218 CALL daxpy( n, -one, b( 1, j ), 1, work, 1 )
219 bnorm = dasum( n, work, 1 )
220 xnorm = dasum( n, x( 1, j ), 1 )
221 IF( xnorm.LE.zero ) THEN
222 resid = one / eps
223 ELSE
224 resid = max( resid, ( ( bnorm / anorm ) / xnorm ) / eps )
225 END IF
226 10 CONTINUE
227*
228 RETURN
229*
230* End of DTRT02
231*
232 END
subroutine dtrt02(uplo, trans, diag, n, nrhs, a, lda, x, ldx, b, ldb, work, resid)
DTRT02
Definition dtrt02.f:150
subroutine daxpy(n, da, dx, incx, dy, incy)
DAXPY
Definition daxpy.f:89
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
Definition dcopy.f:82
subroutine dtrmv(uplo, trans, diag, n, a, lda, x, incx)
DTRMV
Definition dtrmv.f:147