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_symm_posdef_inv (f01ab)

Purpose

nag_matop_real_symm_posdef_inv (f01ab) calculates the accurate inverse of a real symmetric positive definite matrix, using a Cholesky factorization and iterative refinement.

Syntax

[a, b, ifail] = f01ab(a, 'n', n)
[a, b, ifail] = nag_matop_real_symm_posdef_inv(a, 'n', n)

Description

To compute the inverse XX of a real symmetric positive definite matrix AA, nag_matop_real_symm_posdef_inv (f01ab) first computes a Cholesky factorization of AA as A = LLTA=LLT, where LL is lower triangular. An approximation to XX is found by computing L1L-1 and then the product LTL1L-TL-1. The residual matrix R = IAXR=I-AX is calculated using additional precision, and a correction DD to XX is found by solving LLTD = RLLTD=R. XX is replaced by X + DX+D, and this iterative refinement of the inverse is repeated until full machine accuracy has been obtained.

References

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

Parameters

Compulsory Input Parameters

1:     a(lda,n) – double array
lda, the first dimension of the array, must satisfy the constraint ldan + 1ldan+1.
The upper triangle of the nn by nn positive definite symmetric matrix AA. The elements of the array below the diagonal need not be set.

Optional Input Parameters

1:     n – int64int32nag_int scalar
Default: The second dimension of the array a.
nn, the order of the matrix AA.
Constraint: n1n1.

Input Parameters Omitted from the MATLAB Interface

lda ldb z

Output Parameters

1:     a(lda,n) – double array
ldan + 1ldan+1.
The lower triangle of the inverse matrix XX is stored in the elements of the array below the diagonal, in rows 22 to n + 1n+1; xijxij is stored in a(i + 1,j)ai+1j for ijij. The upper triangle of the original matrix is unchanged.
2:     b(ldb,n) – double array
ldbnldbn.
The lower triangle of the inverse matrix XX, with xijxij stored in b(i,j)bij, for ijij.
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:
  ifail = 1ifail=1
The matrix AA is not positive definite, possibly due to rounding errors.
  ifail = 2ifail=2
The refinement process fails to converge, i.e., the matrix AA is ill-conditioned.
  ifail = 3ifail=3
n < 1 n<1 , or lda < n + 1 lda<n+1 , or ldb < n ldb<n .

Accuracy

The computed inverse should be correct to full machine accuracy. For a detailed error analysis see page 40 of Wilkinson and Reinsch (1971).

Further Comments

The time taken by nag_matop_real_symm_posdef_inv (f01ab) is approximately proportional to n3n3.

Example

function nag_matop_real_symm_posdef_inv_example
a = [5, 7, 6, 5;
     7, 10, 8, 7;
     6, 8, 10, 9;
     5, 7, 9, 10;
     0, -0.01238964498066066, 0, 0];
[aOut, b, ifail] = nag_matop_real_symm_posdef_inv(a)
 

aOut =

     5     7     6     5
    68    10     8     7
   -41    25    10     9
   -17    10     5    10
    10    -6    -3     2


b =

    68     0     0     0
   -41    25     0     0
   -17    10     5     0
    10    -6    -3     2


ifail =

                    0


function f01ab_example
a = [5, 7, 6, 5;
     7, 10, 8, 7;
     6, 8, 10, 9;
     5, 7, 9, 10;
     0, -0.01238964498066066, 0, 0];
[aOut, b, ifail] = f01ab(a)
 

aOut =

     5     7     6     5
    68    10     8     7
   -41    25    10     9
   -17    10     5    10
    10    -6    -3     2


b =

    68     0     0     0
   -41    25     0     0
   -17    10     5     0
    10    -6    -3     2


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