md为了做这题恶补了之前一知半解的很多多项式相关姿势，才总算是写掉了这题
首先考虑一个暴力的矩乘做法，我们考虑设一个矩阵(f)，其中(f_{i,j})即为(i)轮后状态为(i)的人的分数，很显然我们可以找出一个矩阵(B)来进行转移，那么最后要求的就是(B^t imes f_0)
然后我们进行一个转化，令(oplus)为三进制下不进位加法，(ominus)为三进制下不退位减法，它们两个互为逆运算
然后我们根据剪刀石头布的定义可以知道当(kin [0,3^m))时，有(B_{i,j}=B_{ioplus k,joplus k}=B_{iominus k,jominus k})
因此我们有(B_{i,j}=B_{iominus j,0})，然后它也可以用归纳法扩展到(B^n)
所以(f_t=B^t imes f_0)的第(i)项就是：
然后我们发现这个式子就是一个卷积的形式，而且只需要(B)的第一行即可
那么我们再化一下式子：
那么我们发现其实只要把(B)的第一行和自己做(t)次卷积即可，然后再卷上(f_0)就是答案
那么现在我们要找到一个三进制不进位加法的卷积变换，那么它其实是一个循环卷积，它的一般形式是：
对于循环卷积的话其实就是对于它的卷积的矩阵(T)的每一行(x)都有：
因此(T)的每一行都是方程组(x_k imes x_t=x_{koplus t})的一组解
但是我们都知道如果卷积没有逆变换（就叫鸡卷好了）就GG了，因此我们的(T)还必须有逆矩阵
然后这也就说明(T)的行列式不为(0)，因此方程组(x_k imes x_t=x_{koplus t})至少要有(n)组不同的解
所以这个矩阵(T)应该怎么构造呢？然后如果你认真学了FFT的姿势就会知道这个矩阵就是范德蒙德矩阵
因此矩阵长这样：
即(T_{i,j}=omega_n^{ij})
那么它的逆矩阵为：
然后因为这题是三进制不进位加法，因此我们选三次单位根(omega)，然后我们发现所有的复数都可以被表示成(a+bomega)的形式
又由单位根的一个性质(omega^2+omega+1=0)，因此有(omega^2=-omega-1)，因此我们可以将乘出来的结果降幂，因此我们发现现在复数的乘法就是：
然后我们由上面的可知：
同时我们知道：
其中(I)是单位矩阵，因此我们就把后面那个当作逆矩阵来算就好了
然后由于DFT和IDFT的本质都是分治做向量对矩阵的乘法，它每次IDFT的时候都多乘了一个(3)，而一共有(log_3 n)层，因此总共多乘了(3^{log_3 n}=n)次，那么最后只要把所有元素都除以一个(n)就好了！
PS：是不是明白了一般的FFT的转移系数是哪里来的，以及做一般卷积的时候要把(n)扩大一倍（这样做循环卷积不会被膜到前面去），以及IDFT之后为什么要除以(n)这些问题
然后提到除以(n)我们就会想到是不是有逆元的问题，因为我们要求(n=3^m)有逆元，那么也就是(3 ot | p)
假设(3|p)，那么我们令(k=frac{p}{3})，那么有：
与假设矛盾！因此得证(3 ot | p)
所以这题终于是做完了……
#include<cstdio>
#define RI register int
#define CI const int&
using namespace std;
const int N=531441,M=20;
int m,t,mod,n,x,b[M+5][M+5];
inline int sum(CI x,CI y)
{
int t=x+y; return t>=mod?t-mod:t;
}
inline int sub(CI x,CI y)
{
int t=x-y; return t<0?t+mod:t;
}
inline void exgcd(int& x,int& y,CI m,CI n)
{
if (!n) return (void)(x=1,y=0);
exgcd(y,x,n,m%n); y-=m/n*x;
}
inline int inv(CI x)
{
int m,n; exgcd(m,n,x,mod);
return (m%mod+mod)%mod;
}
struct Complex
{
int x,y; //x+y*w
inline Complex(CI X=0,CI Y=0)
{
x=X; y=Y;
}
friend inline Complex operator + (const Complex& A,const Complex& B)
{
return Complex(sum(A.x,B.x),sum(A.y,B.y));
}
friend inline Complex operator - (const Complex& A,const Complex& B)
{
return Complex(sub(A.x,B.x),sub(A.y,B.y));
}
friend inline Complex operator * (const Complex& A,const Complex& B)
{
return Complex(sub(1LL*A.x*B.x%mod,1LL*A.y*B.y%mod),
sub(sum(1LL*A.x*B.y%mod,1LL*A.y*B.x%mod),1LL*A.y*B.y%mod));
}
}f[N],g[N];
inline Complex quick_pow(Complex x,int p,Complex mul=Complex(1,0))
{
for (;p;p>>=1,x=x*x) if (p&1) mul=mul*x;
return mul;
}
inline Complex calc1(const Complex& A) //calc A*w
{
return Complex(mod-A.y,sub(A.x,A.y));
}
inline Complex calc2(const Complex& A) //calc A*w^2
{
return Complex(sub(A.y,A.x),mod-A.x);
}
inline void DFT(Complex *f)
{
RI mid,j,k; for (mid=1;mid<n;mid*=3)
for (j=0;j<n;j+=mid*3) for (k=0;k<mid;++k)
{
Complex x=f[j+k],y=f[j+k+mid],z=f[j+k+(mid<<1)];
f[j+k]=x+y+z; f[j+k+mid]=x+calc1(y)+calc2(z);
f[j+k+(mid<<1)]=x+calc2(y)+calc1(z);
}
}
inline void IDFT(Complex *f)
{
RI mid,j,k; for (mid=1;mid<n;mid*=3)
for (j=0;j<n;j+=mid*3) for (k=0;k<mid;++k)
{
Complex x=f[j+k],y=f[j+k+mid],z=f[j+k+(mid<<1)];
f[j+k]=x+y+z; f[j+k+mid]=x+calc2(y)+calc1(z);
f[j+k+(mid<<1)]=x+calc1(y)+calc2(z);
}
int invn=inv(n); for (RI i=0;i<n;++i) f[i].x=1LL*f[i].x*invn%mod;
}
int main()
{
RI i,j; for (scanf("%d%d%d",&m,&t,&mod),n=i=1;i<=m;++i) n*=3;
for (i=0;i<n;++i) scanf("%d",&x),f[i].x=x;
for (i=0;i<=m;++i) for (j=0;j<=m-i;++j) scanf("%d",&b[i][j]);
for (i=0;i<n;++i)
{
int x=i,cw=0,cl=0; while (x)
{
int nw=x%3; if (nw==1) ++cw;
else if (nw==2) ++cl; x/=3;
}
g[i].x=b[cw][cl];
}
for (DFT(f),DFT(g),i=0;i<n;++i) f[i]=f[i]*quick_pow(g[i],t);
for (IDFT(f),i=0;i<n;++i) printf("%d
",f[i].x);
return 0;
}