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*> \date September 2012
124*
125*> \ingroup doubleGTsolve
126*
127*  =====================================================================
128      SUBROUTINE DGTSV( N, NRHS, DL, D, DU, B, LDB, INFO )
129*
130*  -- LAPACK driver routine (version 3.4.2) --
131*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
132*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
133*     September 2012
134*
135*     .. Scalar Arguments ..
136      INTEGER            INFO, LDB, N, NRHS
137*     ..
138*     .. Array Arguments ..
139      DOUBLE PRECISION   B( LDB, * ), D( * ), DL( * ), DU( * )
140*     ..
141*
142*  =====================================================================
143*
144*     .. Parameters ..
145      DOUBLE PRECISION   ZERO
146      PARAMETER          ( ZERO = 0.0D+0 )
147*     ..
148*     .. Local Scalars ..
149      INTEGER            I, J
150      DOUBLE PRECISION   FACT, TEMP
151*     ..
152*     .. Intrinsic Functions ..
153      INTRINSIC          ABS, MAX
154*     ..
155*     .. External Subroutines ..
156      EXTERNAL           XERBLA
157*     ..
158*     .. Executable Statements ..
159*
160      INFO = 0
161      IF( N.LT.0 ) THEN
162         INFO = -1
163      ELSE IF( NRHS.LT.0 ) THEN
164         INFO = -2
165      ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
166         INFO = -7
167      END IF
168      IF( INFO.NE.0 ) THEN
169         CALL XERBLA( 'DGTSV ', -INFO )
170         RETURN
171      END IF
172*
173      IF( N.EQ.0 )
174     $   RETURN
175*
176      IF( NRHS.EQ.1 ) THEN
177         DO 10 I = 1, N - 2
178            IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
179*
180*              No row interchange required
181*
182               IF( D( I ).NE.ZERO ) THEN
183                  FACT = DL( I ) / D( I )
184                  D( I+1 ) = D( I+1 ) - FACT*DU( I )
185                  B( I+1, 1 ) = B( I+1, 1 ) - FACT*B( I, 1 )
186               ELSE
187                  INFO = I
188                  RETURN
189               END IF
190               DL( I ) = ZERO
191            ELSE
192*
193*              Interchange rows I and I+1
194*
195               FACT = D( I ) / DL( I )
196               D( I ) = DL( I )
197               TEMP = D( I+1 )
198               D( I+1 ) = DU( I ) - FACT*TEMP
199               DL( I ) = DU( I+1 )
200               DU( I+1 ) = -FACT*DL( I )
201               DU( I ) = TEMP
202               TEMP = B( I, 1 )
203               B( I, 1 ) = B( I+1, 1 )
204               B( I+1, 1 ) = TEMP - FACT*B( I+1, 1 )
205            END IF
206   10    CONTINUE
207         IF( N.GT.1 ) THEN
208            I = N - 1
209            IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
210               IF( D( I ).NE.ZERO ) THEN
211                  FACT = DL( I ) / D( I )
212                  D( I+1 ) = D( I+1 ) - FACT*DU( I )
213                  B( I+1, 1 ) = B( I+1, 1 ) - FACT*B( I, 1 )
214               ELSE
215                  INFO = I
216                  RETURN
217               END IF
218            ELSE
219               FACT = D( I ) / DL( I )
220               D( I ) = DL( I )
221               TEMP = D( I+1 )
222               D( I+1 ) = DU( I ) - FACT*TEMP
223               DU( I ) = TEMP
224               TEMP = B( I, 1 )
225               B( I, 1 ) = B( I+1, 1 )
226               B( I+1, 1 ) = TEMP - FACT*B( I+1, 1 )
227            END IF
228         END IF
229         IF( D( N ).EQ.ZERO ) THEN
230            INFO = N
231            RETURN
232         END IF
233      ELSE
234         DO 40 I = 1, N - 2
235            IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
236*
237*              No row interchange required
238*
239               IF( D( I ).NE.ZERO ) THEN
240                  FACT = DL( I ) / D( I )
241                  D( I+1 ) = D( I+1 ) - FACT*DU( I )
242                  DO 20 J = 1, NRHS
243                     B( I+1, J ) = B( I+1, J ) - FACT*B( I, J )
244   20             CONTINUE
245               ELSE
246                  INFO = I
247                  RETURN
248               END IF
249               DL( I ) = ZERO
250            ELSE
251*
252*              Interchange rows I and I+1
253*
254               FACT = D( I ) / DL( I )
255               D( I ) = DL( I )
256               TEMP = D( I+1 )
257               D( I+1 ) = DU( I ) - FACT*TEMP
258               DL( I ) = DU( I+1 )
259               DU( I+1 ) = -FACT*DL( I )
260               DU( I ) = TEMP
261               DO 30 J = 1, NRHS
262                  TEMP = B( I, J )
263                  B( I, J ) = B( I+1, J )
264                  B( I+1, J ) = TEMP - FACT*B( I+1, J )
265   30          CONTINUE
266            END IF
267   40    CONTINUE
268         IF( N.GT.1 ) THEN
269            I = N - 1
270            IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
271               IF( D( I ).NE.ZERO ) THEN
272                  FACT = DL( I ) / D( I )
273                  D( I+1 ) = D( I+1 ) - FACT*DU( I )
274                  DO 50 J = 1, NRHS
275                     B( I+1, J ) = B( I+1, J ) - FACT*B( I, J )
276   50             CONTINUE
277               ELSE
278                  INFO = I
279                  RETURN
280               END IF
281            ELSE
282               FACT = D( I ) / DL( I )
283               D( I ) = DL( I )
284               TEMP = D( I+1 )
285               D( I+1 ) = DU( I ) - FACT*TEMP
286               DU( I ) = TEMP
287               DO 60 J = 1, NRHS
288                  TEMP = B( I, J )
289                  B( I, J ) = B( I+1, J )
290                  B( I+1, J ) = TEMP - FACT*B( I+1, J )
291   60          CONTINUE
292            END IF
293         END IF
294         IF( D( N ).EQ.ZERO ) THEN
295            INFO = N
296            RETURN
297         END IF
298      END IF
299*
300*     Back solve with the matrix U from the factorization.
301*
302      IF( NRHS.LE.2 ) THEN
303         J = 1
304   70    CONTINUE
305         B( N, J ) = B( N, J ) / D( N )
306         IF( N.GT.1 )
307     $      B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) / D( N-1 )
308         DO 80 I = N - 2, 1, -1
309            B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DL( I )*
310     $                  B( I+2, J ) ) / D( I )
311   80    CONTINUE
312         IF( J.LT.NRHS ) THEN
313            J = J + 1
314            GO TO 70
315         END IF
316      ELSE
317         DO 100 J = 1, NRHS
318            B( N, J ) = B( N, J ) / D( N )
319            IF( N.GT.1 )
320     $         B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) /
321     $                       D( N-1 )
322            DO 90 I = N - 2, 1, -1
323               B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DL( I )*
324     $                     B( I+2, J ) ) / D( I )
325   90       CONTINUE
326  100    CONTINUE
327      END IF
328*
329      RETURN
330*
331*     End of DGTSV
332*
333      END
334