Is there any fast method of matrix exponentiation?

Question

Is there any faster method of matrix exponentiation to calculate Mⁿ (where M is a matrix and n is an integer) than the simple divide and conquer algorithm?

Hey i found one link in stackoverflow only check it out http://stackoverflow.com/questions/12268516/matrix-exponentiation-using-fermats-theorem — , Sep 07 '12 at 04:57
Expokit is a well known package for performing matrix exponentiations. http://fortranwiki.org/fortran/show/Expokit — Sayan, Sep 07 '12 at 05:00

score 37 · Accepted Answer · edited Apr 24 '22 at 16:10

37

You could factor the matrix into eigenvalues and eigenvectors. Then you get

M = V * D * V^-1

Where V is the eigenvector matrix and D is a diagonal matrix. To raise this to the Nth power, you get something like:

M^n = (V * D * V^-1) * (V * D * V^-1) * ... * (V * D * V^-1)
    = V * D^n * V^-1

Because all the V and V^-1 terms cancel.

Since D is diagonal, you just have to raise a bunch of (real) numbers to the nth power, rather than full matrices. You can do that in logarithmic time in n.

Calculating eigenvalues and eigenvectors is r^3 (where r is the number of rows/columns of M). Depending on the relative sizes of r and n, this might be faster or not.

edited Apr 24 '22 at 16:10

Marek Fiołka

4,825
1
5
20

answered Sep 07 '12 at 05:04

Michael

617
4
3

As far as i know this method has the same complexity as Exponentiation by Squaring. So is there any faster method ? – Akashdeep Saluja Sep 07 '12 at 05:12
1

@AkashdeepSaluja: this is faster then exponentiation by squaring. This is O(r^3) time, exponentiation by squaring is O(r^3 logn) time. – Keith Randall Sep 07 '12 at 06:43
for a better explanation of the method mentioned above http://www.google.co.in/url?sa=t&rct=j&q=pdf%20nth%20power%20of%20matrix&source=web&cd=1&cad=rja&ved=0CCAQFjAA&url=http%3A%2F%2Fwww.qc.edu.hk%2Fmath%2FTeaching_Learning%2FNth%2520power%2520of%2520a%2520square%2520matrix.pdf&ei=Jf9JULrwFsi8rAejh4C4DQ&usg=AFQjCNE7yqQce5jdtyyVLFpSZmYUnoWyVA – Akashdeep Saluja Sep 07 '12 at 14:08
Very nice solution of the problem. Thanks a lot! I will use it in future! – Gloomcore Sep 11 '12 at 08:39
Isn't it necessary for the matrix to be symmetric/hermitian in order to be decomposed into eigenvectors?http://en.wikipedia.org/wiki/Spectral_theorem – pqnet Jun 16 '14 at 09:05
1

Not necessary, sufficient. – Valentin Waeselynck Aug 12 '14 at 22:31
But this isn't a general solution since spectral decomposition isn't possible for defective matrices. – jodag Sep 12 '15 at 08:20
Unfortunately, solving the characteristic equation itself will take logarithmic time because you'll need to binary search, so speed increase is negligible. – SinByCos Sep 20 '17 at 16:32
1

@SinByCos yes, but isn't it logarithmic in the size of the matrix? Squaring is logarithmic in the exponent, so you can't really compare the two. – Stefan Dragnev Jan 17 '18 at 10:07
1

You can always find the [Jordan normal form](https://en.wikipedia.org/wiki/Generalized_eigenvector#Jordan_normal_form) even for defective matrices. Then, D is not diagonal but the sum of a diagonal and a nilpotent matrix which you can still use very efficiently. – WorldSEnder Oct 27 '18 at 21:02
2

@WorldSEnder: Unfortunately the Jordan normal form is not numerically stable (the normal form is a discontinuous function of the matrix), so small rounding errors in computing the matrix can lead to large errors in the result. – Nate Eldredge Apr 13 '20 at 13:16

Gloomcore · Answer 2 · 2012-09-12T10:34:42.543

It's quite simple to use Euler fast power algorith. Use next algorith.

#define SIZE 10

//It's simple E matrix
// 1 0 ... 0
// 0 1 ... 0
// ....
// 0 0 ... 1
void one(long a[SIZE][SIZE])
{
    for (int i = 0; i < SIZE; i++)
        for (int j = 0; j < SIZE; j++)
            a[i][j] = (i == j);
}

//Multiply matrix a to matrix b and print result into a
void mul(long a[SIZE][SIZE], long b[SIZE][SIZE])
{
    long res[SIZE][SIZE] = {{0}};

    for (int i = 0; i < SIZE; i++)
        for (int j = 0; j < SIZE; j++)
            for (int k = 0; k < SIZE; k++)
            {
                res[i][j] += a[i][k] * b[k][j];
            }

    for (int i = 0; i < SIZE; i++)
        for (int j = 0; j < SIZE; j++)
            a[i][j] = res[i][j];
}

//Caluclate a^n and print result into matrix res
void pow(long a[SIZE][SIZE], long n, long res[SIZE][SIZE])
{
    one(res);

    while (n > 0) {
        if (n % 2 == 0)
        {
            mul(a, a);
            n /= 2;
        }
        else {
            mul(res, a);
            n--;
        }
    }
}

Below please find equivalent for numbers:

long power(long num, long pow)
{
    if (pow == 0) return 1;
    if (pow % 2 == 0)
        return power(num*num, pow / 2);
    else
        return power(num, pow - 1) * num;
}

score 4 · Answer 3 · answered Sep 07 '12 at 04:52

4

Exponentiation by squaring is frequently used to get high powers of matrices.

answered Sep 07 '12 at 04:52

MBo

77,366
5
53
86

i know this method but need to speed it up further. – Akashdeep Saluja Sep 07 '12 at 05:00
1

You'd better add this algorithm name into the question to avoid similar answers :) – MBo Sep 07 '12 at 05:19
Faster algorithm are much more complicated. – Ari Sep 07 '12 at 12:08

score 0 · Answer 4 · answered Sep 07 '12 at 05:09

0

I would recommend approach used to calculate Fibbonacci sequence in matrix form. AFAIK, its efficiency is O(log(n)).

answered Sep 07 '12 at 05:09

Rail Suleymanov

361
1
4
16

3

You have to multiply that by the cost of multiplying matrices. The overall running time is O (n^3 log n). – mrk Jun 13 '14 at 13:23

Is there any fast method of matrix exponentiation?

4 Answers4

Linked