Codechef Oct chanllenge Queries on Matrix-JIIT
Codechef Oct chanllenge Queries on Matrix-JIIT
首先发现矩阵的两个维度显然是互不相干的，假设最后操作后有\(x\)列被操作奇数次，\(y\)行操作奇数次
那么最后为奇数的格子个数就是\(x(m-y)+(n-x)y\)
考虑求出\(q\)操作后有\(x\)个位置被操作奇数次的方案数
考虑一个Naive的\(dp\)，令\(dp_{i,j}\)为操作\(i\)次后有\(j\)位置操作奇数的方案数，显然得到转移为
\(dp_{i,j}\cdot j\rightarrow dp_{i+1,j-1}\)
\(dp_{i,j}\cdot (n-j)\rightarrow dp_{i+1,j+1}\)
直接\(dp\)复杂度为\(O(nq)\)，用矩阵优化复杂度为\(O(n^3\log q)\)
没有考虑过求向量列的现行递推式？说不定暴力求递推然后。。。
\(O(n^3)\)
由于笔者不会数学，所以考虑一个非常暴力非常直观的理解，可以完全抛开组合意义
每次是挑选一个位置异或上1，用形式幂级数可能比较蛋疼，不如直接搞成集合幂级数
即令集合\(S\)包含所有被操作奇数次的位置，用一个多项式\(F=\sum_S a_S x^S\)表示答案
那么转移多项式即为\(F=\sum_{i=1}^{n} x^{\lbrace i \rbrace}\)，转移运算为集合对称差运算(哎就是异或)
那么实际就是要求出\(F^q\)，可以直接用\(\text{FWT}\)优化，先\(\text{FWT}\)，然后求出每一项的\(q\)次幂，然后\(\text{FWT}\)回来
(这样不是\(n2^n\)的吗)
显然的可以发现，\(F^i\)的任意位置系数\([x^S]F^i\)只与\(|S|\)有关，所以实际上只需要存下含有\(0,1,2\cdots,n\)个元素的项的系数即可
考虑在这样的多项式上模拟原先的\(\text{FWT}\)过程
按照快速沃尔什变换的式子，令$G=\text{FWT(F)}=\sum_{S}x^S \sum_{T}(-1)^{|S\cap T|}\cdot [x^T]F $
考虑枚举\(|S|,|T|,|S\cap T|\)，然后组合数计算系数，令\(F_i\)表示\([x^S]F(|S|=i)\)
则有\(\begin{aligned} G_i=\sum_{j}F_j\sum_k C_i^k\cdot C_{n-i}^{j-k}(-1)^{k}\end{aligned}\)，其中\(C_{n-i}^{j-k}\)表示从不相交的部分里选出\(j\)中剩下的元素
按照该式即可完成\(O(n^3)\)模拟\(\text{FWT}\)，注意\(\text{IFWT}\)时需要除掉\(2^n\)
算上快速幂复杂度应为\(O(n^3+n\log q)\)
\(O(n^2\log n)\)
用\(\text{NTT}\)优化上式
\(O(n^2)\)
依然考虑上面\(\text{FWT}\)的转移式，发现\(i\)项得到\(j\)项的贡献为一个常数，考虑直接计算这个常数\(W_{i,j}=\sum_k C_i^k\cdot C_{n-i}^{j-k}(-1)^{k}\)
我们知道组合数就是二项展开的结果，所以发现实际就是\(W_{i,j}=[x^j] (-x+1)^i\cdot (x+1)^{n-i}\)
该式可以\(O(n^2)\)按照\(i\)递推求出，每次乘上一个\(\frac{1-x}{1+x}\)即可
实际复杂度为\(O(n^2+n\log q)\)
实际上应该还可以处理一些稍微复杂点的问题，比如每次可以操作若干个位置
不知道能否优化到\(n^2\)以下
// NTT Version
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
typedef double db;
typedef long double ldb;
typedef pair <int,int> Pii;
#define reg register
#define pb push_back
#define mp make_pair
#define Mod1(x) ((x>=P)&&(x-=P))
#define Mod2(x) ((x<0)&&(x+=P))
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
bool Mbe;
const int N=2050,P=998244353;
int n,m,Z; ll k;
int A[N],B[N];
ll qpow(ll x,ll k=P-2) {
ll res=1;
for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
return res;
}
int C[N][N];
int rev[N],T[N],IT[N];
void NTT(int n,int *a,int f){
rep(i,0,n-1) if(i<rev[i]) swap(a[i],a[rev[i]]);
static int e[N>>1];
for(reg int i=e[0]=1;i<n;i<<=1) {
ll t=f==1?T[i]:IT[i];
//qpow(f==1?3:(P+1)/3,(P-1)/i/2);
for(reg int j=i-2;j>=0;j-=2) e[j+1]=t*(e[j]=e[j>>1])%P;
for(reg int l=0;l<n;l+=i*2) {
for(reg int j=l;j<l+i;++j) {
reg int t=1ll*a[j+i]*e[j-l]%P;
a[j+i]=a[j]-t,Mod2(a[j+i]);
a[j]+=t,Mod1(a[j]);
}
}
}
if(f==-1) {
ll base=qpow(n);
rep(i,0,n-1) a[i]=a[i]*base%P;
}
}
int Init(int n){
int R=1,cc=-1;
while(R<=n) R<<=1,cc++;
rep(i,1,R-1) rev[i]=(rev[i>>1]>>1)|((i&1)<<cc);
return R;
}
int mk1[N],mk2[N];
void FWT(int n,int *A,int *B,int f,int *mk){
static int X[N],Y[N];
rep(i,0,n) {
memset(X,0,sizeof X),memset(Y,0,sizeof Y);
rep(k,0,i) X[i-k]=(k&1)?P-C[i][k]:C[i][k];
rep(j,0,n) Y[j]=A[j];
int R=Init(n+i+1);
NTT(R,X,1),NTT(R,Y,1);
rep(j,0,R-1) X[j]=1ll*X[j]*Y[j]%P;
NTT(R,X,-1);
rep(j,0,n-i) B[i]=(B[i]+1ll*C[n-i][j]*X[i+j])%P;
}
}
void Solve(int n,int *A,int *mk) {
static int X[N],Y[N];
memset(X,0,sizeof X),memset(Y,0,sizeof Y);
X[1]=1,FWT(n,X,Y,1,mk);
rep(i,0,n) Y[i]=qpow(Y[i],k);
memset(X,0,sizeof X);
FWT(n,Y,X,2,mk);
ll base=qpow(qpow(2,n),P-2);
rep(i,0,n) A[i]=X[i]*base%P*C[n][i]%P;
}
bool Med;
int main(){
//fprintf(stderr,"%.2lf\n",(&Med-&Mbe)/1024.0/1024.0);
freopen("clone.in","r",stdin),freopen("clone.out","w",stdout);
rep(i,0,N-1) rep(j,C[i][0]=1,i) C[i][j]=(C[i-1][j-1]+C[i-1][j])%P;
for(int i=1;i<N;i<<=1) {
T[i]=qpow(3,(P-1)/i/2);
IT[i]=qpow((P+1)/3,(P-1)/i/2);
}
scanf("%d%d%lld%d",&n,&m,&k,&Z);
rep(i,0,n) rep(j,0,m) {
if(i*(m-j)+j*(n-i)!=Z) continue;
mk1[i]=1,mk2[j]=1;
}
Solve(n,A,mk1),Solve(m,B,mk2);
int ans=0;
rep(i,0,n) rep(j,0,m) {
if(i*(m-j)+j*(n-i)!=Z) continue;
ans=(ans+1ll*A[i]*B[j])%P;
}
printf("%d\n",ans);
}
// n^2
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
#define Mod1(x) ((x>=P)&&(x-=P))
#define Mod2(x) ((x<0)&&(x+=P))
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
const int N=2050,P=998244353;
int n,m,Z; ll k;
int A[N],B[N];
ll qpow(ll x,ll k=P-2) {
ll res=1;
for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
return res;
}
int C[N][N],W[N][N],mk1[N],mk2[N];
void Solve(int n,int *A,int *mk) {
static int X[N],Y[N];
rep(i,0,n) W[0][i]=C[n][i];
rep(i,1,n) {
rep(j,0,n) W[i][j]=W[i-1][j]-(j?W[i][j-1]:0),Mod2(W[i][j]);
drep(j,n,0) W[i][j]=(j?-W[i][j-1]:0)+W[i][j],Mod2(W[i][j]);
}
memset(Y,0,sizeof Y);
rep(i,0,n) X[i]=qpow(W[i][1],k);
rep(i,0,n) if(mk[i]) rep(j,0,n) Y[i]=(Y[i]+1ll*W[i][j]*X[j])%P;
ll base=qpow(qpow(2,n),P-2);
rep(i,0,n) A[i]=Y[i]*base%P*C[n][i]%P;
}
int main(){
freopen("clone.in","r",stdin),freopen("clone.out","w",stdout);
scanf("%d%d%lld%d",&n,&m,&k,&Z);
rep(i,0,max(n,m)) rep(j,C[i][0]=1,i) C[i][j]=C[i-1][j-1]+C[i-1][j],Mod1(C[i][j]);
rep(i,0,n) rep(j,0,m) {
if(i*(m-j)+j*(n-i)!=Z) continue;
mk1[i]=mk2[j]=1;
}
Solve(n,A,mk1),Solve(m,B,mk2);
int ans=0;
rep(i,0,n) rep(j,0,m) {
if(i*(m-j)+j*(n-i)!=Z) continue;
ans=(ans+1ll*A[i]*B[j])%P;
}
printf("%d\n",ans);
}
