something wrong with my new blog!
I can't type matrixs
so I come back. qwq
题目:https://www.codechef.com/problems/BIKE
题解
是我naive了，二维和一维其实差不多
首先，n很小，t很大，什么算法？矩阵乘法！没跑了
然后矩阵里填什么？一条边是两个值啊，还要一个%n一个%(n - 1)，怎么搞
我们设计一个多项式(x^{a}y^{b})，x指数（也就是a）代表前轮加上一条边的值后取模的值，y指数（也就是b）代表后轮加上一条边的值后取模的值
最后每一项前面的系数就是答案
这样我们矩阵里填的就是多项式了，例如一条边的值是(f = 1,r = 2)那么这个邻接矩阵里填的就是(xy^{2})
好了我们做矩阵乘法并且指数分别取模n和(n - 1)就做完啦！（大雾
啥啥啥，为啥这样就做完了啊，我要在矩阵里填多项式，还要资瓷多项式的加和乘，虽然可行？但是复杂度听起来不可靠（项有(n^2)级别啊，我多项式乘暴力搞要(n^4)我FFT还要(n^2 log n)，而我矩阵乘法要做(n^3)次，excuse me？)而且很玄学（我不会写），有没有什么简单优美的办法啊。
考虑插值，我们熟悉的插值都有什么呢，例如FFT，例如拉格朗日插值，都是找一些点值代进去然后反解出系数
考虑继续用FFT的思路插值反解系数，因为单位根有一个很棒的性质就是(w^{n} = 1)（先记住这一点，不过可能大家都知道）
啥，你说这模数不是NTT？不是哇你不是N只有22么，你会发现(P - 1)是2..22所有的数的倍数，所以求出原根来就可以求单位根啦，不是(2^k)而是正好(N)的(所以我们不能FFT了而要暴力DFT，也就是把每个单位根的值代进去求一遍值，显然你(n^2)和大常数(n log n)在22没有什么太大的区别，况且你也不能FFT）
那么我们就可以尝试把(x)的(N)种取值，(y)的(N - 1)种取值结合起来
例如(xy^{2})就是(w_{N}w_{N - 1}^{2})
好惹现在窝们（奇怪的口音？）有了(N(N - 1))个点值了……我们还要画在三维坐标系里不成？
回忆一下学习FFT时候的矩阵
大概是这个样子的
(egin{bmatrix} (omega_{n}^{0})^{0} & (omega_{n}^{0})^{1} & cdots & (omega_{n}^{0})^{n - 1} \ (omega_{n}^{1})^{0} & (omega_{n}^{1})^{1} & cdots & (omega_{n}^{1})^{n - 1} \ vdots &vdots & ddots & vdots \ (omega_{n}^{n - 1})^{0} & (omega_{n}^{n - 1})^{1} & cdots & (omega_{n}^{n - 1})^{n - 1} end{bmatrix} egin{bmatrix} a_{0}\ a_{1}\ vdots \ a_{n - 1} end{bmatrix} = egin{bmatrix} A(omega_{0}) \ A(omega_{1})\ vdots \ A(omega_{n - 1}) end{bmatrix})
感觉没有问题，那么回忆一下是怎么反解出来的
假如有这么一个矩阵
(egin{bmatrix}
(omega_{n}^{-0})^{0} &(omega_{n}^{-0})^{1} &cdots &(omega_{n}^{-0})^{n - 1} \
(omega_{n}^{-1})^{0} &(omega_{n}^{-1})^{1} &cdots &(omega_{n}^{-1})^{n - 1} \
vdots &vdots &ddots &vdots \
(omega_{n}^{-(n - 1)})^{0} &(omega_{n}^{-(n - 1)})^{1} &cdots &(omega_{n}^{-(n - 1)})^{n - 1}
end{bmatrix})
我们让两个矩阵相乘
(egin{bmatrix} (omega_{n}^{0})^{0} & (omega_{n}^{0})^{1} & cdots & (omega_{n}^{0})^{n - 1} \ (omega_{n}^{1})^{0} & (omega_{n}^{1})^{1} & cdots & (omega_{n}^{1})^{n - 1} \ vdots &vdots & ddots & vdots \ (omega_{n}^{n - 1})^{0} & (omega_{n}^{n - 1})^{1} & cdots & (omega_{n}^{n - 1})^{n - 1} end{bmatrix}egin{bmatrix}
(omega_{n}^{-0})^{0} &(omega_{n}^{-0})^{1} &cdots &(omega_{n}^{-0})^{n - 1} \
(omega_{n}^{-1})^{0} &(omega_{n}^{-1})^{1} &cdots &(omega_{n}^{-1})^{n - 1} \
vdots &vdots &ddots &vdots \
(omega_{n}^{-(n - 1)})^{0} &(omega_{n}^{-(n - 1)})^{1} &cdots &(omega_{n}^{-(n - 1)})^{n - 1}
end{bmatrix} = E)
设
(D = egin{bmatrix} (omega_{n}^{0})^{0} & (omega_{n}^{0})^{1} & cdots & (omega_{n}^{0})^{n - 1} \ (omega_{n}^{1})^{0} & (omega_{n}^{1})^{1} & cdots & (omega_{n}^{1})^{n - 1} \ vdots &vdots & ddots & vdots \ (omega_{n}^{n - 1})^{0} & (omega_{n}^{n - 1})^{1} & cdots & (omega_{n}^{n - 1})^{n - 1} end{bmatrix})
(V = egin{bmatrix}
(omega_{n}^{-0})^{0} &(omega_{n}^{-0})^{1} &cdots &(omega_{n}^{-0})^{n - 1} \
(omega_{n}^{-1})^{0} &(omega_{n}^{-1})^{1} &cdots &(omega_{n}^{-1})^{n - 1} \
vdots &vdots &ddots &vdots \
(omega_{n}^{-(n - 1)})^{0} &(omega_{n}^{-(n - 1)})^{1} &cdots &(omega_{n}^{-(n - 1)})^{n - 1}
end{bmatrix})
(sum_{i = 0}^{n - 1}sum_{j = 0}^{n - 1} E_{ij} = D_{ik}V_{kj})
当(i == j)
(E_{ij} = n)
当(i != j)
(E_{ij} = sum_{k = 0}^{n - 1} D_{ik}V_{kj})
(E_{ij} = sum_{k = 0}^{n - 1} omega_{n}^{ik - kj})
(E_{ij} = sum_{k = 0}^{n - 1} omega_{n}^{k(i - j)})
(E_{ij} = frac{1 - (omega_{n}^{i - j})^n}{1 - omega_{n}^{i - j}})
根据单位根神奇的性质，我们有
((omega_{n}^{i - j})^n = 1)！
这样的话
(E = nI)
(I)是单位矩阵
所以就有了
(frac{1}{n} egin{bmatrix} (omega_{n}^{0})^{0} & (omega_{n}^{0})^{1} & cdots & (omega_{n}^{0})^{n - 1} \ (omega_{n}^{1})^{0} & (omega_{n}^{1})^{1} & cdots & (omega_{n}^{1})^{n - 1} \ vdots &vdots & ddots & vdots \ (omega_{n}^{n - 1})^{0} & (omega_{n}^{n - 1})^{1} & cdots & (omega_{n}^{n - 1})^{n - 1} end{bmatrix} egin{bmatrix} A(omega_{0}) \ A(omega_{1})\ vdots \ A(omega_{n - 1}) end{bmatrix} = egin{bmatrix} a_{0}\ a_{1}\ vdots \ a_{n - 1} end{bmatrix})
妙趣横生，我们只需要把原来的FFT的单位根改成反着旋转，最后除一下数组长度即可
然而？？？你说你会，因为wys在WC2018讲过了。。。（大雾
其实二维的毕克在WC2017也讲了，显然这道题是BI KE出的。。。= =
我们考虑一个(n(n - 1))的矩阵！
(egin{bmatrix} (omega_{n}^{0})^{0}(omega_{n - 1}^{0})^{0} &(omega_{n}^{0})^{0}(omega_{n - 1}^{0})^{1} &cdots &(omega_{n}^{0})^{n - 1}(omega_{n - 1}^{0})^{n - 3} &(omega_{n}^{0})^{n - 1}(omega_{n - 1}^{0})^{n - 2} \ (omega_{n}^{0})^{0}(omega_{n - 1}^{1})^{0}&(omega_{n}^{0})^{0}(omega_{n - 1}^{1})^{1} & cdots &(omega_{n}^{0})^{n - 1}(omega_{n - 1}^{1})^{n - 3} & (omega_{n}^{0})^{n - 1}(omega_{n - 1}^{1})^{n - 2}\ vdots &vdots & ddots &vdots & vdots \ (omega_{n}^{n - 1})^{0}(omega_{n - 1}^{n - 3})^{0} & (omega_{n}^{n - 1})^{0}(omega_{n - 1}^{n - 3})^{1} & cdots & vdots & vdots \ (omega_{n}^{n - 1})^{0}(omega_{n - 1}^{n - 2})^{0} & (omega_{n}^{n - 1})^{0}(omega_{n - 1}^{n - 2})^{1} & cdots & (omega_{n}^{n - 1})^{n - 1}(omega_{n - 1}^{n - 2})^{n - 3} & (omega_{n}^{n - 1})^{n - 1}(omega_{n - 1}^{n - 2})^{n - 2} end{bmatrix}egin{bmatrix} a_{0,0}\ a_{0,1}\ vdots \ a_{n-1,n-3}\ a_{n-1,n-2} end{bmatrix} = egin{bmatrix} A(omega_{n}^{0}omega_{n - 1}^{0})\ A(omega_{n}^{0}omega_{n - 1}^{1})\ vdots \ A(omega_{n}^{n - 1}omega_{n - 1}^{n - 3})\ A(omega_{n}^{n - 1}omega_{n - 1}^{n - 2}) end{bmatrix})
按照一样的方法把号都取反
(D = egin{bmatrix}
(omega_{n}^{0})^{0}(omega_{n - 1}^{0})^{0} &(omega_{n}^{0})^{0}(omega_{n - 1}^{0})^{1} &cdots &(omega_{n}^{0})^{n - 1}(omega_{n - 1}^{0})^{n - 3} &(omega_{n}^{0})^{n - 1}(omega_{n - 1}^{0})^{n - 2} \
(omega_{n}^{0})^{0}(omega_{n - 1}^{1})^{0}&(omega_{n}^{0})^{0}(omega_{n - 1}^{1})^{1} & cdots &(omega_{n}^{0})^{n - 1}(omega_{n - 1}^{1})^{n - 3} & (omega_{n}^{0})^{n - 1}(omega_{n - 1}^{1})^{n - 2}\
vdots &vdots & ddots &vdots & vdots \
(omega_{n}^{n - 1})^{0}(omega_{n - 1}^{n - 3})^{0} & (omega_{n}^{n - 1})^{0}(omega_{n - 1}^{n - 3})^{1} & cdots & vdots & vdots \
(omega_{n}^{n - 1})^{0}(omega_{n - 1}^{n - 2})^{0} & (omega_{n}^{n - 1})^{0}(omega_{n - 1}^{n - 2})^{1} & cdots & (omega_{n}^{n - 1})^{n - 1}(omega_{n - 1}^{n - 2})^{n - 3} & (omega_{n}^{n - 1})^{n - 1}(omega_{n - 1}^{n - 2})^{n - 2}
end{bmatrix})
(V = egin{bmatrix}
(omega_{n}^{-0})^{0}(omega_{n - 1}^{-0})^{0} &(omega_{n}^{-0})^{0}(omega_{n - 1}^{-0})^{1} &cdots &(omega_{n}^{-0})^{n - 1}(omega_{n - 1}^{-0})^{n - 3} &(omega_{n}^{-0})^{n - 1}(omega_{n - 1}^{-0})^{n - 2} \
(omega_{n}^{-0})^{0}(omega_{n - 1}^{-1})^{0}&(omega_{n}^{-0})^{0}(omega_{n - 1}^{-1})^{1} & cdots &(omega_{n}^{-0})^{n - 1}(omega_{n - 1}^{-1})^{n - 3} & (omega_{n}^{-0})^{n - 1}(omega_{n - 1}^{-1})^{n - 2}\
vdots &vdots & ddots &vdots & vdots \
(omega_{n}^{-(n - 1)})^{0}(omega_{n - 1}^{-(n - 3)})^{0} & (omega_{n}^{-(n - 1)})^{0}(omega_{n - 1}^{-(n - 3)})^{1} & cdots & vdots & vdots \
(omega_{n}^{-(n - 1)})^{0}(omega_{n - 1}^{-(n - 2)})^{0} & (omega_{n}^{-(n - 1)})^{0}(omega_{n - 1}^{-(n - 2)})^{1} & cdots & (omega_{n}^{-(n - 1)})^{n - 1}(omega_{n - 1}^{-(n - 2)})^{n - 3} & (omega_{n}^{-(n - 1)})^{n - 1}(omega_{n - 1}^{-(n - 2)})^{n - 2}
end{bmatrix})
(E = DV)
(sum_{i = 0}^{n(n - 1) - 1}sum_{j = 0}^{n(n - 1) - 1} E_{i,j} = D_{i,k}V_{k,j})
考虑把标号拆成一个点，类似一些的矩阵标号方法？
设数值为i拆成的点是(i_{x},i_{y})
当(i == j)时
(E_{i,j} = (omega_{n}^{i_{x}})^{k_x}(omega_{n - 1}^{i_y})^{k_y} (omega_{n}^{-k_{x}})^{j_x}(omega_{n - 1}^{-k_y})^{j_y})
(E_{i,j} = n * (n - 1))
而当(i != j)时
(E_{i,j} = (omega_{n}^{i_{x}})^{k_x}(omega_{n - 1}^{i_y})^{k_y} (omega_{n}^{-k_{x}})^{j_x}(omega_{n - 1}^{-k_y})^{j_y})
(E_{i,j} = omega_{n}^{k_x(i_{x} - j_{x})}omega_{n - 1}^{k_y(i_{y} - j_{y})})
嗯？迷一般的。。。眼熟？
可是这样就不是等比数列了？…………等等！
(E_{i,j} = sum_{k_{x} = 0}^{n - 1}omega_{n}^{k_x(i_{x} - j_{x})}sum_{k_{y} = 0}^{n - 2}omega_{n - 1}^{k_y(i_{y} - j_{y})})
(E_{i,j} = frac{1 - omega_{n}^{n(i_{x} - j_{x})}}{1 - omega_{n}^{i_x - j_x}}frac{1 - omega_{n - 1}^{(n - 1)(i_{y} - j_{y})}}{1 - omega_{n - 1}^{i_{y} - j_{y}}})
(E_{i,j} = 0)
也就是！我们把等比数列的公比q当成了(omega_{n}^{i_x - j_x})
后面就不写了，大家都知道了，最后模长除个(n(n - 1))就可以了！
说到这，我想我已经想到了一个绝妙的(n^2 log n)的二维FFT（在界是2^k的情况下）
我猜你们也想到了，大概二进制平摊反转这一个小技巧就不能用了，可以递归
代码
#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#define MAXM 100005
//#define ivorysi
using namespace std;
typedef long long int64;
const int MOD = 1163962801,G = 46;
int N,M,T;
int px[105],py[105];
struct node {int s,t,x,y;} E[MAXM];
int v[25][25][25][25],w[25][25];
struct Matrix {
int f[24][24];
void clear() {
memset(f,0,sizeof(f));
}
friend Matrix operator * (const Matrix &a,const Matrix &b) {
Matrix c;c.clear();
for(int i = 1 ; i <= N ; ++i) {
for(int j = 1 ; j <= N ; ++j) {
for(int k = 1 ; k <= N ; ++k) {
c.f[i][j] = (c.f[i][j] + 1LL * a.f[i][k] * b.f[k][j] % MOD) % MOD;
}
}
}
return c;
}
void set_Matrix(int x,int y) {
clear();
for(int i = 1 ; i <= M ; ++i) {
f[E[i].s][E[i].t] = (f[E[i].s][E[i].t] + 1LL * px[E[i].x * x % N] * py[E[i].y * y % (N - 1)] % MOD) % MOD;
}
}
}A,Ans,I;
int fpow(int x,int c) {
int res = 1,t = x;
while(c) {
if(c & 1) res = 1LL * res * t % MOD;
t = 1LL * t * t % MOD;
c >>= 1;
}
return res;
}
void mul(int c) {
Ans = I;Matrix t = A;
while(c) {
if(c & 1) Ans = Ans * t;
t = t * t;
c >>= 1;
}
}
void Process(int x,int y) {
for(int i = 0 ; i < N ; ++i) {
for(int j = 0 ; j < N - 1 ; ++j) {
w[i][j] = 0;
for(int k = 0 ; k < N ; ++k) {
for(int l = 0 ; l < N - 1 ; ++l) {
w[i][j] = (w[i][j] + 1LL * v[x][y][k][l] * px[N - k * i % N] % MOD
* py[N - 1 - l * j % (N - 1)] % MOD) % MOD;
}
}
w[i][j] = 1LL * w[i][j] * fpow(N * (N - 1),MOD - 2) % MOD;
}
}
}
void Init() {
scanf("%d%d%d",&N,&M,&T);
for(int i = 1 ; i <= M ; ++i) {
scanf("%d%d%d%d",&E[i].s,&E[i].t,&E[i].x,&E[i].y);
E[i].x %= N;
E[i].y %= N - 1;
}
px[0] = 1;py[0] = 1;
px[1] = fpow(G,(MOD - 1) / N);
py[1] = fpow(G,(MOD - 1) / (N - 1));
for(int i = 2 ;i <= N ; ++i) {
px[i] = 1LL * px[i - 1] * px[1] % MOD;
py[i] = 1LL * py[i - 1] * py[1] % MOD;
}
}
void Solve() {
Init();
I.clear();
for(int i = 1 ; i <= N ; ++i) I.f[i][i] = 1;
for(int i = 0 ; i < N ; ++i) {
for(int j = 0 ; j < N - 1 ; ++j) {
A.set_Matrix(i,j);
mul(T);
for(int k = 1 ; k <= N ; ++k) {
for(int l = 1 ; l <= N ; ++l) {
v[k][l][i][j] = Ans.f[k][l];
}
}
}
}
for(int i = 1 ; i <= N ; ++i) {
Process(i,i);
for(int k = 0 ; k < N ; ++k) {
for(int l = 0 ; l < N - 1 ; ++l) {
printf("%d%c",w[k][l],"
"[l == N - 2]);
}
}
}
}
int main() {
#ifdef ivorysi
freopen("f1.in","r",stdin);
#endif
Solve();
}