LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
double precision function dla_gbrcond ( character  TRANS,
integer  N,
integer  KL,
integer  KU,
double precision, dimension( ldab, * )  AB,
integer  LDAB,
double precision, dimension( ldafb, * )  AFB,
integer  LDAFB,
integer, dimension( * )  IPIV,
integer  CMODE,
double precision, dimension( * )  C,
integer  INFO,
double precision, dimension( * )  WORK,
integer, dimension( * )  IWORK 
)

DLA_GBRCOND estimates the Skeel condition number for a general banded matrix.

Download DLA_GBRCOND + dependencies [TGZ] [ZIP] [TXT]

Purpose: 
    DLA_GBRCOND Estimates the Skeel condition number of  op(A) * op2(C)
    where op2 is determined by CMODE as follows
    CMODE =  1    op2(C) = C
    CMODE =  0    op2(C) = I
    CMODE = -1    op2(C) = inv(C)
    The Skeel condition number  cond(A) = norminf( |inv(A)||A| )
    is computed by computing scaling factors R such that
    diag(R)*A*op2(C) is row equilibrated and computing the standard
    infinity-norm condition number.
Parameters
[in]TRANS
          TRANS is CHARACTER*1
     Specifies the form of the system of equations:
       = 'N':  A * X = B     (No transpose)
       = 'T':  A**T * X = B  (Transpose)
       = 'C':  A**H * X = B  (Conjugate Transpose = Transpose)
[in]N
          N is INTEGER
     The number of linear equations, i.e., the order 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]AB
          AB is DOUBLE PRECISION array, dimension (LDAB,N)
     On entry, the matrix A in band storage, in rows 1 to KL+KU+1.
     The j-th column of A is stored in the j-th column of the
     array AB as follows:
     AB(KU+1+i-j,j) = A(i,j) for max(1,j-KU)<=i<=min(N,j+kl)
[in]LDAB
          LDAB is INTEGER
     The leading dimension of the array AB.  LDAB >= KL+KU+1.
[in]AFB
          AFB is DOUBLE PRECISION array, dimension (LDAFB,N)
     Details of the LU factorization of the band matrix A, 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.
[in]LDAFB
          LDAFB is INTEGER
     The leading dimension of the array AFB.  LDAFB >= 2*KL+KU+1.
[in]IPIV
          IPIV is INTEGER array, dimension (N)
     The pivot indices from the factorization A = P*L*U
     as computed by DGBTRF; row i of the matrix was interchanged
     with row IPIV(i).
[in]CMODE
          CMODE is INTEGER
     Determines op2(C) in the formula op(A) * op2(C) as follows:
     CMODE =  1    op2(C) = C
     CMODE =  0    op2(C) = I
     CMODE = -1    op2(C) = inv(C)
[in]C
          C is DOUBLE PRECISION array, dimension (N)
     The vector C in the formula op(A) * op2(C).
[out]INFO
          INFO is INTEGER
       = 0:  Successful exit.
     i > 0:  The ith argument is invalid.
[in]WORK
          WORK is DOUBLE PRECISION array, dimension (5*N).
     Workspace.
[in]IWORK
          IWORK is INTEGER array, dimension (N).
     Workspace.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 172 of file dla_gbrcond.f.

172 *
173 * -- LAPACK computational routine (version 3.4.2) --
174 * -- LAPACK is a software package provided by Univ. of Tennessee, --
175 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
176 * September 2012
177 *
178 * .. Scalar Arguments ..
179  CHARACTER trans
180  INTEGER n, ldab, ldafb, info, kl, ku, cmode
181 * ..
182 * .. Array Arguments ..
183  INTEGER iwork( * ), ipiv( * )
184  DOUBLE PRECISION ab( ldab, * ), afb( ldafb, * ), work( * ),
185  $ c( * )
186 * ..
187 *
188 * =====================================================================
189 *
190 * .. Local Scalars ..
191  LOGICAL notrans
192  INTEGER kase, i, j, kd, ke
193  DOUBLE PRECISION ainvnm, tmp
194 * ..
195 * .. Local Arrays ..
196  INTEGER isave( 3 )
197 * ..
198 * .. External Functions ..
199  LOGICAL lsame
200  EXTERNAL lsame
201 * ..
202 * .. External Subroutines ..
203  EXTERNAL dlacn2, dgbtrs, xerbla
204 * ..
205 * .. Intrinsic Functions ..
206  INTRINSIC abs, max
207 * ..
208 * .. Executable Statements ..
209 *
210  dla_gbrcond = 0.0d+0
211 *
212  info = 0
213  notrans = lsame( trans, 'N' )
214  IF ( .NOT. notrans .AND. .NOT. lsame(trans, 'T')
215  $ .AND. .NOT. lsame(trans, 'C') ) THEN
216  info = -1
217  ELSE IF( n.LT.0 ) THEN
218  info = -2
219  ELSE IF( kl.LT.0 .OR. kl.GT.n-1 ) THEN
220  info = -3
221  ELSE IF( ku.LT.0 .OR. ku.GT.n-1 ) THEN
222  info = -4
223  ELSE IF( ldab.LT.kl+ku+1 ) THEN
224  info = -6
225  ELSE IF( ldafb.LT.2*kl+ku+1 ) THEN
226  info = -8
227  END IF
228  IF( info.NE.0 ) THEN
229  CALL xerbla( 'DLA_GBRCOND', -info )
230  RETURN
231  END IF
232  IF( n.EQ.0 ) THEN
233  dla_gbrcond = 1.0d+0
234  RETURN
235  END IF
236 *
237 * Compute the equilibration matrix R such that
238 * inv(R)*A*C has unit 1-norm.
239 *
240  kd = ku + 1
241  ke = kl + 1
242  IF ( notrans ) THEN
243  DO i = 1, n
244  tmp = 0.0d+0
245  IF ( cmode .EQ. 1 ) THEN
246  DO j = max( i-kl, 1 ), min( i+ku, n )
247  tmp = tmp + abs( ab( kd+i-j, j ) * c( j ) )
248  END DO
249  ELSE IF ( cmode .EQ. 0 ) THEN
250  DO j = max( i-kl, 1 ), min( i+ku, n )
251  tmp = tmp + abs( ab( kd+i-j, j ) )
252  END DO
253  ELSE
254  DO j = max( i-kl, 1 ), min( i+ku, n )
255  tmp = tmp + abs( ab( kd+i-j, j ) / c( j ) )
256  END DO
257  END IF
258  work( 2*n+i ) = tmp
259  END DO
260  ELSE
261  DO i = 1, n
262  tmp = 0.0d+0
263  IF ( cmode .EQ. 1 ) THEN
264  DO j = max( i-kl, 1 ), min( i+ku, n )
265  tmp = tmp + abs( ab( ke-i+j, i ) * c( j ) )
266  END DO
267  ELSE IF ( cmode .EQ. 0 ) THEN
268  DO j = max( i-kl, 1 ), min( i+ku, n )
269  tmp = tmp + abs( ab( ke-i+j, i ) )
270  END DO
271  ELSE
272  DO j = max( i-kl, 1 ), min( i+ku, n )
273  tmp = tmp + abs( ab( ke-i+j, i ) / c( j ) )
274  END DO
275  END IF
276  work( 2*n+i ) = tmp
277  END DO
278  END IF
279 *
280 * Estimate the norm of inv(op(A)).
281 *
282  ainvnm = 0.0d+0
283 
284  kase = 0
285  10 CONTINUE
286  CALL dlacn2( n, work( n+1 ), work, iwork, ainvnm, kase, isave )
287  IF( kase.NE.0 ) THEN
288  IF( kase.EQ.2 ) THEN
289 *
290 * Multiply by R.
291 *
292  DO i = 1, n
293  work( i ) = work( i ) * work( 2*n+i )
294  END DO
295 
296  IF ( notrans ) THEN
297  CALL dgbtrs( 'No transpose', n, kl, ku, 1, afb, ldafb,
298  $ ipiv, work, n, info )
299  ELSE
300  CALL dgbtrs( 'Transpose', n, kl, ku, 1, afb, ldafb, ipiv,
301  $ work, n, info )
302  END IF
303 *
304 * Multiply by inv(C).
305 *
306  IF ( cmode .EQ. 1 ) THEN
307  DO i = 1, n
308  work( i ) = work( i ) / c( i )
309  END DO
310  ELSE IF ( cmode .EQ. -1 ) THEN
311  DO i = 1, n
312  work( i ) = work( i ) * c( i )
313  END DO
314  END IF
315  ELSE
316 *
317 * Multiply by inv(C**T).
318 *
319  IF ( cmode .EQ. 1 ) THEN
320  DO i = 1, n
321  work( i ) = work( i ) / c( i )
322  END DO
323  ELSE IF ( cmode .EQ. -1 ) THEN
324  DO i = 1, n
325  work( i ) = work( i ) * c( i )
326  END DO
327  END IF
328 
329  IF ( notrans ) THEN
330  CALL dgbtrs( 'Transpose', n, kl, ku, 1, afb, ldafb, ipiv,
331  $ work, n, info )
332  ELSE
333  CALL dgbtrs( 'No transpose', n, kl, ku, 1, afb, ldafb,
334  $ ipiv, work, n, info )
335  END IF
336 *
337 * Multiply by R.
338 *
339  DO i = 1, n
340  work( i ) = work( i ) * work( 2*n+i )
341  END DO
342  END IF
343  GO TO 10
344  END IF
345 *
346 * Compute the estimate of the reciprocal condition number.
347 *
348  IF( ainvnm .NE. 0.0d+0 )
349  $ dla_gbrcond = ( 1.0d+0 / ainvnm )
350 *
351  RETURN
352 *
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
dla_gbrcond
double precision function dla_gbrcond(TRANS, N, KL, KU, AB, LDAB, AFB, LDAFB, IPIV, CMODE, C, INFO, WORK, IWORK)
DLA_GBRCOND estimates the Skeel condition number for a general banded matrix.
Definition: dla_gbrcond.f:172
lsame
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
dgbtrs
subroutine dgbtrs(TRANS, N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB, INFO)
DGBTRS
Definition: dgbtrs.f:140
dlacn2
subroutine dlacn2(N, V, X, ISGN, EST, KASE, ISAVE)
DLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vec...
Definition: dlacn2.f:138

Here is the call graph for this function:

Here is the caller graph for this function: