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
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dgtsv.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack