白兔之舞 题解

题意：白兔初始在 \((0,x)\) 的位置，每次跳跃要求第一维单增，第二维从 \(u\) 到 \(v\) 共 \(w_{u,v}\) 条边。第一维不能跳超过 \(L\)，第二维取值在 \(1\) 到 \(n\) 之间。
设步数为 \(m\)，询问 \(m\bmod k=t\) 时的跳跃方案对 \(p\) 取模的结果。

题解：
单位根反演的式子为：\(\frac{1}{n}\sum\limits_{i=0}^{n-1}\omega_{n}^{ik}=[n|k]\)
当 \(n=1\) 时，设 \(W\) 为 \(w_{1,1}\) 的结果，则答案为：

\[\begin{aligned} \sum\limits_{i=0}^L [k|(i-t)]W^i\binom{L}{i}&=\sum\limits_{i=0}^L\frac{1}{k}\sum\limits_{j=0}^{k-1}\omega_k^{j(i-t)}W^i\binom{L}{i}\\ &=\frac{1}{k}\sum\limits_{i=0}^LW^i\binom{L}{i}\sum\limits_{j=0}^{k-1}\omega_k^{-jt+ij}\\ &=\frac{1}{k}\sum\limits_{i=0}^{k-1}\omega_k^{-it}\sum\limits_{j=0}^LW^j\binom{L}{j}\omega_k^{ij}\\ &=\frac{1}{k}\sum\limits_{i=0}^{k-1}\omega^{-it}_k(W\omega _k^i+1)^L \end{aligned} \]

设 \(C_i=(W\omega_k^i+1)^L\)
则原式 \(=\frac{1}{k}\sum\limits_{i=0}^{k-1}\omega_k^{\binom{i}{2}+\binom{t}{2}-\binom{i+t}{2}}C_i=\frac{1}{k}\omega_k^{\binom{t}{2}}\sum\limits_{i=0}^{k-1}\omega_k^{\binom{i}{2}}C_i\omega_k^{-\binom{i+t}{2}}\)
发现为卷积的形式。可以用任意模数 NTT 解决。
当 \(n\gt 1\) 时，设 \(Begin\) 为初始矩阵，（只有 \(\left(1,x\right)\) 位置为 \(1\)），\(W\)为转移矩阵，则 \(C_i\) 为 \(Begin\times\left(W\omega_k^i+1\right)^L\) 这个矩阵的 \(\left(1,y\right)\) 位置的值。

代码：

#include <bits/stdc++.h>
const int mod1 = 998244353;
const int mod2 = 1004535809;
const int mod3 = 469762049;
const int M = 1000005;
const int g = 3;
using namespace std;
int max_L,max_len,rev[M],C[M],W_k[M],n,k,L,x,y,p;
int read(){
	char c=getchar();int ans=0;bool flag=1;
	while (c<'0'||c>'9') flag&=(c!='-'),c=getchar();
	while (c>='0'&&c<='9') ans=ans*10+c-'0',c=getchar();
	return flag?ans:-ans;
}
void Write(long long x){
	if (x<10) putchar(x^48);
	else Write(x/10),putchar((x%10)^48);
	return;
}
int power(int x,int y,int mod){
	long long ans=1,now=x;
	for (int i=y;i;i>>=1,now=now*now%mod)
		if (i&1) ans=ans*now%mod;
	return ans;
}
int add(int u,int v,int mod){u+=v-mod;return u+((u>>31)&mod);}
int sub(int u,int v,int mod){u-=v;return u+((u>>31)&mod);}
int invs(int x,int mod){return x==1?1:(mod-mod/x+0ll)*invs(mod%x,mod)%mod;}
const long long inv1=invs(mod1%mod2,mod2);
const long long mod1_2=mod1*(mod2+0ll);
const long long inv2=invs(mod1_2%mod3,mod3);
int get_root(int p){
	static int S[105];int tmp=p-1;S[0]=0;
	for (register int i=2;i*i<=tmp;i++){
		if (tmp%i) continue;S[++S[0]]=i;
		while (tmp%i==0) tmp/=i;
	}
	if (tmp!=1) S[++S[0]]=tmp;
	for (register int i=2;;i++){
		bool flag=1;
		for (register int j=1;j<=S[0];j++)
			if (power(i,(p-1)/S[j],p)==1){flag=0;break;}
		if (flag) return i;
	}
	return 0;
}
struct number{
	int num1,num2,num3;
	number operator+ (const number b) const{
		return (number){add(num1,b.num1,mod1),add(num2,b.num2,mod2),add(num3,b.num3,mod3)};
	}
	number operator+= (const number &b){return *this=*this+b;}
	number operator* (const number b) const{
		return (number){(int)((num1+0ll)*b.num1%mod1),(int)((num2+0ll)*b.num2%mod2),(int)((num3+0ll)*b.num3%mod3)};
	}
	number operator*= (const number &b){return *this=*this*b;}
	number operator- (const number b) const{
		return (number){sub(num1,b.num1,mod1),sub(num2,b.num2,mod2),sub(num3,b.num3,mod3)};
	}
	number operator-= (const number &b){return *this=*this-b;}
	number operator^ (const int b) const{
		number ans=(number){1,1,1},now=*this;
		for (register int i=b;i;i>>=1,now*=now)
			if (i&1) ans*=now;
		return ans;
	}
}W[M],inv[M],A[M],B[M];
number make_number (int x){return (number){x%mod1,x%mod2,x%mod3};}
long long get_ans(number now){
	long long x4=(now.num2-now.num1+mod2)*inv1%mod2*mod1+now.num1;
	return (x4+mod1_2%p*((now.num3+mod3-x4%mod3)*inv2%mod3))%p;
}
void prepare(){
	max_len=1<<19,max_L=19;
	for (register int i=0;i<max_len;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<max_L-1);
	number wn=(number){power(g,(mod1-1)/max_len,mod1),power(g,(mod2-1)/max_len,mod2),power(g,(mod3-1)/max_len,mod3)};
	W[max_len>>1]=(number){1,1,1};inv[1]=(number){1,1,1};
	for (register int i=2;i<=max_len;i++){
		int X=(mod1-mod1/i+0ll)*inv[mod1%i].num1%mod1;
		int Y=(mod2-mod2/i+0ll)*inv[mod2%i].num2%mod2;
		int Z=(mod3-mod3/i+0ll)*inv[mod3%i].num3%mod3;
		inv[i]=(number){X,Y,Z};
	}
	for (register int i=(max_len>>1)+1;i<max_len;i++) W[i]=W[i-1]*wn;
	for (register int i=(max_len>>1)-1;i;i--) W[i]=W[i<<1];
	return;
}
void NTT(number *A,int len,int L,int flag){
	static number F[M];
	for (register int i=0,s=max_L-L;i<len;i++) F[rev[i]>>s]=A[i];
	for (register int l=1;l<len;l<<=1)
		for (register int i=0;i<len;i+=(l<<1))
			for (register int j=0;j<l;j++){
				number u=F[i|j],v=F[i|j|l]*W[j|l];
				F[i|j]=u+v,F[i|j|l]=u-v;
			}
	if (flag){
		A[0]=F[0]*inv[len];
		for (register int i=1;i<len;i++) A[i]=F[len-i]*inv[len];
	}
	else for (register int i=0;i<len;i++) A[i]=F[i];
	return;
}
struct Matrix{
	int a[3][3];
	Matrix operator* (const Matrix b) const{
		Matrix c;
		for (register int i=0;i<3;i++)
			for (register int j=0;j<3;j++)
				c.a[i][j]=(a[i][0]*(b.a[0][j]+0ll)+a[i][1]*(b.a[1][j]+0ll)+a[i][2]*(b.a[2][j]+0ll))%p;
		return c;
	}
	Matrix operator*= (const Matrix &b){return *this=*this*b;}
	Matrix operator* (const int b) const{
		Matrix c=*this;
		for (register int i=0;i<3;i++)
			for (register int j=0;j<3;j++)
				c.a[i][j]=c.a[i][j]*(b+0ll)%p;
		return c;
	}
	Matrix operator*= (const int &b){return *this=*this*b;}
	Matrix operator+ (const Matrix b) const{
		Matrix c;
		for (register int i=0;i<3;i++)
			for (register int j=0;j<3;j++) c.a[i][j]=add(a[i][j],b.a[i][j],p);
		return c;
	}
	Matrix operator+= (const Matrix &b){return *this=*this+b;}
	Matrix operator^ (const int b) const{
		Matrix ans=*this,now=*this;
		if (!b){
			ans.a[0][0]=ans.a[1][1]=ans.a[2][2]=1;
			ans.a[0][1]=ans.a[1][2]=ans.a[2][0]=0;
			ans.a[0][2]=ans.a[1][0]=ans.a[2][1]=0;
			return ans;
		}
		for (register int i=b-1;i;i>>=1,now*=now)
			if (i&1) ans*=now;
		return ans;
	}
}G,Begin,I;
void mul(number *A,number *B,int len){
	int rlen=1,L=0;while (rlen<=(len<<1)) rlen<<=1,++L;
	NTT(A,rlen,L,0),NTT(B,rlen,L,0);
	for (register int i=0;i<rlen;i++) A[i]*=B[i];
	NTT(A,rlen,L,1);
	for (register int i=len;i<rlen;i++) A[i]=make_number(0);
	return;
}
int main(){
	n=read(),k=read(),L=read(),x=read()-1,y=read()-1,p=read();
	for (register int i=0;i<n;i++)
		for (register int j=0;j<n;j++) G.a[i][j]=read();
	I.a[0][0]=I.a[1][1]=I.a[2][2]=1;
	I.a[0][1]=I.a[1][2]=I.a[2][0]=0;
	I.a[0][2]=I.a[1][0]=I.a[2][1]=0;
	Begin.a[0][x]=1;W_k[0]=1;W_k[1]=power(get_root(p),(p-1)/k,p);prepare();
	for (register int i=2;i<=k;i++) W_k[i]=W_k[i-1]*(W_k[1]+0ll)%p;
	for (register int i=0;i<k;i++) C[i]=(Begin*((G*W_k[i]+I)^L)).a[0][y];
	for (register int i=0;i<k;i++) A[i]=make_number(W_k[i*(i-1ll)/2%k]*(C[i]+0ll)%p);
	for (register int i=0;i<k+k;i++) B[k+k-1-i]=make_number(W_k[k-i*(i-1ll)/2%k]%p);
	mul(A,B,k<<1);int invk=invs(k,p);
	for (register int i=0;i<k;i++)
		Write(get_ans(A[k+k-i-1])*invk%p*W_k[i*(i-1ll)/2%k]%p),putchar('\n');
	return 0;
}