LAPACK 3.3.0

dgbsv.f

Go to the documentation of this file.
00001       SUBROUTINE DGBSV( N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB, INFO )
00002 *
00003 *  -- LAPACK driver routine (version 3.2) --
00004 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00005 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00006 *     November 2006
00007 *
00008 *     .. Scalar Arguments ..
00009       INTEGER            INFO, KL, KU, LDAB, LDB, N, NRHS
00010 *     ..
00011 *     .. Array Arguments ..
00012       INTEGER            IPIV( * )
00013       DOUBLE PRECISION   AB( LDAB, * ), B( LDB, * )
00014 *     ..
00015 *
00016 *  Purpose
00017 *  =======
00018 *
00019 *  DGBSV computes the solution to a real system of linear equations
00020 *  A * X = B, where A is a band matrix of order N with KL subdiagonals
00021 *  and KU superdiagonals, and X and B are N-by-NRHS matrices.
00022 *
00023 *  The LU decomposition with partial pivoting and row interchanges is
00024 *  used to factor A as A = L * U, where L is a product of permutation
00025 *  and unit lower triangular matrices with KL subdiagonals, and U is
00026 *  upper triangular with KL+KU superdiagonals.  The factored form of A
00027 *  is then used to solve the system of equations A * X = B.
00028 *
00029 *  Arguments
00030 *  =========
00031 *
00032 *  N       (input) INTEGER
00033 *          The number of linear equations, i.e., the order of the
00034 *          matrix A.  N >= 0.
00035 *
00036 *  KL      (input) INTEGER
00037 *          The number of subdiagonals within the band of A.  KL >= 0.
00038 *
00039 *  KU      (input) INTEGER
00040 *          The number of superdiagonals within the band of A.  KU >= 0.
00041 *
00042 *  NRHS    (input) INTEGER
00043 *          The number of right hand sides, i.e., the number of columns
00044 *          of the matrix B.  NRHS >= 0.
00045 *
00046 *  AB      (input/output) DOUBLE PRECISION array, dimension (LDAB,N)
00047 *          On entry, the matrix A in band storage, in rows KL+1 to
00048 *          2*KL+KU+1; rows 1 to KL of the array need not be set.
00049 *          The j-th column of A is stored in the j-th column of the
00050 *          array AB as follows:
00051 *          AB(KL+KU+1+i-j,j) = A(i,j) for max(1,j-KU)<=i<=min(N,j+KL)
00052 *          On exit, details of the factorization: U is stored as an
00053 *          upper triangular band matrix with KL+KU superdiagonals in
00054 *          rows 1 to KL+KU+1, and the multipliers used during the
00055 *          factorization are stored in rows KL+KU+2 to 2*KL+KU+1.
00056 *          See below for further details.
00057 *
00058 *  LDAB    (input) INTEGER
00059 *          The leading dimension of the array AB.  LDAB >= 2*KL+KU+1.
00060 *
00061 *  IPIV    (output) INTEGER array, dimension (N)
00062 *          The pivot indices that define the permutation matrix P;
00063 *          row i of the matrix was interchanged with row IPIV(i).
00064 *
00065 *  B       (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS)
00066 *          On entry, the N-by-NRHS right hand side matrix B.
00067 *          On exit, if INFO = 0, the N-by-NRHS solution matrix X.
00068 *
00069 *  LDB     (input) INTEGER
00070 *          The leading dimension of the array B.  LDB >= max(1,N).
00071 *
00072 *  INFO    (output) INTEGER
00073 *          = 0:  successful exit
00074 *          < 0:  if INFO = -i, the i-th argument had an illegal value
00075 *          > 0:  if INFO = i, U(i,i) is exactly zero.  The factorization
00076 *                has been completed, but the factor U is exactly
00077 *                singular, and the solution has not been computed.
00078 *
00079 *  Further Details
00080 *  ===============
00081 *
00082 *  The band storage scheme is illustrated by the following example, when
00083 *  M = N = 6, KL = 2, KU = 1:
00084 *
00085 *  On entry:                       On exit:
00086 *
00087 *      *    *    *    +    +    +       *    *    *   u14  u25  u36
00088 *      *    *    +    +    +    +       *    *   u13  u24  u35  u46
00089 *      *   a12  a23  a34  a45  a56      *   u12  u23  u34  u45  u56
00090 *     a11  a22  a33  a44  a55  a66     u11  u22  u33  u44  u55  u66
00091 *     a21  a32  a43  a54  a65   *      m21  m32  m43  m54  m65   *
00092 *     a31  a42  a53  a64   *    *      m31  m42  m53  m64   *    *
00093 *
00094 *  Array elements marked * are not used by the routine; elements marked
00095 *  + need not be set on entry, but are required by the routine to store
00096 *  elements of U because of fill-in resulting from the row interchanges.
00097 *
00098 *  =====================================================================
00099 *
00100 *     .. External Subroutines ..
00101       EXTERNAL           DGBTRF, DGBTRS, XERBLA
00102 *     ..
00103 *     .. Intrinsic Functions ..
00104       INTRINSIC          MAX
00105 *     ..
00106 *     .. Executable Statements ..
00107 *
00108 *     Test the input parameters.
00109 *
00110       INFO = 0
00111       IF( N.LT.0 ) THEN
00112          INFO = -1
00113       ELSE IF( KL.LT.0 ) THEN
00114          INFO = -2
00115       ELSE IF( KU.LT.0 ) THEN
00116          INFO = -3
00117       ELSE IF( NRHS.LT.0 ) THEN
00118          INFO = -4
00119       ELSE IF( LDAB.LT.2*KL+KU+1 ) THEN
00120          INFO = -6
00121       ELSE IF( LDB.LT.MAX( N, 1 ) ) THEN
00122          INFO = -9
00123       END IF
00124       IF( INFO.NE.0 ) THEN
00125          CALL XERBLA( 'DGBSV ', -INFO )
00126          RETURN
00127       END IF
00128 *
00129 *     Compute the LU factorization of the band matrix A.
00130 *
00131       CALL DGBTRF( N, N, KL, KU, AB, LDAB, IPIV, INFO )
00132       IF( INFO.EQ.0 ) THEN
00133 *
00134 *        Solve the system A*X = B, overwriting B with X.
00135 *
00136          CALL DGBTRS( 'No transpose', N, KL, KU, NRHS, AB, LDAB, IPIV,
00137      $                B, LDB, INFO )
00138       END IF
00139       RETURN
00140 *
00141 *     End of DGBSV
00142 *
00143       END
 All Files Functions