F04 Chapter Contents
F04 Chapter Introduction
NAG Library Manual

# NAG Library Routine DocumentF04YDF

Note:  before using this routine, please read the Users' Note for your implementation to check the interpretation of bold italicised terms and other implementation-dependent details.

## 1  Purpose

F04YDF estimates the $1$-norm of a real rectangular matrix without accessing the matrix explicitly. It uses reverse communication for evaluating matrix products. The routine may be used for estimating condition numbers of square matrices.

## 2  Specification

 SUBROUTINE F04YDF ( IREVCM, M, N, X, LDX, Y, LDY, ESTNRM, T, SEED, WORK, IWORK, IFAIL)
 INTEGER IREVCM, M, N, LDX, LDY, T, SEED, IWORK(2*N+5*T+20), IFAIL REAL (KIND=nag_wp) X(LDX,*), Y(LDY,*), ESTNRM, WORK(M*T)

## 3  Description

F04YDF computes an estimate (a lower bound) for the $1$-norm
 $A1 = max 1≤j≤n ∑ i=1 m aij$ (1)
of an $m$ by $n$ real matrix $A=\left({a}_{ij}\right)$. The routine regards the matrix $A$ as being defined by a user-supplied ‘Black Box’ which, given an $n×t$ matrix $X$ (with $t\ll n$) or an $m×t$ matrix $Y$, can return $AX$ or ${A}^{\mathrm{T}}Y$. A reverse communication interface is used; thus control is returned to the calling program whenever a matrix product is required.
Note:  this routine 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 formula (1) above.
The main use of the routine is for estimating ${‖{B}^{-1}‖}_{1}$ for a square matrix, $B$, and hence the condition number ${\kappa }_{1}\left(B\right)={‖B‖}_{1}{‖{B}^{-1}‖}_{1}$, without forming ${B}^{-1}$ explicitly ($A={B}^{-1}$ above).
If, for example, an $LU$ factorization of $B$ is available, the matrix products ${B}^{-1}X$ and ${B}^{-\mathrm{T}}Y$ required by F04YDF may be computed by back- and forward-substitutions, without computing ${B}^{-1}$.
The routine can also be used to estimate $1$-norms of matrix products such as ${A}^{-1}B$ and $ABC$, without forming the products explicitly. Further applications are described by Higham (1988).
Since ${‖A‖}_{\infty }={‖{A}^{\mathrm{T}}‖}_{1}$, F04YDF can be used to estimate the $\infty$-norm of $A$ by working with ${A}^{\mathrm{T}}$ instead of $A$.
The algorithm used is described in Higham and Tisseur (2000).

## 4  References

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
Higham N J and Tisseur F (2000) A block algorithm for matrix $1$-norm estimation, with an application to $1$-norm pseudospectra SIAM J. Matrix. Anal. Appl. 21 1185–1201

## 5  Parameters

Note: this routine uses reverse communication. Its use involves an initial entry, intermediate exits and re-entries, and a final exit, as indicated by the parameter IREVCM. Between intermediate exits and re-entries, all parameters other than X and Y must remain unchanged.
1:     IREVCM – INTEGERInput/Output
On initial entry: must be set to $0$.
On intermediate exit: ${\mathbf{IREVCM}}=1$ or $2$, and X and Y contain the elements of $n×t$ and $m×t$ matrices respectively. The calling program must
 (a) if ${\mathbf{IREVCM}}=1$, evaluate $AX$ and store the result in Y or if ${\mathbf{IREVCM}}=2$, evaluate ${A}^{\mathrm{T}}Y$ and store the result in X, (b) call F04YDF once again, with all the other parameters unchanged.
On intermediate re-entry: IREVCM must be unchanged.
On final exit: ${\mathbf{IREVCM}}=0$.
2:     M – INTEGERInput
On entry: the number of rows of the matrix $A$.
Constraint: ${\mathbf{M}}\ge 0$.
3:     N – INTEGERInput
On entry: $n$, the number of columns of the matrix $A$.
Constraint: ${\mathbf{N}}\ge 0$.
4:     X(LDX,$*$) – REAL (KIND=nag_wp) arrayInput/Output
Note: the second dimension of the array X must be at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{T}}\right)$.
On initial entry: need not be set.
On intermediate exit: if ${\mathbf{IREVCM}}=1$, contains the current matrix $X$.
On intermediate re-entry: if ${\mathbf{IREVCM}}=2$, must contain ${A}^{\mathrm{T}}Y$.
On final exit: the array is undefined.
5:     LDX – INTEGERInput
On initial entry: the leading dimension of the array X as declared in the (sub)program from which F04YDF is called.
Constraint: ${\mathbf{LDX}}\ge {\mathbf{N}}$.
6:     Y(LDY,$*$) – REAL (KIND=nag_wp) arrayInput/Output
Note: the second dimension of the array Y must be at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{T}}\right)$.
On initial entry: need not be set.
On intermediate exit: if ${\mathbf{IREVCM}}=2$, contains the current matrix $Y$.
On intermediate re-entry: if ${\mathbf{IREVCM}}=1$, must contain $AX$.
On final exit: the array is undefined.
7:     LDY – INTEGERInput
On initial entry: the leading dimension of the array Y as declared in the (sub)program from which F04YDF is called.
Constraint: ${\mathbf{LDY}}\ge {\mathbf{M}}$.
8:     ESTNRM – REAL (KIND=nag_wp)Input/Output
On initial entry: need not be set.
On intermediate re-entry: must not be changed.
On final exit: an estimate (a lower bound) for ${‖A‖}_{1}$.
9:     T – INTEGERInput
On entry: the number of columns $t$ of the matrices $X$ and $Y$. This is a parameter that can be used to control the accuracy and reliability of the estimate and corresponds roughly to the number of columns of $A$ that are visited during each iteration of the algorithm.
If ${\mathbf{T}}\ge 2$ then a partly random starting matrix is used in the algorithm.
Suggested value: ${\mathbf{T}}=2$.
Constraint: $1\le {\mathbf{T}}\le {\mathbf{M}}$.
10:   SEED – INTEGERInput
On entry: the seed used for random number generation.
If ${\mathbf{T}}=1$, SEED is not used.
11:   WORK(${\mathbf{M}}×{\mathbf{T}}$) – REAL (KIND=nag_wp) arrayCommunication Array
12:   IWORK($2×{\mathbf{N}}+5×{\mathbf{T}}+20$) – INTEGER arrayCommunication Array
On initial entry: need not be set.
On intermediate re-entry: must not be changed.
13:   IFAIL – INTEGERInput/Output
On initial entry: IFAIL must be set to $0$, $-1\text{​ or ​}1$. If you are unfamiliar with this parameter you should refer to Section 3.3 in the Essential Introduction for details.
For environments where it might be inappropriate to halt program execution when an error is detected, the value $-1\text{​ or ​}1$ is recommended. If the output of error messages is undesirable, then the value $1$ is recommended. Otherwise, because for this routine the values of the output parameters may be useful even if ${\mathbf{IFAIL}}\ne {\mathbf{0}}$ on exit, the recommended value is $-1$. When the value $-\mathbf{1}\text{​ or ​}1$ is used it is essential to test the value of IFAIL on exit.
On final exit: ${\mathbf{IFAIL}}={\mathbf{0}}$ unless the routine detects an error or a warning has been flagged (see Section 6).

## 6  Error Indicators and Warnings

If on entry ${\mathbf{IFAIL}}={\mathbf{0}}$ or $-{\mathbf{1}}$, explanatory error messages are output on the current error message unit (as defined by X04AAF).
Errors or warnings detected by the routine:
${\mathbf{IFAIL}}=1$
${\mathbf{IFAIL}}=-1$
On entry, ${\mathbf{IREVCM}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{IREVCM}}=0$, $1$ or $2$.
On initial entry, ${\mathbf{IREVCM}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{IREVCM}}=0$.
${\mathbf{IFAIL}}=-2$
On entry, ${\mathbf{M}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{M}}\ge 0$.
${\mathbf{IFAIL}}=-3$
On entry, ${\mathbf{N}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{N}}\ge 0$.
${\mathbf{IFAIL}}=-5$
On entry, ${\mathbf{LDX}}=⟨\mathit{\text{value}}⟩$ and ${\mathbf{N}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{LDX}}\ge {\mathbf{N}}$.
${\mathbf{IFAIL}}=-7$
On entry, ${\mathbf{LDY}}=⟨\mathit{\text{value}}⟩$ and ${\mathbf{M}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{LDY}}\ge {\mathbf{M}}$.
${\mathbf{IFAIL}}=-9$
On entry, ${\mathbf{M}}=⟨\mathit{\text{value}}⟩$ and ${\mathbf{T}}=⟨\mathit{\text{value}}⟩$.
Constraint: $1\le {\mathbf{T}}\le {\mathbf{M}}$.

## 7  Accuracy

In extensive tests on random matrices of size up to $m=n=450$ the estimate ESTNRM has been found always to be within a factor two of ${‖A‖}_{1}$; often the estimate has many correct figures. However, matrices exist for which the estimate is smaller than ${‖A‖}_{1}$ by an arbitrary factor; such matrices are very unlikely to arise in practice. See Higham and Tisseur (2000) for further details.

### 8.1  Timing

For most problems the time taken during calls to F04YDF will be negligible compared with the time spent evaluating matrix products between calls to F04YDF.
The number of matrix products required depends on the matrix $A$. At most six products of the form $Y=AX$ and five products of the form $X={A}^{\mathrm{T}}Y$ will be required. The number of iterations is independent of the choice of $t$.

### 8.2  Overflow

It is your responsibility to guard against potential overflows during evaluation of the matrix products. In particular, when estimating ${‖{B}^{-1}‖}_{1}$ using a triangular factorization of $B$, F04YDF should not be called if one of the factors is exactly singular – otherwise division by zero may occur in the substitutions.

### 8.3  Choice of $t$

The parameter $t$ controls the accuracy and reliability of the estimate. For $t=1$, the algorithm behaves similarly to the LAPACK estimator xLACON. Increasing $t$ typically improves the estimate, without increasing the number of iterations required.
For $t\ge 2$, random matrices are used in the algorithm, so for repeatable results the value of SEED should be kept constant.
A value of $t=2$ is recommended for new users.

### 8.4  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) THEN
IREVCM = 0
10       CALL F04YDF (IREVCM,M,N,X,LDX,Y,LDY,ESTNRM,T,SEED,WORK,   &
IWORK,IFAIL)
IF (IREVCM.NE.0) THEN
IF (IREVCM.EQ.1) THEN
...  code to compute Y=inv(A)X ...
ELSE
...  code to compute X=inv(transpose(A))Y ...
END IF
GO TO 10
END IF
END IF
```
To compute ${A}^{-1}X$ or ${A}^{-\mathrm{T}}Y$, solve the equation $AY=X$ or ${A}^{\mathrm{T}}X=Y$, storing the result in Y or X respectively. The code will vary, depending on the type of the matrix $A$, and the NAG routine used to factorize $A$.
The factorization will normally have been performed by a suitable routine from Chapters F01, F03 or F07. Note also that many of the ‘Black Box’ routines in Chapter F04 for solving systems of equations also return a factorization of the matrix. The example program in Section 9 illustrates how F04YDF can be used in conjunction with NAG Library routines for $LU$ factorization of a real matrix F07ADF (DGETRF).
It is straightforward to use F04YDF for the following other types of matrix, using the named routines for factorization and solution:
For upper or lower triangular matrices, no factorization routine is needed: $Y={A}^{-1}X$ and $X={A}^{-\mathrm{T}}Y$ may be computed by calls to F06PJF (DTRSV) (or F06PKF (DTBSV) if the matrix is banded, or F06PLF (DTPSV) if the matrix is stored in packed form).

## 9  Example

For this routine two examples are provided. There is a single example program for F04YDF, with a main program and the code to solve the two example problems is given in Example 1 (EX1) and Example 2 (EX2).
Example 1 (EX1)
This example estimates the condition number ${‖A‖}_{1}{‖{A}^{-1}‖}_{1}$ of the matrix $A$ given by
 $A= 0.7 -0.2 1.0 0.0 2.0 0.1 0.3 0.7 0.0 1.0 0.9 0.2 0.0 0.0 0.2 0.7 0.0 -1.1 0.0 3.4 -0.7 0.2 0.1 0.1 0.0 -4.0 0.0 1.0 9.0 0.0 0.4 1.2 4.3 0.0 6.2 5.9 .$
Example 2 (EX2)
This example estimates the condition number of the sparse matrix $A$ (stored in symmetric coordinate storage format) given by
 $A = 0.0 0.0 0.0 1.0 0.0 3.0 1.0 0.0 0.0 0.0 0.0 2.0 0.0 2.0 0.0 2.0 0.0 4.0 0.0 5.0 0.0 1.0 2.0 0.0 0.0 .$

### 9.1  Program Text

Program Text (f04ydfe.f90)

### 9.2  Program Data

Program Data (f04ydfe.d)

### 9.3  Program Results

Program Results (f04ydfe.r)