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_matop_complex_gen_matrix_cond_exp (f01kg)

## Purpose

nag_matop_complex_gen_matrix_cond_exp (f01kg) computes an estimate of the relative condition number ${\kappa }_{\mathrm{exp}}\left(A\right)$ of the exponential of a complex $n$ by $n$ matrix $A$, in the $1$-norm. The matrix exponential ${e}^{A}$ is also returned.

## Syntax

[a, condea, ifail] = f01kg(a, 'n', n)
[a, condea, ifail] = nag_matop_complex_gen_matrix_cond_exp(a, 'n', n)

## Description

The Fréchet derivative of the matrix exponential of $A$ is the unique linear mapping $E⟼L\left(A,E\right)$ such that for any matrix $E$
 $eA+E - e A - LA,E = oE .$
The derivative describes the first-order effect of perturbations in $A$ on the exponential ${e}^{A}$.
The relative condition number of the matrix exponential can be defined by
 $κexpA = LA A expA ,$
where $‖L\left(A\right)‖$ is the norm of the Fréchet derivative of the matrix exponential at $A$.
To obtain the estimate of ${\kappa }_{\mathrm{exp}}\left(A\right)$, nag_matop_complex_gen_matrix_cond_exp (f01kg) first estimates $‖L\left(A\right)‖$ by computing an estimate $\gamma$ of a quantity $K\in \left[{n}^{-1}{‖L\left(A\right)‖}_{1},n{‖L\left(A\right)‖}_{1}\right]$, such that $\gamma \le K$.
The algorithms used to compute ${\kappa }_{\mathrm{exp}}\left(A\right)$ are detailed in the Al–Mohy and Higham (2009a) and Al–Mohy and Higham (2009b).
The matrix exponential ${e}^{A}$ is computed using a Padé approximant and the scaling and squaring method. The Padé approximant is differentiated to obtain the Fréchet derivatives $L\left(A,E\right)$ which are used to estimate the condition number.

## References

Al–Mohy A H and Higham N J (2009a) A new scaling and squaring algorithm for the matrix exponential SIAM J. Matrix Anal. 31(3) 970–989
Al–Mohy A H and Higham N J (2009b) Computing the Fréchet derivative of the matrix exponential, with an application to condition number estimation SIAM J. Matrix Anal. Appl. 30(4) 1639–1657
Higham N J (2008) Functions of Matrices: Theory and Computation SIAM, Philadelphia, PA, USA
Moler C B and Van Loan C F (2003) Nineteen dubious ways to compute the exponential of a matrix, twenty-five years later SIAM Rev. 45 3–49

## Parameters

### Compulsory Input Parameters

1:     $\mathrm{a}\left(\mathit{lda},:\right)$ – complex array
The first dimension of the array a must be at least ${\mathbf{n}}$.
The second dimension of the array a must be at least ${\mathbf{n}}$.
The $n$ by $n$ matrix $A$.

### Optional Input Parameters

1:     $\mathrm{n}$int64int32nag_int scalar
Default: the first dimension of the array a.
$n$, the order of the matrix $A$.
Constraint: ${\mathbf{n}}\ge 0$.

### Output Parameters

1:     $\mathrm{a}\left(\mathit{lda},:\right)$ – complex array
The first dimension of the array a will be ${\mathbf{n}}$.
The second dimension of the array a will be ${\mathbf{n}}$.
The $n$ by $n$ matrix exponential ${e}^{A}$.
2:     $\mathrm{condea}$ – double scalar
An estimate of the relative condition number of the matrix exponential ${\kappa }_{\mathrm{exp}}\left(A\right)$.
3:     $\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$
The linear equations to be solved for the Padé approximant are singular; it is likely that this function has been called incorrectly.
${\mathbf{ifail}}=2$
${e}^{A}$ has been computed using an IEEE double precision Padé approximant, although the arithmetic precision is higher than IEEE double precision.
${\mathbf{ifail}}=3$
${\mathbf{ifail}}=-1$
Constraint: ${\mathbf{n}}\ge 0$.
${\mathbf{ifail}}=-3$
Constraint: $\mathit{lda}\ge {\mathbf{n}}$.
${\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

nag_matop_complex_gen_matrix_cond_exp (f01kg) uses the norm estimation function nag_linsys_complex_gen_norm_rcomm (f04zd) to produce an estimate $\gamma$ of a quantity $K\in \left[{n}^{-1}{‖L\left(A\right)‖}_{1},n{‖L\left(A\right)‖}_{1}\right]$, such that $\gamma \le K$. For further details on the accuracy of norm estimation, see the documentation for nag_linsys_complex_gen_norm_rcomm (f04zd).
For a normal matrix $A$ (for which ${A}^{\mathrm{H}}A=A{A}^{\mathrm{H}}$) the computed matrix, ${e}^{A}$, is guaranteed to be close to the exact matrix, that is, the method is forward stable. No such guarantee can be given for non-normal matrices. See Section 10.3 of Higham (2008) for details and further discussion.
For further discussion of the condition of the matrix exponential see Section 10.2 of Higham (2008).

nag_matop_complex_gen_matrix_cond_std (f01ka) uses a similar algorithm to nag_matop_complex_gen_matrix_cond_exp (f01kg) to compute an estimate of the absolute condition number (which is related to the relative condition number by a factor of $‖A‖/‖\mathrm{exp}\left(A\right)‖$). However, the required Fréchet derivatives are computed in a more efficient and stable manner by nag_matop_complex_gen_matrix_cond_exp (f01kg) and so its use is recommended over nag_matop_complex_gen_matrix_cond_std (f01ka).
The cost of the algorithm is $O\left({n}^{3}\right)$ and the complex allocatable memory required is approximately $15{n}^{2}$; see Al–Mohy and Higham (2009a) and Al–Mohy and Higham (2009b) for further details.
If the matrix exponential alone is required, without an estimate of the condition number, then nag_matop_complex_gen_matrix_exp (f01fc) should be used. If the Fréchet derivative of the matrix exponential is required then nag_matop_complex_gen_matrix_frcht_exp (f01kh) should be used.
As well as the excellent book Higham (2008), the classic reference for the computation of the matrix exponential is Moler and Van Loan (2003).

## Example

This example estimates the relative condition number of the matrix exponential ${e}^{A}$, where
 $A = 1+0i 2+0i 2+0i 2+i 3+2i 1i+0 1i+0 2+i 3+2i 2+0i 1i+0 2+i 3+2i 3+2i 3+2i 1+i .$
```function f01kg_example

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

% Exponential and conditioning of matrix A
a =  [1+ i, 2+ i, 2+ i, 2+i;
3+2i, 1,    1,    2+i;
3+2i, 2+ i, 1,    2+i;
3+2i, 3+2i, 3+2i, 1+i];

% Compute exp(a)
[expa, condea, ifail] = f01kg(a);

% Display results
disp('exp(A):');
disp(expa);

fprintf('Estimated condition number is: %6.2f\n', condea);

```
```f01kg example results

exp(A):
1.0e+03 *

-0.1579 - 0.7544i  -0.1947 - 0.5551i  -0.1866 - 0.4755i  -0.1558 - 0.5202i
-0.2069 - 0.6947i  -0.2255 - 0.5054i  -0.2124 - 0.4311i  -0.1866 - 0.4755i
-0.2087 - 0.8082i  -0.2385 - 0.5908i  -0.2255 - 0.5054i  -0.1947 - 0.5551i
-0.1334 - 1.0855i  -0.2087 - 0.8082i  -0.2069 - 0.6947i  -0.1579 - 0.7544i

Estimated condition number is:  15.29
```