LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ dgbt01()

 subroutine dgbt01 ( integer m, integer n, integer kl, integer ku, double precision, dimension( lda, * ) a, integer lda, double precision, dimension( ldafac, * ) afac, integer ldafac, integer, dimension( * ) ipiv, double precision, dimension( * ) work, double precision resid )

DGBT01

Purpose:
``` DGBT01 reconstructs a band matrix A from its L*U factorization and
computes the residual:
norm(L*U - A) / ( N * norm(A) * EPS ),
where EPS is the machine epsilon.

The expression L*U - A is computed one column at a time, so A and
AFAC are not modified.```
Parameters
 [in] M ``` M is INTEGER The number of rows of the matrix A. M >= 0.``` [in] N ``` N is INTEGER The number of columns of the matrix A. N >= 0.``` [in] KL ``` KL is INTEGER The number of subdiagonals within the band of A. KL >= 0.``` [in] KU ``` KU is INTEGER The number of superdiagonals within the band of A. KU >= 0.``` [in,out] A ``` A is DOUBLE PRECISION array, dimension (LDA,N) The original matrix A in band storage, stored in rows 1 to KL+KU+1.``` [in] LDA ``` LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,KL+KU+1).``` [in] AFAC ``` AFAC is DOUBLE PRECISION array, dimension (LDAFAC,N) The factored form of the matrix A. AFAC contains the banded factors L and U from the L*U factorization, as computed by DGBTRF. U is stored as an upper triangular band matrix with KL+KU superdiagonals in rows 1 to KL+KU+1, and the multipliers used during the factorization are stored in rows KL+KU+2 to 2*KL+KU+1. See DGBTRF for further details.``` [in] LDAFAC ``` LDAFAC is INTEGER The leading dimension of the array AFAC. LDAFAC >= max(1,2*KL*KU+1).``` [in] IPIV ``` IPIV is INTEGER array, dimension (min(M,N)) The pivot indices from DGBTRF.``` [out] WORK ` WORK is DOUBLE PRECISION array, dimension (2*KL+KU+1)` [out] RESID ``` RESID is DOUBLE PRECISION norm(L*U - A) / ( N * norm(A) * EPS )```

Definition at line 124 of file dgbt01.f.

126*
127* -- LAPACK test routine --
128* -- LAPACK is a software package provided by Univ. of Tennessee, --
129* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
130*
131* .. Scalar Arguments ..
132 INTEGER KL, KU, LDA, LDAFAC, M, N
133 DOUBLE PRECISION RESID
134* ..
135* .. Array Arguments ..
136 INTEGER IPIV( * )
137 DOUBLE PRECISION A( LDA, * ), AFAC( LDAFAC, * ), WORK( * )
138* ..
139*
140* =====================================================================
141*
142* .. Parameters ..
143 DOUBLE PRECISION ZERO, ONE
144 parameter( zero = 0.0d+0, one = 1.0d+0 )
145* ..
146* .. Local Scalars ..
147 INTEGER I, I1, I2, IL, IP, IW, J, JL, JU, JUA, KD, LENJ
148 DOUBLE PRECISION ANORM, EPS, T
149* ..
150* .. External Functions ..
151 DOUBLE PRECISION DASUM, DLAMCH
152 EXTERNAL dasum, dlamch
153* ..
154* .. External Subroutines ..
155 EXTERNAL daxpy, dcopy
156* ..
157* .. Intrinsic Functions ..
158 INTRINSIC dble, max, min
159* ..
160* .. Executable Statements ..
161*
162* Quick exit if M = 0 or N = 0.
163*
164 resid = zero
165 IF( m.LE.0 .OR. n.LE.0 )
166 \$ RETURN
167*
168* Determine EPS and the norm of A.
169*
170 eps = dlamch( 'Epsilon' )
171 kd = ku + 1
172 anorm = zero
173 DO 10 j = 1, n
174 i1 = max( kd+1-j, 1 )
175 i2 = min( kd+m-j, kl+kd )
176 IF( i2.GE.i1 )
177 \$ anorm = max( anorm, dasum( i2-i1+1, a( i1, j ), 1 ) )
178 10 CONTINUE
179*
180* Compute one column at a time of L*U - A.
181*
182 kd = kl + ku + 1
183 DO 40 j = 1, n
184*
185* Copy the J-th column of U to WORK.
186*
187 ju = min( kl+ku, j-1 )
188 jl = min( kl, m-j )
189 lenj = min( m, j ) - j + ju + 1
190 IF( lenj.GT.0 ) THEN
191 CALL dcopy( lenj, afac( kd-ju, j ), 1, work, 1 )
192 DO 20 i = lenj + 1, ju + jl + 1
193 work( i ) = zero
194 20 CONTINUE
195*
196* Multiply by the unit lower triangular matrix L. Note that L
197* is stored as a product of transformations and permutations.
198*
199 DO 30 i = min( m-1, j ), j - ju, -1
200 il = min( kl, m-i )
201 IF( il.GT.0 ) THEN
202 iw = i - j + ju + 1
203 t = work( iw )
204 CALL daxpy( il, t, afac( kd+1, i ), 1, work( iw+1 ),
205 \$ 1 )
206 ip = ipiv( i )
207 IF( i.NE.ip ) THEN
208 ip = ip - j + ju + 1
209 work( iw ) = work( ip )
210 work( ip ) = t
211 END IF
212 END IF
213 30 CONTINUE
214*
215* Subtract the corresponding column of A.
216*
217 jua = min( ju, ku )
218 IF( jua+jl+1.GT.0 )
219 \$ CALL daxpy( jua+jl+1, -one, a( ku+1-jua, j ), 1,
220 \$ work( ju+1-jua ), 1 )
221*
222* Compute the 1-norm of the column.
223*
224 resid = max( resid, dasum( ju+jl+1, work, 1 ) )
225 END IF
226 40 CONTINUE
227*
228* Compute norm(L*U - A) / ( N * norm(A) * EPS )
229*
230 IF( anorm.LE.zero ) THEN
231 IF( resid.NE.zero )
232 \$ resid = one / eps
233 ELSE
234 resid = ( ( resid / dble( n ) ) / anorm ) / eps
235 END IF
236*
237 RETURN
238*
239* End of DGBT01
240*
double precision function dasum(n, dx, incx)
DASUM
Definition dasum.f:71
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
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
Here is the call graph for this function:
Here is the caller graph for this function: