您的位置:首页 > 其它

sgu 197 Nice Patterns Strike Back 矩阵快速幂

2013-12-03 14:03 351 查看
    这题如果n不大的话显然是个状压dp,但n的范围是10^100...但m的范围很小,最大只有5,状态最多也就是2^5,所以可以根据m暴力构造出2^m * 2^m 的状态转移矩阵A,然后A^(n-1),统计结果矩阵中每个数的和就是答案了,但是n太大了,所以还要写个高精/2,和高精--...太烦了- -...另外当m==1的时候,直接求2^n就是答案...

     #include <iostream>
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <string>
#include <cstring>
using namespace std;
const int maxn=10000+10;
const int N=11000;
int n,mod;
int m,p,q,tt;
int mat[50][50];
int a[10],b[10];
char num[1200];
int phi[11000];

struct Mat
{
int mat[50][50];
};
struct bnum
{
int dt[120];
void init(char* s,int len)
{

memset(dt,0,sizeof dt);
for (int i=len-1,j=1; i>=0; i--,j++)
dt[j]=s[i]-'0';
dt[0]=len;

dt[1]--;

for (int i=1; i<=dt[0]; i++)
if (dt[i]<0)
{
dt[i+1]--;
dt[i]+=10;
}

}
bool half()
{
int i,ad=0;
bool flag=(dt[1] & 1);
for (int i=dt[0]; i>=1; i--)
{
ad=ad*10+dt[i];
dt[i]=ad/2;
ad%=2;
}
return flag;
}
bool zero()
{
for (int i=1; i<=dt[0]; i++)
if (dt[i]!=0) return false;
return true;
}
}bx;
void Unit(Mat &a)
{
for(int i = 0; i < n; i++)
for(int j = 0; j < n; j++)
a.mat[i][j] = (i == j);
}

Mat operator*(Mat a, Mat b)
{
Mat c;
memset(c.mat, 0, sizeof c.mat);
for(int i = 0; i < n; i++)
{
for(int j = 0; j < n; j++)
{
for(int k = 0; k < n; k++)
{
if(a.mat[i][k] && b.mat[k][j])
{
c.mat[i][j] = (c.mat[i][j] + a.mat[i][k] * b.mat[k][j]);
if(c.mat[i][j] >= mod)
c.mat[i][j] %= mod;
}
}
}
}
return c;
}

Mat pw(Mat A)
{
Mat c;
Unit(c); // 初始化为单位矩阵
while(!bx.zero())
{
if(bx.half()) c = c * A;
A = A * A;
}
return c;
}
void euler_phi()
{
for(int i=1;i<N;i++) phi[i]=i;
for(int i=2;i<N;i++) if(phi[i]==i)
for(int j=i;j<N;j+=i)
phi[j]=phi[j]/i*(i-1);
}

bool check(int p,int q)
{
int tp=p,tq=q;
int cnt1=0,cnt2=0;
memset(a,0,sizeof a);
memset(b,0,sizeof b);
while (tp)
{
a[cnt1++]=(tp & 1);
tp>>=1;
}
while (tq)
{
b[cnt2++]=(tq & 1);
tq>>=1;
}
for (int i=0; i<tt-1; i++)
{
if (a[i]==a[i+1] && a[i]==b[i] && b[i]==b[i+1])
return false;
}

return true;
}
int pow_mod(int a,int b)
{
int ans=1;
while(b)
{
if(b%2==1) ans=ans*a%mod;
a=a*a%mod;
b/=2;
}
return ans % mod;
}

int main()
{
// freopen("in.txt","r",stdin);
// freopen("out.txt","w",stdout);
cin>>num>>tt>>mod;
int len=strlen(num);
int kk=0;
euler_phi();
int ph=phi[mod];
bool flag=false;

if (tt==1)
{
for (int i=0; i<len; i++)
{
kk*=10;
kk+=(num[i]-'0');
if (kk>=ph) flag=true;
kk %=ph;
}
cout<<pow_mod(2,kk)<<endl;
return 0;
}

bx.init(num,len);

memset(mat,0,sizeof mat);
int tot=(1<<tt);
n=tot;
for (int i=0; i<tot; i++)
for (int j=i; j<tot; j++)
if (check(i,j)) mat[i][j]=mat[j][i]=1;
else mat[i][j]=mat[j][i]=0;
Mat a,b;
for (int i=0; i<tot; i++)
for (int j=0; j<tot; j++)
a.mat[i][j]=mat[i][j];

// for (int i=0; i<tot; i++)
// {
// for (int j=0; j<tot; j++)
// cout<<mat[i][j]<<" ";
// cout<<endl;
// }
// cout<<endl;

b=pw(a);
int ans=0;
for (int i=0; i<n; i++)
for (int j=0; j<n; j++)
ans+=b.mat[i][j], ans%=mod;
cout<<ans<<endl;
return 0;
}


    
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: