Higham N J and Lin L (2013) An improved Schur–Padé algorithm for fractional powers of a matrix and their Fréchet derivatives MIMS Eprint 2013.1 Manchester Institute for Mathematical Sciences, School of Mathematics, University of Manchester http://eprints.ma.man.ac.uk/

Parameters

Compulsory Input Parameters

1: $a (lda :)$ – complex array: The first dimension of the array a must be at least $n$ .
The second dimension of the array a must be at least $n$ .

The $n$ by $n$ matrix $A$ .
2: $p$ – double scalar: The required power of $A$ .

Optional Input Parameters

1: $n$ – int64int32nag_int scalar: Default: the first dimension of the array a.
$n$ , the order of the matrix $A$ .

Constraint: $n \geq 0$ .

Output Parameters

1: $a (lda :)$ – complex array: The first dimension of the array a will be $n$ .
The second dimension of the array a will be $n$ .

If $ifail = 0$ , the $n$ by $n$ matrix $p$ th power, $A^{p}$ . Alternatively, if $ifail = 1$ , contains an $n$ by $n$ non-principal power of $A$ .
2: $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:

$ifail = 1$: $A$ has eigenvalues on the negative real line. The principal $p$ th power is not defined so a non-principal power is returned.

$ifail = 2$: $A$ is singular so the $p$ th power cannot be computed.

$ifail = 3$: $A^{p}$ has been computed using an IEEE double precision Padé approximant, although the arithmetic precision is higher than IEEE double precision.

$ifail = 4$: An unexpected internal error occurred. This failure should not occur and suggests that the function has been called incorrectly.

$ifail = - 1$: Constraint: $n \geq 0$ .

$ifail = - 3$: Constraint: $lda \geq n$ .

$ifail = - 99$: An unexpected error has been triggered by this routine. Please contact NAG.

$ifail = - 399$: Your licence key may have expired or may not have been installed correctly.

$ifail = - 999$: Dynamic memory allocation failed.

Accuracy

For positive integer

p

, the algorithm reduces to a sequence of matrix multiplications. For negative integer

p

, the algorithm consists of a combination of matrix inversion and matrix multiplications.

For a normal matrix

A

(for which

A^{H} A = A A^{H}

) and non-integer

p

, the Schur decomposition is diagonal and the algorithm reduces to evaluating powers of the eigenvalues of

A

and then constructing

A^{p}

using the Schur vectors. This should give a very accurate result. In general however, no error bounds are available for the algorithm.

Further Comments

The cost of the algorithm is

O (n^{3})

. The exact cost depends on the matrix

A

but if

p \in (- 1, 1)

then the cost is independent of

p

O (4 \times n^{2})

complex allocatable memory is required by the function.

If estimates of the condition number of

A^{p}

are required then nag_matop_complex_gen_matrix_cond_pow (f01ke) should be used.

Example

This example finds

A^{p}

where

p = 0.2

and

A = (\begin{array}{r} 2 & 3 & 2 & 1 + 3 i \\ 2 + i & 1 & 1 & 2 + 2 i \\ 2 + i & 2 + 2 i & 2 i & 2 + 4 i \\ 3 & 2 + 2 i & 3 & 1 \end{array}) .

Open in the MATLAB editor: f01fq_example

function f01fq_example


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

% Principal power p of matrix A

a = [ 2    3     2     1+3i;
      2+i  1     1     2+2i;
      2+i  2+2i    2i  2+4i;
      3    2+2i  3     1];

p = 0.2;

[pa, ifail] = f01fq(a,p);

disp('A^p:');
disp(pa);

f01fq example results

A^p:
   1.1766 - 0.0758i   0.1375 + 0.2241i   0.2742 - 0.2223i  -0.1435 + 0.0816i
   0.2074 - 0.0998i   1.1118 - 0.0039i  -0.1343 - 0.0404i   0.1794 + 0.3590i
  -0.0859 - 0.0824i   0.5224 - 0.0530i   1.0616 + 0.3921i   0.2308 + 0.1856i
   0.3313 + 0.1303i  -0.1507 + 0.0982i   0.2178 - 0.0061i   1.1710 - 0.2136i

PDF version (NAG web site, 64-bit version, 64-bit version)

Chapter Contents

Chapter Introduction

NAG Toolbox