LP3321 [SDOI2015]序列统计 [NTT/DP]
qwq...
Problem
https://www.luogu.org/problemnew/show/P3321
Solution
由于 m 是质数所以把 $ab \equiv x \pmod {m}$ 换成 $g^{i+j} \equiv g^k \pmod {m}$ , $g$ 是模 $m$ 意义下的原根。
那么考虑 dp, $f_{k, i}=\sum f_{k-1, j}~d_{i-j}$ 。 $f_{k,i}$ 是前 $k$ 局得到 $g_i$ 的方案数,而 $d_i$ 是 取 1 个数 能得到的 $g_i$ 的方案个数。
根据费马小定理做多项式乘法并将前后两半合并,就可以按局数进行多项式快速幂来转移。
记得长度 2 倍。
记得质数不要筛错。
记得该上 1ll 的地方不要忘记上。
Code
#include<bits/stdc++.h>
#define MOD (1004535809)
using namespace std;
typedef long long ll;
int qpow(int x, int y, int mo) {
x%=mo;
int ret=1;
while(y) {
if(y&1) ret=1ll*ret*x%mo;
x=1ll*x*x%mo, y>>=1;
}
return ret;
}
int g(int x) {
x--;
int pri[30], p=0, rx=x;
for(int i=2;i*i<=x;i++)
if(x%i==0) {
while(x%i==0) x/=i;
pri[p++]=i;
}
if(x>1) pri[p++]=x;
char ok;
for(int g=2; g<=rx; g++) {
ok=1;
for(int i=0; i<p; i++)
if(qpow(g, rx/pri[i], rx+1)==1) { ok=false; break; }
if(ok) return g;
}
return -1;
}
const int N=4e4+10;
int r[N], inv, T, M, X, S;
int n, m, id[N], a[N], b[N], x[N], y[N];
void fft(int x[], int y) {
for(int i=1; i<n; i++) if(i<r[i]) swap(x[i], x[r[i]]);
int t, w, o;
for(int i=1; i<n; i<<=1) {
o=qpow(3, y*(MOD-1)/(i<<1)+MOD-1, MOD);
for(int j=0; j<n; j+=(i<<1)) {
w=1;
for(int k=0; k<i; k++, w=1ll*w*o%MOD)
t=1ll*x[i+j+k]*w%MOD, x[i+j+k]=(x[j+k]-t+MOD)%MOD, (x[j+k]+=t)%=MOD;
}
}
}
void work(int a[], int b[]) {
fill(x, x+n, 0), fill(y, y+n, 0);
for(int i=0;i<=M-2;i++) x[i]=a[i];
for(int i=0;i<=M-2;i++) y[i]=b[i];
fft(x, 1), fft(y, 1); for(int i=0; i<n; i++) x[i]=1ll*x[i]*y[i]%MOD;
fft(x, -1); for(int i=0; i<=M-2; i++) a[i]=1ll*(x[i]+x[i+M-1])*inv%MOD;
}
int main() {
ios::sync_with_stdio(false), cin.tie(0);
cin>>T>>M>>X>>S;
int G=g(M), na;
memset(id, -1, sizeof(id));
for(int i=0, j=1;i<M-1;i++, j=j*G%M) id[j]=i;
for(int i=1; i<=S; i++) {
cin>>na;
if(na) b[id[na]]=1;
}
a[0]=1;
m=(M-2)<<1; int l=0; for(n=1; n<=m; n<<=1) l++;
for(int i=1; i<n; i++) r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
inv=qpow(n, MOD-2, MOD);
while(T) {
if(T&1) work(a, b);
work(b, b), T>>=1;
}
cout<<a[id[X]]<<endl;
return 0;
}
本文链接:https://pst.iorinn.moe/archives/sol-luogu-3321.html
许可: https://pst.iorinn.moe/license.html若无特别说明,博客内的文章默认将采用 CC BY 4.0 许可协议 进行许可☆
%%%
话说为什么是矩阵快速幂啊....
已改emm
不要在意细节(雾