hide long namesshow long names
hide short namesshow short names
Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

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 variable-bandwidth matrix.

Syntax

[al, d, ifail] = f01mc(a, nrow, 'n', n, 'lal', lal)
[al, d, ifail] = nag_matop_real_vband_posdef_fac(a, nrow, 'n', n, 'lal', lal)

Description

nag_matop_real_vband_posdef_fac (f01mc) determines the unit lower triangular matrix LL and the diagonal matrix DD in the Cholesky factorization A = LDLTA=LDLT of a symmetric positive definite variable-bandwidth matrix AA of order nn. (Such a matrix is sometimes called a ‘sky-line’ matrix.)
The matrix AA 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 [Example] for an example). The width nrow(i)nrowi of the iith 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 nn (also see Section [Further Comments]).
The method is based on the property that during Cholesky factorization there is no fill-in outside the envelope.
The determination of LL and DD is normally the first of two steps in the solution of the system of equations Ax = bAx=b. The remaining step, viz. the solution of LDLTx = bLDLTx=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:     a(lal) – double array
lal, the dimension of the array, must satisfy the constraint lalnrow(1) + nrow(2) + + nrow(n)lalnrow1+nrow2++nrown.
The elements within the envelope of the lower triangle of the positive definite symmetric matrix AA, 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 = i-nrow(i)+1:i
     k = k + 1;
     a(k) = matrix(i,j);
   end
 end
2:     nrow(n) – int64int32nag_int array
n, the dimension of the array, must satisfy the constraint n1n1.
nrow(i)nrowi must contain the width of row ii of the matrix AA, i.e., the number of elements between the first (leftmost) nonzero element and the element on the diagonal, inclusive.
Constraint: 1nrow(i)i1nrowii, for i = 1,2,,ni=1,2,,n.

Optional Input Parameters

1:     n – int64int32nag_int scalar
Default: The dimension of the array nrow.
nn, the order of the matrix AA.
Constraint: n1n1.
2:     lal – int64int32nag_int scalar
Default: The dimension of the array a.
The dimension of the arrays a and al as declared in the (sub)program from which nag_matop_real_vband_posdef_fac (f01mc) is called.
Constraint: lalnrow(1) + nrow(2) + + nrow(n)lalnrow1+nrow2++nrown.

Input Parameters Omitted from the MATLAB Interface

None.

Output Parameters

1:     al(lal) – double array
The elements within the envelope of the lower triangular matrix LL, taken in row by row order. The envelope of LL is identical to that of the lower triangle of AA. The unit diagonal elements of LL are stored explicitly. See also Section [Further Comments].
2:     d(n) – double array
The diagonal elements of the diagonal matrix DD. Note that the determinant of AA 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:     ifail – int64int32nag_int scalar
ifail = 0ifail=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.

  ifail = 1ifail=1
On entry,n < 1n<1,
orfor some ii, nrow(i) < 1nrowi<1 or nrow(i) > inrowi>i,
orlal < nrow(1) + nrow(2) + + nrow(n)lal<nrow1+nrow2++nrown.
  ifail = 2ifail=2
AA is not positive definite, or this property has been destroyed by rounding errors. The factorization has not been completed.
W ifail = 3ifail=3
AA is not positive definite, or this property has been destroyed by rounding errors. The factorization has not been completed.

Accuracy

If ifail = 0ifail=0 on exit, then the computed LL and DD satisfy the relation LDLT = A + FLDLT=A+F, where
F2km2ε × max aii
i
F2 km2 ε × maxi aii
and
F2 km2 ε × A2 ,
F2 km2 ε × A2 ,
where kk is a constant of order unity, mm is the largest value of nrow(i)nrowi, and εε 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 nrow(i)nrowi.
The distribution of row widths may be very non-uniform 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

function nag_matop_real_vband_posdef_fac_example
a = [1;
     2;
     5;
     3;
     13;
     16;
     5;
     14;
     18;
     8;
     55;
     24;
     17;
     77];
nrow = [int64(1);2;2;1;5;3];
[al, d, ifail] = nag_matop_real_vband_posdef_fac(a, nrow)
 

al =

    1.0000
    2.0000
    1.0000
    3.0000
    1.0000
    1.0000
    5.0000
    4.0000
    1.5000
    0.5000
    1.0000
    1.5000
    5.0000
    1.0000


d =

     1
     1
     4
    16
     1
    16


ifail =

                    0


function f01mc_example
a = [1;
     2;
     5;
     3;
     13;
     16;
     5;
     14;
     18;
     8;
     55;
     24;
     17;
     77];
nrow = [int64(1);2;2;1;5;3];
[al, d, ifail] = f01mc(a, nrow)
 

al =

    1.0000
    2.0000
    1.0000
    3.0000
    1.0000
    1.0000
    5.0000
    4.0000
    1.5000
    0.5000
    1.0000
    1.5000
    5.0000
    1.0000


d =

     1
     1
     4
    16
     1
    16


ifail =

                    0



PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2013