Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

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 X$X$ of a real symmetric positive definite matrix A$A$, nag_matop_real_symm_posdef_inv (f01ab) first computes a Cholesky factorization of A$A$ as A = LLT$A=L{L}^{\mathrm{T}}$, where L$L$ is lower triangular. An approximation to X$X$ is found by computing L1${L}^{-1}$ and then the product LTL1${L}^{-\mathrm{T}}{L}^{-1}$. The residual matrix R = IAX$R=I-AX$ is calculated using additional precision, and a correction D$D$ to X$X$ is found by solving LLTD = R$L{L}^{\mathrm{T}}D=R$. X$X$ is replaced by X + D$X+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 + 1$\mathit{lda}\ge {\mathbf{n}}+1$.
The upper triangle of the n$n$ by n$n$ positive definite symmetric matrix A$A$. 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.
n$n$, the order of the matrix A$A$.
Constraint: n1${\mathbf{n}}\ge 1$.

lda ldb z

### Output Parameters

1:     a(lda,n) – double array
ldan + 1$\mathit{lda}\ge {\mathbf{n}}+1$.
The lower triangle of the inverse matrix X$X$ is stored in the elements of the array below the diagonal, in rows 2$2$ to n + 1$n+1$; xij${x}_{ij}$ is stored in a(i + 1,j)${\mathbf{a}}\left(i+1,j\right)$ for ij$i\ge j$. The upper triangle of the original matrix is unchanged.
2:     b(ldb,n) – double array
ldbn$\mathit{ldb}\ge {\mathbf{n}}$.
The lower triangle of the inverse matrix X$X$, with xij${x}_{ij}$ stored in b(i,j)${\mathbf{b}}\left(i,j\right)$, for ij$i\ge j$.
3:     ifail – int64int32nag_int scalar
${\mathrm{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:
ifail = 1${\mathbf{ifail}}=1$
The matrix A$A$ is not positive definite, possibly due to rounding errors.
ifail = 2${\mathbf{ifail}}=2$
The refinement process fails to converge, i.e., the matrix A$A$ is ill-conditioned.
ifail = 3${\mathbf{ifail}}=3$
n < 1 ${\mathbf{n}}<1$, or lda < n + 1 $\mathit{lda}<{\mathbf{n}}+1$, or ldb < n $\mathit{ldb}<{\mathbf{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).

The time taken by nag_matop_real_symm_posdef_inv (f01ab) is approximately proportional to n3${n}^{3}$.

## 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

```