NAG FL Interface
f01lef (real_gen_tridiag_lu)
1
Purpose
f01lef computes an $LU$ factorization of a real tridiagonal matrix, using Gaussian elimination with partial pivoting.
2
Specification
Fortran Interface
Integer, Intent (In) 
:: 
n 
Integer, Intent (Inout) 
:: 
ifail 
Integer, Intent (Out) 
:: 
ipiv(n) 
Real (Kind=nag_wp), Intent (In) 
:: 
lambda, tol 
Real (Kind=nag_wp), Intent (Inout) 
:: 
a(n), b(n), c(n) 
Real (Kind=nag_wp), Intent (Out) 
:: 
d(n) 

C Header Interface
#include <nag.h>
void 
f01lef_ (const Integer *n, double a[], const double *lambda, double b[], double c[], const double *tol, double d[], Integer ipiv[], Integer *ifail) 

C++ Header Interface
#include <nag.h> extern "C" {
void 
f01lef_ (const Integer &n, double a[], const double &lambda, double b[], double c[], const double &tol, double d[], Integer ipiv[], Integer &ifail) 
}

The routine may be called by the names f01lef or nagf_matop_real_gen_tridiag_lu.
3
Description
The matrix
$T\lambda I$, where
$T$ is a real
$n$ by
$n$ tridiagonal matrix, is factorized as
where
$P$ is a permutation matrix,
$L$ is a unit lower triangular matrix with at most one nonzero subdiagonal element per column, and
$U$ is an upper triangular matrix with at most two nonzero superdiagonal elements per column.
The factorization is obtained by Gaussian elimination with partial pivoting and implicit row scaling.
An indication of whether or not the matrix
$T\lambda I$ is nearly singular is returned in the
$n$th element of the array
ipiv. If it is important that
$T\lambda I$ is nonsingular, as is usually the case when solving a system of tridiagonal equations, then it is strongly recommended that
${\mathbf{ipiv}}\left(n\right)$ is inspected on return from
f01lef. (See the argument
ipiv and
Section 9 for further details.)
The argument
$\lambda $ is included in the routine so that
f01lef may be used, in conjunction with
f04lef, to obtain eigenvectors of
$T$ by inverse iteration.
4
References
Wilkinson J H (1965) The Algebraic Eigenvalue Problem Oxford University Press, Oxford
Wilkinson J H and Reinsch C (1971) Handbook for Automatic Computation II, Linear Algebra Springer–Verlag
5
Arguments

1:
$\mathbf{n}$ – Integer
Input

On entry: $n$, the order of the matrix $T$.
Constraint:
${\mathbf{n}}\ge 1$.

2:
$\mathbf{a}\left({\mathbf{n}}\right)$ – Real (Kind=nag_wp) array
Input/Output

On entry: the diagonal elements of $T$.
On exit: the diagonal elements of the upper triangular matrix $U$.

3:
$\mathbf{lambda}$ – Real (Kind=nag_wp)
Input

On entry: the scalar $\lambda $. f01lef factorizes $T\lambda I$.

4:
$\mathbf{b}\left({\mathbf{n}}\right)$ – Real (Kind=nag_wp) array
Input/Output

On entry: the superdiagonal elements of $T$, stored in ${\mathbf{b}}\left(2\right)$ to ${\mathbf{b}}\left(n\right)$; ${\mathbf{b}}\left(1\right)$ is not used.
On exit: the elements of the first superdiagonal of $U$, stored in ${\mathbf{b}}\left(2\right)$ to ${\mathbf{b}}\left(n\right)$.

5:
$\mathbf{c}\left({\mathbf{n}}\right)$ – Real (Kind=nag_wp) array
Input/Output

On entry: the subdiagonal elements of $T$, stored in ${\mathbf{c}}\left(2\right)$ to ${\mathbf{c}}\left(n\right)$; ${\mathbf{c}}\left(1\right)$ is not used.
On exit: the subdiagonal elements of $L$, stored in ${\mathbf{c}}\left(2\right)$ to ${\mathbf{c}}\left(n\right)$.

6:
$\mathbf{tol}$ – Real (Kind=nag_wp)
Input

On entry: a relative tolerance used to indicate whether or not the matrix (
$T\lambda I$) is nearly singular.
tol should normally be chosen as approximately the largest relative error in the elements of
$T$. For example, if the elements of
$T$ are correct to about
$4$ significant figures, then
tol should be set to about
$5\times {10}^{4}$. See
Section 9 for further details on how
tol is used. If
tol is supplied as less than
$\epsilon $, where
$\epsilon $ is the
machine precision, then the value
$\epsilon $ is used in place of
tol.

7:
$\mathbf{d}\left({\mathbf{n}}\right)$ – Real (Kind=nag_wp) array
Output

On exit: the elements of the second superdiagonal of $U$, stored in ${\mathbf{d}}\left(3\right)$ to ${\mathbf{d}}\left(n\right)$; ${\mathbf{d}}\left(1\right)$ and ${\mathbf{d}}\left(2\right)$ are not used.

8:
$\mathbf{ipiv}\left({\mathbf{n}}\right)$ – Integer array
Output

On exit: details of the permutation matrix
$P$. If an interchange occurred at the
$k$th step of the elimination, then
${\mathbf{ipiv}}\left(k\right)=1$, otherwise
${\mathbf{ipiv}}\left(k\right)=0$. If a diagonal element of
$U$ is small, indicating that
$\left(T\lambda I\right)$ is nearly singular, then the element
${\mathbf{ipiv}}\left(n\right)$ is returned as positive. Otherwise
${\mathbf{ipiv}}\left(n\right)$ is returned as
$0$. See
Section 9 for further details. If the application is such that it is important that
$\left(T\lambda I\right)$ is not nearly singular, then it is strongly recommended that
${\mathbf{ipiv}}\left(n\right)$ is inspected on return.

9:
$\mathbf{ifail}$ – Integer
Input/Output

On entry:
ifail must be set to
$0$,
$1$ or
$1$ to set behaviour on detection of an error; these values have no effect when no error is detected.
A value of $0$ causes the printing of an error message and program execution will be halted; otherwise program execution continues. A value of $1$ means that an error message is printed while a value of $1$ means that it is not.
If halting is not appropriate, the value
$1$ or
$1$ is recommended. If message printing is undesirable, then the value
$1$ is recommended. Otherwise, the value
$0$ is recommended.
When the value $\mathbf{1}$ or $\mathbf{1}$ is used it is essential to test the value of ifail on exit.
On exit:
${\mathbf{ifail}}={\mathbf{0}}$ unless the routine detects an error or a warning has been flagged (see
Section 6).
6
Error Indicators and Warnings
If on entry
${\mathbf{ifail}}=0$ or
$1$, explanatory error messages are output on the current error message unit (as defined by
x04aaf).
Errors or warnings detected by the routine:
 ${\mathbf{ifail}}=1$

On entry, ${\mathbf{n}}=\u2329\mathit{\text{value}}\u232a$.
Constraint: ${\mathbf{n}}\ge 1$.
 ${\mathbf{ifail}}=99$
An unexpected error has been triggered by this routine. Please
contact
NAG.
See
Section 7 in the Introduction to the NAG Library FL Interface for further information.
 ${\mathbf{ifail}}=399$
Your licence key may have expired or may not have been installed correctly.
See
Section 8 in the Introduction to the NAG Library FL Interface for further information.
 ${\mathbf{ifail}}=999$
Dynamic memory allocation failed.
See
Section 9 in the Introduction to the NAG Library FL Interface for further information.
7
Accuracy
The computed factorization will satisfy the equation
where
where
$\epsilon $ is the
machine precision.
8
Parallelism and Performance
f01lef makes calls to BLAS and/or LAPACK routines, which may be threaded within the vendor library used by this implementation. Consult the documentation for the vendor library for further information.
Please consult the
X06 Chapter Introduction for information on how to control and interrogate the OpenMP environment used within this routine. Please also consult the
Users' Note for your implementation for any additional implementationspecific information.
The time taken by f01lef is approximately proportional to $n$.
The factorization of a tridiagonal matrix proceeds in $\left(n1\right)$ steps, each step eliminating one subdiagonal element of the tridiagonal matrix. In order to avoid small pivot elements and to prevent growth in the size of the elements of $L$, rows $k$ and ($k+1$) will, if necessary, be interchanged at the $k$th step prior to the elimination.
The element
${\mathbf{ipiv}}\left(n\right)$ returns the smallest integer,
$j$, for which
where
${\Vert {\left(T\lambda I\right)}_{j}\Vert}_{1}$ denotes the sum of the absolute values of the
$j$th row of the matrix (
$T\lambda I$). If no such
$j$ exists, then
${\mathbf{ipiv}}\left(n\right)$ is returned as zero. If such a
$j$ exists, then
$\left{u}_{jj}\right$ is small and hence (
$T\lambda I$) is singular or nearly singular.
This routine may be followed by
f04lef to solve systems of tridiagonal equations. If you wish to solve single systems of tridiagonal equations you should be aware of
f07caf, which solves tridiagonal systems with a single call.
f07caf requires less storage and will generally be faster than the combination of
f01lef and
f04lef, but no test for near singularity is included in
f07caf and so it should only be used when the equations are known to be nonsingular.
10
Example
This example factorizes the tridiagonal matrix
$T$ where
and then to solve the equations
$Tx=y$, where
by a call to
f04lef. The example program sets
${\mathbf{tol}}=5\times {10}^{5}$ and, of course, sets
${\mathbf{lambda}}=0$.
10.1
Program Text
10.2
Program Data
10.3
Program Results