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_sparse_direct_real_gen_solve (f11mf)

## Purpose

nag_sparse_direct_real_gen_solve (f11mf) solves a real sparse system of linear equations with multiple right-hand sides given an LU $LU$ factorization of the sparse matrix computed by nag_sparse_direct_real_gen_lu (f11me).

## Syntax

[b, ifail] = f11mf(trans, iprm, il, lval, iu, uval, b, 'n', n, 'nrhs_p', nrhs_p)
[b, ifail] = nag_sparse_direct_real_gen_solve(trans, iprm, il, lval, iu, uval, b, 'n', n, 'nrhs_p', nrhs_p)

## Description

nag_sparse_direct_real_gen_solve (f11mf) solves a real system of linear equations with multiple right-hand sides AX = B$AX=B$ or ATX = B${A}^{\mathrm{T}}X=B$, according to the value of the parameter trans, where the matrix factorization Pr A Pc = LU ${P}_{r}A{P}_{c}=LU$ corresponds to an LU $LU$ decomposition of a sparse matrix stored in compressed column (Harwell–Boeing) format, as computed by nag_sparse_direct_real_gen_lu (f11me).
In the above decomposition L$L$ is a lower triangular sparse matrix with unit diagonal elements and U$U$ is an upper triangular sparse matrix; Pr${P}_{r}$ and Pc${P}_{c}$ are permutation matrices.

None.

## Parameters

### Compulsory Input Parameters

1:     trans – string (length ≥ 1)
Specifies whether AX = B$AX=B$ or ATX = B${A}^{\mathrm{T}}X=B$ is solved.
trans = 'N'${\mathbf{trans}}=\text{'N'}$
AX = B$AX=B$ is solved.
trans = 'T'${\mathbf{trans}}=\text{'T'}$
ATX = B${A}^{\mathrm{T}}X=B$ is solved.
Constraint: trans = 'N'${\mathbf{trans}}=\text{'N'}$ or 'T'$\text{'T'}$.
2:     iprm(7 × n$7×{\mathbf{n}}$) – int64int32nag_int array
The column permutation which defines Pc${P}_{c}$, the row permutation which defines Pr${P}_{r}$, plus associated data structures as computed by nag_sparse_direct_real_gen_lu (f11me).
3:     il( : $:$) – int64int32nag_int array
Note: the dimension of the array il must be at least as large as the dimension of the array of the same name in nag_sparse_direct_real_gen_lu (f11me).
Records the sparsity pattern of matrix L$L$ as computed by nag_sparse_direct_real_gen_lu (f11me).
4:     lval( : $:$) – double array
Note: the dimension of the array lval must be at least as large as the dimension of the array of the same name in nag_sparse_direct_real_gen_lu (f11me).
Records the nonzero values of matrix L$L$ and some nonzero values of matrix U$U$ as computed by nag_sparse_direct_real_gen_lu (f11me).
5:     iu( : $:$) – int64int32nag_int array
Note: the dimension of the array iu must be at least as large as the dimension of the array of the same name in nag_sparse_direct_real_gen_lu (f11me).
Records the sparsity pattern of matrix U$U$ as computed by nag_sparse_direct_real_gen_lu (f11me).
6:     uval( : $:$) – double array
Note: the dimension of the array uval must be at least as large as the dimension of the array of the same name in nag_sparse_direct_real_gen_lu (f11me).
Records some nonzero values of matrix U$U$ as computed by nag_sparse_direct_real_gen_lu (f11me).
7:     b(ldb, : $:$) – double array
The first dimension of the array b must be at least max (1,n)$\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$
The second dimension of the array must be at least max (1,nrhs)$\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{nrhs}}\right)$
The n${\mathbf{n}}$ by nrhs${\mathbf{nrhs}}$ right-hand side matrix B$B$.

### Optional Input Parameters

1:     n – int64int32nag_int scalar
Default: The first dimension of the array b.
n$n$, the order of the matrix A$A$.
Constraint: n0${\mathbf{n}}\ge 0$.
2:     nrhs_p – int64int32nag_int scalar
Default: The second dimension of the array b.
nrhs$\mathit{nrhs}$, the number of right-hand sides in B$B$.
Constraint: nrhs0${\mathbf{nrhs}}\ge 0$.

ldb

### Output Parameters

1:     b(ldb, : $:$) – double array
The first dimension of the array b will be max (1,n)$\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$
The second dimension of the array will be max (1,nrhs)$\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{nrhs}}\right)$
ldbmax (1,n)$\mathit{ldb}\ge \mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$.
The n${\mathbf{n}}$ by nrhs${\mathbf{nrhs}}$ solution matrix X$X$.
2:     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$
 On entry, trans ≠ 'N'${\mathbf{trans}}\ne \text{'N'}$ or 'T'$\text{'T'}$, or n < 0${\mathbf{n}}<0$, or nrhs < 0${\mathbf{nrhs}}<0$, or ldb < max (1,n)$\mathit{ldb}<\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{n}}\right)$.
ifail = 2${\mathbf{ifail}}=2$
Ill-defined row permutation in array iprm${\mathbf{iprm}}$. Internal checks have revealed that the iprm${\mathbf{iprm}}$ array is corrupted.
ifail = 3${\mathbf{ifail}}=3$
Ill-defined column permutations in array iprm${\mathbf{iprm}}$. Internal checks have revealed that the iprm${\mathbf{iprm}}$ array is corrupted.
ifail = 301${\mathbf{ifail}}=301$
Unable to allocate required internal workspace.

## Accuracy

For each right-hand side vector b$b$, the computed solution x$x$ is the exact solution of a perturbed system of equations (A + E)x = b$\left(A+E\right)x=b$, where
 |E| ≤ c(n)ε|L||U|, $|E|≤c(n)ε|L||U|,$
c(n)$c\left(n\right)$ is a modest linear function of n$n$, and ε$\epsilon$ is the machine precision, when partial pivoting is used.
If $\stackrel{^}{x}$ is the true solution, then the computed solution x$x$ satisfies a forward error bound of the form
 ( ‖x − x̂‖∞ )/( ‖x‖∞ ) ≤ c(n)cond(A,x)ε $‖x-x^‖∞ ‖x‖∞ ≤c(n)cond(A,x)ε$
where cond(A,x) = |A1||A||x| / xcond(A) = |A1||A|κ (A)$\mathrm{cond}\left(A,x\right)={‖|{A}^{-1}||A||x|‖}_{\infty }/{‖x‖}_{\infty }\le \mathrm{cond}\left(A\right)={‖|{A}^{-1}||A|‖}_{\infty }\le {\kappa }_{\infty }\left(A\right)$. Note that cond(A,x)$\mathrm{cond}\left(A,x\right)$ can be much smaller than cond(A)$\mathrm{cond}\left(A\right)$, and cond(AT)$\mathrm{cond}\left({A}^{\mathrm{T}}\right)$ can be much larger (or smaller) than cond(A)$\mathrm{cond}\left(A\right)$.
Forward and backward error bounds can be computed by calling nag_sparse_direct_real_gen_refine (f11mh), and an estimate for κ(A)${\kappa }_{\infty }\left(A\right)$ can be obtained by calling nag_sparse_direct_real_gen_cond (f11mg).

nag_sparse_direct_real_gen_solve (f11mf) may be followed by a call to nag_sparse_direct_real_gen_refine (f11mh) to refine the solution and return an error estimate.

## Example

```function nag_sparse_direct_real_gen_solve_example
n = int64(5);
nz = int64(11);
icolzp = [int64(1); 3; 5; 7; 9; 12];
irowix = [int64(1); 3; 1; 5; 2; 3; 2; 4; 3; 4; 5];
a = [2; 4; 1; -2; 1; 1; -1; 1; 1; 2; 3];
b = [  1.56,  3.12;
-0.25, -0.50;
3.60,  7.20;
1.33,  2.66;
0.52,  1.04];
iprm = zeros(1, 7*n, 'int64');
spec = 'M';
thresh = 1;
nzlmx = int64(8*nz);
nzlumx = int64(8*nz);
nzumx = int64(8*nz);
trans = 'N';

% Calculate COLAMD permutation
[iprm, ifail] = nag_sparse_direct_real_gen_setup(spec, n, icolzp, irowix, iprm);

% Factorise
[iprm, nzlumx, il, lval, iu, uval, nnzl, nnzu, flop, ifail] = ...
nag_sparse_direct_real_gen_lu(n, irowix, a, iprm, thresh, nzlmx, nzlumx, nzumx);

% Solve
[b, ifail] = nag_sparse_direct_real_gen_solve(trans, iprm, il, lval, iu, uval, b)
```
```

b =

0.7000    1.4000
0.1600    0.3200
0.5200    1.0400
0.7700    1.5400
0.2800    0.5600

ifail =

0

```
```function f11mf_example
n = int64(5);
nz = int64(11);
icolzp = [int64(1); 3; 5; 7; 9; 12];
irowix = [int64(1); 3; 1; 5; 2; 3; 2; 4; 3; 4; 5];
a = [2; 4; 1; -2; 1; 1; -1; 1; 1; 2; 3];
b = [  1.56,  3.12;
-0.25, -0.50;
3.60,  7.20;
1.33,  2.66;
0.52,  1.04];
iprm = zeros(1, 7*n, 'int64');
spec = 'M';
thresh = 1;
nzlmx = int64(8*nz);
nzlumx = int64(8*nz);
nzumx = int64(8*nz);
trans = 'N';

% Calculate COLAMD permutation
[iprm, ifail] = f11md(spec, n, icolzp, irowix, iprm);

% Factorise
[iprm, nzlumx, il, lval, iu, uval, nnzl, nnzu, flop, ifail] = ...
f11me(n, irowix, a, iprm, thresh, nzlmx, nzlumx, nzumx);

% Solve
[b, ifail] = f11mf(trans, iprm, il, lval, iu, uval, b)
```
```

b =

0.7000    1.4000
0.1600    0.3200
0.5200    1.0400
0.7700    1.5400
0.2800    0.5600

ifail =

0

```