NAG CL Interface
f01mcc (real_vband_posdef_fac)
1
Purpose
f01mcc computes the Cholesky factorization of a real symmetric positive definite variablebandwidth matrix.
2
Specification
void 
f01mcc (Integer n,
const double a[],
Integer lal,
Integer row[],
double al[],
double d[],
NagError *fail) 

The function may be called by the names: f01mcc, nag_matop_real_vband_posdef_fac or nag_real_cholesky_skyline.
3
Description
f01mcc determines the unit lower triangular matrix $L$ and the diagonal matrix $D$ in the Cholesky factorization $A=LD{L}^{\mathrm{T}}$ of a symmetric positive definite variablebandwidth matrix $A$ of order $n$. (Such a matrix is sometimes called a ‘skyline’ matrix.)
The matrix
$A$ is represented by the elements lying within the
envelope of its lower triangular part, that is, between the first nonzero of each row and the diagonal (see
Section 10 for an example). The width
${\mathbf{row}}\left[i\right]$ of the
$i$th row is the number of elements between the first nonzero element and the element on the diagonal, inclusive. Although, of course, any matrix possesses an envelope as defined, this function is primarily intended for the factorization of symmetric positive definite matrices with an
average bandwidth which is small compared with
$n$ (also see
Section 9).
The method is based on the property that during Cholesky factorization there is no fillin outside the envelope.
The determination of
$L$ and
$D$ is normally the first of two steps in the solution of the system of equations
$Ax=b$. The remaining step, viz. the solution of
${LDL}^{\mathrm{T}}x=b$ may be carried out using
f04mcc.
4
References
Jennings A (1966) A compact storage scheme for the solution of symmetric linear simultaneous equations Comput. J. 9 281–285
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 $A$.
Constraint:
${\mathbf{n}}\ge 1$.

2:
$\mathbf{a}\left[{\mathbf{lal}}\right]$ – const double
Input

On entry: the elements within the envelope of the lower triangle of the positive definite symmetric matrix
$A$, taken in row by row order. The following code assigns the matrix elements within the envelope to the correct elements of the array
k=0;
for(i=0; i<n; ++i)
for(j=irow[i]+1; j<=i; ++j)
a[k++]=matrix[i][j];
See also
Section 9

3:
$\mathbf{lal}$ – Integer
Input

On entry: the smaller of the dimensions of the arrays
a and
al as declared in the function from which
f01mcc is called.
Constraint:
${\mathbf{lal}}\ge {\mathbf{row}}\left[0\right]+{\mathbf{row}}\left[1\right]+\cdots +{\mathbf{row}}\left[n1\right]$.

4:
$\mathbf{row}\left[{\mathbf{n}}\right]$ – Integer
Input

On entry: ${\mathbf{row}}\left[i\right]$ must contain the width of row $i$ of the matrix $A$, i.e., the number of elements between the first (leftmost) nonzero element and the element on the diagonal, inclusive.
Constraint:
$1\le {\mathbf{row}}\left[\mathit{i}\right]\le \mathit{i}+1$, for $\mathit{i}=0,1,\dots ,n1$.

5:
$\mathbf{al}\left[{\mathbf{lal}}\right]$ – double
Output

On exit: the elements within the envelope of the lower triangular matrix
$L$, taken in row by row order. The envelope of
$L$ is identical to that of the lower triangle of
$A$. The unit diagonal elements of
$L$ are stored explicitly. See also
Section 9

6:
$\mathbf{d}\left[{\mathbf{n}}\right]$ – double
Output

On exit: the diagonal elements of the diagonal matrix $D$. Note that the determinant of $A$ is equal to the product of these diagonal elements. If the value of the determinant is required it should not be determined by forming the product explicitly, because of the possibility of overflow or underflow. The logarithm of the determinant may safely be formed from the sum of the logarithms of the diagonal elements.

7:
$\mathbf{fail}$ – NagError *
Input/Output

The NAG error argument (see
Section 7 in the Introduction to the NAG Library CL Interface).
6
Error Indicators and Warnings
 NE_2_INT_ARG_GT

On entry, ${\mathbf{row}}\left[\u2329\mathit{\text{value}}\u232a\right]=\u2329\mathit{\text{value}}\u232a$ while $i=\u2329\mathit{\text{value}}\u232a$. These arguments must satisfy ${\mathbf{row}}\left[i\right]\le i+1$.
 NE_2_INT_ARG_LT

On entry, ${\mathbf{lal}}=\u2329\mathit{\text{value}}\u232a$ while ${\mathbf{row}}\left[0\right]+\cdots +{\mathbf{row}}\left[n1\right]=\u2329\mathit{\text{value}}\u232a$. These arguments must satisfy ${\mathbf{lal}}\ge {\mathbf{row}}\left[0\right]+\cdots +{\mathbf{row}}\left[n1\right]$.
 NE_INT_ARG_LT

On entry, ${\mathbf{n}}=\u2329\mathit{\text{value}}\u232a$.
Constraint: ${\mathbf{n}}\ge 1$.
On entry, ${\mathbf{row}}\left[\u2329\mathit{\text{value}}\u232a\right]$ must not be less than 1: ${\mathbf{row}}\left[\u2329\mathit{\text{value}}\u232a\right]=\u2329\mathit{\text{value}}\u232a$.
 NE_NOT_POS_DEF

The matrix is not positive definite, possibly due to rounding errors.
 NE_NOT_POS_DEF_FACT

The matrix is not positive definite, possibly due to rounding errors. The factorization has been completed but may be very inaccurate.
7
Accuracy
On successful exit then the
computed $L$ and
$D$ satisfy the relation
${LDL}^{\mathrm{T}}=A+F$, where
and
where
$k$ is a constant of order unity,
$m$ is the largest value of
${\mathbf{row}}\left[i\right]$, and
$\epsilon $ is the
machine precision. See pages 25–27 and 54–55 or
Wilkinson and Reinsch (1971). If the error
NE_NOT_POS_DEF_FACT is reported then the factorization has been completed although the matrix was not positive definite. However the factorization may be very inaccurate and should be used only with great caution. For instance, if it is used to solve a set of equations
$Ax=b$ using
f04mcc, the residual vector
$bAx$ should be checked.
8
Parallelism and Performance
f01mcc is not threaded in any implementation.
The time taken by f01mcc is approximately proportional to the sum of squares of the values of ${\mathbf{row}}\left[i\right]$.
The distribution of row widths may be very nonuniform without undue loss of efficiency. Moreover, the function has been designed to be as competitive as possible in speed with functions designed for full or uniformly banded matrices, when applied to such matrices.
The function may be called with the same actual array supplied for arguments
a and
al, in which case
$L$ overwrites the lower triangle of
$A$.
10
Example
To obtain the Cholesky factorization of the symmetric matrix, whose lower triangle is
For this matrix, the elements of
row must be set to
$1$,
$2$,
$2$,
$1$,
$5$,
$3$, and the elements within the envelope must be supplied in row order as
10.1
Program Text
10.2
Program Data
10.3
Program Results