dgbcon_8f_source.html

*> \brief \b DGBCON

*

*  =========== DOCUMENTATION ===========

*

* Online html documentation available at

*            http://www.netlib.org/lapack/explore-html/

*

*> \htmlonly

*> Download DGBCON + dependencies

*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgbcon.f">

*> [TGZ]</a>

*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgbcon.f">

*> [ZIP]</a>

*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgbcon.f">

*> [TXT]</a>

*> \endhtmlonly

*

*  Definition:

*  ===========

*

*       SUBROUTINE DGBCON( NORM, N, KL, KU, AB, LDAB, IPIV, ANORM, RCOND,

*                          WORK, IWORK, INFO )

*

*       .. Scalar Arguments ..

*       CHARACTER          NORM

*       INTEGER            INFO, KL, KU, LDAB, N

*       DOUBLE PRECISION   ANORM, RCOND

*       ..

*       .. Array Arguments ..

*       INTEGER            IPIV( * ), IWORK( * )

*       DOUBLE PRECISION   AB( LDAB, * ), WORK( * )

*       ..

*

*

*> \par Purpose:

*  =============

*>

*> \verbatim

*>

*> DGBCON estimates the reciprocal of the condition number of a real

*> general band matrix A, in either the 1-norm or the infinity-norm,

*> using the LU factorization computed by DGBTRF.

*>

*> An estimate is obtained for norm(inv(A)), and the reciprocal of the

*> condition number is computed as

*>    RCOND = 1 / ( norm(A) * norm(inv(A)) ).

*> \endverbatim

*

*  Arguments:

*  ==========

*

*> \param[in] NORM

*> \verbatim

*>          NORM is CHARACTER*1

*>          Specifies whether the 1-norm condition number or the

*>          infinity-norm condition number is required:

*>          = '1' or 'O':  1-norm;

*>          = 'I':         Infinity-norm.

*> \endverbatim

*>

*> \param[in] N

*> \verbatim

*>          N is INTEGER

*>          The order of the matrix A.  N >= 0.

*> \endverbatim

*>

*> \param[in] KL

*> \verbatim

*>          KL is INTEGER

*>          The number of subdiagonals within the band of A.  KL >= 0.

*> \endverbatim

*>

*> \param[in] KU

*> \verbatim

*>          KU is INTEGER

*>          The number of superdiagonals within the band of A.  KU >= 0.

*> \endverbatim

*>

*> \param[in] AB

*> \verbatim

*>          AB is DOUBLE PRECISION array, dimension (LDAB,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.

*> \endverbatim

*>

*> \param[in] LDAB

*> \verbatim

*>          LDAB is INTEGER

*>          The leading dimension of the array AB.  LDAB >= 2*KL+KU+1.

*> \endverbatim

*>

*> \param[in] IPIV

*> \verbatim

*>          IPIV is INTEGER array, dimension (N)

*>          The pivot indices; for 1 <= i <= N, row i of the matrix was

*>          interchanged with row IPIV(i).

*> \endverbatim

*>

*> \param[in] ANORM

*> \verbatim

*>          ANORM is DOUBLE PRECISION

*>          If NORM = '1' or 'O', the 1-norm of the original matrix A.

*>          If NORM = 'I', the infinity-norm of the original matrix A.

*> \endverbatim

*>

*> \param[out] RCOND

*> \verbatim

*>          RCOND is DOUBLE PRECISION

*>          The reciprocal of the condition number of the matrix A,

*>          computed as RCOND = 1/(norm(A) * norm(inv(A))).

*> \endverbatim

*>

*> \param[out] WORK

*> \verbatim

*>          WORK is DOUBLE PRECISION array, dimension (3*N)

*> \endverbatim

*>

*> \param[out] IWORK

*> \verbatim

*>          IWORK is INTEGER array, dimension (N)

*> \endverbatim

*>

*> \param[out] INFO

*> \verbatim

*>          INFO is INTEGER

*>          = 0:  successful exit

*>          < 0: if INFO = -i, the i-th argument had an illegal value

*> \endverbatim

*

*  Authors:

*  ========

*

*> \author Univ. of Tennessee

*> \author Univ. of California Berkeley

*> \author Univ. of Colorado Denver

*> \author NAG Ltd.

*

*> \ingroup doubleGBcomputational

*

*  =====================================================================


      SUBROUTINE dgbcon( NORM, N, KL, KU, AB, LDAB, IPIV, ANORM, RCOND,

     $                   WORK, IWORK, INFO )

*

*  -- LAPACK computational routine --

*  -- LAPACK is a software package provided by Univ. of Tennessee,    --

*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--

*

*     .. Scalar Arguments ..

      CHARACTER          NORM

      INTEGER            INFO, KL, KU, LDAB, N

      DOUBLE PRECISION   ANORM, RCOND

*     ..

*     .. Array Arguments ..

      INTEGER            IPIV( * ), IWORK( * )

      DOUBLE PRECISION   AB( LDAB, * ), WORK( * )

*     ..

*

*  =====================================================================

*

*     .. Parameters ..

      DOUBLE PRECISION   ONE, ZERO

      parameter( one = 1.0d+0, zero = 0.0d+0 )

*     ..

*     .. Local Scalars ..

      LOGICAL            LNOTI, ONENRM

      CHARACTER          NORMIN

      INTEGER            IX, J, JP, KASE, KASE1, KD, LM

      DOUBLE PRECISION   AINVNM, SCALE, SMLNUM, T

*     ..

*     .. Local Arrays ..

      INTEGER            ISAVE( 3 )

*     ..

*     .. External Functions ..

      LOGICAL            LSAME

      INTEGER            IDAMAX

      DOUBLE PRECISION   DDOT, DLAMCH

      EXTERNAL           lsame, idamax, ddot, dlamch

*     ..

*     .. External Subroutines ..

      EXTERNAL           daxpy, dlacn2, dlatbs, drscl, xerbla

*     ..

*     .. Intrinsic Functions ..

      INTRINSIC          abs, min

*     ..

*     .. Executable Statements ..

*

*     Test the input parameters.

*

      info = 0

      onenrm = norm.EQ.'1' .OR. lsame( norm, 'o' )

.NOT..AND..NOT.      IF( ONENRM  LSAME( NORM, 'i' ) ) THEN

         INFO = -1

.LT.      ELSE IF( N0 ) THEN

         INFO = -2

.LT.      ELSE IF( KL0 ) THEN

         INFO = -3

.LT.      ELSE IF( KU0 ) THEN

         INFO = -4

.LT.      ELSE IF( LDAB2*KL+KU+1 ) THEN

         INFO = -6

.LT.      ELSE IF( ANORMZERO ) THEN

         INFO = -8

      END IF

.NE.      IF( INFO0 ) THEN

         CALL XERBLA( 'dgbcon', -INFO )

         RETURN

      END IF

*

*     Quick return if possible

*

      RCOND = ZERO

.EQ.      IF( N0 ) THEN

         RCOND = ONE

         RETURN

.EQ.      ELSE IF( ANORMZERO ) THEN

         RETURN

      END IF

*

      SMLNUM = DLAMCH( 'safe minimum' )

*

*     Estimate the norm of inv(A).

*

      AINVNM = ZERO

      NORMIN = 'n'

      IF( ONENRM ) THEN

         KASE1 = 1

      ELSE

         KASE1 = 2

      END IF

      KD = KL + KU + 1

.GT.      LNOTI = KL0

      KASE = 0

   10 CONTINUE

      CALL DLACN2( N, WORK( N+1 ), WORK, IWORK, AINVNM, KASE, ISAVE )

.NE.      IF( KASE0 ) THEN

.EQ.         IF( KASEKASE1 ) THEN

*

*           Multiply by inv(L).

*

            IF( LNOTI ) THEN

               DO 20 J = 1, N - 1

                  LM = MIN( KL, N-J )

                  JP = IPIV( J )

                  T = WORK( JP )

.NE.                  IF( JPJ ) THEN

                     WORK( JP ) = WORK( J )

                     WORK( J ) = T

                  END IF

                  CALL DAXPY( LM, -T, AB( KD+1, J ), 1, WORK( J+1 ), 1 )

   20          CONTINUE

            END IF

*

*           Multiply by inv(U).

*

            CALL DLATBS( 'upper', 'no transpose', 'non-unit', NORMIN, N,

     $                   KL+KU, AB, LDAB, WORK, SCALE, WORK( 2*N+1 ),

     $                   INFO )

         ELSE

*

*           Multiply by inv(U**T).

*

            CALL DLATBS( 'upper', 'transpose', 'non-unit', NORMIN, N,

     $                   KL+KU, AB, LDAB, WORK, SCALE, WORK( 2*N+1 ),

     $                   INFO )

*

*           Multiply by inv(L**T).

*

            IF( LNOTI ) THEN

               DO 30 J = N - 1, 1, -1

                  LM = MIN( KL, N-J )

                  WORK( J ) = WORK( J ) - DDOT( LM, AB( KD+1, J ), 1,

     $                        WORK( J+1 ), 1 )

                  JP = IPIV( J )

.NE.                  IF( JPJ ) THEN

                     T = WORK( JP )

                     WORK( JP ) = WORK( J )

                     WORK( J ) = T

                  END IF

   30          CONTINUE

            END IF

         END IF

*

*        Divide X by 1/SCALE if doing so will not cause overflow.

*

         NORMIN = 'y'

.NE.         IF( SCALEONE ) THEN

            IX = IDAMAX( N, WORK, 1 )

.LT..OR..EQ.            IF( SCALEABS( WORK( IX ) )*SMLNUM  SCALEZERO )

     $         GO TO 40

            CALL DRSCL( N, SCALE, WORK, 1 )

         END IF

         GO TO 10

      END IF

*

*     Compute the estimate of the reciprocal condition number.

*

.NE.      IF( AINVNMZERO )

     $   RCOND = ( ONE / AINVNM ) / ANORM

*

   40 CONTINUE

      RETURN

*

*     End of DGBCON

*


      END

xerbla
subroutine xerbla(srname, info)
XERBLA
Definition xerbla.f:60

dgbcon
subroutine dgbcon(norm, n, kl, ku, ab, ldab, ipiv, anorm, rcond, work, iwork, info)
DGBCON
Definition dgbcon.f:146

drscl
subroutine drscl(n, sa, sx, incx)
DRSCL multiplies a vector by the reciprocal of a real scalar.
Definition drscl.f:84

dlatbs
subroutine dlatbs(uplo, trans, diag, normin, n, kd, ab, ldab, x, scale, cnorm, info)
DLATBS solves a triangular banded system of equations.
Definition dlatbs.f:242

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:136

daxpy
subroutine daxpy(n, da, dx, incx, dy, incy)
DAXPY
Definition daxpy.f:89

min
#define min(a, b)
Definition macros.h:20