UOJ #272. 【清华集训2016】石家庄的工人阶级队伍比较坚强
md为了做这题恶补了之前一知半解的很多多项式相关姿势，才总算是写掉了这题
首先考虑一个暴力的矩乘做法，我们考虑设一个矩阵\(f\)，其中\(f_{i,j}\)即为\(i\)轮后状态为\(i\)的人的分数，很显然我们可以找出一个矩阵\(B\)来进行转移，那么最后要求的就是\(B^t\times f_0\)
然后我们进行一个转化，令\(\oplus\)为三进制下不进位加法，\(\ominus\)为三进制下不退位减法，它们两个互为逆运算
然后我们根据剪刀石头布的定义可以知道当\(k\in [0,3^m)\)时，有\(B_{i,j}=B_{i\oplus k,j\oplus k}=B_{i\ominus k,j\ominus k}\)
因此我们有\(B_{i,j}=B_{i\ominus j,0}\)，然后它也可以用归纳法扩展到\(B^n\)
所以\(f_t=B^t\times f_0\)的第\(i\)项就是：
然后我们发现这个式子就是一个卷积的形式，而且只需要\(B\)的第一行即可
那么我们再化一下式子：
那么我们发现其实只要把\(B\)的第一行和自己做\(t\)次卷积即可，然后再卷上\(f_0\)就是答案
那么现在我们要找到一个三进制不进位加法的卷积变换，那么它其实是一个循环卷积，它的一般形式是：
对于循环卷积的话其实就是对于它的卷积的矩阵\(T\)的每一行\(x\)都有：
因此\(T\)的每一行都是方程组\(x_k\times x_t=x_{k\oplus t}\)的一组解
但是我们都知道如果卷积没有逆变换（就叫鸡卷好了）就GG了，因此我们的\(T\)还必须有逆矩阵
然后这也就说明\(T\)的行列式不为\(0\)，因此方程组\(x_k\times x_t=x_{k\oplus t}\)至少要有\(n\)组不同的解
所以这个矩阵\(T\)应该怎么构造呢？然后如果你认真学了FFT的姿势就会知道这个矩阵就是范德蒙德矩阵
因此矩阵长这样：
即\(T_{i,j}=\omega_n^{ij}\)
那么它的逆矩阵为：
然后因为这题是三进制不进位加法，因此我们选三次单位根\(\omega\)，然后我们发现所有的复数都可以被表示成\(a+b\omega\)的形式
又由单位根的一个性质\(\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\not | \ p\)
假设\(3|p\)，那么我们令\(k=\frac{p}{3}\)，那么有：
与假设矛盾！因此得证\(3\not |\ 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\n",f[i].x);
return 0;
}
