对于式子的值,我们若直接算出等比数列的和再模上M,那么中间结果可能会溢出。对于这个问题,我们可以用二分求和来做,这样一来,它的中间结果每次都在模M,任何时候都不会有溢出的危险。对上述式子进行如下化简:
根据,以上化简的结果,我们可以得出递归的代码:
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 #include <cstdio> 2 #include <cstring> 3 using namespace std; 4 const int MOD = 10000; 5 int Pow(int a, int n) 6 { 7 int r=1; 8 while(n) 9 { 10 if(n&1) 11 r = (r*a)%MOD; 12 a = (a*a)%MOD; 13 n >>= 1; 14 } 15 return r; 16 } 17 // 计算Sn = (a^1+a^2+...+a^n)%MOD 18 int Sum(int a, int n) 19 { 20 if(n==1) return a; 21 int t = Sum(a, n/2); 22 if(n&1) 23 { 24 int c = Pow(a, n/2+1); 25 t = (t%MOD+(t%MOD*c%MOD)%MOD+c%MOD)%MOD; 26 } 27 else 28 { 29 int c = Pow(a, n/2); 30 t = (t + t*c) % MOD; 31 } 32 return t; 33 } 34 int main() 35 { 36 int a, n; 37 while(~scanf("%d%d",&a, &n)) 38 { 39 printf("%d ",Sum(a, n)); 40 } 41 return 0; 42 }
相应的可以得出矩阵二分求和的代码:
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 #include <cstdio> 2 #include <cstring> 3 #include <algorithm> 4 using namespace std; 5 #define INF 0x3f3f3f3f 6 #define N 35 7 int MOD; 8 9 struct Matrix 10 { 11 int map[N][N]; 12 int n; 13 Matrix() 14 { 15 n=0; 16 memset(map, 0, sizeof map); 17 } 18 Matrix(int k) 19 { 20 n = 0; 21 memset(map, 0, sizeof map); 22 for(int i=0; i<N; i++) map[i][i] = k; 23 } 24 Matrix operator * (const Matrix &a) const 25 { 26 Matrix e; 27 e.n = n; 28 for(int i=0; i<n; i++) 29 for(int j=0; j<n; j++) 30 for(int k=0; k<n; k++) 31 e.map[i][j] = (e.map[i][j]+a.map[i][k]*map[k][j])%MOD; 32 return e; 33 } 34 Matrix operator + (const Matrix &a) const 35 { 36 Matrix e; 37 e.n = n; 38 for(int i=0; i<n; i++) 39 for(int j=0; j<n; j++) 40 e.map[i][j] = (a.map[i][j] + map[i][j])%MOD; 41 return e; 42 } 43 void print() const 44 { 45 for(int i=0; i<n; i++) 46 { 47 for(int j=0; j<n-1; j++) 48 printf("%d ",map[i][j]); 49 printf("%d ",map[i][n-1]); 50 } 51 } 52 }; 53 54 Matrix Pow(Matrix a, int b) 55 { 56 Matrix e(1); 57 e.n = a.n; 58 while(b) 59 { 60 if(b&1) 61 e = a * e; 62 a = a * a; 63 b >>= 1; 64 } 65 return e; 66 } 67 68 // 计算 Sn = (a+a^2+a^3+...+a^n)%MOD; 其中a为矩阵. 69 Matrix Sum(Matrix a, int n) 70 { 71 if(n==1)return a; 72 Matrix t = Sum(a, n/2); 73 if(n&1) 74 { 75 Matrix c = Pow(a, n/2+1); 76 t = t + (t * c); 77 t = t + c; 78 } 79 else 80 { 81 Matrix c = Pow(a, n/2); 82 t = t + (t * c); 83 } 84 return t; 85 } 86 87 int main() 88 { 89 int n, m, k; 90 while(~scanf("%d%d%d",&n, &k, &m) && n+m) 91 { 92 MOD = m; 93 Matrix A, B; 94 A.n = n; 95 for(int i=0; i<n; i++) 96 for(int j=0; j<n; j++) 97 scanf("%d",&A.map[i][j]); 98 B = Sum(A, k); 99 B.print(); 100 } 101 return 0; 102 }