LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
dgtsv.f
Go to the documentation of this file.
00001 *> \brief \b DGTSV
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *> \htmlonly
00009 *> Download DGTSV + dependencies 
00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgtsv.f"> 
00011 *> [TGZ]</a> 
00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgtsv.f"> 
00013 *> [ZIP]</a> 
00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgtsv.f"> 
00015 *> [TXT]</a>
00016 *> \endhtmlonly 
00017 *
00018 *  Definition:
00019 *  ===========
00020 *
00021 *       SUBROUTINE DGTSV( N, NRHS, DL, D, DU, B, LDB, INFO )
00022 * 
00023 *       .. Scalar Arguments ..
00024 *       INTEGER            INFO, LDB, N, NRHS
00025 *       ..
00026 *       .. Array Arguments ..
00027 *       DOUBLE PRECISION   B( LDB, * ), D( * ), DL( * ), DU( * )
00028 *       ..
00029 *  
00030 *
00031 *> \par Purpose:
00032 *  =============
00033 *>
00034 *> \verbatim
00035 *>
00036 *> DGTSV  solves the equation
00037 *>
00038 *>    A*X = B,
00039 *>
00040 *> where A is an n by n tridiagonal matrix, by Gaussian elimination with
00041 *> partial pivoting.
00042 *>
00043 *> Note that the equation  A**T*X = B  may be solved by interchanging the
00044 *> order of the arguments DU and DL.
00045 *> \endverbatim
00046 *
00047 *  Arguments:
00048 *  ==========
00049 *
00050 *> \param[in] N
00051 *> \verbatim
00052 *>          N is INTEGER
00053 *>          The order of the matrix A.  N >= 0.
00054 *> \endverbatim
00055 *>
00056 *> \param[in] NRHS
00057 *> \verbatim
00058 *>          NRHS is INTEGER
00059 *>          The number of right hand sides, i.e., the number of columns
00060 *>          of the matrix B.  NRHS >= 0.
00061 *> \endverbatim
00062 *>
00063 *> \param[in,out] DL
00064 *> \verbatim
00065 *>          DL is DOUBLE PRECISION array, dimension (N-1)
00066 *>          On entry, DL must contain the (n-1) sub-diagonal elements of
00067 *>          A.
00068 *>
00069 *>          On exit, DL is overwritten by the (n-2) elements of the
00070 *>          second super-diagonal of the upper triangular matrix U from
00071 *>          the LU factorization of A, in DL(1), ..., DL(n-2).
00072 *> \endverbatim
00073 *>
00074 *> \param[in,out] D
00075 *> \verbatim
00076 *>          D is DOUBLE PRECISION array, dimension (N)
00077 *>          On entry, D must contain the diagonal elements of A.
00078 *>
00079 *>          On exit, D is overwritten by the n diagonal elements of U.
00080 *> \endverbatim
00081 *>
00082 *> \param[in,out] DU
00083 *> \verbatim
00084 *>          DU is DOUBLE PRECISION array, dimension (N-1)
00085 *>          On entry, DU must contain the (n-1) super-diagonal elements
00086 *>          of A.
00087 *>
00088 *>          On exit, DU is overwritten by the (n-1) elements of the first
00089 *>          super-diagonal of U.
00090 *> \endverbatim
00091 *>
00092 *> \param[in,out] B
00093 *> \verbatim
00094 *>          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
00095 *>          On entry, the N by NRHS matrix of right hand side matrix B.
00096 *>          On exit, if INFO = 0, the N by NRHS solution matrix X.
00097 *> \endverbatim
00098 *>
00099 *> \param[in] LDB
00100 *> \verbatim
00101 *>          LDB is INTEGER
00102 *>          The leading dimension of the array B.  LDB >= max(1,N).
00103 *> \endverbatim
00104 *>
00105 *> \param[out] INFO
00106 *> \verbatim
00107 *>          INFO is INTEGER
00108 *>          = 0: successful exit
00109 *>          < 0: if INFO = -i, the i-th argument had an illegal value
00110 *>          > 0: if INFO = i, U(i,i) is exactly zero, and the solution
00111 *>               has not been computed.  The factorization has not been
00112 *>               completed unless i = N.
00113 *> \endverbatim
00114 *
00115 *  Authors:
00116 *  ========
00117 *
00118 *> \author Univ. of Tennessee 
00119 *> \author Univ. of California Berkeley 
00120 *> \author Univ. of Colorado Denver 
00121 *> \author NAG Ltd. 
00122 *
00123 *> \date November 2011
00124 *
00125 *> \ingroup doubleOTHERcomputational
00126 *
00127 *  =====================================================================
00128       SUBROUTINE DGTSV( N, NRHS, DL, D, DU, B, LDB, INFO )
00129 *
00130 *  -- LAPACK computational routine (version 3.4.0) --
00131 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00132 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00133 *     November 2011
00134 *
00135 *     .. Scalar Arguments ..
00136       INTEGER            INFO, LDB, N, NRHS
00137 *     ..
00138 *     .. Array Arguments ..
00139       DOUBLE PRECISION   B( LDB, * ), D( * ), DL( * ), DU( * )
00140 *     ..
00141 *
00142 *  =====================================================================
00143 *
00144 *     .. Parameters ..
00145       DOUBLE PRECISION   ZERO
00146       PARAMETER          ( ZERO = 0.0D+0 )
00147 *     ..
00148 *     .. Local Scalars ..
00149       INTEGER            I, J
00150       DOUBLE PRECISION   FACT, TEMP
00151 *     ..
00152 *     .. Intrinsic Functions ..
00153       INTRINSIC          ABS, MAX
00154 *     ..
00155 *     .. External Subroutines ..
00156       EXTERNAL           XERBLA
00157 *     ..
00158 *     .. Executable Statements ..
00159 *
00160       INFO = 0
00161       IF( N.LT.0 ) THEN
00162          INFO = -1
00163       ELSE IF( NRHS.LT.0 ) THEN
00164          INFO = -2
00165       ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
00166          INFO = -7
00167       END IF
00168       IF( INFO.NE.0 ) THEN
00169          CALL XERBLA( 'DGTSV ', -INFO )
00170          RETURN
00171       END IF
00172 *
00173       IF( N.EQ.0 )
00174      $   RETURN
00175 *
00176       IF( NRHS.EQ.1 ) THEN
00177          DO 10 I = 1, N - 2
00178             IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
00179 *
00180 *              No row interchange required
00181 *
00182                IF( D( I ).NE.ZERO ) THEN
00183                   FACT = DL( I ) / D( I )
00184                   D( I+1 ) = D( I+1 ) - FACT*DU( I )
00185                   B( I+1, 1 ) = B( I+1, 1 ) - FACT*B( I, 1 )
00186                ELSE
00187                   INFO = I
00188                   RETURN
00189                END IF
00190                DL( I ) = ZERO
00191             ELSE
00192 *
00193 *              Interchange rows I and I+1
00194 *
00195                FACT = D( I ) / DL( I )
00196                D( I ) = DL( I )
00197                TEMP = D( I+1 )
00198                D( I+1 ) = DU( I ) - FACT*TEMP
00199                DL( I ) = DU( I+1 )
00200                DU( I+1 ) = -FACT*DL( I )
00201                DU( I ) = TEMP
00202                TEMP = B( I, 1 )
00203                B( I, 1 ) = B( I+1, 1 )
00204                B( I+1, 1 ) = TEMP - FACT*B( I+1, 1 )
00205             END IF
00206    10    CONTINUE
00207          IF( N.GT.1 ) THEN
00208             I = N - 1
00209             IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
00210                IF( D( I ).NE.ZERO ) THEN
00211                   FACT = DL( I ) / D( I )
00212                   D( I+1 ) = D( I+1 ) - FACT*DU( I )
00213                   B( I+1, 1 ) = B( I+1, 1 ) - FACT*B( I, 1 )
00214                ELSE
00215                   INFO = I
00216                   RETURN
00217                END IF
00218             ELSE
00219                FACT = D( I ) / DL( I )
00220                D( I ) = DL( I )
00221                TEMP = D( I+1 )
00222                D( I+1 ) = DU( I ) - FACT*TEMP
00223                DU( I ) = TEMP
00224                TEMP = B( I, 1 )
00225                B( I, 1 ) = B( I+1, 1 )
00226                B( I+1, 1 ) = TEMP - FACT*B( I+1, 1 )
00227             END IF
00228          END IF
00229          IF( D( N ).EQ.ZERO ) THEN
00230             INFO = N
00231             RETURN
00232          END IF
00233       ELSE
00234          DO 40 I = 1, N - 2
00235             IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
00236 *
00237 *              No row interchange required
00238 *
00239                IF( D( I ).NE.ZERO ) THEN
00240                   FACT = DL( I ) / D( I )
00241                   D( I+1 ) = D( I+1 ) - FACT*DU( I )
00242                   DO 20 J = 1, NRHS
00243                      B( I+1, J ) = B( I+1, J ) - FACT*B( I, J )
00244    20             CONTINUE
00245                ELSE
00246                   INFO = I
00247                   RETURN
00248                END IF
00249                DL( I ) = ZERO
00250             ELSE
00251 *
00252 *              Interchange rows I and I+1
00253 *
00254                FACT = D( I ) / DL( I )
00255                D( I ) = DL( I )
00256                TEMP = D( I+1 )
00257                D( I+1 ) = DU( I ) - FACT*TEMP
00258                DL( I ) = DU( I+1 )
00259                DU( I+1 ) = -FACT*DL( I )
00260                DU( I ) = TEMP
00261                DO 30 J = 1, NRHS
00262                   TEMP = B( I, J )
00263                   B( I, J ) = B( I+1, J )
00264                   B( I+1, J ) = TEMP - FACT*B( I+1, J )
00265    30          CONTINUE
00266             END IF
00267    40    CONTINUE
00268          IF( N.GT.1 ) THEN
00269             I = N - 1
00270             IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
00271                IF( D( I ).NE.ZERO ) THEN
00272                   FACT = DL( I ) / D( I )
00273                   D( I+1 ) = D( I+1 ) - FACT*DU( I )
00274                   DO 50 J = 1, NRHS
00275                      B( I+1, J ) = B( I+1, J ) - FACT*B( I, J )
00276    50             CONTINUE
00277                ELSE
00278                   INFO = I
00279                   RETURN
00280                END IF
00281             ELSE
00282                FACT = D( I ) / DL( I )
00283                D( I ) = DL( I )
00284                TEMP = D( I+1 )
00285                D( I+1 ) = DU( I ) - FACT*TEMP
00286                DU( I ) = TEMP
00287                DO 60 J = 1, NRHS
00288                   TEMP = B( I, J )
00289                   B( I, J ) = B( I+1, J )
00290                   B( I+1, J ) = TEMP - FACT*B( I+1, J )
00291    60          CONTINUE
00292             END IF
00293          END IF
00294          IF( D( N ).EQ.ZERO ) THEN
00295             INFO = N
00296             RETURN
00297          END IF
00298       END IF
00299 *
00300 *     Back solve with the matrix U from the factorization.
00301 *
00302       IF( NRHS.LE.2 ) THEN
00303          J = 1
00304    70    CONTINUE
00305          B( N, J ) = B( N, J ) / D( N )
00306          IF( N.GT.1 )
00307      $      B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) / D( N-1 )
00308          DO 80 I = N - 2, 1, -1
00309             B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DL( I )*
00310      $                  B( I+2, J ) ) / D( I )
00311    80    CONTINUE
00312          IF( J.LT.NRHS ) THEN
00313             J = J + 1
00314             GO TO 70
00315          END IF
00316       ELSE
00317          DO 100 J = 1, NRHS
00318             B( N, J ) = B( N, J ) / D( N )
00319             IF( N.GT.1 )
00320      $         B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) /
00321      $                       D( N-1 )
00322             DO 90 I = N - 2, 1, -1
00323                B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DL( I )*
00324      $                     B( I+2, J ) ) / D( I )
00325    90       CONTINUE
00326   100    CONTINUE
00327       END IF
00328 *
00329       RETURN
00330 *
00331 *     End of DGTSV
00332 *
00333       END
 All Files Functions