SRC/dgbtrs.f

   1 *> \brief \b DGBTRS
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download DGBTRS + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgbtrs.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgbtrs.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgbtrs.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE DGBTRS( TRANS, N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB,
  22 *                          INFO )
  23 *
  24 *       .. Scalar Arguments ..
  25 *       CHARACTER          TRANS
  26 *       INTEGER            INFO, KL, KU, LDAB, LDB, N, NRHS
  27 *       ..
  28 *       .. Array Arguments ..
  29 *       INTEGER            IPIV( * )
  30 *       DOUBLE PRECISION   AB( LDAB, * ), B( LDB, * )
  31 *       ..
  32 *
  33 *
  34 *> \par Purpose:
  35 *  =============
  36 *>
  37 *> \verbatim
  38 *>
  39 *> DGBTRS solves a system of linear equations
  40 *>    A * X = B  or  A**T * X = B
  41 *> with a general band matrix A using the LU factorization computed
  42 *> by DGBTRF.
  43 *> \endverbatim
  44 *
  45 *  Arguments:
  46 *  ==========
  47 *
  48 *> \param[in] TRANS
  49 *> \verbatim
  50 *>          TRANS is CHARACTER*1
  51 *>          Specifies the form of the system of equations.
  52 *>          = 'N':  A * X = B  (No transpose)
  53 *>          = 'T':  A**T* X = B  (Transpose)
  54 *>          = 'C':  A**T* X = B  (Conjugate transpose = Transpose)
  55 *> \endverbatim
  56 *>
  57 *> \param[in] N
  58 *> \verbatim
  59 *>          N is INTEGER
  60 *>          The order of the matrix A.  N >= 0.
  61 *> \endverbatim
  62 *>
  63 *> \param[in] KL
  64 *> \verbatim
  65 *>          KL is INTEGER
  66 *>          The number of subdiagonals within the band of A.  KL >= 0.
  67 *> \endverbatim
  68 *>
  69 *> \param[in] KU
  70 *> \verbatim
  71 *>          KU is INTEGER
  72 *>          The number of superdiagonals within the band of A.  KU >= 0.
  73 *> \endverbatim
  74 *>
  75 *> \param[in] NRHS
  76 *> \verbatim
  77 *>          NRHS is INTEGER
  78 *>          The number of right hand sides, i.e., the number of columns
  79 *>          of the matrix B.  NRHS >= 0.
  80 *> \endverbatim
  81 *>
  82 *> \param[in] AB
  83 *> \verbatim
  84 *>          AB is DOUBLE PRECISION array, dimension (LDAB,N)
  85 *>          Details of the LU factorization of the band matrix A, as
  86 *>          computed by DGBTRF.  U is stored as an upper triangular band
  87 *>          matrix with KL+KU superdiagonals in rows 1 to KL+KU+1, and
  88 *>          the multipliers used during the factorization are stored in
  89 *>          rows KL+KU+2 to 2*KL+KU+1.
  90 *> \endverbatim
  91 *>
  92 *> \param[in] LDAB
  93 *> \verbatim
  94 *>          LDAB is INTEGER
  95 *>          The leading dimension of the array AB.  LDAB >= 2*KL+KU+1.
  96 *> \endverbatim
  97 *>
  98 *> \param[in] IPIV
  99 *> \verbatim
 100 *>          IPIV is INTEGER array, dimension (N)
 101 *>          The pivot indices; for 1 <= i <= N, row i of the matrix was
 102 *>          interchanged with row IPIV(i).
 103 *> \endverbatim
 104 *>
 105 *> \param[in,out] B
 106 *> \verbatim
 107 *>          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
 108 *>          On entry, the right hand side matrix B.
 109 *>          On exit, the solution matrix X.
 110 *> \endverbatim
 111 *>
 112 *> \param[in] LDB
 113 *> \verbatim
 114 *>          LDB is INTEGER
 115 *>          The leading dimension of the array B.  LDB >= max(1,N).
 116 *> \endverbatim
 117 *>
 118 *> \param[out] INFO
 119 *> \verbatim
 120 *>          INFO is INTEGER
 121 *>          = 0:  successful exit
 122 *>          < 0: if INFO = -i, the i-th argument had an illegal value
 123 *> \endverbatim
 124 *
 125 *  Authors:
 126 *  ========
 127 *
 128 *> \author Univ. of Tennessee
 129 *> \author Univ. of California Berkeley
 130 *> \author Univ. of Colorado Denver
 131 *> \author NAG Ltd.
 132 *
 133 *> \date November 2011
 134 *
 135 *> \ingroup doubleGBcomputational
 136 *
 137 *  =====================================================================
 138       SUBROUTINE DGBTRS( TRANS, N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB,
 139      $                   INFO )
 140 *
 141 *  -- LAPACK computational routine (version 3.4.0) --
 142 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 143 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 144 *     November 2011
 145 *
 146 *     .. Scalar Arguments ..
 147       CHARACTER          TRANS
 148       INTEGER            INFO, KL, KU, LDAB, LDB, N, NRHS
 149 *     ..
 150 *     .. Array Arguments ..
 151       INTEGER            IPIV( * )
 152       DOUBLE PRECISION   AB( LDAB, * ), B( LDB, * )
 153 *     ..
 154 *
 155 *  =====================================================================
 156 *
 157 *     .. Parameters ..
 158       DOUBLE PRECISION   ONE
 159       PARAMETER          ( ONE = 1.0D+0 )
 160 *     ..
 161 *     .. Local Scalars ..
 162       LOGICAL            LNOTI, NOTRAN
 163       INTEGER            I, J, KD, L, LM
 164 *     ..
 165 *     .. External Functions ..
 166       LOGICAL            LSAME
 167       EXTERNAL           LSAME
 168 *     ..
 169 *     .. External Subroutines ..
 170       EXTERNAL           DGEMV, DGER, DSWAP, DTBSV, XERBLA
 171 *     ..
 172 *     .. Intrinsic Functions ..
 173       INTRINSIC          MAX, MIN
 174 *     ..
 175 *     .. Executable Statements ..
 176 *
 177 *     Test the input parameters.
 178 *
 179       INFO = 0
 180       NOTRAN = LSAME( TRANS, 'N' )
 181       IF( .NOT.NOTRAN .AND. .NOT.LSAME( TRANS, 'T' ) .AND. .NOT.
 182      $    LSAME( TRANS, 'C' ) ) THEN
 183          INFO = -1
 184       ELSE IF( N.LT.0 ) THEN
 185          INFO = -2
 186       ELSE IF( KL.LT.0 ) THEN
 187          INFO = -3
 188       ELSE IF( KU.LT.0 ) THEN
 189          INFO = -4
 190       ELSE IF( NRHS.LT.0 ) THEN
 191          INFO = -5
 192       ELSE IF( LDAB.LT.( 2*KL+KU+1 ) ) THEN
 193          INFO = -7
 194       ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
 195          INFO = -10
 196       END IF
 197       IF( INFO.NE.0 ) THEN
 198          CALL XERBLA( 'DGBTRS', -INFO )
 199          RETURN
 200       END IF
 201 *
 202 *     Quick return if possible
 203 *
 204       IF( N.EQ.0 .OR. NRHS.EQ.0 )
 205      $   RETURN
 206 *
 207       KD = KU + KL + 1
 208       LNOTI = KL.GT.0
 209 *
 210       IF( NOTRAN ) THEN
 211 *
 212 *        Solve  A*X = B.
 213 *
 214 *        Solve L*X = B, overwriting B with X.
 215 *
 216 *        L is represented as a product of permutations and unit lower
 217 *        triangular matrices L = P(1) * L(1) * ... * P(n-1) * L(n-1),
 218 *        where each transformation L(i) is a rank-one modification of
 219 *        the identity matrix.
 220 *
 221          IF( LNOTI ) THEN
 222             DO 10 J = 1, N - 1
 223                LM = MIN( KL, N-J )
 224                L = IPIV( J )
 225                IF( L.NE.J )
 226      $            CALL DSWAP( NRHS, B( L, 1 ), LDB, B( J, 1 ), LDB )
 227                CALL DGER( LM, NRHS, -ONE, AB( KD+1, J ), 1, B( J, 1 ),
 228      $                    LDB, B( J+1, 1 ), LDB )
 229    10       CONTINUE
 230          END IF
 231 *
 232          DO 20 I = 1, NRHS
 233 *
 234 *           Solve U*X = B, overwriting B with X.
 235 *
 236             CALL DTBSV( 'Upper', 'No transpose', 'Non-unit', N, KL+KU,
 237      $                  AB, LDAB, B( 1, I ), 1 )
 238    20    CONTINUE
 239 *
 240       ELSE
 241 *
 242 *        Solve A**T*X = B.
 243 *
 244          DO 30 I = 1, NRHS
 245 *
 246 *           Solve U**T*X = B, overwriting B with X.
 247 *
 248             CALL DTBSV( 'Upper', 'Transpose', 'Non-unit', N, KL+KU, AB,
 249      $                  LDAB, B( 1, I ), 1 )
 250    30    CONTINUE
 251 *
 252 *        Solve L**T*X = B, overwriting B with X.
 253 *
 254          IF( LNOTI ) THEN
 255             DO 40 J = N - 1, 1, -1
 256                LM = MIN( KL, N-J )
 257                CALL DGEMV( 'Transpose', LM, NRHS, -ONE, B( J+1, 1 ),
 258      $                     LDB, AB( KD+1, J ), 1, ONE, B( J, 1 ), LDB )
 259                L = IPIV( J )
 260                IF( L.NE.J )
 261      $            CALL DSWAP( NRHS, B( L, 1 ), LDB, B( J, 1 ), LDB )
 262    40       CONTINUE
 263          END IF
 264       END IF
 265       RETURN
 266 *
 267 *     End of DGBTRS
 268 *
 269       END