1*> \brief <b> DGTSV computes the solution to system of linear equations A * X = B for GT matrices </b>
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download DGTSV + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgtsv.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgtsv.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgtsv.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18*  Definition:
19*  ===========
20*
21*       SUBROUTINE DGTSV( N, NRHS, DL, D, DU, B, LDB, INFO )
22*
23*       .. Scalar Arguments ..
24*       INTEGER            INFO, LDB, N, NRHS
25*       ..
26*       .. Array Arguments ..
27*       DOUBLE PRECISION   B( LDB, * ), D( * ), DL( * ), DU( * )
28*       ..
29*
30*
31*> \par Purpose:
32*  =============
33*>
34*> \verbatim
35*>
36*> DGTSV  solves the equation
37*>
38*>    A*X = B,
39*>
40*> where A is an n by n tridiagonal matrix, by Gaussian elimination with
41*> partial pivoting.
42*>
43*> Note that the equation  A**T*X = B  may be solved by interchanging the
44*> order of the arguments DU and DL.
45*> \endverbatim
46*
47*  Arguments:
48*  ==========
49*
50*> \param[in] N
51*> \verbatim
52*>          N is INTEGER
53*>          The order of the matrix A.  N >= 0.
54*> \endverbatim
55*>
56*> \param[in] NRHS
57*> \verbatim
58*>          NRHS is INTEGER
59*>          The number of right hand sides, i.e., the number of columns
60*>          of the matrix B.  NRHS >= 0.
61*> \endverbatim
62*>
63*> \param[in,out] DL
64*> \verbatim
65*>          DL is DOUBLE PRECISION array, dimension (N-1)
66*>          On entry, DL must contain the (n-1) sub-diagonal elements of
67*>          A.
68*>
69*>          On exit, DL is overwritten by the (n-2) elements of the
70*>          second super-diagonal of the upper triangular matrix U from
71*>          the LU factorization of A, in DL(1), ..., DL(n-2).
72*> \endverbatim
73*>
74*> \param[in,out] D
75*> \verbatim
76*>          D is DOUBLE PRECISION array, dimension (N)
77*>          On entry, D must contain the diagonal elements of A.
78*>
79*>          On exit, D is overwritten by the n diagonal elements of U.
80*> \endverbatim
81*>
82*> \param[in,out] DU
83*> \verbatim
84*>          DU is DOUBLE PRECISION array, dimension (N-1)
85*>          On entry, DU must contain the (n-1) super-diagonal elements
86*>          of A.
87*>
88*>          On exit, DU is overwritten by the (n-1) elements of the first
89*>          super-diagonal of U.
90*> \endverbatim
91*>
92*> \param[in,out] B
93*> \verbatim
94*>          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
95*>          On entry, the N by NRHS matrix of right hand side matrix B.
96*>          On exit, if INFO = 0, the N by NRHS solution matrix X.
97*> \endverbatim
98*>
99*> \param[in] LDB
100*> \verbatim
101*>          LDB is INTEGER
102*>          The leading dimension of the array B.  LDB >= max(1,N).
103*> \endverbatim
104*>
105*> \param[out] INFO
106*> \verbatim
107*>          INFO is INTEGER
108*>          = 0: successful exit
109*>          < 0: if INFO = -i, the i-th argument had an illegal value
110*>          > 0: if INFO = i, U(i,i) is exactly zero, and the solution
111*>               has not been computed.  The factorization has not been
112*>               completed unless i = N.
113*> \endverbatim
114*
115*  Authors:
116*  ========
117*
118*> \author Univ. of Tennessee
119*> \author Univ. of California Berkeley
120*> \author Univ. of Colorado Denver
121*> \author NAG Ltd.
122*
123*> \ingroup doubleGTsolve
124*
125*  =====================================================================
126      SUBROUTINE DGTSV( N, NRHS, DL, D, DU, B, LDB, INFO )
127*
128*  -- LAPACK driver routine --
129*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
130*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
131*
132*     .. Scalar Arguments ..
133      INTEGER            INFO, LDB, N, NRHS
134*     ..
135*     .. Array Arguments ..
136      DOUBLE PRECISION   B( LDB, * ), D( * ), DL( * ), DU( * )
137*     ..
138*
139*  =====================================================================
140*
141*     .. Parameters ..
142      DOUBLE PRECISION   ZERO
143      PARAMETER          ( ZERO = 0.0D+0 )
144*     ..
145*     .. Local Scalars ..
146      INTEGER            I, J
147      DOUBLE PRECISION   FACT, TEMP
148*     ..
149*     .. Intrinsic Functions ..
150      INTRINSIC          ABS, MAX
151*     ..
152*     .. External Subroutines ..
153      EXTERNAL           XERBLA
154*     ..
155*     .. Executable Statements ..
156*
157      INFO = 0
158      IF( N.LT.0 ) THEN
159         INFO = -1
160      ELSE IF( NRHS.LT.0 ) THEN
161         INFO = -2
162      ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
163         INFO = -7
164      END IF
165      IF( INFO.NE.0 ) THEN
166         CALL XERBLA( 'DGTSV ', -INFO )
167         RETURN
168      END IF
169*
170      IF( N.EQ.0 )
171     $   RETURN
172*
173      IF( NRHS.EQ.1 ) THEN
174         DO 10 I = 1, N - 2
175            IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
176*
177*              No row interchange required
178*
179               IF( D( I ).NE.ZERO ) THEN
180                  FACT = DL( I ) / D( I )
181                  D( I+1 ) = D( I+1 ) - FACT*DU( I )
182                  B( I+1, 1 ) = B( I+1, 1 ) - FACT*B( I, 1 )
183               ELSE
184                  INFO = I
185                  RETURN
186               END IF
187               DL( I ) = ZERO
188            ELSE
189*
190*              Interchange rows I and I+1
191*
192               FACT = D( I ) / DL( I )
193               D( I ) = DL( I )
194               TEMP = D( I+1 )
195               D( I+1 ) = DU( I ) - FACT*TEMP
196               DL( I ) = DU( I+1 )
197               DU( I+1 ) = -FACT*DL( I )
198               DU( I ) = TEMP
199               TEMP = B( I, 1 )
200               B( I, 1 ) = B( I+1, 1 )
201               B( I+1, 1 ) = TEMP - FACT*B( I+1, 1 )
202            END IF
203   10    CONTINUE
204         IF( N.GT.1 ) THEN
205            I = N - 1
206            IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
207               IF( D( I ).NE.ZERO ) THEN
208                  FACT = DL( I ) / D( I )
209                  D( I+1 ) = D( I+1 ) - FACT*DU( I )
210                  B( I+1, 1 ) = B( I+1, 1 ) - FACT*B( I, 1 )
211               ELSE
212                  INFO = I
213                  RETURN
214               END IF
215            ELSE
216               FACT = D( I ) / DL( I )
217               D( I ) = DL( I )
218               TEMP = D( I+1 )
219               D( I+1 ) = DU( I ) - FACT*TEMP
220               DU( I ) = TEMP
221               TEMP = B( I, 1 )
222               B( I, 1 ) = B( I+1, 1 )
223               B( I+1, 1 ) = TEMP - FACT*B( I+1, 1 )
224            END IF
225         END IF
226         IF( D( N ).EQ.ZERO ) THEN
227            INFO = N
228            RETURN
229         END IF
230      ELSE
231         DO 40 I = 1, N - 2
232            IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
233*
234*              No row interchange required
235*
236               IF( D( I ).NE.ZERO ) THEN
237                  FACT = DL( I ) / D( I )
238                  D( I+1 ) = D( I+1 ) - FACT*DU( I )
239                  DO 20 J = 1, NRHS
240                     B( I+1, J ) = B( I+1, J ) - FACT*B( I, J )
241   20             CONTINUE
242               ELSE
243                  INFO = I
244                  RETURN
245               END IF
246               DL( I ) = ZERO
247            ELSE
248*
249*              Interchange rows I and I+1
250*
251               FACT = D( I ) / DL( I )
252               D( I ) = DL( I )
253               TEMP = D( I+1 )
254               D( I+1 ) = DU( I ) - FACT*TEMP
255               DL( I ) = DU( I+1 )
256               DU( I+1 ) = -FACT*DL( I )
257               DU( I ) = TEMP
258               DO 30 J = 1, NRHS
259                  TEMP = B( I, J )
260                  B( I, J ) = B( I+1, J )
261                  B( I+1, J ) = TEMP - FACT*B( I+1, J )
262   30          CONTINUE
263            END IF
264   40    CONTINUE
265         IF( N.GT.1 ) THEN
266            I = N - 1
267            IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
268               IF( D( I ).NE.ZERO ) THEN
269                  FACT = DL( I ) / D( I )
270                  D( I+1 ) = D( I+1 ) - FACT*DU( I )
271                  DO 50 J = 1, NRHS
272                     B( I+1, J ) = B( I+1, J ) - FACT*B( I, J )
273   50             CONTINUE
274               ELSE
275                  INFO = I
276                  RETURN
277               END IF
278            ELSE
279               FACT = D( I ) / DL( I )
280               D( I ) = DL( I )
281               TEMP = D( I+1 )
282               D( I+1 ) = DU( I ) - FACT*TEMP
283               DU( I ) = TEMP
284               DO 60 J = 1, NRHS
285                  TEMP = B( I, J )
286                  B( I, J ) = B( I+1, J )
287                  B( I+1, J ) = TEMP - FACT*B( I+1, J )
288   60          CONTINUE
289            END IF
290         END IF
291         IF( D( N ).EQ.ZERO ) THEN
292            INFO = N
293            RETURN
294         END IF
295      END IF
296*
297*     Back solve with the matrix U from the factorization.
298*
299      IF( NRHS.LE.2 ) THEN
300         J = 1
301   70    CONTINUE
302         B( N, J ) = B( N, J ) / D( N )
303         IF( N.GT.1 )
304     $      B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) / D( N-1 )
305         DO 80 I = N - 2, 1, -1
306            B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DL( I )*
307     $                  B( I+2, J ) ) / D( I )
308   80    CONTINUE
309         IF( J.LT.NRHS ) THEN
310            J = J + 1
311            GO TO 70
312         END IF
313      ELSE
314         DO 100 J = 1, NRHS
315            B( N, J ) = B( N, J ) / D( N )
316            IF( N.GT.1 )
317     $         B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) /
318     $                       D( N-1 )
319            DO 90 I = N - 2, 1, -1
320               B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DL( I )*
321     $                     B( I+2, J ) ) / D( I )
322   90       CONTINUE
323  100    CONTINUE
324      END IF
325*
326      RETURN
327*
328*     End of DGTSV
329*
330      END
331