/*矩阵乘法模板题给定矩阵A,请快速计算出A^n(n个A相乘)的结果,输出的每个数都mod p。 由于矩阵乘法具有结合律,因此A^4 = A * A * A * A = (A*A) * (A*A) = A^2 * A^2。我们可以得到这样的结论:当n为偶数时,A^n = A^(n/2) * A^(n/2);当n为奇数时,A^n = A^(n/2) * A^(n/2) * A (其中n/2取整)。这就告诉我们,计算A^n也可以使用二分快速求幂的方法。例如,为了算出A^25的值,我们只需要递归地计算出A^12、A^6、A^3的值即可。根据这里的一些结果,我们可以在计算过程中不断取模,避免高精度运算。递归实现POW函数Matrix POW( Matrix t,int k ){ if( k == 1 ) return t; Matrix t1 = POW( t, k/2 ); t1 = t1*t1; if( k & 1 ) return t1 * t; else return t1;}递归的容易理解,但时间花费较多。*/#include <stdio.h> #include <stdlib.h> #include <string.h> typedef struct Node { int m[11][11]; }Matrix; Matrix init, unit; //初始化输入矩阵,单位矩阵如果用递归写Pow函数可以不用单位矩阵 int n, K; void Init( ) //初始化 { scanf( "%d%d", &n, &K ); for( int i=0; i<n; ++ i ) for( int j=0; j<n; ++ j ) { scanf( "%d", &init.m[i][j] ); unit.m[i][j]=( i == j ); } } Matrix Mul( Matrix a, Matrix b ) //矩阵乘法 { Matrix c; for( int i=0; i<n; ++ i ) { for( int j=0; j<n; ++ j ) { c.m[i][j]=0; //特别注意 for( int k=0; k<n; ++ k ) { c.m[i][j] += a.m[i][k]*b.m[k][j]; } c.m[i][j] %= 9973; } } return c; } Matrix Pow( Matrix a, Matrix b, int k ) { while( k>1 ) { if( k&1 ) // k为奇数时 { k --; b=Mul( a, b ); } else // k为偶数 { k >>= 1; a=Mul( a, a ); } } a=Mul( a, b ); return a; } int main( ) { int T; scanf( "%d", &T ); while( T -- ) { Matrix x; Init( ); x=Pow( init, unit, K ); int sum=0, i=0; n--; while( n >= 0 ) { sum += x.m[n][n]; sum%=9973; n --; } printf( "%d\n", sum%9973 ); } }
A Simple Math Problem
/*If x < 10 f(x) = x.If x >= 10 f(x) = a0 * f(x-1) + a1 * f(x-2) + a2 * f(x-3) + …… + a9 * f(x-10);|f(10) | |a0 a1 a2 ...a8 a9| |f(9)|| f(9) | | 1 0 0 ... 0 0| |f(8)|| .....| = |.. ... ... ... ..| | .. || f(2) | | 0 0 0 ... 0 0| |f(1)|| f(1) | | 0 0 0 ... 1 0| |f(0)|另A举证为10*10的举证,如上图。可以推出:(f(n),f(n-1),...,f(n-9))^(-1) = A^(n-9)*(f(9),f(8),...,f(0))^(-1)*/#include<iostream>using namespace std;__int64 k,m;struct mat{ int mar[10][10];}a,b,tmp;mat matrixmul(mat a,mat b){ int i,j,k; for(i = 0;i < 10;i++) for(j = 0;j < 10;j++) { tmp.mar[i][j] = 0; for(k = 0;k < 10;k++) tmp.mar[i][j] += (a.mar[i][k] * b.mar[k][j]) % m; tmp.mar[i][j] %= m; } return tmp;}void matrix_binary(){ while(k) { if(k & 1) b = matrixmul(b,a); a = matrixmul(a,a); k = k >> 1; }}int main(){ int i; while (scanf("%I64d%I64d",&k,&m) != EOF) { memset(a.mar,0,sizeof(a.mar)); for(i = 1;i < 10;i++) a.mar[i][i-1] = 1; memset(b.mar,0,sizeof(b.mar)); for(i = 0;i < 10;i++) b.mar[i][i] = 1; for(i = 0;i < 10;i++) scanf("%d",&a.mar[0][i]); if(k < 10) { printf("%d\n", k % m); continue; } k -= 9; matrix_binary(); int res = 0; for(i = 0;i < 10;i++) res += (b.mar[0][i] * (9-i)) % m; printf("%d\n",res%m); } return 0;}
Pendant
/*本题从表面上看是排列组合题,但要推出公式还是有相当难度。所以想到用DP来做。方程可通过枚举几种情况来推出。 以F[i][j]表示长度为i的pendant,用了j种珍珠,所构成的方案数,则F[i][j]=F[i-1][j]*j+F[i-1][j-1]*(k-j+1)。结果就是F[1][k]+…+F[n][k]。但注意到N的范围很大,申请那么大的数组会MLE。 注意到当前状态只与上一个状态有关,那么可以使用循环数组来,并累加上每个F[I][K]的方法来做。但这样的复杂度为O(NK),会TLE。 优化的方法是使用矩阵来做。将F[i-1]到F[i]的转移用矩阵来描述,相当于一个k*k的线性变换矩阵。因此F[i]=A*F[i-1],这里A是转移矩阵,即F[i]=Ai-1*F[1],所以F[1]+…+F[n]=A0*F[1]+…+An-1*F[1]=(E+A+A2+…+An-1)*F[1]。 F[i-1]这个矩阵是F[i-1][0] F[i-1][1] F[i-1][2] ...... F[i-1][k]0 0 0 …… 0 . . . . .. . . . .0 0 0 0 0(这是个K+1阶的阶阵) A矩阵是0 k 0 0 0 00 1 k-1 0 0 00 0 2 0 0 0. . . . . .. . . . . .0 0 0 0 k-1 10 0 0 0 0 k 这两个矩阵为什么是这样的,用矩阵的乘法来乘一下,将结果与DP方程对比下就知道了。接下来要解决的问题就是F[1]要怎样填写。每次再另写一个程序来填F[1]?大可不必。因为F[1]=F[0]*A,那么只要填F[0]就行了。将F[0]变为E(单位阵)就可以了。最后的问题就是求矩阵的平方和了。这里的方法很多,到网上一找也很容易找到,就不多说了。下面用的是二分求和,二分求幂。*/#include<iostream> void MIN(int& a,int b) { if( a>b ) a=b; }void MAX(int& a,int b) { if( a<b ) a=b; }#define CLR(NAME,VALUE) memset(NAME,VALUE,sizeof(NAME)) using namespace std; const int M=1234567891; int n,k;int a[31][31],b[31][31],c[31][31];int d[31][31],z[31][31]; void Mul(int x[31][31],int y[31][31]) { //矩阵乘法,结果放在X中 int i,j,jj; for(i=0;i<=k;++i) { for(j=0;j<=k;++j) { z[i][j]=0; for(jj=0;jj<=k;++jj) { z[i][j]=((__int64)z[i][j]+((__int64)x[i][jj]*y[jj][j])%M)%M; } } } for(i=0;i<=k;++i) { for(j=0;j<=k;++j) { x[i][j]=z[i][j]; } }} void Cal(int p) { //递归二分算和式:A^1+A^2+...+A^N if( p==0 ) return; Cal(p/2); int i,j; for(i=0;i<=k;++i) { for(j=0;j<=k;++j) { d[i][j]=b[i][j]; } } for(i=0;i<=k;++i) { //D=B+E d[i][i]=(d[i][i]+1)%M; } Mul(c,d); //C=C*D=C*B+C; Mul(b,b); //B=B*B; if( p&1 ) { Mul(b,a); //B=B*A for(i=0;i<=k;++i) { //C=C+B for(j=0;j<=k;++j) { c[i][j]=((__int64)c[i][j]+b[i][j])%M; } } }} int main() { int ca,i; scanf("%d",&ca); while( ca-- ) { scanf("%d%d",&n,&k); CLR(a,0); CLR(b,0); CLR(c,0); //将B矩阵变为单位阵 for(i=0;i<=k;++i) b[i][i]=1; //F[i][j]=F[i-1][j]*j+F[i-1][j-1]*(k-j+1) for(i=1;i<=k;++i) { //转移矩阵A a[i-1][i]=k-i+1; a[i][i]=i; } Cal(n); printf("%d\n",c[0][k]); } return 0;}
奥运
#include <iostream>#include <cstring>#include <cstdio>#define MAX 35using namespace std;struct matrix{ int num[MAX][MAX]; matrix() { memset(num,0,sizeof(num)); }};int n;long long map[MAX][MAX];matrix A[10005];long long city[35];matrix operator*(matrix &a,matrix &b){ matrix t; int i,j,k; for(i=0;i<n;i++) for(j=0;j<n;j++) { t.num[i][j]=0; for(k=0;k<n;k++) t.num[i][j]+=(a.num[i][k]*b.num[k][j]); t.num[i][j]%=2008; } return t;}void init(){ int i,j; for(i=0;i<n;i++) { for(j=0;j<n;j++) { A[0].num[i][j]=map[i][j]; } }}int main(){ int nn; int Find(int u); while(scanf("%d",&nn)!=EOF) { int i; memset(map,0,sizeof(map)); memset(city,0,sizeof(city)); n=0; for(i=0;i<nn;i++) { int u,v; scanf("%d%d",&u,&v); int lo1=Find(u); int lo2=Find(v); map[lo1][lo2]++; } int k; init(); for(i=1;i<10005;i++) A[i]=A[i-1]*A[0]; scanf("%d",&k); while(k--) { int v1,v2,t1,t2; scanf("%d%d%d%d",&v1,&v2,&t1,&t2); if(t1>t2) { int t=t1; t1=t2; t2=t; } int l=n; int l1=Find(v1); int l2=Find(v2); if(l1>=l||l2>=l) { printf("0\n"); n=l; continue; } int sum=0; for(i=t1-1;i<t2;i++) sum=(sum%2008+A[i].num[l1][l2]%2008)%2008; printf("%d\n",sum); } } return 0;}int Find(int u){ int i; for(i=0;i<n;i++) { if(city[i]==u) break; } if(i>=n) { n++; city[i]=u; return i; } else { return i; }}
Kiki & Little Kiki 2
/*题目大意:题目大意给定一系列灯的初始状态,0代表暗,1代表亮,每一秒所有的灯都有可能发生状态切换,切换规则:当前灯的左边第一个灯是亮的,则当前的灯切换状态,如果当前灯的左边第一盏灯是暗的,则当前灯的状态无需变化!注意:最左边的参考左右边那栈灯。题目分析;首先有两种情况:左边第一盏灯亮的: 当前灯的动作: 1->0; 0->1;左边第一盏灯暗的: 当前灯的动作: 1->1; 0->0;我们可以看到的是, 可以用一个方程来模拟这个过程: F[i] = ( f [ i] +f [i+n-2]%n+1 )%2;所以我们只要计算这个值就OK啦。然后由于数据过大,开数组肯定会爆掉~这里我们要想到的是 矩阵是一个模拟递归的高效算法这里我们要构造一个 可以计算如上的方程的矩阵: 1 0 0...0 1 1 1 0...0 0 0 1 1..0 0 0 0 1..0 0 . . . 0 .... . . .0..... 0 0 0..0 1 */#include <iostream> #include <cstdio> #include <cstring> using namespace std; const int N = 105; const int mod = 2; struct Mat { int num[N][N]; Mat() { for(int i = 0; i < N; i++) for(int j = 0; j < N; j++) num[i][j] = 0; } }; Mat mul(Mat a, Mat b, int n) { Mat r; for(int i = 0; i < n; i++) for(int k = 0; k < n; k++) { if(a.num[i][k] == 0) continue; for(int j = 0; j < n; j++) { if(b.num[k][j] == 0) continue; r.num[i][j] = (r.num[i][j] + a.num[i][k] * b.num[k][j]) % mod; } } return r; } Mat mal(Mat init, Mat unit, int m, int n) { while(m) { if(m & 1) { unit = mul(init, unit, n); m--; } else { init = mul(init, init, n); m >>= 1; } } return unit; } int main() { int m; while(cin >> m) { char str[N]; cin >> str; Mat init, unit; int len = strlen(str); for(int i = 0; i < len; i++) { unit.num[i][0] = str[i] - '0'; } init.num[0][0] = 1; init.num[0][len - 1] = 1; int k = 0; for(int i = 1; i < len; i++) { for(int j = k; j < k + 2; j++) init.num[i][j] = 1; k++; if(k == len - 1) break; } Mat ans; ans = mal(init, unit, m, len); for(int i = 0 ; i < len; i++) cout << ans.num[i][0]; cout << endl; } }
Tower
/*由 a(n) = 2*a(2)*a(n-1) - a(n-2)let p = 2*a(2);==>a(n) =p*a(n-1) - a(n-2)==>a(n)^2 = p^2*a(n-1)^2 - 2*p*a(n-1)*a(n-2) + a(n-2)^2................(1)let s(n) =sum( a(i) );s(n) = s(n-1) + a(n)^2 ................................................................... (2)a(n)*a(n-1) = p*a(n-1)^2 - a(n-1)*a(n-2)...................(3)so we have:a(n-1)^2 0 , 1 , 0 , 0 a(n-2)^2a(n)^2 = 1 , p^2, -2*p, 0 * a(n-1)^2a(n)*a(n-1) 0 , p , -1 , 0 a(n-1)*a(n-2)s(n-1) 0 , 1 , 0 , 1 s(n-2)use matrix to do this.Problem is sovle now.*/#include<cstdio> #define FOR(i,a,b) for(int (i)=(a);(i)<=(b);(i)++) #define nMax 10 #define LL long long LL mod; template<typename T> struct Matrix{ T mtix[nMax][nMax]; int n; Matrix(int n=1):n(n) { FOR(i,0,n-1) FOR(j,0,n-1) mtix[i][j] = 0; FOR(i,0,n-1) mtix[i][i] = 1; } friend Matrix operator * (const Matrix& a,const Matrix& b) { Matrix c(a.n); FOR(i,0,a.n-1) FOR(j,0,a.n-1) { c.mtix[i][j] = 0; FOR(k,0,a.n-1) if(a.mtix[i][k] && b.mtix[k][j]){ c.mtix[i][j] += a.mtix[i][k] * b.mtix[k][j] ; c.mtix[i][j]%=mod; } } return c; } friend Matrix Exp(Matrix a,int n){ Matrix c(a.n); Matrix b = a; while(n){ if(n&1) c=c*b; n >>= 1; b = b*b; } return c; } void init(T p) { p%=mod; FOR(i,0,n-1) FOR(j,0,n-1) mtix[i][j] = 0 ; mtix[1][1] = p*p%mod; mtix[0][1] = mtix[1][0] = mtix[3][1] = mtix[3][3] = 1; mtix[1][2] = ((-2*p) % mod + mod) % mod; mtix[2][1] = p; mtix[2][2] = ((-1 ) % mod + mod) % mod; } }; LL b[5]; int main(){ //freopen("input.txt","r",stdin); int n; LL a; int t; scanf("%d",&t); while(t--){ scanf("%I64d%d%I64d",&a,&n,&mod); a%=mod; if(mod == 1){ printf("0\n"); continue; } Matrix<LL> s(4); s.init(a*2); b[1]=a*a%mod,b[0]=1LL,b[2]=a,b[3]=1; s = Exp( s ,n-1 ); LL ans = 0LL; FOR(i,0,3) ans = (ans+s.mtix[3][i]*b[i]) % mod; printf("%I64d\n",ans); } return 0; }
Lucky Coins Sequence
//参考文献http://blog.sina.com.cn/s/blog_83ccc39d0101jv6d.html//代码参考http://blog.csdn.net/acm_davidcn/article/details/5819193//矩阵快速幂#include<iostream>#define MOD 10007using namespace std;int f[4][4], res[4][4];void A(int a[4][4], int b[4][4]){int t[4][4];memset(t, 0, sizeof(t));int i, j, k;for( i=1; i < 4; i++ )for( j=1; j < 4; j++ )for( k=1; k < 4; k++ )t[i][j] += a[i][k] * b[k][j] % MOD;for( i=1; i < 4; i++ )for( j=1; j < 4; j++ )a[i][j] = t[i][j];}int main(){int n;int a[5] = {0, 0, 0, 2, 6};while(scanf("%d", &n) != EOF ){if(n <= 4) { printf("%d\n", a[n]); continue;}memset(f, 0, sizeof(f));memset(res, 0, sizeof(res));f[1][1] = f[1][2] = f[1][3] = f[2][1] = 1;f[3][3] = 2;res[1][1] = res[2][2] = res[3][3] = 1;n -= 4;while(n > 0){if(n % 2 == 0)n /=2, A(f, f);elsen--, A(res, f);}printf("%d\n", (res[1][1]*6 + res[1][2]*2 + res[1][3]*8) % MOD);}return 0;}
Buge's Fibonacci Number Problem
#include<cstdio>#include<iostream>#include<cstring>using namespace std;const int mm=55;__int64 x[mm][mm],y[mm][mm],s[mm][mm],A[mm],B[mm],F1[mm],F2[mm],tmp;int f1,f2,a,b,K,n,m;int i,j,k,t;void mul(__int64 x[mm][mm],__int64 y[mm][mm]){ int i,j,k; for(i=0;i<K+2;++i) for(j=0;j<K+2;++j) for(s[i][j]=k=0;k<K+2;++k) { s[i][j]+=x[i][k]*y[k][j]; if(s[i][j]>m)s[i][j]%=m; } for(i=0;i<K+2;++i) for(j=0;j<K+2;++j) x[i][j]=s[i][j];}int main(){ scanf("%d",&t); while(t--) { scanf("%d%d%d%d%d%d%d",&f1,&f2,&a,&b,&K,&n,&m); A[0]=B[0]=F1[0]=F2[0]=1; for(i=1;i<K+2;++i) { A[i]=A[i-1]*a; B[i]=B[i-1]*b; F1[i]=F1[i-1]*f1; F2[i]=F2[i-1]*f2; if(A[i]>=m)A[i]%=m; if(B[i]>=m)B[i]%=m; if(F1[i]>=m)F1[i]%=m; if(F2[i]>=m)F2[i]%=m; } memset(x,0,sizeof(x)); memset(y,0,sizeof(y)); for(i=0;i<K+2;++i)x[i][i]=1; y[0][0]=y[0][1]=1; for(i=1;i<K+2;++i) for(j=1,tmp=1,k=K-i+2;k--;++j) { y[i][j]=(tmp%m)*A[k]; if(y[i][j]>=m)y[i][j]%=m; y[i][j]*=B[j-1]; if(y[i][j]>=m)y[i][j]%=m; tmp*=k; tmp/=j; } --n; while(n) { if(n&1)mul(x,y); mul(y,y); n>>=1; } tmp=F1[K]*x[0][0]; if(tmp>=m)tmp%=m; for(i=1;i<K+2;++i) { tmp+=((F2[K-i+1]*F1[i-1])%m)*x[0][i]; if(tmp>=m)tmp%=m; } printf("%I64d\n",tmp); } return 0;}