H: Magical Balls
总Time Limit: 1000ms Memory Limit: 262144kB
Description
Wenwen has a magical ball. When put on an infinite plane, it will keep duplicating itself forever.
Initially, Wenwen puts the ball on the location (x0, y0) of the plane. Then the ball starts to duplicate itself right away. For every unit of time, each existing ball on the plane will duplicate itself, and the new balls will be put on the adjacent locations. The duplication rule of these balls is, during the i-th unit of time, a ball, which locates at (x, y), will duplicate ui balls to (x, y+1), di balls to (x, y-1), li balls to (x-1, y) and ri balls to (x+1, y).
The duplication rule has a period of M. In other words, ui=ui-M, di=di-M, li=li-M, ri=ri-M for i=M+1, M+2, …
Wenwen is very happy because she will get many balls. It is easy to calculate how many balls she will get after N units of time. However, she wants to know the sum of x-coordinates and y-coordinates of all balls after N units of time. This is a bit difficult for her. Could you help her? Since the sum might be very large, you should give the sum modulo 1,000,000,007 to her.
Input
The first line contains an integer T (1 ≤ T ≤ 25), indicating the number of test cases.
For each test case:
The first line contains four integers N (1 ≤ N ≤ 10^18), M (1 ≤ M ≤ 20,000), x0 and y0 (-10^18 ≤ x0, y0 ≤ 10^18);
Then follows M lines, the i-th line contains four integers: ui, di, li and ri (0 ≤ ui, di, li, ri ≤ 10,000).
Output
For each test case, output one integer on a single line, indicating the sum of x-coordinates and y-coordinates of all balls after N units of time, modulo 1,000,000,007.
Sample Input
1
2 2 1 1
2 0 0 0
0 0 0 1
Sample Output
19
Hint
In the Sample Input:
Initially, there is 1 ball on (1,1).
After 1 unit of time, there is 1 ball on (1,1) and 2 balls on (1,2);
After 2 units of time, there is 1 ball on (1,1), 2 balls on (1,2), 1 ball on (2,1) and 2 balls on (2,2).
Therefore, after 2 units of time, the sum of x-coordinates and y-coordinates of all balls is
(1+1)*1+(1+2)*2+(2+1)*1+(2+2)*2=19.
Solution
首先我们来考察题目里蕴含的数学关系,用 num[i] 表示时间 i 时平面上球的总数,用 sum[i] 表示时间 i 时平面上所有球横纵坐标的总和,那么根据平面上每一个球都会分别向四个方向复制,可以得到以下递推式
$\qquad \mathrm{num}[i] = \mathrm{num}[i - 1] * (1 + u[i] + d[i] + r[i] + l[i])$
$\qquad \mathrm{sum}[i] = \mathrm{sum}[i - 1] * (1 + u[i] + d[i] + r[i] + l[i]) + \mathrm{num}[i - 1] * (u[i] - d[i] + r[i] - l[i])$
为简便我们记
$\qquad A[i] = 1 + u[i] + d[i] + r[i] + l[i]$
$\qquad B[i] = u[i] - d[i] + r[i] - l[i]$
那么递推式就可以写成
$\qquad \mathrm{num}[i] = \mathrm{num}[i - 1] * A[i]$
$\qquad \mathrm{sum}[i] = \mathrm{sum}[i - 1] * A[i] + \mathrm{num}[i - 1] * B[i]$
边界条件为 num[0] = 1, sum[0] = x0 + y0。于是可以进行线性递推,复杂度为O(n),考虑到 n<=10^18 且操作具有周期性,我们考虑矩阵快速幂优化,将递推式写成
$\qquad\displaystyle\left(\begin{array}{c}\mathrm{num}[i]\\ \mathrm{sum}[i]\end{array}\right)=\left(\begin{array}{cc}A[i] & 0 \\ B[i] & A[i]\end{array}\right)\left(\begin{array}{c}\mathrm{num}[i-1]\\ \mathrm{sum}[i-1]\end{array}\right)$
于是可以先计算出一个周期(m个单位时间)内的矩阵,然后利用快速幂计算其 [n/m] 次方,最后把剩下的 n%m 次操作乘上去,这样就把递推过程的复杂度降为O(log(n/m))。
务必要注意,初始坐标(x, y)可能为负,所以取模的时候要再加一次再取模。
Code
1 /* Problem: pkuCampus2016H: Magical Balls 2 Category: Recursion with Matrix multiplication 3 Author: Niwatori 4 Date: 2016/07/18 */ 5 6 #include <stdio.h> 7 #define MOD 1000000007 8 #define MAXN 20005 9 10 struct Matrix{ 11 long long c[2][2]; 12 Matrix(long long a1 = 1, long long a2 = 0, long long a3 = 0, long long a4 = 1) 13 {c[0][0] = a1; c[0][1] = a2; c[1][0] = a3; c[1][1] = a4;} 14 long long & operator()(int i, int j) {return c[i][j];} 15 }; 16 17 Matrix mul(Matrix a, Matrix b) 18 { 19 Matrix tmp(0, 0, 0, 0); 20 for (int i = 0; i < 2; ++i) 21 for (int j = 0; j < 2; ++j) 22 for (int k = 0; k < 2; ++k) 23 tmp(i, j) = (tmp(i, j) + a(i, k) * b(k, j)) % MOD; 24 return tmp; 25 } 26 27 Matrix power(Matrix x, long long n) // 快速幂 28 { 29 Matrix ans, tmp = x; 30 while (n) 31 { 32 if (n & 1) ans = mul(ans, tmp); 33 tmp = mul(tmp, tmp); 34 n = n >> 1; 35 } 36 return ans; 37 } 38 39 int main() 40 { 41 int t; scanf("%d", &t); 42 while (t--) 43 { 44 long long n, m, x, y, u, d, l, r, A[MAXN], B[MAXN]; 45 scanf("%lld%lld%lld%lld", &n, &m, &x, &y); 46 x = (x % MOD + MOD) % MOD; // 注意加模再取模 47 y = (y % MOD + MOD) % MOD; 48 for (int i = 0; i < m; ++i) 49 { 50 scanf("%lld%lld%lld%lld", &u, &d, &l, &r); 51 A[i] = 1 + u + d + l + r; 52 B[i] = u - d + r - l; 53 } 54 55 Matrix base; 56 for (int i = 0; i < m; ++i) 57 base = mul(Matrix(A[i], 0, B[i], A[i]), base); 58 base = power(base, n / m); 59 for (int i = 0; i < n % m; ++i) 60 base = mul(Matrix(A[i], 0, B[i], A[i]), base); 61 62 printf("%lld\n", (base(1, 0) + base(1, 1) * (x + y)) % MOD); 63 } 64 return 0; 65 }