LAPACK  3.6.0
LAPACK: Linear Algebra PACKage
Collaboration diagram for double:

Functions

subroutine dgtcon (NORM, N, DL, D, DU, DU2, IPIV, ANORM, RCOND, WORK, IWORK, INFO)
 DGTCON More...
 
subroutine dgtrfs (TRANS, N, NRHS, DL, D, DU, DLF, DF, DUF, DU2, IPIV, B, LDB, X, LDX, FERR, BERR, WORK, IWORK, INFO)
 DGTRFS More...
 
subroutine dgttrf (N, DL, D, DU, DU2, IPIV, INFO)
 DGTTRF More...
 
subroutine dgttrs (TRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB, INFO)
 DGTTRS More...
 
subroutine dgtts2 (ITRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB)
 DGTTS2 solves a system of linear equations with a tridiagonal matrix using the LU factorization computed by sgttrf. More...
 

Detailed Description

This is the group of double computational functions for GT matrices

Function Documentation

subroutine dgtcon ( character  NORM,
integer  N,
double precision, dimension( * )  DL,
double precision, dimension( * )  D,
double precision, dimension( * )  DU,
double precision, dimension( * )  DU2,
integer, dimension( * )  IPIV,
double precision  ANORM,
double precision  RCOND,
double precision, dimension( * )  WORK,
integer, dimension( * )  IWORK,
integer  INFO 
)

DGTCON

Download DGTCON + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DGTCON estimates the reciprocal of the condition number of a real
 tridiagonal matrix A using the LU factorization as computed by
 DGTTRF.

 An estimate is obtained for norm(inv(A)), and the reciprocal of the
 condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))).
Parameters
[in]NORM
          NORM is CHARACTER*1
          Specifies whether the 1-norm condition number or the
          infinity-norm condition number is required:
          = '1' or 'O':  1-norm;
          = 'I':         Infinity-norm.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in]DL
          DL is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) multipliers that define the matrix L from the
          LU factorization of A as computed by DGTTRF.
[in]D
          D is DOUBLE PRECISION array, dimension (N)
          The n diagonal elements of the upper triangular matrix U from
          the LU factorization of A.
[in]DU
          DU is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) elements of the first superdiagonal of U.
[in]DU2
          DU2 is DOUBLE PRECISION array, dimension (N-2)
          The (n-2) elements of the second superdiagonal of U.
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, row i of the matrix was
          interchanged with row IPIV(i).  IPIV(i) will always be either
          i or i+1; IPIV(i) = i indicates a row interchange was not
          required.
[in]ANORM
          ANORM is DOUBLE PRECISION
          If NORM = '1' or 'O', the 1-norm of the original matrix A.
          If NORM = 'I', the infinity-norm of the original matrix A.
[out]RCOND
          RCOND is DOUBLE PRECISION
          The reciprocal of the condition number of the matrix A,
          computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an
          estimate of the 1-norm of inv(A) computed in this routine.
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (2*N)
[out]IWORK
          IWORK is INTEGER array, dimension (N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 148 of file dgtcon.f.

148 *
149 * -- LAPACK computational routine (version 3.4.2) --
150 * -- LAPACK is a software package provided by Univ. of Tennessee, --
151 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
152 * September 2012
153 *
154 * .. Scalar Arguments ..
155  CHARACTER norm
156  INTEGER info, n
157  DOUBLE PRECISION anorm, rcond
158 * ..
159 * .. Array Arguments ..
160  INTEGER ipiv( * ), iwork( * )
161  DOUBLE PRECISION d( * ), dl( * ), du( * ), du2( * ), work( * )
162 * ..
163 *
164 * =====================================================================
165 *
166 * .. Parameters ..
167  DOUBLE PRECISION one, zero
168  parameter( one = 1.0d+0, zero = 0.0d+0 )
169 * ..
170 * .. Local Scalars ..
171  LOGICAL onenrm
172  INTEGER i, kase, kase1
173  DOUBLE PRECISION ainvnm
174 * ..
175 * .. Local Arrays ..
176  INTEGER isave( 3 )
177 * ..
178 * .. External Functions ..
179  LOGICAL lsame
180  EXTERNAL lsame
181 * ..
182 * .. External Subroutines ..
183  EXTERNAL dgttrs, dlacn2, xerbla
184 * ..
185 * .. Executable Statements ..
186 *
187 * Test the input arguments.
188 *
189  info = 0
190  onenrm = norm.EQ.'1' .OR. lsame( norm, 'O' )
191  IF( .NOT.onenrm .AND. .NOT.lsame( norm, 'I' ) ) THEN
192  info = -1
193  ELSE IF( n.LT.0 ) THEN
194  info = -2
195  ELSE IF( anorm.LT.zero ) THEN
196  info = -8
197  END IF
198  IF( info.NE.0 ) THEN
199  CALL xerbla( 'DGTCON', -info )
200  RETURN
201  END IF
202 *
203 * Quick return if possible
204 *
205  rcond = zero
206  IF( n.EQ.0 ) THEN
207  rcond = one
208  RETURN
209  ELSE IF( anorm.EQ.zero ) THEN
210  RETURN
211  END IF
212 *
213 * Check that D(1:N) is non-zero.
214 *
215  DO 10 i = 1, n
216  IF( d( i ).EQ.zero )
217  $ RETURN
218  10 CONTINUE
219 *
220  ainvnm = zero
221  IF( onenrm ) THEN
222  kase1 = 1
223  ELSE
224  kase1 = 2
225  END IF
226  kase = 0
227  20 CONTINUE
228  CALL dlacn2( n, work( n+1 ), work, iwork, ainvnm, kase, isave )
229  IF( kase.NE.0 ) THEN
230  IF( kase.EQ.kase1 ) THEN
231 *
232 * Multiply by inv(U)*inv(L).
233 *
234  CALL dgttrs( 'No transpose', n, 1, dl, d, du, du2, ipiv,
235  $ work, n, info )
236  ELSE
237 *
238 * Multiply by inv(L**T)*inv(U**T).
239 *
240  CALL dgttrs( 'Transpose', n, 1, dl, d, du, du2, ipiv, work,
241  $ n, info )
242  END IF
243  GO TO 20
244  END IF
245 *
246 * Compute the estimate of the reciprocal condition number.
247 *
248  IF( ainvnm.NE.zero )
249  $ rcond = ( one / ainvnm ) / anorm
250 *
251  RETURN
252 *
253 * End of DGTCON
254 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine dgttrs(TRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB, INFO)
DGTTRS
Definition: dgttrs.f:140
subroutine dlacn2(N, V, X, ISGN, EST, KASE, ISAVE)
DLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vec...
Definition: dlacn2.f:138

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine dgtrfs ( character  TRANS,
integer  N,
integer  NRHS,
double precision, dimension( * )  DL,
double precision, dimension( * )  D,
double precision, dimension( * )  DU,
double precision, dimension( * )  DLF,
double precision, dimension( * )  DF,
double precision, dimension( * )  DUF,
double precision, dimension( * )  DU2,
integer, dimension( * )  IPIV,
double precision, dimension( ldb, * )  B,
integer  LDB,
double precision, dimension( ldx, * )  X,
integer  LDX,
double precision, dimension( * )  FERR,
double precision, dimension( * )  BERR,
double precision, dimension( * )  WORK,
integer, dimension( * )  IWORK,
integer  INFO 
)

DGTRFS

Download DGTRFS + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DGTRFS improves the computed solution to a system of linear
 equations when the coefficient matrix is tridiagonal, and provides
 error bounds and backward error estimates for the solution.
Parameters
[in]TRANS
          TRANS is CHARACTER*1
          Specifies the form of the system of equations:
          = 'N':  A * X = B     (No transpose)
          = 'T':  A**T * X = B  (Transpose)
          = 'C':  A**H * X = B  (Conjugate transpose = Transpose)
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in]NRHS
          NRHS is INTEGER
          The number of right hand sides, i.e., the number of columns
          of the matrix B.  NRHS >= 0.
[in]DL
          DL is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) subdiagonal elements of A.
[in]D
          D is DOUBLE PRECISION array, dimension (N)
          The diagonal elements of A.
[in]DU
          DU is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) superdiagonal elements of A.
[in]DLF
          DLF is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) multipliers that define the matrix L from the
          LU factorization of A as computed by DGTTRF.
[in]DF
          DF is DOUBLE PRECISION array, dimension (N)
          The n diagonal elements of the upper triangular matrix U from
          the LU factorization of A.
[in]DUF
          DUF is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) elements of the first superdiagonal of U.
[in]DU2
          DU2 is DOUBLE PRECISION array, dimension (N-2)
          The (n-2) elements of the second superdiagonal of U.
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, row i of the matrix was
          interchanged with row IPIV(i).  IPIV(i) will always be either
          i or i+1; IPIV(i) = i indicates a row interchange was not
          required.
[in]B
          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
          The right hand side matrix B.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
[in,out]X
          X is DOUBLE PRECISION array, dimension (LDX,NRHS)
          On entry, the solution matrix X, as computed by DGTTRS.
          On exit, the improved solution matrix X.
[in]LDX
          LDX is INTEGER
          The leading dimension of the array X.  LDX >= max(1,N).
[out]FERR
          FERR is DOUBLE PRECISION array, dimension (NRHS)
          The estimated forward error bound for each solution vector
          X(j) (the j-th column of the solution matrix X).
          If XTRUE is the true solution corresponding to X(j), FERR(j)
          is an estimated upper bound for the magnitude of the largest
          element in (X(j) - XTRUE) divided by the magnitude of the
          largest element in X(j).  The estimate is as reliable as
          the estimate for RCOND, and is almost always a slight
          overestimate of the true error.
[out]BERR
          BERR is DOUBLE PRECISION array, dimension (NRHS)
          The componentwise relative backward error of each solution
          vector X(j) (i.e., the smallest relative change in
          any element of A or B that makes X(j) an exact solution).
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (3*N)
[out]IWORK
          IWORK is INTEGER array, dimension (N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
Internal Parameters:
  ITMAX is the maximum number of steps of iterative refinement.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 211 of file dgtrfs.f.

211 *
212 * -- LAPACK computational routine (version 3.4.2) --
213 * -- LAPACK is a software package provided by Univ. of Tennessee, --
214 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
215 * September 2012
216 *
217 * .. Scalar Arguments ..
218  CHARACTER trans
219  INTEGER info, ldb, ldx, n, nrhs
220 * ..
221 * .. Array Arguments ..
222  INTEGER ipiv( * ), iwork( * )
223  DOUBLE PRECISION b( ldb, * ), berr( * ), d( * ), df( * ),
224  $ dl( * ), dlf( * ), du( * ), du2( * ), duf( * ),
225  $ ferr( * ), work( * ), x( ldx, * )
226 * ..
227 *
228 * =====================================================================
229 *
230 * .. Parameters ..
231  INTEGER itmax
232  parameter( itmax = 5 )
233  DOUBLE PRECISION zero, one
234  parameter( zero = 0.0d+0, one = 1.0d+0 )
235  DOUBLE PRECISION two
236  parameter( two = 2.0d+0 )
237  DOUBLE PRECISION three
238  parameter( three = 3.0d+0 )
239 * ..
240 * .. Local Scalars ..
241  LOGICAL notran
242  CHARACTER transn, transt
243  INTEGER count, i, j, kase, nz
244  DOUBLE PRECISION eps, lstres, s, safe1, safe2, safmin
245 * ..
246 * .. Local Arrays ..
247  INTEGER isave( 3 )
248 * ..
249 * .. External Subroutines ..
250  EXTERNAL daxpy, dcopy, dgttrs, dlacn2, dlagtm, xerbla
251 * ..
252 * .. Intrinsic Functions ..
253  INTRINSIC abs, max
254 * ..
255 * .. External Functions ..
256  LOGICAL lsame
257  DOUBLE PRECISION dlamch
258  EXTERNAL lsame, dlamch
259 * ..
260 * .. Executable Statements ..
261 *
262 * Test the input parameters.
263 *
264  info = 0
265  notran = lsame( trans, 'N' )
266  IF( .NOT.notran .AND. .NOT.lsame( trans, 'T' ) .AND. .NOT.
267  $ lsame( trans, 'C' ) ) THEN
268  info = -1
269  ELSE IF( n.LT.0 ) THEN
270  info = -2
271  ELSE IF( nrhs.LT.0 ) THEN
272  info = -3
273  ELSE IF( ldb.LT.max( 1, n ) ) THEN
274  info = -13
275  ELSE IF( ldx.LT.max( 1, n ) ) THEN
276  info = -15
277  END IF
278  IF( info.NE.0 ) THEN
279  CALL xerbla( 'DGTRFS', -info )
280  RETURN
281  END IF
282 *
283 * Quick return if possible
284 *
285  IF( n.EQ.0 .OR. nrhs.EQ.0 ) THEN
286  DO 10 j = 1, nrhs
287  ferr( j ) = zero
288  berr( j ) = zero
289  10 CONTINUE
290  RETURN
291  END IF
292 *
293  IF( notran ) THEN
294  transn = 'N'
295  transt = 'T'
296  ELSE
297  transn = 'T'
298  transt = 'N'
299  END IF
300 *
301 * NZ = maximum number of nonzero elements in each row of A, plus 1
302 *
303  nz = 4
304  eps = dlamch( 'Epsilon' )
305  safmin = dlamch( 'Safe minimum' )
306  safe1 = nz*safmin
307  safe2 = safe1 / eps
308 *
309 * Do for each right hand side
310 *
311  DO 110 j = 1, nrhs
312 *
313  count = 1
314  lstres = three
315  20 CONTINUE
316 *
317 * Loop until stopping criterion is satisfied.
318 *
319 * Compute residual R = B - op(A) * X,
320 * where op(A) = A, A**T, or A**H, depending on TRANS.
321 *
322  CALL dcopy( n, b( 1, j ), 1, work( n+1 ), 1 )
323  CALL dlagtm( trans, n, 1, -one, dl, d, du, x( 1, j ), ldx, one,
324  $ work( n+1 ), n )
325 *
326 * Compute abs(op(A))*abs(x) + abs(b) for use in the backward
327 * error bound.
328 *
329  IF( notran ) THEN
330  IF( n.EQ.1 ) THEN
331  work( 1 ) = abs( b( 1, j ) ) + abs( d( 1 )*x( 1, j ) )
332  ELSE
333  work( 1 ) = abs( b( 1, j ) ) + abs( d( 1 )*x( 1, j ) ) +
334  $ abs( du( 1 )*x( 2, j ) )
335  DO 30 i = 2, n - 1
336  work( i ) = abs( b( i, j ) ) +
337  $ abs( dl( i-1 )*x( i-1, j ) ) +
338  $ abs( d( i )*x( i, j ) ) +
339  $ abs( du( i )*x( i+1, j ) )
340  30 CONTINUE
341  work( n ) = abs( b( n, j ) ) +
342  $ abs( dl( n-1 )*x( n-1, j ) ) +
343  $ abs( d( n )*x( n, j ) )
344  END IF
345  ELSE
346  IF( n.EQ.1 ) THEN
347  work( 1 ) = abs( b( 1, j ) ) + abs( d( 1 )*x( 1, j ) )
348  ELSE
349  work( 1 ) = abs( b( 1, j ) ) + abs( d( 1 )*x( 1, j ) ) +
350  $ abs( dl( 1 )*x( 2, j ) )
351  DO 40 i = 2, n - 1
352  work( i ) = abs( b( i, j ) ) +
353  $ abs( du( i-1 )*x( i-1, j ) ) +
354  $ abs( d( i )*x( i, j ) ) +
355  $ abs( dl( i )*x( i+1, j ) )
356  40 CONTINUE
357  work( n ) = abs( b( n, j ) ) +
358  $ abs( du( n-1 )*x( n-1, j ) ) +
359  $ abs( d( n )*x( n, j ) )
360  END IF
361  END IF
362 *
363 * Compute componentwise relative backward error from formula
364 *
365 * max(i) ( abs(R(i)) / ( abs(op(A))*abs(X) + abs(B) )(i) )
366 *
367 * where abs(Z) is the componentwise absolute value of the matrix
368 * or vector Z. If the i-th component of the denominator is less
369 * than SAFE2, then SAFE1 is added to the i-th components of the
370 * numerator and denominator before dividing.
371 *
372  s = zero
373  DO 50 i = 1, n
374  IF( work( i ).GT.safe2 ) THEN
375  s = max( s, abs( work( n+i ) ) / work( i ) )
376  ELSE
377  s = max( s, ( abs( work( n+i ) )+safe1 ) /
378  $ ( work( i )+safe1 ) )
379  END IF
380  50 CONTINUE
381  berr( j ) = s
382 *
383 * Test stopping criterion. Continue iterating if
384 * 1) The residual BERR(J) is larger than machine epsilon, and
385 * 2) BERR(J) decreased by at least a factor of 2 during the
386 * last iteration, and
387 * 3) At most ITMAX iterations tried.
388 *
389  IF( berr( j ).GT.eps .AND. two*berr( j ).LE.lstres .AND.
390  $ count.LE.itmax ) THEN
391 *
392 * Update solution and try again.
393 *
394  CALL dgttrs( trans, n, 1, dlf, df, duf, du2, ipiv,
395  $ work( n+1 ), n, info )
396  CALL daxpy( n, one, work( n+1 ), 1, x( 1, j ), 1 )
397  lstres = berr( j )
398  count = count + 1
399  GO TO 20
400  END IF
401 *
402 * Bound error from formula
403 *
404 * norm(X - XTRUE) / norm(X) .le. FERR =
405 * norm( abs(inv(op(A)))*
406 * ( abs(R) + NZ*EPS*( abs(op(A))*abs(X)+abs(B) ))) / norm(X)
407 *
408 * where
409 * norm(Z) is the magnitude of the largest component of Z
410 * inv(op(A)) is the inverse of op(A)
411 * abs(Z) is the componentwise absolute value of the matrix or
412 * vector Z
413 * NZ is the maximum number of nonzeros in any row of A, plus 1
414 * EPS is machine epsilon
415 *
416 * The i-th component of abs(R)+NZ*EPS*(abs(op(A))*abs(X)+abs(B))
417 * is incremented by SAFE1 if the i-th component of
418 * abs(op(A))*abs(X) + abs(B) is less than SAFE2.
419 *
420 * Use DLACN2 to estimate the infinity-norm of the matrix
421 * inv(op(A)) * diag(W),
422 * where W = abs(R) + NZ*EPS*( abs(op(A))*abs(X)+abs(B) )))
423 *
424  DO 60 i = 1, n
425  IF( work( i ).GT.safe2 ) THEN
426  work( i ) = abs( work( n+i ) ) + nz*eps*work( i )
427  ELSE
428  work( i ) = abs( work( n+i ) ) + nz*eps*work( i ) + safe1
429  END IF
430  60 CONTINUE
431 *
432  kase = 0
433  70 CONTINUE
434  CALL dlacn2( n, work( 2*n+1 ), work( n+1 ), iwork, ferr( j ),
435  $ kase, isave )
436  IF( kase.NE.0 ) THEN
437  IF( kase.EQ.1 ) THEN
438 *
439 * Multiply by diag(W)*inv(op(A)**T).
440 *
441  CALL dgttrs( transt, n, 1, dlf, df, duf, du2, ipiv,
442  $ work( n+1 ), n, info )
443  DO 80 i = 1, n
444  work( n+i ) = work( i )*work( n+i )
445  80 CONTINUE
446  ELSE
447 *
448 * Multiply by inv(op(A))*diag(W).
449 *
450  DO 90 i = 1, n
451  work( n+i ) = work( i )*work( n+i )
452  90 CONTINUE
453  CALL dgttrs( transn, n, 1, dlf, df, duf, du2, ipiv,
454  $ work( n+1 ), n, info )
455  END IF
456  GO TO 70
457  END IF
458 *
459 * Normalize error.
460 *
461  lstres = zero
462  DO 100 i = 1, n
463  lstres = max( lstres, abs( x( i, j ) ) )
464  100 CONTINUE
465  IF( lstres.NE.zero )
466  $ ferr( j ) = ferr( j ) / lstres
467 *
468  110 CONTINUE
469 *
470  RETURN
471 *
472 * End of DGTRFS
473 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine dlagtm(TRANS, N, NRHS, ALPHA, DL, D, DU, X, LDX, BETA, B, LDB)
DLAGTM performs a matrix-matrix product of the form C = αAB+βC, where A is a tridiagonal matrix...
Definition: dlagtm.f:147
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine dgttrs(TRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB, INFO)
DGTTRS
Definition: dgttrs.f:140
subroutine dlacn2(N, V, X, ISGN, EST, KASE, ISAVE)
DLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vec...
Definition: dlacn2.f:138
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:53
subroutine daxpy(N, DA, DX, INCX, DY, INCY)
DAXPY
Definition: daxpy.f:54

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine dgttrf ( integer  N,
double precision, dimension( * )  DL,
double precision, dimension( * )  D,
double precision, dimension( * )  DU,
double precision, dimension( * )  DU2,
integer, dimension( * )  IPIV,
integer  INFO 
)

DGTTRF

Download DGTTRF + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DGTTRF computes an LU factorization of a real tridiagonal matrix A
 using elimination with partial pivoting and row interchanges.

 The factorization has the form
    A = L * U
 where L is a product of permutation and unit lower bidiagonal
 matrices and U is upper triangular with nonzeros in only the main
 diagonal and first two superdiagonals.
Parameters
[in]N
          N is INTEGER
          The order of the matrix A.
[in,out]DL
          DL is DOUBLE PRECISION array, dimension (N-1)
          On entry, DL must contain the (n-1) sub-diagonal elements of
          A.

          On exit, DL is overwritten by the (n-1) multipliers that
          define the matrix L from the LU factorization of A.
[in,out]D
          D is DOUBLE PRECISION array, dimension (N)
          On entry, D must contain the diagonal elements of A.

          On exit, D is overwritten by the n diagonal elements of the
          upper triangular matrix U from the LU factorization of A.
[in,out]DU
          DU is DOUBLE PRECISION array, dimension (N-1)
          On entry, DU must contain the (n-1) super-diagonal elements
          of A.

          On exit, DU is overwritten by the (n-1) elements of the first
          super-diagonal of U.
[out]DU2
          DU2 is DOUBLE PRECISION array, dimension (N-2)
          On exit, DU2 is overwritten by the (n-2) elements of the
          second super-diagonal of U.
[out]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, row i of the matrix was
          interchanged with row IPIV(i).  IPIV(i) will always be either
          i or i+1; IPIV(i) = i indicates a row interchange was not
          required.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -k, the k-th argument had an illegal value
          > 0:  if INFO = k, U(k,k) is exactly zero. The factorization
                has been completed, but the factor U is exactly
                singular, and division by zero will occur if it is used
                to solve a system of equations.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 126 of file dgttrf.f.

126 *
127 * -- LAPACK computational routine (version 3.4.2) --
128 * -- LAPACK is a software package provided by Univ. of Tennessee, --
129 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
130 * September 2012
131 *
132 * .. Scalar Arguments ..
133  INTEGER info, n
134 * ..
135 * .. Array Arguments ..
136  INTEGER ipiv( * )
137  DOUBLE PRECISION d( * ), dl( * ), du( * ), du2( * )
138 * ..
139 *
140 * =====================================================================
141 *
142 * .. Parameters ..
143  DOUBLE PRECISION zero
144  parameter( zero = 0.0d+0 )
145 * ..
146 * .. Local Scalars ..
147  INTEGER i
148  DOUBLE PRECISION fact, temp
149 * ..
150 * .. Intrinsic Functions ..
151  INTRINSIC abs
152 * ..
153 * .. External Subroutines ..
154  EXTERNAL xerbla
155 * ..
156 * .. Executable Statements ..
157 *
158  info = 0
159  IF( n.LT.0 ) THEN
160  info = -1
161  CALL xerbla( 'DGTTRF', -info )
162  RETURN
163  END IF
164 *
165 * Quick return if possible
166 *
167  IF( n.EQ.0 )
168  $ RETURN
169 *
170 * Initialize IPIV(i) = i and DU2(I) = 0
171 *
172  DO 10 i = 1, n
173  ipiv( i ) = i
174  10 CONTINUE
175  DO 20 i = 1, n - 2
176  du2( i ) = zero
177  20 CONTINUE
178 *
179  DO 30 i = 1, n - 2
180  IF( abs( d( i ) ).GE.abs( dl( i ) ) ) THEN
181 *
182 * No row interchange required, eliminate DL(I)
183 *
184  IF( d( i ).NE.zero ) THEN
185  fact = dl( i ) / d( i )
186  dl( i ) = fact
187  d( i+1 ) = d( i+1 ) - fact*du( i )
188  END IF
189  ELSE
190 *
191 * Interchange rows I and I+1, eliminate DL(I)
192 *
193  fact = d( i ) / dl( i )
194  d( i ) = dl( i )
195  dl( i ) = fact
196  temp = du( i )
197  du( i ) = d( i+1 )
198  d( i+1 ) = temp - fact*d( i+1 )
199  du2( i ) = du( i+1 )
200  du( i+1 ) = -fact*du( i+1 )
201  ipiv( i ) = i + 1
202  END IF
203  30 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  dl( i ) = fact
210  d( i+1 ) = d( i+1 ) - fact*du( i )
211  END IF
212  ELSE
213  fact = d( i ) / dl( i )
214  d( i ) = dl( i )
215  dl( i ) = fact
216  temp = du( i )
217  du( i ) = d( i+1 )
218  d( i+1 ) = temp - fact*d( i+1 )
219  ipiv( i ) = i + 1
220  END IF
221  END IF
222 *
223 * Check for a zero on the diagonal of U.
224 *
225  DO 40 i = 1, n
226  IF( d( i ).EQ.zero ) THEN
227  info = i
228  GO TO 50
229  END IF
230  40 CONTINUE
231  50 CONTINUE
232 *
233  RETURN
234 *
235 * End of DGTTRF
236 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine dgttrs ( character  TRANS,
integer  N,
integer  NRHS,
double precision, dimension( * )  DL,
double precision, dimension( * )  D,
double precision, dimension( * )  DU,
double precision, dimension( * )  DU2,
integer, dimension( * )  IPIV,
double precision, dimension( ldb, * )  B,
integer  LDB,
integer  INFO 
)

DGTTRS

Download DGTTRS + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DGTTRS solves one of the systems of equations
    A*X = B  or  A**T*X = B,
 with a tridiagonal matrix A using the LU factorization computed
 by DGTTRF.
Parameters
[in]TRANS
          TRANS is CHARACTER*1
          Specifies the form of the system of equations.
          = 'N':  A * X = B  (No transpose)
          = 'T':  A**T* X = B  (Transpose)
          = 'C':  A**T* X = B  (Conjugate transpose = Transpose)
[in]N
          N is INTEGER
          The order of the matrix A.
[in]NRHS
          NRHS is INTEGER
          The number of right hand sides, i.e., the number of columns
          of the matrix B.  NRHS >= 0.
[in]DL
          DL is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) multipliers that define the matrix L from the
          LU factorization of A.
[in]D
          D is DOUBLE PRECISION array, dimension (N)
          The n diagonal elements of the upper triangular matrix U from
          the LU factorization of A.
[in]DU
          DU is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) elements of the first super-diagonal of U.
[in]DU2
          DU2 is DOUBLE PRECISION array, dimension (N-2)
          The (n-2) elements of the second super-diagonal of U.
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, row i of the matrix was
          interchanged with row IPIV(i).  IPIV(i) will always be either
          i or i+1; IPIV(i) = i indicates a row interchange was not
          required.
[in,out]B
          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
          On entry, the matrix of right hand side vectors B.
          On exit, B is overwritten by the solution vectors X.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 140 of file dgttrs.f.

140 *
141 * -- LAPACK computational routine (version 3.4.2) --
142 * -- LAPACK is a software package provided by Univ. of Tennessee, --
143 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
144 * September 2012
145 *
146 * .. Scalar Arguments ..
147  CHARACTER trans
148  INTEGER info, ldb, n, nrhs
149 * ..
150 * .. Array Arguments ..
151  INTEGER ipiv( * )
152  DOUBLE PRECISION b( ldb, * ), d( * ), dl( * ), du( * ), du2( * )
153 * ..
154 *
155 * =====================================================================
156 *
157 * .. Local Scalars ..
158  LOGICAL notran
159  INTEGER itrans, j, jb, nb
160 * ..
161 * .. External Functions ..
162  INTEGER ilaenv
163  EXTERNAL ilaenv
164 * ..
165 * .. External Subroutines ..
166  EXTERNAL dgtts2, xerbla
167 * ..
168 * .. Intrinsic Functions ..
169  INTRINSIC max, min
170 * ..
171 * .. Executable Statements ..
172 *
173  info = 0
174  notran = ( trans.EQ.'N' .OR. trans.EQ.'n' )
175  IF( .NOT.notran .AND. .NOT.( trans.EQ.'T' .OR. trans.EQ.
176  $ 't' ) .AND. .NOT.( trans.EQ.'C' .OR. trans.EQ.'c' ) ) THEN
177  info = -1
178  ELSE IF( n.LT.0 ) THEN
179  info = -2
180  ELSE IF( nrhs.LT.0 ) THEN
181  info = -3
182  ELSE IF( ldb.LT.max( n, 1 ) ) THEN
183  info = -10
184  END IF
185  IF( info.NE.0 ) THEN
186  CALL xerbla( 'DGTTRS', -info )
187  RETURN
188  END IF
189 *
190 * Quick return if possible
191 *
192  IF( n.EQ.0 .OR. nrhs.EQ.0 )
193  $ RETURN
194 *
195 * Decode TRANS
196 *
197  IF( notran ) THEN
198  itrans = 0
199  ELSE
200  itrans = 1
201  END IF
202 *
203 * Determine the number of right-hand sides to solve at a time.
204 *
205  IF( nrhs.EQ.1 ) THEN
206  nb = 1
207  ELSE
208  nb = max( 1, ilaenv( 1, 'DGTTRS', trans, n, nrhs, -1, -1 ) )
209  END IF
210 *
211  IF( nb.GE.nrhs ) THEN
212  CALL dgtts2( itrans, n, nrhs, dl, d, du, du2, ipiv, b, ldb )
213  ELSE
214  DO 10 j = 1, nrhs, nb
215  jb = min( nrhs-j+1, nb )
216  CALL dgtts2( itrans, n, jb, dl, d, du, du2, ipiv, b( 1, j ),
217  $ ldb )
218  10 CONTINUE
219  END IF
220 *
221 * End of DGTTRS
222 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine dgtts2(ITRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB)
DGTTS2 solves a system of linear equations with a tridiagonal matrix using the LU factorization compu...
Definition: dgtts2.f:130
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine dgtts2 ( integer  ITRANS,
integer  N,
integer  NRHS,
double precision, dimension( * )  DL,
double precision, dimension( * )  D,
double precision, dimension( * )  DU,
double precision, dimension( * )  DU2,
integer, dimension( * )  IPIV,
double precision, dimension( ldb, * )  B,
integer  LDB 
)

DGTTS2 solves a system of linear equations with a tridiagonal matrix using the LU factorization computed by sgttrf.

Download DGTTS2 + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DGTTS2 solves one of the systems of equations
    A*X = B  or  A**T*X = B,
 with a tridiagonal matrix A using the LU factorization computed
 by DGTTRF.
Parameters
[in]ITRANS
          ITRANS is INTEGER
          Specifies the form of the system of equations.
          = 0:  A * X = B  (No transpose)
          = 1:  A**T* X = B  (Transpose)
          = 2:  A**T* X = B  (Conjugate transpose = Transpose)
[in]N
          N is INTEGER
          The order of the matrix A.
[in]NRHS
          NRHS is INTEGER
          The number of right hand sides, i.e., the number of columns
          of the matrix B.  NRHS >= 0.
[in]DL
          DL is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) multipliers that define the matrix L from the
          LU factorization of A.
[in]D
          D is DOUBLE PRECISION array, dimension (N)
          The n diagonal elements of the upper triangular matrix U from
          the LU factorization of A.
[in]DU
          DU is DOUBLE PRECISION array, dimension (N-1)
          The (n-1) elements of the first super-diagonal of U.
[in]DU2
          DU2 is DOUBLE PRECISION array, dimension (N-2)
          The (n-2) elements of the second super-diagonal of U.
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices; for 1 <= i <= n, row i of the matrix was
          interchanged with row IPIV(i).  IPIV(i) will always be either
          i or i+1; IPIV(i) = i indicates a row interchange was not
          required.
[in,out]B
          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
          On entry, the matrix of right hand side vectors B.
          On exit, B is overwritten by the solution vectors X.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 130 of file dgtts2.f.

130 *
131 * -- LAPACK computational routine (version 3.4.2) --
132 * -- LAPACK is a software package provided by Univ. of Tennessee, --
133 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
134 * September 2012
135 *
136 * .. Scalar Arguments ..
137  INTEGER itrans, ldb, n, nrhs
138 * ..
139 * .. Array Arguments ..
140  INTEGER ipiv( * )
141  DOUBLE PRECISION b( ldb, * ), d( * ), dl( * ), du( * ), du2( * )
142 * ..
143 *
144 * =====================================================================
145 *
146 * .. Local Scalars ..
147  INTEGER i, ip, j
148  DOUBLE PRECISION temp
149 * ..
150 * .. Executable Statements ..
151 *
152 * Quick return if possible
153 *
154  IF( n.EQ.0 .OR. nrhs.EQ.0 )
155  $ RETURN
156 *
157  IF( itrans.EQ.0 ) THEN
158 *
159 * Solve A*X = B using the LU factorization of A,
160 * overwriting each right hand side vector with its solution.
161 *
162  IF( nrhs.LE.1 ) THEN
163  j = 1
164  10 CONTINUE
165 *
166 * Solve L*x = b.
167 *
168  DO 20 i = 1, n - 1
169  ip = ipiv( i )
170  temp = b( i+1-ip+i, j ) - dl( i )*b( ip, j )
171  b( i, j ) = b( ip, j )
172  b( i+1, j ) = temp
173  20 CONTINUE
174 *
175 * Solve U*x = b.
176 *
177  b( n, j ) = b( n, j ) / d( n )
178  IF( n.GT.1 )
179  $ b( n-1, j ) = ( b( n-1, j )-du( n-1 )*b( n, j ) ) /
180  $ d( n-1 )
181  DO 30 i = n - 2, 1, -1
182  b( i, j ) = ( b( i, j )-du( i )*b( i+1, j )-du2( i )*
183  $ b( i+2, j ) ) / d( i )
184  30 CONTINUE
185  IF( j.LT.nrhs ) THEN
186  j = j + 1
187  GO TO 10
188  END IF
189  ELSE
190  DO 60 j = 1, nrhs
191 *
192 * Solve L*x = b.
193 *
194  DO 40 i = 1, n - 1
195  IF( ipiv( i ).EQ.i ) THEN
196  b( i+1, j ) = b( i+1, j ) - dl( i )*b( i, j )
197  ELSE
198  temp = b( i, j )
199  b( i, j ) = b( i+1, j )
200  b( i+1, j ) = temp - dl( i )*b( i, j )
201  END IF
202  40 CONTINUE
203 *
204 * Solve U*x = b.
205 *
206  b( n, j ) = b( n, j ) / d( n )
207  IF( n.GT.1 )
208  $ b( n-1, j ) = ( b( n-1, j )-du( n-1 )*b( n, j ) ) /
209  $ d( n-1 )
210  DO 50 i = n - 2, 1, -1
211  b( i, j ) = ( b( i, j )-du( i )*b( i+1, j )-du2( i )*
212  $ b( i+2, j ) ) / d( i )
213  50 CONTINUE
214  60 CONTINUE
215  END IF
216  ELSE
217 *
218 * Solve A**T * X = B.
219 *
220  IF( nrhs.LE.1 ) THEN
221 *
222 * Solve U**T*x = b.
223 *
224  j = 1
225  70 CONTINUE
226  b( 1, j ) = b( 1, j ) / d( 1 )
227  IF( n.GT.1 )
228  $ b( 2, j ) = ( b( 2, j )-du( 1 )*b( 1, j ) ) / d( 2 )
229  DO 80 i = 3, n
230  b( i, j ) = ( b( i, j )-du( i-1 )*b( i-1, j )-du2( i-2 )*
231  $ b( i-2, j ) ) / d( i )
232  80 CONTINUE
233 *
234 * Solve L**T*x = b.
235 *
236  DO 90 i = n - 1, 1, -1
237  ip = ipiv( i )
238  temp = b( i, j ) - dl( i )*b( i+1, j )
239  b( i, j ) = b( ip, j )
240  b( ip, j ) = temp
241  90 CONTINUE
242  IF( j.LT.nrhs ) THEN
243  j = j + 1
244  GO TO 70
245  END IF
246 *
247  ELSE
248  DO 120 j = 1, nrhs
249 *
250 * Solve U**T*x = b.
251 *
252  b( 1, j ) = b( 1, j ) / d( 1 )
253  IF( n.GT.1 )
254  $ b( 2, j ) = ( b( 2, j )-du( 1 )*b( 1, j ) ) / d( 2 )
255  DO 100 i = 3, n
256  b( i, j ) = ( b( i, j )-du( i-1 )*b( i-1, j )-
257  $ du2( i-2 )*b( i-2, j ) ) / d( i )
258  100 CONTINUE
259  DO 110 i = n - 1, 1, -1
260  IF( ipiv( i ).EQ.i ) THEN
261  b( i, j ) = b( i, j ) - dl( i )*b( i+1, j )
262  ELSE
263  temp = b( i+1, j )
264  b( i+1, j ) = b( i, j ) - dl( i )*temp
265  b( i, j ) = temp
266  END IF
267  110 CONTINUE
268  120 CONTINUE
269  END IF
270  END IF
271 *
272 * End of DGTTS2
273 *

Here is the caller graph for this function: