NAG Toolbox: nag_sparse_complex_gen_basic_solver (f11bs)

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

1:     $\mathrm{irevcm}$ – int64int32nag_int scalar

On initial entry: $\mathbf{irevcm}=0$, otherwise an error condition will be raised.

On intermediate re-entry: must either be unchanged from its previous exit value, or can have one of the following values.

$\mathbf{irevcm}=5$

Tidy termination: the computation will terminate at the end of the current iteration. Further reverse communication exits may occur depending on when the termination request is issued. nag_sparse_complex_gen_basic_solver (f11bs) will then return with the termination code $\mathbf{irevcm}=4$. Note that before calling nag_sparse_complex_gen_basic_solver (f11bs) with $\mathbf{irevcm}=5$ the calling program must have performed the tasks required by the value of irevcm returned by the previous call to nag_sparse_complex_gen_basic_solver (f11bs), otherwise subsequently returned values may be invalid.

$\mathbf{irevcm}=6$

Immediate termination: nag_sparse_complex_gen_basic_solver (f11bs) will return immediately with termination code $\mathbf{irevcm}=4$ and with any useful information available. This includes the last iterate of the solution.
Immediate termination may be useful, for example, when errors are detected during matrix-vector multiplication or during the solution of the preconditioning equation.

Changing irevcm to any other value between calls will result in an error.

Constraint: on initial entry, $\mathbf{irevcm}=0$; on re-entry, either irevcm must remain unchanged or be reset to $\mathbf{irevcm}=5$ or $6$.

2:     $\mathrm{u}\left(:\right)$ – complex array

The dimension of the array u must be at least $\mathit{n}$

On initial entry: an initial estimate, $x_0$, of the solution of the system of equations $Ax=b$.

On intermediate re-entry: must remain unchanged.

3:     $\mathrm{v}\left(:\right)$ – complex array

The dimension of the array v must be at least $\mathit{n}$

On initial entry: the right-hand side $b$ of the system of equations $Ax=b$.

On intermediate re-entry: the returned value of irevcm determines the contents of v as follows.

If $\mathbf{irevcm}=-1$, $1$ or $2$, v must store the vector $v$, the result of the operation specified by the value of irevcm returned by the previous call to nag_sparse_complex_gen_basic_solver (f11bs).

If $\mathbf{irevcm}=3$, v must remain unchanged.

4:     $\mathrm{wgt}\left(:\right)$ – double array

The dimension of the array wgt must be at least $\max \left(1,\mathit{n}\right)$

The user-supplied weights, if these are to be used in the computation of the vector norms in the termination criterion (see Description and Arguments in nag_sparse_complex_gen_basic_setup (f11br)).

Constraint: if weights are to be used, at least one element of wgt must be nonzero.

5:     $\mathrm{work}\left(\mathbf{lwork}\right)$ – complex array

On initial entry: the array work as returned by nag_sparse_complex_gen_basic_setup (f11br) (see also Arguments in nag_sparse_complex_gen_basic_setup (f11br)).

On intermediate re-entry: must remain unchanged.

Optional Input Parameters

1:     $\mathrm{lwork}$ – int64int32nag_int scalar

Default: the dimension of the array work.

On initial entry: the dimension of the array work (see also Description and Arguments in nag_sparse_complex_gen_basic_setup (f11br)). the required amount of workspace is as follows:

Method	Requirements
RGMRES	$\mathbf{lwork}=120+\mathit{n}\left(m+3\right)+m\left(m+5\right)+1$, where $m$ is the dimension of the basis
CGS	$\mathbf{lwork}=120+7\mathit{n}$
Bi-CGSTAB($\ell$)	$\mathbf{lwork}=120+\left(2\mathit{n}+\ell \right)\left(\ell +2\right)+p$, where $\ell$ is the order of the method
TFQMR	$\mathbf{lwork}=120+10\mathit{n}$,

where

• $p=2\mathit{n}$, if $\ell >1$ and $\mathbf{iterm}=2$ was supplied;
• $p=\mathit{n}$, if $\ell >1$ and a preconditioner is used or $\mathbf{iterm}=2$ was supplied;
• $p=0$, otherwise.

Constraint: $\mathbf{lwork}\ge \mathbf{lwreq}$, where lwreq is returned by nag_sparse_complex_gen_basic_setup (f11br).

Output Parameters

1:     $\mathrm{irevcm}$ – int64int32nag_int scalar

On intermediate exit: has the following meanings.

$\mathbf{irevcm}=-1$

The calling program must compute the matrix-vector product $v=A^{\mathrm{H}}u$, where $u$ and $v$ are stored in u and v, respectively; RGMRES, CGS and Bi-CGSTAB($\ell$) methods return $\mathbf{irevcm}=-1$ only if the matrix norm ${‖A‖}_1$ or ${‖A‖}_{\infty }$ is estimated internally using Higham's method. This can only happen if $\mathbf{iterm}=1$ in nag_sparse_complex_gen_basic_setup (f11br).

$\mathbf{irevcm}=1$

The calling program must compute the matrix-vector product $v=Au$, where $u$ and $v$ are stored in u and v, respectively.

$\mathbf{irevcm}=2$

The calling program must solve the preconditioning equation $Mv=u$, where $u$ and $v$ are stored in u and v, respectively.

$\mathbf{irevcm}=3$

Monitoring step: the solution and residual at the current iteration are returned in the arrays u and v, respectively. No action by the calling program is required. nag_sparse_complex_gen_basic_diag (f11bt) can be called at this step to return additional information.

On final exit: $\mathbf{irevcm}=4$: nag_sparse_complex_gen_basic_solver (f11bs) has completed its tasks. The value of ifail determines whether the iteration has been successfully completed, errors have been detected or the calling program has requested termination.

2:     $\mathrm{u}\left(:\right)$ – complex array

The dimension of the array u will be $\mathit{n}$

On intermediate exit: the returned value of irevcm determines the contents of u as follows.

If $\mathbf{irevcm}=-1$, $1$ or $2$, u holds the vector $u$ on which the operation specified by irevcm is to be carried out.

If $\mathbf{irevcm}=3$, u holds the current iterate of the solution vector.

On final exit: if $\mathbf{ifail}=\mathbf{3}$ or $\mathbf{ifail}<\mathbf{0}$, the array u is unchanged from the initial entry to nag_sparse_complex_gen_basic_solver (f11bs).

If $\mathbf{ifail}=\mathbf{1}$, the array u is unchanged from the last entry to nag_sparse_complex_gen_basic_solver (f11bs).

Otherwise, u holds the last available iterate of the solution of the system of equations, for all returned values of ifail.

3:     $\mathrm{v}\left(:\right)$ – complex array

The dimension of the array v will be $\mathit{n}$

On intermediate exit: if $\mathbf{irevcm}=3$, v holds the current iterate of the residual vector. Note that this is an approximation to the true residual vector. Otherwise, it does not contain any useful information.

On final exit: if $\mathbf{ifail}=\mathbf{3}$ or $\mathbf{ifail}<\mathbf{0}$, the array v is unchanged from the initial entry to nag_sparse_complex_gen_basic_solver (f11bs).

If $\mathbf{ifail}=\mathbf{1}$, the array v is unchanged from the last entry to nag_sparse_complex_gen_basic_solver (f11bs).

If $\mathbf{ifail}=\mathbf{0}$ or $\mathbf{2}$, the array v contains the true residual vector of the system of equations (see also Error Indicators and Warnings).

Otherwise, v stores the last available iterate of the residual vector unless $\mathbf{ifail}=\mathbf{8}$ is returned on last entry, in which case v is set to $0.0$.

4:     $\mathrm{work}\left(\mathbf{lwork}\right)$ – complex array

5:     $\mathrm{ifail}$ – int64int32nag_int scalar

On final exit: $\mathbf{ifail}=\mathbf{0}$ unless the function detects an error (see Error Indicators and Warnings).

Error Indicators and Warnings

Cases prefixed with W are classified as warnings and do not generate an error of type NAG:error_n. See nag_issue_warnings.

   $\mathbf{ifail}=-i$

If $\mathbf{ifail}=-i$, parameter $i$ had an illegal value on entry. The parameters are numbered as follows:

1: irevcm, 2: u, 3: v, 4: wgt, 5: work, 6: lwork, 7: ifail.

W  $\mathbf{ifail}=1$

nag_sparse_complex_gen_basic_solver (f11bs) has been called again after returning the termination code $\mathbf{irevcm}=4$. No further computation has been carried out and all input data and data stored for access by nag_sparse_complex_gen_basic_diag (f11bt) have remained unchanged.

W  $\mathbf{ifail}=2$

The required accuracy could not be obtained. However, nag_sparse_complex_gen_basic_solver (f11bs) has terminated with reasonable accuracy: the last iterate of the residual satisfied the termination criterion but the exact residual $r=b-Ax$, did not. After the first occurrence of this situation, the iteration was restarted once, but nag_sparse_complex_gen_basic_solver (f11bs) could not improve on the accuracy. This error code usually implies that your problem has been fully and satisfactorily solved to within or close to the accuracy available on your system. Further iterations are unlikely to improve on this situation. You should call nag_sparse_complex_gen_basic_diag (f11bt) to check the values of the left- and right-hand sides of the termination condition.

   $\mathbf{ifail}=3$

nag_sparse_complex_gen_basic_setup (f11br) was either not called before calling nag_sparse_complex_gen_basic_solver (f11bs) or it returned an error. The arguments u and v remain unchanged.

W  $\mathbf{ifail}=4$

The calling program requested a tidy termination before the solution had converged. The arrays u and v return the last iterates available of the solution and of the residual vector, respectively.

W  $\mathbf{ifail}=5$

The solution did not converge within the maximum number of iterations allowed. The arrays u and v return the last iterates available of the solution and of the residual vector, respectively.

   $\mathbf{ifail}=6$

Algorithmic breakdown. The last available iterates of the solution and residuals are returned, although it is possible that they are completely inaccurate.

W  $\mathbf{ifail}=8$

The calling program requested an immediate termination. However, the array u returns the last iterate of the solution, the array v returns the last iterate of the residual vector, for the CGS and TFQMR methods only.

   $\mathbf{ifail}=10$

User-supplied weights are to be used, but all elements of the array wgt are zero.

   $\mathbf{ifail}=-99$

An unexpected error has been triggered by this routine. Please contact NAG.

   $\mathbf{ifail}=-399$

Your licence key may have expired or may not have been installed correctly.

   $\mathbf{ifail}=-999$

Dynamic memory allocation failed.

Accuracy

Further Comments

Example

Open in the MATLAB editor: f11bs_example

function f11bs_example


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

% Solve sparse system Ax = b iteratively using ILU preconditioner

% Define sparse matrix A and RHS B
nz   = int64(24);
n    = int64(8);
a    = complex(zeros(3*nz,1));
irow = zeros(3*nz,1,'int64');
icol = zeros(3*nz,1,'int64');
a(1:nz) = [ 2 + 1i, -1 + 1i,  1 - 3i,  4 + 7i, -3 + 0i,  2 + 4i, ...
           -7 - 5i,  2 + 1i,  3 + 2i, -4 + 2i,  0 + 1i,  5 - 3i, ...
           -1 + 2i,  8 + 6i, -3 - 4i, -6 - 2i,  5 - 2i,  2 + 0i, ...
            0 - 5i, -1 + 5i,  6 + 2i, -1 + 4i,  2 + 0i,  3 + 3i];
b       = [ 7 + 11i;
            1 + 24i;
          -13 - 18i;
          -10 +  3i;
           23 + 14i;
           17 -  7i;
           15 -  3i;
           -3 + 20i];
irow(1:nz) = int64(...
             [1; 1; 1; 2; 2; 2; 3; 3; 4; 4; 4; 4;
              5; 5; 5; 6; 6; 6; 7; 7; 7; 8; 8; 8]);
icol(1:nz) = int64(...
             [1; 4; 8; 1; 2; 5; 3; 6; 1; 3; 4; 7;
              2; 5; 7; 1; 3; 6; 3; 5; 7; 2; 6; 8]);

% ILU preconditioner
%     Input parameters
lfill = int64(0);
dtol  = 0;
milu  = 'No modification';
ipivp = zeros(n, 1, 'int64');
ipivq = zeros(n, 1, 'int64');
% Compute preconditioner
[a, irow, icol, ipivp, ipivq, istr, idiag, nnzc, npivm, ifail] = ...
    f11dn(...
          nz, a, irow, icol, lfill, dtol, milu, ipivp, ipivq);


% Iterative Solution Setup
%      Input parameters
method = 'TFQMR   ';
precon = 'P';
lpoly  = int64(1);
tol    = sqrt(x02aj);
maxitn = int64(20);
anorm  = 0;
sigmax = 0;
monit  = int64(2);
lwork  = int64(6000);

[lwreq, work, ifail] = ...
     f11br(...
           method, precon, n, lpoly, tol, maxitn, anorm, sigmax, ...
           monit, lwork, 'norm_p', '1');

irevcm = int64(0);
wgt    = zeros(n,1);
u      = complex(zeros(n,1));
v      = b;

while (irevcm ~= 4)
  [irevcm, u, v, work, ifail] = f11bs(...
                                      irevcm, u, v, wgt, work);

  if (irevcm == -1)
    % v = A^Hu
    [v, ifail] = f11xn(...
                       'T', a(1:nz), irow(1:nz), icol(1:nz), 'N', u);
  elseif (irevcm == 1)
    % v = Au
    [v, ifail] = f11xn(...
                       'N', a(1:nz), irow(1:nz), icol(1:nz), 'N', u);
  elseif (irevcm == 2)
    % Solve ILU v = u  
    [v, ifail] = f11dp(...
                       'N', a, irow, icol, ipivp, ipivq, istr, idiag, 'N', u);
  elseif (irevcm == 3)
    % Monitoring stage
    [itn, stplhs, stprhs, anorm, sigmax, work, ifail] = ...
       f11bt(work);
    fprintf('\nMonitoring at iteration number %2d\n',itn);
    fprintf('residual norm:              %14.4e\n', stplhs);
    fprintf('\n   Solution Vector\n');
    disp(u);
    fprintf('\n   Residual Vector\n');
    disp(v);
  end
end

% Get information about the computation
[itn, stplhs, stprhs, anorm, sigmax, work, ifail] = ...
      f11bt(work);

fprintf('\nNumber of iterations for convergence:     %4d\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\n');
disp(u);
fprintf('\n   Residual Vector\n');
disp(v);

f11bs example results


Monitoring at iteration number  2
residual norm:                  8.2345e+01

   Solution Vector
   0.6905 + 1.4236i
   0.0739 - 1.1880i
   1.4778 + 0.4785i
   5.6572 - 3.0786i
   1.4243 - 1.1246i
   0.1037 + 1.9740i
   0.4498 - 1.2715i
   2.5704 + 1.7578i


   Residual Vector
   1.7772 + 4.6797i
   1.0774 + 6.4600i
  -3.2812 -11.3135i
  -3.8698 - 1.6438i
   8.9912 +11.1004i
   9.7428 - 0.4622i
   3.1668 + 2.8721i
 -10.3231 + 1.5837i


Number of iterations for convergence:        4
Residual norm:                               1.3396e-11
Right-hand side of termination criteria:     9.3882e-06
i-norm of matrix a:                          2.7000e+01

   Solution Vector
   1.0000 + 1.0000i
   2.0000 - 1.0000i
   3.0000 + 1.0000i
   4.0000 - 1.0000i
   3.0000 - 1.0000i
   2.0000 + 1.0000i
   1.0000 - 1.0000i
  -0.0000 + 3.0000i


   Residual Vector
   1.0e-11 *

  -0.0060 - 0.0782i
   0.1551 - 0.1364i
   0.0822 + 0.0703i
   0.0988 - 0.0218i
  -0.1254 - 0.0417i
  -0.0735 + 0.0527i
   0.0284 + 0.0162i
   0.1112 + 0.2416i