![[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|
Return to dgtsv.f CVS log | Up to [local] / rpl / lapack / lapack |
1.1 bertrand 1: SUBROUTINE DGTSV( N, NRHS, DL, D, DU, B, LDB, INFO )
2: *
3: * -- LAPACK routine (version 3.2) --
4: * -- LAPACK is a software package provided by Univ. of Tennessee, --
5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
6: * November 2006
7: *
8: * .. Scalar Arguments ..
9: INTEGER INFO, LDB, N, NRHS
10: * ..
11: * .. Array Arguments ..
12: DOUBLE PRECISION B( LDB, * ), D( * ), DL( * ), DU( * )
13: * ..
14: *
15: * Purpose
16: * =======
17: *
18: * DGTSV solves the equation
19: *
20: * A*X = B,
21: *
22: * where A is an n by n tridiagonal matrix, by Gaussian elimination with
23: * partial pivoting.
24: *
25: * Note that the equation A'*X = B may be solved by interchanging the
26: * order of the arguments DU and DL.
27: *
28: * Arguments
29: * =========
30: *
31: * N (input) INTEGER
32: * The order of the matrix A. N >= 0.
33: *
34: * NRHS (input) INTEGER
35: * The number of right hand sides, i.e., the number of columns
36: * of the matrix B. NRHS >= 0.
37: *
38: * DL (input/output) DOUBLE PRECISION array, dimension (N-1)
39: * On entry, DL must contain the (n-1) sub-diagonal elements of
40: * A.
41: *
42: * On exit, DL is overwritten by the (n-2) elements of the
43: * second super-diagonal of the upper triangular matrix U from
44: * the LU factorization of A, in DL(1), ..., DL(n-2).
45: *
46: * D (input/output) DOUBLE PRECISION array, dimension (N)
47: * On entry, D must contain the diagonal elements of A.
48: *
49: * On exit, D is overwritten by the n diagonal elements of U.
50: *
51: * DU (input/output) DOUBLE PRECISION array, dimension (N-1)
52: * On entry, DU must contain the (n-1) super-diagonal elements
53: * of A.
54: *
55: * On exit, DU is overwritten by the (n-1) elements of the first
56: * super-diagonal of U.
57: *
58: * B (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS)
59: * On entry, the N by NRHS matrix of right hand side matrix B.
60: * On exit, if INFO = 0, the N by NRHS solution matrix X.
61: *
62: * LDB (input) INTEGER
63: * The leading dimension of the array B. LDB >= max(1,N).
64: *
65: * INFO (output) INTEGER
66: * = 0: successful exit
67: * < 0: if INFO = -i, the i-th argument had an illegal value
68: * > 0: if INFO = i, U(i,i) is exactly zero, and the solution
69: * has not been computed. The factorization has not been
70: * completed unless i = N.
71: *
72: * =====================================================================
73: *
74: * .. Parameters ..
75: DOUBLE PRECISION ZERO
76: PARAMETER ( ZERO = 0.0D+0 )
77: * ..
78: * .. Local Scalars ..
79: INTEGER I, J
80: DOUBLE PRECISION FACT, TEMP
81: * ..
82: * .. Intrinsic Functions ..
83: INTRINSIC ABS, MAX
84: * ..
85: * .. External Subroutines ..
86: EXTERNAL XERBLA
87: * ..
88: * .. Executable Statements ..
89: *
90: INFO = 0
91: IF( N.LT.0 ) THEN
92: INFO = -1
93: ELSE IF( NRHS.LT.0 ) THEN
94: INFO = -2
95: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
96: INFO = -7
97: END IF
98: IF( INFO.NE.0 ) THEN
99: CALL XERBLA( 'DGTSV ', -INFO )
100: RETURN
101: END IF
102: *
103: IF( N.EQ.0 )
104: $ RETURN
105: *
106: IF( NRHS.EQ.1 ) THEN
107: DO 10 I = 1, N - 2
108: IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
109: *
110: * No row interchange required
111: *
112: IF( D( I ).NE.ZERO ) THEN
113: FACT = DL( I ) / D( I )
114: D( I+1 ) = D( I+1 ) - FACT*DU( I )
115: B( I+1, 1 ) = B( I+1, 1 ) - FACT*B( I, 1 )
116: ELSE
117: INFO = I
118: RETURN
119: END IF
120: DL( I ) = ZERO
121: ELSE
122: *
123: * Interchange rows I and I+1
124: *
125: FACT = D( I ) / DL( I )
126: D( I ) = DL( I )
127: TEMP = D( I+1 )
128: D( I+1 ) = DU( I ) - FACT*TEMP
129: DL( I ) = DU( I+1 )
130: DU( I+1 ) = -FACT*DL( I )
131: DU( I ) = TEMP
132: TEMP = B( I, 1 )
133: B( I, 1 ) = B( I+1, 1 )
134: B( I+1, 1 ) = TEMP - FACT*B( I+1, 1 )
135: END IF
136: 10 CONTINUE
137: IF( N.GT.1 ) THEN
138: I = N - 1
139: IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
140: IF( D( I ).NE.ZERO ) THEN
141: FACT = DL( I ) / D( I )
142: D( I+1 ) = D( I+1 ) - FACT*DU( I )
143: B( I+1, 1 ) = B( I+1, 1 ) - FACT*B( I, 1 )
144: ELSE
145: INFO = I
146: RETURN
147: END IF
148: ELSE
149: FACT = D( I ) / DL( I )
150: D( I ) = DL( I )
151: TEMP = D( I+1 )
152: D( I+1 ) = DU( I ) - FACT*TEMP
153: DU( I ) = TEMP
154: TEMP = B( I, 1 )
155: B( I, 1 ) = B( I+1, 1 )
156: B( I+1, 1 ) = TEMP - FACT*B( I+1, 1 )
157: END IF
158: END IF
159: IF( D( N ).EQ.ZERO ) THEN
160: INFO = N
161: RETURN
162: END IF
163: ELSE
164: DO 40 I = 1, N - 2
165: IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
166: *
167: * No row interchange required
168: *
169: IF( D( I ).NE.ZERO ) THEN
170: FACT = DL( I ) / D( I )
171: D( I+1 ) = D( I+1 ) - FACT*DU( I )
172: DO 20 J = 1, NRHS
173: B( I+1, J ) = B( I+1, J ) - FACT*B( I, J )
174: 20 CONTINUE
175: ELSE
176: INFO = I
177: RETURN
178: END IF
179: DL( I ) = ZERO
180: ELSE
181: *
182: * Interchange rows I and I+1
183: *
184: FACT = D( I ) / DL( I )
185: D( I ) = DL( I )
186: TEMP = D( I+1 )
187: D( I+1 ) = DU( I ) - FACT*TEMP
188: DL( I ) = DU( I+1 )
189: DU( I+1 ) = -FACT*DL( I )
190: DU( I ) = TEMP
191: DO 30 J = 1, NRHS
192: TEMP = B( I, J )
193: B( I, J ) = B( I+1, J )
194: B( I+1, J ) = TEMP - FACT*B( I+1, J )
195: 30 CONTINUE
196: END IF
197: 40 CONTINUE
198: IF( N.GT.1 ) THEN
199: I = N - 1
200: IF( ABS( D( I ) ).GE.ABS( DL( I ) ) ) THEN
201: IF( D( I ).NE.ZERO ) THEN
202: FACT = DL( I ) / D( I )
203: D( I+1 ) = D( I+1 ) - FACT*DU( I )
204: DO 50 J = 1, NRHS
205: B( I+1, J ) = B( I+1, J ) - FACT*B( I, J )
206: 50 CONTINUE
207: ELSE
208: INFO = I
209: RETURN
210: END IF
211: ELSE
212: FACT = D( I ) / DL( I )
213: D( I ) = DL( I )
214: TEMP = D( I+1 )
215: D( I+1 ) = DU( I ) - FACT*TEMP
216: DU( I ) = TEMP
217: DO 60 J = 1, NRHS
218: TEMP = B( I, J )
219: B( I, J ) = B( I+1, J )
220: B( I+1, J ) = TEMP - FACT*B( I+1, J )
221: 60 CONTINUE
222: END IF
223: END IF
224: IF( D( N ).EQ.ZERO ) THEN
225: INFO = N
226: RETURN
227: END IF
228: END IF
229: *
230: * Back solve with the matrix U from the factorization.
231: *
232: IF( NRHS.LE.2 ) THEN
233: J = 1
234: 70 CONTINUE
235: B( N, J ) = B( N, J ) / D( N )
236: IF( N.GT.1 )
237: $ B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) / D( N-1 )
238: DO 80 I = N - 2, 1, -1
239: B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DL( I )*
240: $ B( I+2, J ) ) / D( I )
241: 80 CONTINUE
242: IF( J.LT.NRHS ) THEN
243: J = J + 1
244: GO TO 70
245: END IF
246: ELSE
247: DO 100 J = 1, NRHS
248: B( N, J ) = B( N, J ) / D( N )
249: IF( N.GT.1 )
250: $ B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) /
251: $ D( N-1 )
252: DO 90 I = N - 2, 1, -1
253: B( I, J ) = ( B( I, J )-DU( I )*B( I+1, J )-DL( I )*
254: $ B( I+2, J ) ) / D( I )
255: 90 CONTINUE
256: 100 CONTINUE
257: END IF
258: *
259: RETURN
260: *
261: * End of DGTSV
262: *
263: END