-
-
Notifications
You must be signed in to change notification settings - Fork 7.3k
/
matrix_exponentiation.cpp
148 lines (132 loc) · 3.85 KB
/
matrix_exponentiation.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
/**
@file
@brief Matrix Exponentiation.
The problem can be solved with DP but constraints are high.
<br/>\f$a_i = b_i\f$ (for \f$i <= k\f$)
<br/>\f$a_i = c_1 a_{i-1} + c_2 a_{i-2} + ... + c_k a_{i-k}\f$ (for \f$i > k\f$)
<br/>Taking the example of Fibonacci series, \f$k=2\f$
<br/>\f$b_1 = 1,\; b_2=1\f$
<br/>\f$c_1 = 1,\; c_2=1\f$
<br/>\f$a = \begin{bmatrix}0& 1& 1& 2& \ldots\end{bmatrix}\f$
<br/>This way you can find the \f$10^{18}\f$ fibonacci number%MOD.
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 \f$T^{n-1}\f$ [transpose matrix to the power n-1]
4. Multiply with F to get the last matrix of size (1\f$\times\f$k).
The first element of this matrix is the required result.
*/
#include <iostream>
#include <vector>
using std::cin;
using std::cout;
using std::vector;
/*! shorthand definition for `int64_t` */
#define ll int64_t
/*! shorthand definition for `std::endl` */
#define endl std::endl
/*! shorthand definition for `int64_t` */
#define pb push_back
#define MOD 1000000007
/** global variable mat_size
* @todo @stepfencurryxiao add documetnation
*/
ll mat_size;
/** global vector variables used in the ::ans function.
* @todo @stepfencurryxiao add documetnation
*/
vector<ll> fib_b, fib_c;
/** To multiply 2 matrices
* \param [in] A matrix 1 of size (m\f$\times\f$n)
* \param [in] B \p matrix 2 of size (p\f$\times\f$q)\n\note \f$p=n\f$
* \result matrix of dimension (m\f$\times\f$q)
*/
vector<vector<ll>> multiply(const vector<vector<ll>> &A,
const vector<vector<ll>> &B) {
vector<vector<ll>> C(mat_size + 1, vector<ll>(mat_size + 1));
for (ll i = 1; i <= mat_size; i++) {
for (ll j = 1; j <= mat_size; j++) {
for (ll z = 1; z <= mat_size; z++) {
C[i][j] = (C[i][j] + (A[i][z] * B[z][j]) % MOD) % MOD;
}
}
}
return C;
}
/** computing integer power of a matrix using recursive multiplication.
* @note A must be a square matrix for this algorithm.
* \param [in] A base matrix
* \param [in] p exponent
* \return matrix of same dimension as A
*/
vector<vector<ll>> power(const vector<vector<ll>> &A, ll p) {
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);
}
}
/*! Wrapper for Fibonacci
* \param[in] n \f$n^\text{th}\f$ Fibonacci number
* \return \f$n^\text{th}\f$ Fibonacci number
*/
ll ans(ll n) {
if (n == 0)
return 0;
if (n <= mat_size)
return fib_b[n - 1];
// F1
vector<ll> F1(mat_size + 1);
for (ll i = 1; i <= mat_size; i++) F1[i] = fib_b[i - 1];
// Transpose matrix
vector<vector<ll>> T(mat_size + 1, vector<ll>(mat_size + 1));
for (ll i = 1; i <= mat_size; i++) {
for (ll j = 1; j <= mat_size; j++) {
if (i < mat_size) {
if (j == i + 1)
T[i][j] = 1;
else
T[i][j] = 0;
continue;
}
T[i][j] = fib_c[mat_size - j];
}
}
// T^n-1
T = power(T, n - 1);
// T*F1
ll res = 0;
for (ll i = 1; i <= mat_size; i++) {
res = (res + (T[1][i] * F1[i]) % MOD) % MOD;
}
return res;
}
/** Main function */
int main() {
cin.tie(0);
cout.tie(0);
ll t;
cin >> t;
ll i, j, x;
while (t--) {
cin >> mat_size;
for (i = 0; i < mat_size; i++) {
cin >> x;
fib_b.pb(x);
}
for (i = 0; i < mat_size; i++) {
cin >> x;
fib_c.pb(x);
}
cin >> x;
cout << ans(x) << endl;
fib_b.clear();
fib_c.clear();
}
return 0;
}