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_linsys_real_blkdiag_fac_solve (f04lh)

Purpose

nag_linsys_real_blkdiag_fac_solve (f04lh) calculates the approximate solution of a set of real linear equations with multiple right-hand sides, $AX=B$ or ${A}^{\mathrm{T}}X=B$, where $A$ is an almost block-diagonal matrix which has been factorized by nag_matop_real_gen_blkdiag_lu (f01lh).

Syntax

[b, ifail] = f04lh(trans, blkstr, a, pivot, b, 'n', n, 'nbloks', nbloks, 'lena', lena, 'ir', ir)
[b, ifail] = nag_linsys_real_blkdiag_fac_solve(trans, blkstr, a, pivot, b, 'n', n, 'nbloks', nbloks, 'lena', lena, 'ir', ir)

Description

nag_linsys_real_blkdiag_fac_solve (f04lh) solves a set of real linear equations $AX=B$ or ${A}^{\mathrm{T}}X=B$, where $A$ is almost block-diagonal. $A$ must first be factorized by nag_matop_real_gen_blkdiag_lu (f01lh). nag_linsys_real_blkdiag_fac_solve (f04lh) then computes $X$ by forward and backward substitution over the blocks.

References

Diaz J C, Fairweather G and Keast P (1983) Fortran packages for solving certain almost block diagonal linear systems by modified alternate row and column elimination ACM Trans. Math. Software 9 358–375

Parameters

Compulsory Input Parameters

1:     $\mathrm{trans}$ – string (length ≥ 1)
Specifies the equations to be solved.
${\mathbf{trans}}=\text{'N'}$
Solve $AX=B$.
${\mathbf{trans}}=\text{'T'}$
Solve ${A}^{\mathrm{T}}X=B$.
Constraint: ${\mathbf{trans}}=\text{'N'}$ or $\text{'T'}$.
2:     $\mathrm{blkstr}\left(3,{\mathbf{nbloks}}\right)$int64int32nag_int array
Information which describes the block structure of $A$, as supplied to nag_linsys_real_blkdiag_fac_solve (f04lh).
3:     $\mathrm{a}\left({\mathbf{lena}}\right)$ – double array
The elements in the factorization of $A$, as returned by nag_linsys_real_blkdiag_fac_solve (f04lh).
4:     $\mathrm{pivot}\left({\mathbf{n}}\right)$int64int32nag_int array
Details of the interchanges in the factorization, as returned by nag_linsys_real_blkdiag_fac_solve (f04lh).
5:     $\mathrm{b}\left(\mathit{ldb},{\mathbf{ir}}\right)$ – double array
ldb, the first dimension of the array, must satisfy the constraint $\mathit{ldb}\ge {\mathbf{n}}$.
The $n$ by $r$ right-hand side matrix $B$.

Optional Input Parameters

1:     $\mathrm{n}$int64int32nag_int scalar
Default: the dimension of the array pivot and the first dimension of the array b. (An error is raised if these dimensions are not equal.)
$n$, the order of the matrix $A$.
Constraint: ${\mathbf{n}}>0$.
2:     $\mathrm{nbloks}$int64int32nag_int scalar
Default: the dimension of the array blkstr.
The total number of blocks of the matrix $A$, as supplied to nag_linsys_real_blkdiag_fac_solve (f04lh).
Constraint: $0<{\mathbf{nbloks}}\le {\mathbf{n}}$.
3:     $\mathrm{lena}$int64int32nag_int scalar
Default: the dimension of the array a.
The dimension of the array a.
Constraint: ${\mathbf{lena}}\ge \sum _{k=1}^{{\mathbf{nbloks}}}{\mathbf{blkstr}}\left(1,k\right)×{\mathbf{blkstr}}\left(2,k\right)$.
4:     $\mathrm{ir}$int64int32nag_int scalar
Default: the second dimension of the array b.
$r$, the number of right-hand sides.
Constraint: ${\mathbf{ir}}>0$.

Output Parameters

1:     $\mathrm{b}\left(\mathit{ldb},{\mathbf{ir}}\right)$ – double array
b stores the $n$ by $r$ solution matrix $X$.
2:     $\mathrm{ifail}$int64int32nag_int scalar
${\mathbf{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:
${\mathbf{ifail}}=1$
 On entry, ${\mathbf{n}}<1$, or ${\mathbf{nbloks}}<1$, or ${\mathbf{ir}}<1$, or $\mathit{ldb}<{\mathbf{n}}$, or ${\mathbf{n}}<{\mathbf{nbloks}}$, or lena is too small, or illegal values detected in blkstr, or ${\mathbf{trans}}\ne \text{'N'}$ or $\text{'T'}$.
${\mathbf{ifail}}=-99$
${\mathbf{ifail}}=-399$
Your licence key may have expired or may not have been installed correctly.
${\mathbf{ifail}}=-999$
Dynamic memory allocation failed.

Accuracy

The accuracy of the computed solution depends on the conditioning of the original matrix $A$.

None.

Example

This example solves the set of linear equations $Ax=b$ where
 $A= -1.00 -0.98 -0.79 -0.15 -1.00 -0.25 -0.87 0.35 0.78 0.31 -0.85 0.89 -0.69 -0.98 -0.76 -0.82 0.12 -0.01 0.75 0.32 -1.00 -0.53 -0.83 -0.98 -0.58 0.04 0.87 0.38 -1.00 -0.21 -0.93 -0.84 0.37 -0.94 -0.96 -1.00 -0.99 -0.91 -0.28 0.90 0.78 -0.93 -0.76 0.48 -0.87 -0.14 -1.00 -0.59 -0.99 0.21 -0.73 -0.48 -0.93 -0.91 0.10 -0.89 -0.68 -0.09 -0.58 -0.21 0.85 -0.39 0.79 -0.71 0.39 -0.99 -0.12 -0.75 0.17 -1.37 1.29 -1.59 1.10 -1.63 -1.01 -0.27 0.08 0.61 0.54 -0.41 0.16 -0.46 -0.67 0.56 -0.99 0.16 -0.16 0.98 -0.24 -0.41 0.40 -0.93 0.70 0.43 0.71 -0.97 -0.60 -0.30 0.18 -0.47 -0.98 -0.73 0.07 0.04 -0.25 -0.92 -0.52 -0.46 -0.58 0.89 -0.94 -0.54 -1.00 -0.36$
and
 $b= -2.92 -1.17 -1.30 -1.17 -2.10 -4.51 -1.71 -4.59 -4.19 -0.93 -3.31 0.52 -0.12 -0.05 -0.98 -2.07 -2.73 -1.95$
The exact solution is
 $x=1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1T.$
```function f04lh_example

fprintf('f04lh example results\n\n');

% Block structure of A
n = int64(18);
blkstr = [int64(2),4,5,3,4;
4, 7,8,6,5;
3, 4,2,3,0];
a1 = [-1.00  -0.98  -0.79  -0.15;
-1.00   0.25  -0.87   0.35];
a2 = [ 0.78   0.31  -0.85   0.89  -0.69  -0.98 -0.76;
-0.82   0.12  -0.01   0.75   0.32  -1.00 -0.53;
-0.83  -0.98  -0.58   0.04   0.87   0.38 -1.00;
-0.21  -0.93  -0.84   0.37  -0.94  -0.96 -1.00];
a3 = [-0.99  -0.91  -0.28   0.90   0.78  -0.93  -0.76   0.48;
-0.87  -0.14  -1.00  -0.59  -0.99   0.21  -0.73  -0.48;
-0.93  -0.91   0.10  -0.89  -0.68  -0.09  -0.58  -0.21;
0.85  -0.39   0.79  -0.71   0.39  -0.99  -0.12  -0.75;
0.17  -1.37   1.29  -1.59   1.10  -1.63  -1.01  -0.27];
a4 = [ 0.08   0.61   0.54  -0.41   0.16  -0.46;
-0.67   0.56  -0.99   0.16  -0.16   0.98;
-0.24  -0.41   0.40  -0.93   0.70   0.43];
a5 = [ 0.71  -0.97  -0.60  -0.30   0.18;
-0.47  -0.98  -0.73   0.07   0.04;
-0.25  -0.92  -0.52  -0.46  -0.58;
0.89  -0.94  -0.54  -1.00  -0.36];
% Flatten A
a = [reshape(a1,[ 8,1]);
reshape(a2,[28,1]);
reshape(a3,[40,1]);
reshape(a4,[18,1]);
reshape(a5,[20,1])];

% Right hand side
b  = [-2.92;  -1.27;  -1.30;  -1.17;  -2.10;  -4.51;  -1.71;  -4.59;
-4.19;  -0.93;  -3.31;   0.52;  -0.12;  -0.05;  -0.98;  -2.07;
-2.73;  -1.95];

% Factorize A
tol = 0;
[AF, pivot, tol, index, ifail] = ...
f01lh(n, blkstr, a, tol);

% Solve system
trans = 'N';
[x, ifail] = f04lh( ...
trans, blkstr, AF, pivot, b);
disp('Component solution');
disp(x);

```
```f04lh example results

Component solution
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000

```