思路
很容易想到\(O(3^{3m}\log T)\)的暴力大矩乘，显然过不了。
我们分析一下每次转移的性质。题目给的转移方程是填表法，我们试着改成刷表法看看……
发现好像没啥用。
注意到游戏的规则是1吃0，2吃1，0吃2，也就是在\(x-y=1\pmod 3\)的时候\(x\)吃\(y\)。
我们枚举\(j\)，然后再枚举\(i\ominus j\)（这里减法是每一位不退位减法），根据\(i\ominus j\)的状态来更新\(f_i\)。
换句话说，枚举\(j,k\)，然后用\(f_j\times b_{cnt1[k],cnt2[k]}\)来更新\(f_{j\oplus k}\)。
发现这就是一个高维循环卷积的形式，用FWT来做。
由于模数的特殊性质，可以发现模数一定和3互质，所以最后是可以得到答案的。
（不太懂网上那么多题解为什么一定要强行分析转移矩阵的性质啊qwq）
代码
我一开始写得太丑了卡不过去，学习了minamoto神仙的代码之后才过。
#include<bits/stdc++.h>
clock_t t=clock();
namespace my_std{
using namespace std;
#define pii pair<int,int>
#define fir first
#define sec second
#define MP make_pair
#define rep(i,x,y) for (int i=(x);i<=(y);i++)
#define drep(i,x,y) for (int i=(x);i>=(y);i--)
#define go(x) for (int i=head[x];i;i=edge[i].nxt)
#define templ template<typename T>
#define sz 550000
typedef long long ll;
typedef double db;
mt19937 rng(chrono::steady_clock::now().time_since_epoch().count());
templ inline T rnd(T l,T r) {return uniform_int_distribution<T>(l,r)(rng);}
templ inline bool chkmax(T &x,T y){return x<y?x=y,1:0;}
templ inline bool chkmin(T &x,T y){return x>y?x=y,1:0;}
templ inline void read(T& t)
{
t=0;char f=0,ch=getchar();double d=0.1;
while(ch>'9'||ch<'0') f|=(ch=='-'),ch=getchar();
while(ch<='9'&&ch>='0') t=t*10+ch-48,ch=getchar();
if(ch=='.'){ch=getchar();while(ch<='9'&&ch>='0') t+=d*(ch^48),d*=0.1,ch=getchar();}
t=(f?-t:t);
}
template<typename T,typename... Args>inline void read(T& t,Args&... args){read(t); read(args...);}
char __sr[1<<21],__z[20];int __C=-1,__zz=0;
inline void Ot(){fwrite(__sr,1,__C+1,stdout),__C=-1;}
inline void print(register int x)
{
if(__C>1<<20)Ot();if(x<0)__sr[++__C]='-',x=-x;
while(__z[++__zz]=x%10+48,x/=10);
while(__sr[++__C]=__z[__zz],--__zz);__sr[++__C]='\n';
}
void file()
{
#ifdef NTFOrz
freopen("a.in","r",stdin);
#endif
}
inline void chktime()
{
#ifndef ONLINE_JUDGE
cout<<(clock()-t)/1000.0<<'\n';
#endif
}
// inline ll mul(ll a,ll b){ll d=(ll)(a*(double)b/mod+0.5);ll ret=a*b-d*mod;if (ret<0) ret+=mod;return ret;}
}
using namespace my_std;
int n,m,T,mod;
ll b[15][15];
int add(int x,int y)
{
int cur=1,ret=0;
rep(k,1,m) ret+=(x+y)%3*cur,cur*=3,x/=3,y/=3;
return ret;
}
struct hh{ll x,y;hh(ll X=0,ll Y=0){x=X,y=Y;}}f[sz],w[sz]; // x + y \times \omega
inline ll M1(ll x){return x-((mod-x-1)>>31&mod);}
inline ll M2(ll x){return x+(x>>31&mod);}
hh operator * (hh a,hh b) { return hh(M2(a.x*b.x%mod-a.y*b.y%mod),M2(M1(a.x*b.y%mod+a.y*b.x%mod)-a.y*b.y%mod)); }
hh operator + (hh a,hh b) { return hh(M1(a.x+b.x),M1(a.y+b.y)); }
bool operator < (hh a,hh b){return a.x==b.x?a.y<b.y:a.x<b.x;}
ll ksm(ll x,int y){ll ret=1;for (;y;y>>=1,x=x*x%mod) if (y&1) ret=ret*x%mod;return ret;}
ll Ksm(ll x,int y){ll ret=1;for (;y;y>>=1,x=x*x) if (y&1) ret=ret*x;return ret;}
hh kksm(hh x,int y){hh ret=hh(1,0);for (;y;y>>=1,x=x*x) if (y&1) ret=ret*x;return ret;}
void exgcd(ll a,ll b,ll &x,ll &y)
{
if (a==1&&!b) return (void)(x=1,y=0);
exgcd(b,a%b,y,x);y=M2(y-a/b*x%mod);
}
ll inv(ll x){ll a,b;exgcd(x,mod,a,b);return (a+mod)%mod;}
hh calc1(hh a){return hh(M2(-a.y),M2(a.x-a.y));}
hh calc2(hh a){return hh(M2(a.y-a.x),M2(-a.x));}
void FWT(hh *a,int type)
{
for (int len=1,i=0;i<m;i++,len*=3) for (int j=0;j<n;j+=len*3) rep(k,0,len-1)
{
hh x=a[j+k],y=a[j+k+len],z=a[j+k+len*2],Y1=calc1(y),Y2=calc2(y),Z1=calc1(z),Z2=calc2(z);
if (type==-1) swap(Y1,Y2),swap(Z1,Z2);
a[j+k]=x+y+z;
a[j+k+len]=x+Y1+Z2;
a[j+k+len+len]=x+Y2+Z1;
}
}
int main()
{
file();
read(m,T,mod);n=Ksm(3,m);
if (mod==1) { rep(i,0,n-1) puts("0"); return 0; }
rep(i,0,n-1) read(f[i].x);
rep(i,1,m+1) rep(j,1,m+2-i) read(b[i-1][j-1]);
rep(i,0,n-1)
{
int x=i,y=i,c1=0,c2=0;
rep(k,1,m) c1+=(x%3==1),x/=3;
rep(k,1,m) c2+=(y%3==2),y/=3;
w[i].x=b[c1][c2];
}
FWT(w,1);FWT(f,1);
map<hh,hh>M;
rep(i,0,n-1) w[i]=(M.count(w[i])?M[w[i]]:M[w[i]]=kksm(w[i],T))*f[i];
FWT(w,-1);
ll I=inv(n);
rep(i,0,n-1) printf("%lld\n",w[i].x*I%mod);
return 0;
}原文地址：https://www.cnblogs.com/p-b-p-b/p/11404257.html
时间： 2024-10-25 13:53:44