LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
dlagtf.f
Go to the documentation of this file.
1 *> \brief \b DLAGTF computes an LU factorization of a matrix T-λI, where T is a general tridiagonal matrix, and λ a scalar, using partial pivoting with row interchanges.
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download DLAGTF + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlagtf.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlagtf.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlagtf.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE DLAGTF( N, A, LAMBDA, B, C, TOL, D, IN, INFO )
22 *
23 * .. Scalar Arguments ..
24 * INTEGER INFO, N
25 * DOUBLE PRECISION LAMBDA, TOL
26 * ..
27 * .. Array Arguments ..
28 * INTEGER IN( * )
29 * DOUBLE PRECISION A( * ), B( * ), C( * ), D( * )
30 * ..
31 *
32 *
33 *> \par Purpose:
34 * =============
35 *>
36 *> \verbatim
37 *>
38 *> DLAGTF factorizes the matrix (T - lambda*I), where T is an n by n
39 *> tridiagonal matrix and lambda is a scalar, as
40 *>
41 *> T - lambda*I = PLU,
42 *>
43 *> where P is a permutation matrix, L is a unit lower tridiagonal matrix
44 *> with at most one non-zero sub-diagonal elements per column and U is
45 *> an upper triangular matrix with at most two non-zero super-diagonal
46 *> elements per column.
47 *>
48 *> The factorization is obtained by Gaussian elimination with partial
49 *> pivoting and implicit row scaling.
50 *>
51 *> The parameter LAMBDA is included in the routine so that DLAGTF may
52 *> be used, in conjunction with DLAGTS, to obtain eigenvectors of T by
53 *> inverse iteration.
54 *> \endverbatim
55 *
56 * Arguments:
57 * ==========
58 *
59 *> \param[in] N
60 *> \verbatim
61 *> N is INTEGER
62 *> The order of the matrix T.
63 *> \endverbatim
64 *>
65 *> \param[in,out] A
66 *> \verbatim
67 *> A is DOUBLE PRECISION array, dimension (N)
68 *> On entry, A must contain the diagonal elements of T.
69 *>
70 *> On exit, A is overwritten by the n diagonal elements of the
71 *> upper triangular matrix U of the factorization of T.
72 *> \endverbatim
73 *>
74 *> \param[in] LAMBDA
75 *> \verbatim
76 *> LAMBDA is DOUBLE PRECISION
77 *> On entry, the scalar lambda.
78 *> \endverbatim
79 *>
80 *> \param[in,out] B
81 *> \verbatim
82 *> B is DOUBLE PRECISION array, dimension (N-1)
83 *> On entry, B must contain the (n-1) super-diagonal elements of
84 *> T.
85 *>
86 *> On exit, B is overwritten by the (n-1) super-diagonal
87 *> elements of the matrix U of the factorization of T.
88 *> \endverbatim
89 *>
90 *> \param[in,out] C
91 *> \verbatim
92 *> C is DOUBLE PRECISION array, dimension (N-1)
93 *> On entry, C must contain the (n-1) sub-diagonal elements of
94 *> T.
95 *>
96 *> On exit, C is overwritten by the (n-1) sub-diagonal elements
97 *> of the matrix L of the factorization of T.
98 *> \endverbatim
99 *>
100 *> \param[in] TOL
101 *> \verbatim
102 *> TOL is DOUBLE PRECISION
103 *> On entry, a relative tolerance used to indicate whether or
104 *> not the matrix (T - lambda*I) is nearly singular. TOL should
105 *> normally be chose as approximately the largest relative error
106 *> in the elements of T. For example, if the elements of T are
107 *> correct to about 4 significant figures, then TOL should be
108 *> set to about 5*10**(-4). If TOL is supplied as less than eps,
109 *> where eps is the relative machine precision, then the value
110 *> eps is used in place of TOL.
111 *> \endverbatim
112 *>
113 *> \param[out] D
114 *> \verbatim
115 *> D is DOUBLE PRECISION array, dimension (N-2)
116 *> On exit, D is overwritten by the (n-2) second super-diagonal
117 *> elements of the matrix U of the factorization of T.
118 *> \endverbatim
119 *>
120 *> \param[out] IN
121 *> \verbatim
122 *> IN is INTEGER array, dimension (N)
123 *> On exit, IN contains details of the permutation matrix P. If
124 *> an interchange occurred at the kth step of the elimination,
125 *> then IN(k) = 1, otherwise IN(k) = 0. The element IN(n)
126 *> returns the smallest positive integer j such that
127 *>
128 *> abs( u(j,j) ) <= norm( (T - lambda*I)(j) )*TOL,
129 *>
130 *> where norm( A(j) ) denotes the sum of the absolute values of
131 *> the jth row of the matrix A. If no such j exists then IN(n)
132 *> is returned as zero. If IN(n) is returned as positive, then a
133 *> diagonal element of U is small, indicating that
134 *> (T - lambda*I) is singular or nearly singular,
135 *> \endverbatim
136 *>
137 *> \param[out] INFO
138 *> \verbatim
139 *> INFO is INTEGER
140 *> = 0: successful exit
141 *> < 0: if INFO = -k, the kth argument had an illegal value
142 *> \endverbatim
143 *
144 * Authors:
145 * ========
146 *
147 *> \author Univ. of Tennessee
148 *> \author Univ. of California Berkeley
149 *> \author Univ. of Colorado Denver
150 *> \author NAG Ltd.
151 *
152 *> \ingroup auxOTHERcomputational
153 *
154 * =====================================================================
155  SUBROUTINE dlagtf( N, A, LAMBDA, B, C, TOL, D, IN, INFO )
156 *
157 * -- LAPACK computational routine --
158 * -- LAPACK is a software package provided by Univ. of Tennessee, --
159 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
160 *
161 * .. Scalar Arguments ..
162  INTEGER INFO, N
163  DOUBLE PRECISION LAMBDA, TOL
164 * ..
165 * .. Array Arguments ..
166  INTEGER IN( * )
167  DOUBLE PRECISION A( * ), B( * ), C( * ), D( * )
168 * ..
169 *
170 * =====================================================================
171 *
172 * .. Parameters ..
173  DOUBLE PRECISION ZERO
174  parameter( zero = 0.0d+0 )
175 * ..
176 * .. Local Scalars ..
177  INTEGER K
178  DOUBLE PRECISION EPS, MULT, PIV1, PIV2, SCALE1, SCALE2, TEMP, TL
179 * ..
180 * .. Intrinsic Functions ..
181  INTRINSIC abs, max
182 * ..
183 * .. External Functions ..
184  DOUBLE PRECISION DLAMCH
185  EXTERNAL dlamch
186 * ..
187 * .. External Subroutines ..
188  EXTERNAL xerbla
189 * ..
190 * .. Executable Statements ..
191 *
192  info = 0
193  IF( n.LT.0 ) THEN
194  info = -1
195  CALL xerbla( 'DLAGTF', -info )
196  RETURN
197  END IF
198 *
199  IF( n.EQ.0 )
200  $ RETURN
201 *
202  a( 1 ) = a( 1 ) - lambda
203  in( n ) = 0
204  IF( n.EQ.1 ) THEN
205  IF( a( 1 ).EQ.zero )
206  $ in( 1 ) = 1
207  RETURN
208  END IF
209 *
210  eps = dlamch( 'Epsilon' )
211 *
212  tl = max( tol, eps )
213  scale1 = abs( a( 1 ) ) + abs( b( 1 ) )
214  DO 10 k = 1, n - 1
215  a( k+1 ) = a( k+1 ) - lambda
216  scale2 = abs( c( k ) ) + abs( a( k+1 ) )
217  IF( k.LT.( n-1 ) )
218  $ scale2 = scale2 + abs( b( k+1 ) )
219  IF( a( k ).EQ.zero ) THEN
220  piv1 = zero
221  ELSE
222  piv1 = abs( a( k ) ) / scale1
223  END IF
224  IF( c( k ).EQ.zero ) THEN
225  in( k ) = 0
226  piv2 = zero
227  scale1 = scale2
228  IF( k.LT.( n-1 ) )
229  $ d( k ) = zero
230  ELSE
231  piv2 = abs( c( k ) ) / scale2
232  IF( piv2.LE.piv1 ) THEN
233  in( k ) = 0
234  scale1 = scale2
235  c( k ) = c( k ) / a( k )
236  a( k+1 ) = a( k+1 ) - c( k )*b( k )
237  IF( k.LT.( n-1 ) )
238  $ d( k ) = zero
239  ELSE
240  in( k ) = 1
241  mult = a( k ) / c( k )
242  a( k ) = c( k )
243  temp = a( k+1 )
244  a( k+1 ) = b( k ) - mult*temp
245  IF( k.LT.( n-1 ) ) THEN
246  d( k ) = b( k+1 )
247  b( k+1 ) = -mult*d( k )
248  END IF
249  b( k ) = temp
250  c( k ) = mult
251  END IF
252  END IF
253  IF( ( max( piv1, piv2 ).LE.tl ) .AND. ( in( n ).EQ.0 ) )
254  $ in( n ) = k
255  10 CONTINUE
256  IF( ( abs( a( n ) ).LE.scale1*tl ) .AND. ( in( n ).EQ.0 ) )
257  $ in( n ) = n
258 *
259  RETURN
260 *
261 * End of DLAGTF
262 *
263  END
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine dlagtf(N, A, LAMBDA, B, C, TOL, D, IN, INFO)
DLAGTF computes an LU factorization of a matrix T-λI, where T is a general tridiagonal matrix,...
Definition: dlagtf.f:156