dgbtf2()

subroutine dgbtf2	(	integer	m,
		integer	n,
		integer	kl,
		integer	ku,
		double precision, dimension( ldab, * )	ab,
		integer	ldab,
		integer, dimension( * )	ipiv,
		integer	info
	)

DGBTF2 computes the LU factorization of a general band matrix using the unblocked version of the algorithm.

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

Purpose:

 DGBTF2 computes an LU factorization of a real m-by-n band matrix A
 using partial pivoting with row interchanges.

 This is the unblocked version of the algorithm, calling Level 2 BLAS.

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]	AB	AB is DOUBLE PRECISION array, dimension (LDAB,N) On entry, the matrix A in band storage, in rows KL+1 to 2*KL+KU+1; rows 1 to KL of the array need not be set. The j-th column of A is stored in the j-th column of the array AB as follows: AB(kl+ku+1+i-j,j) = A(i,j) for max(1,j-ku)<=i<=min(m,j+kl) On exit, details of the factorization: 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 below for further details.
[in]	LDAB	LDAB is INTEGER The leading dimension of the array AB. LDAB >= 2*KL+KU+1.
[out]	IPIV	IPIV is INTEGER array, dimension (min(M,N)) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i).
[out]	INFO	INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = +i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, and division by zero will occur if it is used to solve a system of equations.

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Further Details:

  The band storage scheme is illustrated by the following example, when
  M = N = 6, KL = 2, KU = 1:

  On entry:                       On exit:

      *    *    *    +    +    +       *    *    *   u14  u25  u36
      *    *    +    +    +    +       *    *   u13  u24  u35  u46
      *   a12  a23  a34  a45  a56      *   u12  u23  u34  u45  u56
     a11  a22  a33  a44  a55  a66     u11  u22  u33  u44  u55  u66
     a21  a32  a43  a54  a65   *      m21  m32  m43  m54  m65   *
     a31  a42  a53  a64   *    *      m31  m42  m53  m64   *    *

  Array elements marked * are not used by the routine; elements marked
  + need not be set on entry, but are required by the routine to store
  elements of U, because of fill-in resulting from the row
  interchanges.

Definition at line 144 of file dgbtf2.f.

*
*  -- 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 ..
      INTEGER            INFO, KL, KU, LDAB, M, N
*     ..
*     .. Array Arguments ..
      INTEGER            IPIV( * )
      DOUBLE PRECISION   AB( LDAB, * )
*     ..
*
*  =====================================================================
*
*     .. Parameters ..
      DOUBLE PRECISION   ONE, ZERO
      parameter( one = 1.0d+0, zero = 0.0d+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I, J, JP, JU, KM, KV
*     ..
*     .. External Functions ..
      INTEGER            IDAMAX
      EXTERNAL           idamax
*     ..
*     .. External Subroutines ..
      EXTERNAL           dger, dscal, dswap, xerbla
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          max, min
*     ..
*     .. Executable Statements ..
*
*     KV is the number of superdiagonals in the factor U, allowing for
*     fill-in.
*
      kv = ku + kl
*
*     Test the input parameters.
*
      info = 0
      IF( m.LT.0 ) THEN
         info = -1
      ELSE IF( n.LT.0 ) THEN
         info = -2
      ELSE IF( kl.LT.0 ) THEN
         info = -3
      ELSE IF( ku.LT.0 ) THEN
         info = -4
      ELSE IF( ldab.LT.kl+kv+1 ) THEN
         info = -6
      END IF
      IF( info.NE.0 ) THEN
         CALL xerbla( 'DGBTF2', -info )
         RETURN
      END IF
*
*     Quick return if possible
*
      IF( m.EQ.0 .OR. n.EQ.0 )
     $   RETURN
*
*     Gaussian elimination with partial pivoting
*
*     Set fill-in elements in columns KU+2 to KV to zero.
*
      DO 20 j = ku + 2, min( kv, n )
         DO 10 i = kv - j + 2, kl
            ab( i, j ) = zero
   10    CONTINUE
   20 CONTINUE
*
*     JU is the index of the last column affected by the current stage
*     of the factorization.
*
      ju = 1
*
      DO 40 j = 1, min( m, n )
*
*        Set fill-in elements in column J+KV to zero.
*
         IF( j+kv.LE.n ) THEN
            DO 30 i = 1, kl
               ab( i, j+kv ) = zero
   30       CONTINUE
         END IF
*
*        Find pivot and test for singularity. KM is the number of
*        subdiagonal elements in the current column.
*
         km = min( kl, m-j )
         jp = idamax( km+1, ab( kv+1, j ), 1 )
         ipiv( j ) = jp + j - 1
         IF( ab( kv+jp, j ).NE.zero ) THEN
            ju = max( ju, min( j+ku+jp-1, n ) )
*
*           Apply interchange to columns J to JU.
*
            IF( jp.NE.1 )
     $         CALL dswap( ju-j+1, ab( kv+jp, j ), ldab-1,
     $                     ab( kv+1, j ), ldab-1 )
*
            IF( km.GT.0 ) THEN
*
*              Compute multipliers.
*
               CALL dscal( km, one / ab( kv+1, j ), ab( kv+2, j ), 1 )
*
*              Update trailing submatrix within the band.
*
               IF( ju.GT.j )
     $            CALL dger( km, ju-j, -one, ab( kv+2, j ), 1,
     $                       ab( kv, j+1 ), ldab-1, ab( kv+1, j+1 ),
     $                       ldab-1 )
            END IF
         ELSE
*
*           If pivot is zero, set INFO to the index of the pivot
*           unless a zero pivot has already been found.
*
            IF( info.EQ.0 )
     $         info = j
         END IF
   40 CONTINUE
      RETURN
*
*     End of DGBTF2
*

Here is the call graph for this function:

Here is the caller graph for this function: