PDF version (NAG web site
, 64bit version, 64bit version)
NAG Toolbox: nag_matop_real_vband_posdef_fac (f01mc)
Purpose
nag_matop_real_vband_posdef_fac (f01mc) computes the Cholesky factorization of a real symmetric positive definite variablebandwidth matrix.
Syntax
Description
nag_matop_real_vband_posdef_fac (f01mc) 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
Example for an example). The
width
${\mathbf{nrow}}\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
Further Comments).
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
$LD{L}^{\mathrm{T}}x=b$, may be carried out using
nag_linsys_real_posdef_vband_solve (f04mc).
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
Parameters
Compulsory Input Parameters
 1:
$\mathrm{a}\left({\mathbf{lal}}\right)$ – double array

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 = 1:n
for j = inrow(i)+1:i
k = k + 1;
a(k) = matrix(i,j);
end
end
 2:
$\mathrm{nrow}\left({\mathbf{n}}\right)$ – int64int32nag_int array

${\mathbf{nrow}}\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{nrow}}\left(\mathit{i}\right)\le \mathit{i}$, for $\mathit{i}=1,2,\dots ,n$.
Optional Input Parameters
 1:
$\mathrm{n}$ – int64int32nag_int scalar

Default:
the dimension of the array
nrow.
$n$, the order of the matrix $A$.
Constraint:
${\mathbf{n}}\ge 1$.
 2:
$\mathrm{lal}$ – int64int32nag_int scalar

Default:
the dimension of the array
a.
The dimension of the arrays
a and
al.
Constraint:
${\mathbf{lal}}\ge {\mathbf{nrow}}\left(1\right)+{\mathbf{nrow}}\left(2\right)+\dots +{\mathbf{nrow}}\left(n\right)$.
Output Parameters
 1:
$\mathrm{al}\left({\mathbf{lal}}\right)$ – double array

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
Further Comments.
 2:
$\mathrm{d}\left({\mathbf{n}}\right)$ – double array

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.
 3:
$\mathrm{ifail}$ – int64int32nag_int scalar
${\mathbf{ifail}}={\mathbf{0}}$ unless the function detects an error (see
Error Indicators and Warnings).
Error Indicators and Warnings
Errors or warnings detected by the function:
Cases prefixed with W are classified as warnings and
do not generate an error of type NAG:error_n. See nag_issue_warnings.
 ${\mathbf{ifail}}=1$

On entry,  ${\mathbf{n}}<1$, 
or  for some $i$, ${\mathbf{nrow}}\left(i\right)<1$ or ${\mathbf{nrow}}\left(i\right)>i$, 
or  ${\mathbf{lal}}<{\mathbf{nrow}}\left(1\right)+{\mathbf{nrow}}\left(2\right)+\dots +{\mathbf{nrow}}\left({\mathbf{n}}\right)$. 
 ${\mathbf{ifail}}=2$

$A$ is not positive definite, or this property has been destroyed by rounding errors. The factorization has not been completed.
 W ${\mathbf{ifail}}=3$

$A$ is not positive definite, or this property has been destroyed by rounding errors. The factorization has not been completed.
 ${\mathbf{ifail}}=99$
An unexpected error has been triggered by this routine. Please
contact
NAG.
 ${\mathbf{ifail}}=399$
Your licence key may have expired or may not have been installed correctly.
 ${\mathbf{ifail}}=999$
Dynamic memory allocation failed.
Accuracy
If
${\mathbf{ifail}}={\mathbf{0}}$ on exit, then the
computed
$L$ and
$D$ satisfy the relation
$LD{L}^{\mathrm{T}}=A+F$, where
and
where
$k$ is a constant of order unity,
$m$ is the largest value of
${\mathbf{nrow}}\left(i\right)$, and
$\epsilon $ is the
machine precision. See pages 25–27 and 54–55 of
Wilkinson and Reinsch (1971).
Further Comments
The time taken by nag_matop_real_vband_posdef_fac (f01mc) is approximately proportional to the sum of squares of the values of ${\mathbf{nrow}}\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.
Example
This example obtains the Cholesky factorization of the symmetric matrix, whose lower triangle is:
For this matrix, the elements of
nrow must be set to
$1$,
$2$,
$2$,
$1$,
$5$,
$3$, and the elements within the envelope must be supplied in row order as:
Open in the MATLAB editor:
f01mc_example
function f01mc_example
fprintf('f01mc example results\n\n');
a = [1;
2; 5;
3; 13;
16;
5; 14; 18; 8; 55;
24; 17; 77];
nrow = [int64(1); 2; 2; 1; 5; 3];
[L, D, ifail] = f01mc(a, nrow);
fprintf(' D L\n');
k = 0;
for i = 1:6
fprintf(' %6.3f ',D(i));
for j = 1:inrow(i)
fprintf('%8s',' ');
end
fprintf('%8.3f',L(k+1:k+nrow(i)));
fprintf('\n');
k = k + nrow(i);
end
f01mc example results
D L
1.000 1.000
1.000 2.000 1.000
4.000 3.000 1.000
16.000 1.000
1.000 5.000 4.000 1.500 0.500 1.000
16.000 1.500 5.000 1.000
PDF version (NAG web site
, 64bit version, 64bit version)
© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2015