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_real_symm_basic_diag (f11gf)

Purpose

nag_sparse_real_symm_basic_diag (f11gf) is the third in a suite of three functions for the iterative solution of a symmetric system of simultaneous linear equations (see Golub and Van Loan (1996)). nag_sparse_real_symm_basic_diag (f11gf) returns information about the computations during an iteration and/or after this has been completed. The first function of the suite, nag_sparse_real_symm_basic_setup (f11gd), is a setup function, the second function, nag_sparse_real_symm_basic_solver (f11ge) is the proper iterative solver.
These three functions are suitable for the solution of large sparse symmetric systems of equations.

Syntax

[itn, stplhs, stprhs, anorm, sigmax, its, sigerr, ifail] = f11gf(work)
[itn, stplhs, stprhs, anorm, sigmax, its, sigerr, ifail] = nag_sparse_real_symm_basic_diag(work)
Note: the interface to this routine has changed since earlier releases of the toolbox:
Mark 22: lwork has been removed from the interface
.

Description

nag_sparse_real_symm_basic_diag (f11gf) returns information about the solution process. It can be called both during a monitoring step of the solver nag_sparse_real_symm_basic_solver (f11ge), or after this solver has completed its tasks. Calling nag_sparse_real_symm_basic_diag (f11gf) at any other time will result in an error condition being raised.
For further information you should read the documentation for nag_sparse_real_symm_basic_setup (f11gd) and nag_sparse_real_symm_basic_solver (f11ge).

References

Golub G H and Van Loan C F (1996) Matrix Computations (3rd Edition) Johns Hopkins University Press, Baltimore

Parameters

Compulsory Input Parameters

1:     work(lwork) – double array
lwork, the dimension of the array, must satisfy the constraint lwork120$\mathit{lwork}\ge 120$.
The array work as returned by nag_sparse_real_symm_basic_solver (f11ge) (see also Section [Description] in (f11ge)).

None.

lwork

Output Parameters

1:     itn – int64int32nag_int scalar
The number of iterations carried out by nag_sparse_real_symm_basic_solver (f11ge).
2:     stplhs – double scalar
The current value of the left-hand side of the termination criterion used by nag_sparse_real_symm_basic_solver (f11ge).
3:     stprhs – double scalar
The current value of the right-hand side of the termination criterion used by nag_sparse_real_symm_basic_solver (f11ge).
4:     anorm – double scalar
For CG and SYMMQ methods, the norm A1 = A${‖A‖}_{1}={‖A‖}_{\infty }$ when either it has been supplied to nag_sparse_real_symm_basic_setup (f11gd) or it has been estimated by nag_sparse_real_symm_basic_solver (f11ge) (see also Sections [Description] and [Parameters] in (f11gd)). Otherwise, anorm = 0.0${\mathbf{anorm}}=0.0$ is returned.
For MINRES method, an estimate of the infinity norm of the preconditioned matrix operator.
5:     sigmax – double scalar
For CG and SYMMQ methods, the current estimate of the largest singular value σ1(A)${\sigma }_{1}\left(\stackrel{-}{A}\right)$ of the preconditioned iteration matrix A = E1AET$\stackrel{-}{A}={E}^{-1}A{E}^{-\mathrm{T}}$, when either it has been supplied to nag_sparse_real_symm_basic_setup (f11gd) or it has been estimated by nag_sparse_real_symm_basic_solver (f11ge) (see also Sections [Description] and [Parameters] in (f11gd)). Note that if ${\mathbf{its}}<{\mathbf{itn}}$ then sigmax contains the final estimate. If, on final exit from nag_sparse_real_symm_basic_solver (f11ge), ${\mathbf{its}}={\mathbf{itn}}$, then the estimation of σ1(A)${\sigma }_{1}\left(\stackrel{-}{A}\right)$ may have not converged; in this case you should look at the value returned in sigerr. Otherwise, sigmax = 0.0${\mathbf{sigmax}}=0.0$ is returned.
For MINRES method, an estimate of the final transformed residual.
6:     its – int64int32nag_int scalar
For CG and SYMMQ methods, the number of iterations employed so far in the computation of the estimate of σ1(A)${\sigma }_{1}\left(\stackrel{-}{A}\right)$, the largest singular value of the preconditioned matrix A = E1AET$\stackrel{-}{A}={E}^{-1}A{E}^{-\mathrm{T}}$, when σ1(A)${\sigma }_{1}\left(\stackrel{-}{A}\right)$ has been estimated by nag_sparse_real_symm_basic_solver (f11ge) using the bisection method (see also Sections [Description], [Parameters] and [Further Comments] in (f11gd)). Otherwise, its = 0${\mathbf{its}}=0$ is returned.
7:     sigerr – double scalar
For CG and SYMMQ methods, if σ1(A)${\sigma }_{1}\left(\stackrel{-}{A}\right)$ has been estimated by nag_sparse_real_symm_basic_solver (f11ge) using bisection,
 sigerr = max ((|σ1(k) − σ1(k − 1)|)/(σ1(k)),(|σ1(k) − σ1(k − 2)|)/(σ1(k))) , $sigerr=max(|σ1(k)-σ1(k-1)|σ1(k),|σ1(k)-σ1(k-2)|σ1(k)) ,$
where k = its$k={\mathbf{its}}$ denotes the iteration number. The estimation has converged if ${\mathbf{sigerr}}\le {\mathbf{sigtol}}$ where sigtol is an input parameter to nag_sparse_real_symm_basic_setup (f11gd). Otherwise, sigerr = 0.0${\mathbf{sigerr}}=0.0$ is returned.
For MINRES method, an estimate of the condition number of the preconditioned matrix.
8:     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 = i${\mathbf{ifail}}=-i$
If ifail = i${\mathbf{ifail}}=-i$, parameter i$i$ had an illegal value on entry. The parameters are numbered as follows:
1: itn, 2: stplhs, 3: stprhs, 4: anorm, 5: sigmax, 6: its, 7: sigerr, 8: work, 9: lwork, 10: ifail.
It is possible that ifail refers to a parameter that is omitted from the MATLAB interface. This usually indicates that an error in one of the other input parameters has caused an incorrect value to be inferred.
ifail = 1${\mathbf{ifail}}=1$
nag_sparse_real_symm_basic_diag (f11gf) has been called out of sequence. For example, the last call to nag_sparse_real_symm_basic_solver (f11ge) did not return irevcm = 3${\mathbf{irevcm}}=3$ or 4$4$.

Not applicable.

None.

Example

```function nag_sparse_real_symm_basic_diag_example
nz = 16;
a = zeros(10000,1);
a(1:16) = [4;
1;
5;
2;
2;
3;
-1;
1;
4;
1;
-2;
3;
2;
-1;
-2;
5];
irow = zeros(10000,1,'int64');
irow(1:16) = [int64(1);2;2;3;4;4;5;5;5;6;6;6;7;7;7;7];
icol = zeros(10000,1,'int64');
icol(1:16) = [int64(1);1;2;3;2;4;1;4;5;2;5;6;1;2;3;7];
lfill = int64(0);
dtol = 0;
mic = 'N';
dscale = 0;
ipiv = zeros(7, 1, 'int64');
method = 'CG    ';
precon = 'P';
n = 7;
tol = 1e-06;
maxitn = int64(20);
anorm = 0;
sigmax = 0;
maxits = int64(7);
monit = int64(2);
irevcm = int64(0);
u = [0; 0; 0; 0; 0; 0; 0];
v = [15; 18; -8; 21; 11; 10; 29];
wgt = [0; 0; 0; 0; 0; 0; 0];

[a, irow, icol, ipiv, istr, nnzc, npivm, ifail] = ...
nag_sparse_real_symm_precon_ichol(int64(nz), a, irow, icol, lfill, ...
dtol, mic, dscale, ipiv);

[lwreq, work, ifail] = ...
nag_sparse_real_symm_basic_setup(method, precon, int64(n), tol, maxitn, anorm, ...
sigmax, maxits, monit, ...
'sigcmp', 'S', 'norm_p', '1');

while (irevcm ~= 4)
[irevcm, u, v, work, ifail] = ...
nag_sparse_real_symm_basic_solver(irevcm, u, v, wgt, work);
if (irevcm == 1)
[v, ifail] = ...
nag_sparse_real_symm_matvec(a(1:nz), irow(1:nz), icol(1:nz), 'N', u);
if (ifail ~= 0)
irecvm = 6;
end
elseif (irevcm == 2)
[v, ifail] = ...
nag_sparse_real_symm_precon_ichol_solve(a, irow, icol, ipiv, istr, 'N', u);
if (ifail ~= 0)
irecvm = 6;
end
elseif (irevcm == 3)
[itn, stplhs, stprhs, anorm, sigmax, its, sigerr, ifail] = ...
nag_sparse_real_symm_basic_diag(work);
if (ifail ~= 0)
irecvm = 6;
end
fprintf('\nMonitoring at iteration number %d\nresidual norm: %14.4e\n', itn, stplhs);
fprintf('\n   Solution Vector  Residual Vector\n');
for i = 1:n
fprintf('%+16.4e %+16.4e\n', u(i), v(i));
end
end
end

% Get information about the computation
[itn, stplhs, stprhs, anorm, sigmax, its, sigerr, ifail] = ...
nag_sparse_real_symm_basic_diag(work);
fprintf('\nNumber of iterations for convergence: %d\n', itn);
fprintf('Residual norm:                           %14.4e\n', stplhs);
fprintf('Right-hand side of termination criteria: %14.4e\n', stprhs);
fprintf('i-norm of matrix a:                      %14.4e\n', anorm);
fprintf('\n   Solution Vector  Residual Vector\n');
for i = 1:n
fprintf('%+16.4e %+16.4e\n', u(i), v(i));
end
```
```

Monitoring at iteration number 2
residual norm:     1.9938e+00

Solution Vector  Residual Vector
+9.6320e-01      -2.2960e-01
+1.9934e+00      +2.2254e-01
+3.0583e+00      +9.5827e-02
+4.1453e+00      -2.5155e-01
+4.8289e+00      -1.7160e-01
+5.6630e+00      +6.7533e-01
+7.1062e+00      -3.4737e-01

Monitoring at iteration number 4
residual norm:     6.6574e-03

Solution Vector  Residual Vector
+9.9940e-01      -1.0551e-03
+2.0011e+00      -2.4675e-03
+3.0008e+00      -1.7116e-05
+3.9996e+00      +4.4929e-05
+4.9991e+00      +2.1359e-03
+5.9993e+00      -8.7482e-04
+7.0007e+00      +6.2045e-05

Number of iterations for convergence: 5
Residual norm:                               7.1054e-15
Right-hand side of termination criteria:     3.9200e-04
i-norm of matrix a:                          1.0000e+01

Solution Vector  Residual Vector
+1.0000e+00      +0.0000e+00
+2.0000e+00      +0.0000e+00
+3.0000e+00      +0.0000e+00
+4.0000e+00      +0.0000e+00
+5.0000e+00      -3.5527e-15
+6.0000e+00      -3.5527e-15
+7.0000e+00      +0.0000e+00

```
```function f11gf_example
nz = 16;
a = zeros(10000,1);
a(1:16) = [4;
1;
5;
2;
2;
3;
-1;
1;
4;
1;
-2;
3;
2;
-1;
-2;
5];
irow = zeros(10000,1,'int64');
irow(1:16) = [int64(1);2;2;3;4;4;5;5;5;6;6;6;7;7;7;7];
icol = zeros(10000,1,'int64');
icol(1:16) = [int64(1);1;2;3;2;4;1;4;5;2;5;6;1;2;3;7];
lfill = int64(0);
dtol = 0;
mic = 'N';
dscale = 0;
ipiv = zeros(7, 1, 'int64');
method = 'CG    ';
precon = 'P';
n = 7;
tol = 1e-06;
maxitn = int64(20);
anorm = 0;
sigmax = 0;
maxits = int64(7);
monit = int64(2);
irevcm = int64(0);
u = [0; 0; 0; 0; 0; 0; 0];
v = [15; 18; -8; 21; 11; 10; 29];
wgt = [0; 0; 0; 0; 0; 0; 0];

[a, irow, icol, ipiv, istr, nnzc, npivm, ifail] = ...
f11ja(int64(nz), a, irow, icol, lfill, dtol, mic, dscale, ipiv);

[lwreq, work, ifail] = ...
f11gd(method, precon, int64(n), tol, maxitn, anorm, sigmax, maxits, monit, ...
'sigcmp', 'S', 'norm_p', '1');

while (irevcm ~= 4)
[irevcm, u, v, work, ifail] = f11ge(irevcm, u, v, wgt, work);
if (irevcm == 1)
[v, ifail] = f11xe(a(1:nz), irow(1:nz), icol(1:nz), 'N', u);
if (ifail ~= 0)
irecvm = 6;
end
elseif (irevcm == 2)
[v, ifail] = f11jb(a, irow, icol, ipiv, istr, 'N', u);
if (ifail ~= 0)
irecvm = 6;
end
elseif (irevcm == 3)
[itn, stplhs, stprhs, anorm, sigmax, its, sigerr, ifail] = f11gf(work);
if (ifail ~= 0)
irecvm = 6;
end
fprintf('\nMonitoring at iteration number %d\nresidual norm: %14.4e\n', ...
itn, stplhs);
fprintf('\n   Solution Vector  Residual Vector\n');
for i = 1:n
fprintf('%+16.4e %+16.4e\n', u(i), v(i));
end
end
end

% Get information about the computation
[itn, stplhs, stprhs, anorm, sigmax, its, sigerr, ifail] = f11gf(work);
fprintf('\nNumber of iterations for convergence: %d\n', itn);
fprintf('Residual norm:                           %14.4e\n', stplhs);
fprintf('Right-hand side of termination criteria: %14.4e\n', stprhs);
fprintf('i-norm of matrix a:                      %14.4e\n', anorm);
fprintf('\n   Solution Vector  Residual Vector\n');
for i = 1:n
fprintf('%+16.4e %+16.4e\n', u(i), v(i));
end
```
```

Monitoring at iteration number 2
residual norm:     1.9938e+00

Solution Vector  Residual Vector
+9.6320e-01      -2.2960e-01
+1.9934e+00      +2.2254e-01
+3.0583e+00      +9.5827e-02
+4.1453e+00      -2.5155e-01
+4.8289e+00      -1.7160e-01
+5.6630e+00      +6.7533e-01
+7.1062e+00      -3.4737e-01

Monitoring at iteration number 4
residual norm:     6.6574e-03

Solution Vector  Residual Vector
+9.9940e-01      -1.0551e-03
+2.0011e+00      -2.4675e-03
+3.0008e+00      -1.7116e-05
+3.9996e+00      +4.4929e-05
+4.9991e+00      +2.1359e-03
+5.9993e+00      -8.7482e-04
+7.0007e+00      +6.2045e-05

Number of iterations for convergence: 5
Residual norm:                               1.3323e-14
Right-hand side of termination criteria:     3.9200e-04
i-norm of matrix a:                          1.0000e+01

Solution Vector  Residual Vector
+1.0000e+00      -1.7764e-15
+2.0000e+00      +0.0000e+00
+3.0000e+00      -8.8818e-16
+4.0000e+00      -3.5527e-15
+5.0000e+00      -5.3291e-15
+6.0000e+00      +1.7764e-15
+7.0000e+00      +0.0000e+00

```