matrix_exponentiation.cpp File Reference

Matrix Exponentiation. More...

#include <iostream>
#include <vector>

Include dependency graph for matrix_exponentiation.cpp:

Macros
#define	ll   int64_t
#define	endl   std::endl
#define	pb   push_back
#define	MOD   1000000007

Functions
ll	ab (ll x)
vector< vector< ll > >	multiply (const vector< vector< ll >> &A, const vector< vector< ll >> &B)
vector< vector< ll > >	power (const vector< vector< ll >> &A, ll p)
ll	ans (ll n)
int	main ()

Variables
ll	k
vector< ll >	a
vector< ll >	b
vector< ll >	c

Detailed Description

Matrix Exponentiation.

The problem can be solved with DP but constraints are high.
\(a_i = b_i\) (for \(i <= k\))
\(a_i = c_1 a_{i-1} + c_2 a_{i-2} + ... + c_k a_{i-k}\) (for \(i > k\))
Taking the example of Fibonacci series, \(k=2\)
\(b_1 = 1,\; b_2=1\)
\(c_1 = 1,\; c_2=1\)
\(a = \begin{bmatrix}0& 1& 1& 2& \ldots\end{bmatrix}\)
This way you can find the \(10^{18}\) fibonacci numberMOD. I have given a general way to use it. The program takes the input of B and C matrix.

Steps for Matrix Expo

1. Create vector F1 : which is the copy of B.
2. Create transpose matrix (Learn more about it on the internet)
3. Perform \(T^{n-1}\) [transpose matrix to the power n-1]
4. Multiply with F to get the last matrix of size (1 \(\times\)k).

The first element of this matrix is the required result.

Macro Definition Documentation

endl

#define endl   std::endl

shorthand definition for std::endl

ll

#define ll   int64_t

shorthand definition for int64_t

pb

#define pb   push_back

shorthand definition for int64_t

Function Documentation

ab()

ll ab ( ll  x )

inline

returns absolute value

{ return x > 0LL ? x : -x; }

ans()

ll ans

(

ll

n

)

Wrapper for Fibonacci

Parameters

[in]

n

\(n^\text{th}\) Fibonacci number

Returns

\(n^\text{th}\) Fibonacci number

             {
    if (n == 0)
        return 0;
    if (n <= k)
        return b[n - 1];
    // F1
    vector<ll> F1(k + 1);
    for (ll i = 1; i <= k; i++) F1[i] = b[i - 1];
 
    // Transpose matrix
    vector<vector<ll>> T(k + 1, vector<ll>(k + 1));
    for (ll i = 1; i <= k; i++) {
        for (ll j = 1; j <= k; j++) {
            if (i < k) {
                if (j == i + 1)
                    T[i][j] = 1;
                else
                    T[i][j] = 0;
                continue;
            }
            T[i][j] = c[k - j];
        }
    }
    // T^n-1
    T = power(T, n - 1);
 
    // T*F1
    ll res = 0;
    for (ll i = 1; i <= k; i++) {
        res = (res + (T[1][i] * F1[i]) % MOD) % MOD;
    }
    return res;
}

main()

int main

(

)

Main function

           {
    cin.tie(0);
    cout.tie(0);
    ll t;
    cin >> t;
    ll i, j, x;
    while (t--) {
        cin >> k;
        for (i = 0; i < k; i++) {
            cin >> x;
            b.pb(x);
        }
        for (i = 0; i < k; i++) {
            cin >> x;
            c.pb(x);
        }
        cin >> x;
        cout << ans(x) << endl;
        b.clear();
        c.clear();
    }
    return 0;
}

multiply()

vector<vector<ll> > multiply	(	const vector< vector< ll >> &	A,
		const vector< vector< ll >> &	B
	)

To multiply 2 matrices

Parameters

[in]	A	matrix 1 of size (m \(\times\)n)
[in]	B	matrix 2 of size (p \(\times\)q)

Note

\(p=n\)

Returns

matrix of dimension (m \(\times\)q)

                                                         {
    vector<vector<ll>> C(k + 1, vector<ll>(k + 1));
    for (ll i = 1; i <= k; i++) {
        for (ll j = 1; j <= k; j++) {
            for (ll z = 1; z <= k; z++) {
                C[i][j] = (C[i][j] + (A[i][z] * B[z][j]) % MOD) % MOD;
            }
        }
    }
    return C;
}

power()

vector<vector<ll> > power	(	const vector< vector< ll >> &	A,
		ll	p
	)

computing integer power of a matrix using recursive multiplication.

Note

A must be a square matrix for this algorithm.

Parameters

[in]	A	base matrix
[in]	p	exponent

Returns

matrix of same dimension as A

                                                            {
    if (p == 1)
        return A;
    if (p % 2 == 1) {
        return multiply(A, power(A, p - 1));
    } else {
        vector<vector<ll>> X = power(A, p / 2);
        return multiply(X, X);
    }
}

Here is the call graph for this function:

Variable Documentation

a

vector<ll> a

global vector variables

Todo:

@stepfencurryxiao add documetnation

k

ll k

global variable k

Todo:

@stepfencurryxiao add documetnation