hide long namesshow long names
hide short namesshow short names
Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

NAG Toolbox: nag_linsys_withdraw_complex_norm_rcomm (f04zc)


    1  Purpose
    2  Syntax
    7  Accuracy
    9  Example


nag_linsys_complex_norm_rcomm (f04zc) estimates the 1-norm of a complex matrix without accessing the matrix explicitly. It uses reverse communication for evaluating matrix-vector products. The function may be used for estimating matrix condition numbers.
Note: this function is scheduled to be withdrawn, please see f04zc in Advice on Replacement Calls for Withdrawn/Superseded Routines..


[icase, x, estnrm, work, ifail] = f04zc(icase, x, estnrm, work, 'n', n)
[icase, x, estnrm, work, ifail] = nag_linsys_withdraw_complex_norm_rcomm(icase, x, estnrm, work, 'n', n)


nag_linsys_complex_norm_rcomm (f04zc) computes an estimate (a lower bound) for the 1-norm
A1 = max 1jn i=1 n aij (1)
of an n by n complex matrix A=aij. The function regards the matrix A as being defined by a user-supplied ‘Black Box’ which, given an input vector x, can return either of the matrix-vector products Ax or AHx, where AH is the complex conjugate transpose. A reverse communication interface is used; thus control is returned to the calling program whenever a matrix-vector product is required.
Note:  this function is not recommended for use when the elements of A are known explicitly; it is then more efficient to compute the 1-norm directly from the formula (1) above.
The main use of the function is for estimating B-11, and hence the condition number κ1B=B1B-11, without forming B-1 explicitly (A=B-1 above).
If, for example, an LU factorization of B is available, the matrix-vector products B-1x and B-Hx required by nag_linsys_complex_norm_rcomm (f04zc) may be computed by back- and forward-substitutions, without computing B-1.
The function can also be used to estimate 1-norms of matrix products such as A-1B and ABC, without forming the products explicitly. Further applications are described in Higham (1988).
Since A=AH1, nag_linsys_complex_norm_rcomm (f04zc) can be used to estimate the -norm of A by working with AH instead of A.
The algorithm used is based on a method given in Hager (1984) and is described in Higham (1988). A comparison of several techniques for condition number estimation is given in Higham (1987).
Note: nag_linsys_complex_gen_norm_rcomm (f04zd) can also be used to estimate the norm of a real matrix. nag_linsys_complex_gen_norm_rcomm (f04zd) uses a more recent algorithm than nag_linsys_complex_norm_rcomm (f04zc) and it is recommended that nag_linsys_complex_gen_norm_rcomm (f04zd) be used in place of nag_linsys_complex_norm_rcomm (f04zc).


Hager W W (1984) Condition estimates SIAM J. Sci. Statist. Comput. 5 311–316
Higham N J (1987) A survey of condition number estimation for triangular matrices SIAM Rev. 29 575–596
Higham N J (1988) FORTRAN codes for estimating the one-norm of a real or complex matrix, with applications to condition estimation ACM Trans. Math. Software 14 381–396


Note:  this function uses reverse communication. Its use involves an initial entry, intermediate exits and re-entries, and a final exit, as indicated by the argument icase. Between intermediate exits and re-entries, all arguments other than x must remain unchanged.

Compulsory Input Parameters

1:     icase int64int32nag_int scalar
On initial entry: must be set to 0.
2:     xn – complex array
On initial entry: need not be set.
On intermediate re-entry: must contain Ax (if icase=1) or AHx (if icase=2).
3:     estnrm – double scalar
On initial entry: need not be set.
4:     workn – complex array
On initial entry: need not be set.

Optional Input Parameters

1:     n int64int32nag_int scalar
Default: the dimension of the arrays x, work. (An error is raised if these dimensions are not equal.)
On initial entry: n, the order of the matrix A.
Constraint: n1.

Output Parameters

1:     icase int64int32nag_int scalar
On intermediate exit: icase=1 or 2, and xi, for i=1,2,,n, contain the elements of a vector x. The calling program must
(a) evaluate Ax (if icase=1) or AHx (if icase=2), where AH is the complex conjugate transpose;
(b) place the result in x; and,
(c) call nag_linsys_complex_norm_rcomm (f04zc) once again, with all the other arguments unchanged.
On final exit: icase=0.
2:     xn – complex array
On intermediate exit: contains the current vector x.
On final exit: the array is undefined.
3:     estnrm – double scalar
On intermediate exit: should not be changed.
On final exit: an estimate (a lower bound) for A1.
4:     workn – complex array
On final exit: contains a vector v such that v=Aw where estnrm=v1/w1 (w is not returned). If A=B-1 and estnrm is large, then v is an approximate null vector for B.
5:     ifail int64int32nag_int scalar
ifail=0 unless the function detects an error (see Error Indicators and Warnings).

Error Indicators and Warnings

Errors or warnings detected by the function:
On entry,n<1.
An unexpected error has been triggered by this routine. Please contact NAG.
Your licence key may have expired or may not have been installed correctly.
Dynamic memory allocation failed.


In extensive tests on random matrices of size up to n=100 the estimate estnrm has been found always to be within a factor eleven of A1; often the estimate has many correct figures. However, matrices exist for which the estimate is smaller than A1 by an arbitrary factor; such matrices are very unlikely to arise in practice. See Higham (1988) for further details.

Further Comments


The total time taken by nag_linsys_complex_norm_rcomm (f04zc) is proportional to n. For most problems the time taken during calls to nag_linsys_complex_norm_rcomm (f04zc) will be negligible compared with the time spent evaluating matrix-vector products between calls to nag_linsys_complex_norm_rcomm (f04zc).
The number of matrix-vector products required varies from 5 to 11 (or is 1 if n=1). In most cases 5 products are required; it is rare for more than 7 to be needed.


It is your responsibility to guard against potential overflows during evaluation of the matrix-vector products. In particular, when estimating B-11 using a triangular factorization of B, nag_linsys_complex_norm_rcomm (f04zc) should not be called if one of the factors is exactly singular – otherwise division by zero may occur in the substitutions.

Use in Conjunction with NAG Library Routines

To estimate the 1-norm of the inverse of a matrix A, the following skeleton code can normally be used:
...  code to factorize A ...
if (A is not singular)
  icase = 0
  [icase, x, estnrm, work, ifail] = f04zc(icase, x, estnrm, work);
  while (icase ~= 0)
    if (icase == 1) 
      ...  code to compute A(-1)x ...
      ...  code to compute (A(-1)(H)) x ...
    [icase, x, estnrm, work, ifail] = f04zc(icase, x, estnrm, work);
To compute A-1x or A-Hx, solve the equation Ay=x or AHy=x for y, overwriting y on x. The code will vary, depending on the type of the matrix A, and the NAG function used to factorize A.
Note that if A is any type of Hermitian matrix, then A=AH, and the if statement after the while can be reduced to:
    ...  code to compute A(-1)x ...
The example program in Example illustrates how nag_linsys_complex_norm_rcomm (f04zc) can be used in conjunction with NAG Toolbox functions for complex band matrices (factorized by nag_lapack_zgbtrf (f07br)).
It is also straightforward to use nag_linsys_complex_norm_rcomm (f04zc) for Hermitian positive definite matrices, using nag_lapack_zpotrf (f07fr) and nag_lapack_zpotrs (f07fs) for factorization and solution.


This example estimates the condition number A1A-11 of the order 5 matrix
A = 1+0i 2+0i 1+2i 0i+0 0i+0 2i 3+5i 1+3i 2+0i 0i+0 0i+0 -2+6i 5+7i 6i 1-0i 0i+0 0i+0 3+9i 4i 4-3i 0i+0 0i+0 0i+0 -1+8i 10-3i  
where A is a band matrix stored in the packed format required by nag_lapack_zgbtrf (f07br) and nag_lapack_zgbtrs (f07bs).
Further examples of the technique for condition number estimation in the case of double matrices can be seen in the example program section of nag_linsys_real_norm_rcomm (f04yc).
function f04zc_example

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

a = [ 1.0 + 1.0i,  2.0 + 1.0i,  1.0 + 2.0i,  0.0 + 0.0i,  0.0 + 0.0i;
      0.0 + 2.0i,  3.0 + 5.0i,  1.0 + 3.0i,  2.0 + 1.0i,  0.0 + 0.0i,
      0.0 + 0.0i, -2.0 + 6.0i,  5.0 + 7.0i,  0.0 + 6.0i,  1.0 - 1.0i;
      0.0 + 0.0i,  0.0 + 0.0i,  3.0 + 9.0i,  0.0 + 4.0i,  4.0 - 3.0i;
      0.0 + 0.0i,  0.0 + 0.0i,  0.0 + 0.0i, -1.0 + 8.0i, 10.0 - 3.0i];

x     = complex(zeros(5, 1));
work  = complex(zeros(5,1));
anorm = norm(a,1);
icase = int64(0);
estnrm = 0;

done = false;
while (~done)
  [icase, x, estnrm, work, ifail] = ...
    f04zc(icase, x, estnrm, work);
  if (icase == 0)
    done = true;
  elseif (icase == 1)
    x = inv(a)*x;
    x = conj(transpose(inv(a)))*x;
fprintf('Computed norm of a              = %6.4g\n', anorm);
fprintf('Estimated norm of inverse(A)    = %6.4g\n', estnrm);
fprintf('Estimated condition number of A = %6.1f\n', estnrm*anorm);

f04zc example results

Computed norm of a              =  23.49
Estimated norm of inverse(A)    =  37.04
Estimated condition number of A =  870.0

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2015