dgtts2.f

Go to the documentation of this file.

*> \brief \b DGTTS2 solves a system of linear equations with a tridiagonal matrix using the LU factorization computed by sgttrf.
*
*  =========== DOCUMENTATION ===========
*
* Online html documentation available at 
*            http://www.netlib.org/lapack/explore-html/ 
*
*> \htmlonly
*> Download DGTTS2 + dependencies 
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgtts2.f"> 
*> [TGZ]</a> 
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgtts2.f"> 
*> [ZIP]</a> 
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgtts2.f"> 
*> [TXT]</a>
*> \endhtmlonly 
*
*  Definition:
*  ===========
*
*       SUBROUTINE DGTTS2( ITRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB )
* 
*       .. Scalar Arguments ..
*       INTEGER            ITRANS, LDB, N, NRHS
*       ..
*       .. Array Arguments ..
*       INTEGER            IPIV( * )
*       DOUBLE PRECISION   B( LDB, * ), D( * ), DL( * ), DU( * ), DU2( * )
*       ..
*  
*
*> \par Purpose:
*  =============
*>
*> \verbatim
*>
*> DGTTS2 solves one of the systems of equations
*>    A*X = B  or  A**T*X = B,
*> with a tridiagonal matrix A using the LU factorization computed
*> by DGTTRF.
*> \endverbatim
*
*  Arguments:
*  ==========
*
*> \param[in] ITRANS
*> \verbatim
*>          ITRANS is INTEGER
*>          Specifies the form of the system of equations.
*>          = 0:  A * X = B  (No transpose)
*>          = 1:  A**T* X = B  (Transpose)
*>          = 2:  A**T* X = B  (Conjugate transpose = Transpose)
*> \endverbatim
*>
*> \param[in] N
*> \verbatim
*>          N is INTEGER
*>          The order of the matrix A.
*> \endverbatim
*>
*> \param[in] NRHS
*> \verbatim
*>          NRHS is INTEGER
*>          The number of right hand sides, i.e., the number of columns
*>          of the matrix B.  NRHS >= 0.
*> \endverbatim
*>
*> \param[in] DL
*> \verbatim
*>          DL is DOUBLE PRECISION array, dimension (N-1)
*>          The (n-1) multipliers that define the matrix L from the
*>          LU factorization of A.
*> \endverbatim
*>
*> \param[in] D
*> \verbatim
*>          D is DOUBLE PRECISION array, dimension (N)
*>          The n diagonal elements of the upper triangular matrix U from
*>          the LU factorization of A.
*> \endverbatim
*>
*> \param[in] DU
*> \verbatim
*>          DU is DOUBLE PRECISION array, dimension (N-1)
*>          The (n-1) elements of the first super-diagonal of U.
*> \endverbatim
*>
*> \param[in] DU2
*> \verbatim
*>          DU2 is DOUBLE PRECISION array, dimension (N-2)
*>          The (n-2) elements of the second super-diagonal of U.
*> \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).  IPIV(i) will always be either
*>          i or i+1; IPIV(i) = i indicates a row interchange was not
*>          required.
*> \endverbatim
*>
*> \param[in,out] B
*> \verbatim
*>          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
*>          On entry, the matrix of right hand side vectors B.
*>          On exit, B is overwritten by the solution vectors X.
*> \endverbatim
*>
*> \param[in] LDB
*> \verbatim
*>          LDB is INTEGER
*>          The leading dimension of the array B.  LDB >= max(1,N).
*> \endverbatim
*
*  Authors:
*  ========
*
*> \author Univ. of Tennessee 
*> \author Univ. of California Berkeley 
*> \author Univ. of Colorado Denver 
*> \author NAG Ltd. 
*
*> \date September 2012
*
*> \ingroup doubleGTcomputational
*
*  =====================================================================
      SUBROUTINE dgtts2( ITRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB )
*
*  -- LAPACK computational routine (version 3.4.2) --
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*     September 2012
*
*     .. Scalar Arguments ..
      INTEGER            itrans, ldb, n, nrhs
*     ..
*     .. Array Arguments ..
      INTEGER            ipiv( * )
      DOUBLE PRECISION   b( ldb, * ), d( * ), dl( * ), du( * ), du2( * )
*     ..
*
*  =====================================================================
*
*     .. Local Scalars ..
      INTEGER            i, ip, j
      DOUBLE PRECISION   temp
*     ..
*     .. Executable Statements ..
*
*     Quick return if possible
*
      IF( n.EQ.0 .OR. nrhs.EQ.0 )
     $   return
*
      IF( itrans.EQ.0 ) THEN
*
*        Solve A*X = B using the LU factorization of A,
*        overwriting each right hand side vector with its solution.
*
         IF( nrhs.LE.1 ) THEN
            j = 1
   10       continue
*
*           Solve L*x = b.
*
            DO 20 i = 1, n - 1
               ip = ipiv( i )
               temp = b( i+1-ip+i, j ) - dl( i )*b( ip, j )
               b( i, j ) = b( ip, j )
               b( i+1, j ) = temp
   20       continue
*
*           Solve U*x = b.
*
            b( n, j ) = b( n, j ) / d( n )
            IF( n.GT.1 )
     $         b( n-1, j ) = ( b( n-1, j )-du( n-1 )*b( n, j ) ) /
     $                       d( n-1 )
            DO 30 i = n - 2, 1, -1
               b( i, j ) = ( b( i, j )-du( i )*b( i+1, j )-du2( i )*
     $                     b( i+2, j ) ) / d( i )
   30       continue
            IF( j.LT.nrhs ) THEN
               j = j + 1
               go to 10
            END IF
         ELSE
            DO 60 j = 1, nrhs
*
*              Solve L*x = b.
*
               DO 40 i = 1, n - 1
                  IF( ipiv( i ).EQ.i ) THEN
                     b( i+1, j ) = b( i+1, j ) - dl( i )*b( i, j )
                  ELSE
                     temp = b( i, j )
                     b( i, j ) = b( i+1, j )
                     b( i+1, j ) = temp - dl( i )*b( i, j )
                  END IF
   40          continue
*
*              Solve U*x = b.
*
               b( n, j ) = b( n, j ) / d( n )
               IF( n.GT.1 )
     $            b( n-1, j ) = ( b( n-1, j )-du( n-1 )*b( n, j ) ) /
     $                          d( n-1 )
               DO 50 i = n - 2, 1, -1
                  b( i, j ) = ( b( i, j )-du( i )*b( i+1, j )-du2( i )*
     $                        b( i+2, j ) ) / d( i )
   50          continue
   60       continue
         END IF
      ELSE
*
*        Solve A**T * X = B.
*
         IF( nrhs.LE.1 ) THEN
*
*           Solve U**T*x = b.
*
            j = 1
   70       continue
            b( 1, j ) = b( 1, j ) / d( 1 )
            IF( n.GT.1 )
     $         b( 2, j ) = ( b( 2, j )-du( 1 )*b( 1, j ) ) / d( 2 )
            DO 80 i = 3, n
               b( i, j ) = ( b( i, j )-du( i-1 )*b( i-1, j )-du2( i-2 )*
     $                     b( i-2, j ) ) / d( i )
   80       continue
*
*           Solve L**T*x = b.
*
            DO 90 i = n - 1, 1, -1
               ip = ipiv( i )
               temp = b( i, j ) - dl( i )*b( i+1, j )
               b( i, j ) = b( ip, j )
               b( ip, j ) = temp
   90       continue
            IF( j.LT.nrhs ) THEN
               j = j + 1
               go to 70
            END IF
*
         ELSE
            DO 120 j = 1, nrhs
*
*              Solve U**T*x = b.
*
               b( 1, j ) = b( 1, j ) / d( 1 )
               IF( n.GT.1 )
     $            b( 2, j ) = ( b( 2, j )-du( 1 )*b( 1, j ) ) / d( 2 )
               DO 100 i = 3, n
                  b( i, j ) = ( b( i, j )-du( i-1 )*b( i-1, j )-
     $                        du2( i-2 )*b( i-2, j ) ) / d( i )
  100          continue
               DO 110 i = n - 1, 1, -1
                  IF( ipiv( i ).EQ.i ) THEN
                     b( i, j ) = b( i, j ) - dl( i )*b( i+1, j )
                  ELSE
                     temp = b( i+1, j )
                     b( i+1, j ) = b( i, j ) - dl( i )*temp
                     b( i, j ) = temp
                  END IF
  110          continue
  120       continue
         END IF
      END IF
*
*     End of DGTTS2
*
      END