NAG Library Function Document

nag_real_cholesky_skyline_solve (f04mcc)

 Contents

    1  Purpose
    7  Accuracy

1
Purpose

nag_real_cholesky_skyline_solve (f04mcc) computes the approximate solution of a system of real linear equations with multiple right-hand sides, AX = B , where A  is a symmetric positive definite variable-bandwidth matrix, which has previously been factorized by nag_real_cholesky_skyline (f01mcc). Related systems may also be solved.

2
Specification

#include <nag.h>
#include <nagf04.h>
void  nag_real_cholesky_skyline_solve (Nag_SolveSystem selct, Integer n, Integer nrhs, const double al[], Integer lal, const double d[], const Integer row[], const double b[], Integer tdb, double x[], Integer tdx, NagError *fail)

3
Description

The normal use of nag_real_cholesky_skyline_solve (f04mcc) is the solution of the systems AX = B , following a call of nag_real_cholesky_skyline (f01mcc) to determine the Cholesky factorization A = LD LT  of the symmetric positive definite variable-bandwidth matrix A .
However, the function may be used to solve any one of the following systems of linear algebraic equations:
LDLT X =B ​ (usual system) (1)
LDX =B ​ (lower triangular system) (2)
DLT X =B ​ (upper triangular system) (3)
LLT X =B (4)
LX =B ​ (unit lower triangular system) (5)
LT X =B ​ (unit upper triangular system) (6)
L  denotes a unit lower triangular variable-bandwidth matrix of order n , D  a diagonal matrix of order n , and B  a set of right-hand sides.
The matrix L  is represented by the elements lying within its envelope, i.e., between the first nonzero of each row and the diagonal (see Section 10 for an example). The width row[i]  of the i th row is the number of elements between the first nonzero element and the element on the diagonal inclusive.

4
References

Wilkinson J H and Reinsch C (1971) Handbook for Automatic Computation II, Linear Algebra Springer–Verlag

5
Arguments

1:     selct Nag_SolveSystemInput
On entry: selct must specify the type of system to be solved, as follows:
  • if selct=Nag_LDLTX: solve LDL TX = B ;
  • if selct=Nag_LDX: solve LDX = B ;
  • if selct=Nag_DLTX: solve DLT X = B ;
  • if selct=Nag_LLTX: solve LLT X = B ;
  • if selct=Nag_LX: solve LX = B ;
  • if selct=Nag_LTX: solve LT X = B .
Constraint: selct=Nag_LDLTX, Nag_LDX, Nag_DLTX, Nag_LLTX, Nag_LX or Nag_LTX.
2:     n IntegerInput
On entry: n , the order of the matrix L .
Constraint: n1 .
3:     nrhs IntegerInput
On entry: r , the number of right-hand sides.
Constraint: nrhs1 .
4:     al[lal] const doubleInput
On entry: the elements within the envelope of the lower triangular matrix L , taken in row by row order, as returned by nag_real_cholesky_skyline (f01mcc). The unit diagonal elements of L  must be stored explicitly.
5:     lal IntegerInput
On entry: the dimension of the array al.
Constraint: lal row[0] + row[1] + + row[n-1] .
6:     d[n] const doubleInput
On entry: the diagonal elements of the diagonal matrix D . d is not referenced if selct=Nag_LLTX, Nag_LX or Nag_LTX
7:     row[n] const IntegerInput
On entry: row[i]  must contain the width of row i  of L , i.e., the number of elements between the first (left-most) nonzero element and the element on the diagonal, inclusive.
Constraint: 1 row[i] i + 1  for i = 0 , 1 , , n - 1 .
8:     b[n×tdb] const doubleInput
Note: the i,jth element of the matrix B is stored in b[i-1×tdb+j-1].
On entry: the n  by r  right-hand side matrix B . See also Section 9.
9:     tdb IntegerInput
On entry: the stride separating matrix column elements in the array b.
Constraint: tdbnrhs .
10:   x[n×tdx] doubleOutput
Note: the i,jth element of the matrix X is stored in x[i-1×tdx+j-1].
On exit: the n  by r  solution matrix X . See also Section 9.
11:   tdx IntegerInput
On entry: the stride separating matrix column elements in the array x.
Constraint: tdxnrhs .
12:   fail NagError *Input/Output
The NAG error argument (see Section 3.7 in How to Use the NAG Library and its Documentation).

6
Error Indicators and Warnings

NE_2_INT_ARG_GT
On entry, row[i] = value while i=value . These arguments must satisfy row[i] i + 1 .
NE_2_INT_ARG_LT
On entry, lal=value  while row[0] + + row[n-1] = value. These arguments must satisfy lal row[0] + + row[n-1] .
On entry, tdb=value  while nrhs=value . These arguments must satisfy tdbnrhs .
On entry, tdx=value  while nrhs=value . These arguments must satisfy tdxnrhs .
NE_BAD_PARAM
On entry, argument selct had an illegal value.
NE_INT_ARG_LT
On entry, n=value.
Constraint: n1.
On entry, nrhs=value.
Constraint: nrhs1.
On entry, row[value]  must not be less than 1: row[value] = value.
NE_NOT_UNIT_DIAG
The lower triangular matrix L  has at least one diagonal element which is not equal to unity. The first non-unit element has been located in the array al[value] .
NE_ZERO_DIAG
The diagonal matrix D  is singular as it has at least one zero element. The first zero element has been located in the array d[value] .

7
Accuracy

The usual backward error analysis of the solution of triangular system applies: each computed solution vector is exact for slightly perturbed matrices L  and D , as appropriate (see pages 25-27 and 54-55 of Wilkinson and Reinsch (1971)).

8
Parallelism and Performance

nag_real_cholesky_skyline_solve (f04mcc) is not threaded in any implementation.

9
Further Comments

The time taken by nag_real_cholesky_skyline_solve (f04mcc) is approximately proportional to pr , where p = row[0] + row[1] + + row[n-1] .
The function may be called with the same actual array supplied for the arguments b and x, in which case the solution matrix will overwrite the right-hand side matrix.

10
Example

To solve the system of equations AX = B , where
A = 1 2 0 0 5 0 2 5 3 0 14 0 0 3 13 0 18 0 0 0 0 16 8 24 5 14 18 8 55 17 0 0 0 24 17 77   and   B = 6 -10 15 -21 11 0-3 00 -24 51 -39 46 -67 .  
Here A  is symmetric and positive definite and must first be factorized by nag_real_cholesky_skyline (f01mcc).

10.1
Program Text

Program Text (f04mcce.c)

10.2
Program Data

Program Data (f04mcce.d)

10.3
Program Results

Program Results (f04mcce.r)

© The Numerical Algorithms Group Ltd, Oxford, UK. 2017