zhangas

Pay more attention

0%

Matrix-Exponentiation

Let’s first consider below simple question.

What is the minimun complexity to find n’th Fibonacci Number?

We can find n’th Fibonacci Number in O(log n) time using $Matrix\quad Exponentiation$. In this post, a general implementation of Matrix Exponentiation is discussed.

For solving the problem, we are assuming a linear recurrence eqution like below:
F(n) = a * F(n - 1) + b * F(n - 2) + c * F(n - 3) for n >= 3  ....Equation (1)
where a, b, c are constants.

For this recurrence relation, it depends on three previous value.
Now we will try to represent Equation(1) in terms of the matrix.
|F(n)|                    |F(n - 1)|
|F(n - 1)| = Matrix 'C' * |F(n - 2)|
|F(n - 2)|                |F(n - 3)|
Now we need to fill the Matrix 'C'.
so according to our equation.
             |a b c|
Matrix 'C' = |1 0 0|
             |0 1 0|
so for n, the equation(1) changes to 
|F(n)|                     |F(2)|
|F(n - 1)| = C ^ (n - 2) * |F(1)|
|F(n - 2)|                 |F(0)|
So we can simply multiply our matrix 'C' n - 2 times and the multiply it with the matrix below to get the result.
|F(2)|
|F(1)|
|F(0)|
Multiplication can be done in O(Log n) time using Divie and Conquer algorithm for fast power.
There is rough code below.


typedef struct Mat{
    int m[105][105];
}Mat;
Mat a;
Mat multiply(Mat x, Mat y){
    // Creating an new matrix to store elements 
    // of the multiplication matrix
    Mat c;
    for(int i = 1; i <= 3; ++ i){
        for(int j = 1; j <= 3; ++ j){
            c.m[i][j] = 0;
            for(int k = 1; k <= 3; ++ k)
                c.m[i][j] += x.m[i][k] * y.m[k][j];
        }
    }
    return c;
    // return the result
}
Mat Matrix_fastpow(Mat a, Mat k){
    //initial values for matrix 'C'
    Mat res;
    res.m[1][1] = a;res.m[1][2] = b;res.m[1][3] = c;
    res.m[2][1] = 1;res.m[2][2] = 0;res.m[2][3] = 0;
    res.m[3][1] = 0;res.m[3][2] = 1;res.m[3][3] = 0;

    while(k){
        if(k & 1){
            res = multiple(a, res);
        }
        k >>= 1;
        a = multiple(a, a);
    }
    return res;
}