LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
dgtsv.f
Go to the documentation of this file.
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