题解：P10890 【烂题杯 Round 1】可持久化糖果树
qczrz6v4nhp6u
·
·
题解
很牛的题！拜谢出题人。
个人感觉 3^m 的做法更加自然，因为这是直接单位根反演的结果。
Solution
考虑单位根反演：
[k|n]=\frac 1k\sum_{i=0}^{k-1}\omega_k^{in}
我们应用它来刻画 \sum_{x=1}^ka_{i,j,x}\equiv 0\pmod 3。一个节点 i 是合法的当且仅当：
\frac 1{3^m}\prod_{j=1}^m(1+\omega_3^{\sum a_{i,j,x}}+\omega_3^{2\sum a_{i,j,x}})=1
则答案即为
\frac 1{3^m}\sum_{i=1}^n\sum_{j=1}^m(1+\omega_3^{\sum a_{i,j,x}}+\omega_3^{2\sum a_{i,j,x}})
将式子展开，得到
&\frac 1{3^m}\sum_{i=1}^n\sum_b\omega_3^{\sum_{j=1}^mb_j\sum_{x=1}^k a_{i,j,x}}\\
=&\frac 1{3^m}\sum_{i=1}^n\sum_{b}\omega_3^{\sum_{x=1}^k\sum_{j=1}^ma_{i,j,x}b_j}
\end{aligned}
其中序列 b 满足 |b|=m,b_i\in\{0,1,2\}。
现在我们已经可以做到 n3^mk 或更快的单次询问了，考虑支持修改。设 a_{i,j,x}\gets a_{i,j,x}\times d_x，继续推式子：
&\frac 1{3^m}\sum_{i=1}^n\sum_{b}\omega_3^{\sum_{x=1}^k \sum_{j=1}^md_xa_{i,j,x}b_j}\\
=&\frac 1{3^m}\sum_{i=1}^n\sum_{b}\omega_3^{\sum_{x=1}^kd_x\sum_{j=1}^ma_{i,j,x}b_j}
\end{aligned}
对于确定的 a_i 和 b，我们可以得到序列 c，满足 c_x=\sum_{j=1}^ma_{i,j,x}b_j。记录 cnt_c 表示上述 c 的出现次数，答案即为：
\frac 1{3^m}\sum_ccnt_c\prod_{x=1}^k\omega_3^{c_xd_x}
总的复杂度即为 $O(n3^mmk+k3^k+qk)$。
Bonus：其实上述的 $k3^k$ DP 即为 3-FWT。
### Code
朴素的实现需要扩域。但我们可以选一个大于等于 $n$ 的质数 $p$ 满足 $3|(p-1)$，这样在 $\mathbb F_p$ 中就存在 $\omega_3$，算答案时直接对 $p$ 取模即可。
```cpp
#include<bits/stdc++.h>
using namespace std;
using ui=unsigned int;
using ll=long long;
using ull=unsigned long long;
using i128=__int128;
using u128=__uint128_t;
using pii=pair<int,int>;
#define fi first
#define se second
constexpr int N=1e6+5,M=5,K=13,mod=1145140831,w1=431040359,w2=714100471;
int n,m,k,X,q,a[M][K],b[K],p3[K],num[N],tmp[K];
ll F[N];
inline ll qpow(ll a,ll b){
ll res=1;
for(;b;b>>=1,a=a*a%mod)
if(b&1)res=res*a%mod;
return res;
}
void dfs(int x){
if(x==m+1){
int cur=0,sum;
for(int i=1;i<=k;i++){
sum=0;
for(int j=1;j<=m;j++)
sum+=b[j]*a[j][i];
cur+=p3[i-1]*(sum%3);
}
++F[cur];
return;
}
for(int i=0;i<3;i++){
b[x]=i;
dfs(x+1);
}
}
void FWT(ll *F,int n){
for(int i=1;i<n;i*=3)
for(int j=0;j<n;j+=i*3)
for(int k=0;k<i;k++){
ll x=F[j+k],y=F[i+j+k],z=F[(i<<1)+j+k];
F[j+k]=(x+y+z)%mod;
F[i+j+k]=(x+w1*y+w2*z)%mod;
F[(i<<1)+j+k]=(x+w2*y+w1*z)%mod;
}
}
int main(){
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
cin>>n>>m>>k>>X;
p3[0]=1;
for(int i=1;i<=k;i++)p3[i]=p3[i-1]*3;
for(int i=1;i<=n;i++){
for(int j=1;j<=m;j++)
for(int x=1;x<=k;x++)
a[j][x]=(X*(i+1ll)+(X^(j*x)))%1000000000%3;
dfs(1);
}
FWT(F,p3[k]);
ll inv=qpow(3,mod-m-1);
for(int i=0;i<p3[k];i++)F[i]=F[i]*inv%mod;
for(int i=1;i<=k;i++)num[0]+=p3[i-1];
int ans=F[num[0]];
cin>>q;
for(int i=1;i<=q;i++){
int x=(X^i)%i,y=(X^i)%k+1,z=(X+(X^i))%999999999%3;
for(int j=1;j<=k;j++)tmp[j]=num[x]/p3[j-1]%3;
tmp[y]=(tmp[y]*z)%3;
for(int j=1;j<=k;j++)num[i]+=tmp[j]*p3[j-1];
ans^=F[num[i]];
}
cout<<ans<<'\n';
return 0;
}
```