题意:
大概是求
∏
i
=
1
n
(
x
+
i
)
\prod_{i=1}^n(x+i)
∏i=1n(x+i)系数模
p
p
p意义下的分布。
题解:
分为
(
∏
i
=
1
p
−
1
(
x
+
i
)
)
⌊
n
p
⌋
(\prod_{i=1}^{p-1}(x+i))^{\lfloor\frac{n}{p}\rfloor}
(∏i=1p−1(x+i))⌊pn⌋,以及
(
∏
i
=
1
n
 
m
o
d
 
p
)
(
x
+
i
)
(\prod_{i=1}^{n \bmod p})(x+i)
(∏i=1nmodp)(x+i)。
前者等于 ( x p − 1 − 1 ) ⌊ n p ⌋ (x^{p-1}-1)^{\lfloor\frac{n}{p}\rfloor} (xp−1−1)⌊pn⌋,后者最高次数小于 p − 1 p-1 p−1(等于把前面的次数加一即可)。
发现两部分系数不影响,可以分开算,后面直接倍增FFT了,前面二项式展开然后Lucas定理合并即可,合并可以用原根+FFT,时间复杂度 O ( p log p + p log n ) O(p \log p + p \log_n) O(plogp+plogn)。
#include <bits/stdc++.h>
using namespace std;
typedef double DB;
typedef long long LL;
const int N=3e5+50, mod=998244353;
LL n,k; int p,G,ind[N],ori[N];
inline int add(int x,int y,int md) {return (x+y>=md) ? (x+y-md) : (x+y);}
inline int dec(int x,int y,int md) {return (x-y<0) ? (x-y+md) : (x-y);}
inline int mul(int x,int y,int md) {return (long long)x*y%md;}
inline int power(int a,int b,int md,int rs=1) {for(;b;b>>=1,a=mul(a,a,md)) if(b&1) rs=mul(rs,a,md); return rs;}
namespace FFT {
const int N=3e6+50, M=sqrt(mod);
const DB PI2=acos(-1.0)*2;
int k,pos[N];
struct CP {
DB r,i;
CP(DB r=0,DB i=0) : r(r),i(i) {}
friend inline CP operator +(const CP &a,const CP &b) {return CP(a.r+b.r,a.i+b.i);}
friend inline CP operator -(const CP &a,const CP &b) {return CP(a.r-b.r,a.i-b.i);}
friend inline CP operator *(const CP &a,const CP &b) {return CP(a.r*b.r-a.i*b.i,a.r*b.i+a.i*b.r);}
} A[N],B[N],C[N],w[N];
inline void clear(CP *a) {for(int i=0;i<k;i++) a[i]=CP(0,0);}
inline void init(int nn) {
for(k=1;k<=nn;k<<=1);
for(int i=0;i<k;i++) pos[i]=(i&1) ? (pos[i>>1]>>1)^(k>>1) : (pos[i>>1]>>1);
}
inline void mul(CP *a,CP *b,CP *c) {
for(int i=0;i<k;i++) c[i]=a[i]*b[i];
}
inline void dft(CP *a,int o=1) {
for(int i=1;i<k;i++)
if(pos[i]>i) swap(a[pos[i]],a[i]);
for(int bl=1;bl<k;bl<<=1) {
int tl=bl<<1;
for(int i=0;i<bl;i++) w[i]=CP(cos(PI2/tl*i),sin(PI2/tl*i*o));
for(int bg=0;bg<k;bg+=tl)
for(int j=0;j<bl;j++) {
CP &t1=a[bg+j], &t2=a[bg+j+bl], t=t2*w[j];
t2=t1-t; t1=t1+t;
}
}
if(!~o) for(int i=0;i<k;i++) a[i].r/=k, a[i].i/=k;
}
}
using FFT::CP;
using FFT::A;
using FFT::B;
using FFT::C;
using FFT::M;
struct poly {
vector <int> a;
poly(int d=0,int t=0) {a.resize(d+1); a[d]=t;}
inline int deg() const {return a.size()-1;}
inline int& operator [](const int &i) {return a[i];}
inline const int& operator[](const int &i) const {return a[i];}
friend inline poly mul(const poly &a,const poly &b,const int &md) {
poly c(a.deg()+b.deg());
FFT::init(c.deg());
FFT::clear(A); FFT::clear(B);
for(int i=0;i<=a.deg();i++) A[i]=CP(a[i]/M,a[i]%M);
FFT::dft(A);
for(int i=0;i<=b.deg();i++) B[i]=CP(b[i]/M,0);
FFT::dft(B);
FFT::mul(A,B,C);
FFT::dft(C,-1);
for(int i=0;i<=c.deg();i++) c[i]=((LL)(C[i].r+0.5)%md*M%md*M%md+(LL)(C[i].i+0.5)%md*M%md)%md;
FFT::clear(B);
for(int i=0;i<=b.deg();i++) B[i]=CP(b[i]%M,0);
FFT::dft(B);
FFT::mul(A,B,C);
FFT::dft(C,-1);
for(int i=0;i<=c.deg();i++)
c[i]=add(c[i],((LL)(C[i].r+0.5)%md*M%md+(LL)(C[i].i+0.5))%md,md);
return c;
}
inline void pt() {
for(int i=0;i<=deg();i++)
cout<<a[i]<<' ';
cout<<'\n';
}
};
struct combin {
int fac[N],ifac[N];
inline void init() {
fac[0]=1;
for(int i=1;i<p;i++) fac[i]=mul(fac[i-1],i,p);
ifac[p-1]=power(fac[p-1],p-2,p);
for(int i=p-2;~i;i--) ifac[i]=mul(ifac[i+1],i+1,p);
}
inline int C(int a,int b) {return mul(fac[a],mul(ifac[b],ifac[a-b],p),p);}
} cb;
namespace ORD {
vector <int> ff;
inline bool check(int v) {
for(auto u:ff)
if(power(v,(p-1)/u,p)==1) return false;
return true;
}
inline void find_G() {
int nn=p-1;
for(int i=2;i*i<=nn;i++)
if(!(nn%i)) {
while(!(nn%i)) nn/=i;
ff.push_back(i);
}
if(nn!=1) ff.push_back(nn);
for(G=1;;G++)
if(check(G)) break;
}
}
inline poly calc_fac(int nn) {
if(!nn) return poly(0,1);
if(nn&1) {
poly t=poly(1,1); t[0]=nn;
return mul(calc_fac(nn-1),t,p);
} else {
int l=nn/2;
poly f=calc_fac(l);
// f.pt();
poly c(l,0),d(l,0);
for(int i=0,pw=1;i<=l;i++) c[i]=mul(f[i],mul(pw,cb.fac[i],p),p), pw=mul(pw,l,p);
for(int i=0;i<=l;i++) d[l-i]=cb.ifac[i];
poly g=mul(c,d,p);
poly h(l,0);
const int inv=power(l,p-2,p);
for(int i=0,iv=1;i<=l;i++) h[i]=mul(iv,mul(g[i+l],cb.ifac[i],p),p), iv=mul(iv,inv,p);
return mul(f,h,p);
}
}
inline poly calc_pow(LL k) {
poly f(0,1);
long long base=1;
while(k) {
int z=k%p;
k/=p;
poly h(p,0);
for(int i=0;i<=z;i++) {
int tag=(base*i%2) ? (p-cb.C(z,i)) : cb.C(z,i);
h[ind[tag]]++;
}
f=mul(f,h,mod);
if(f.deg()>=p-1) {
for(int i=p-1;i<=f.deg();i++)
f[i%(p-1)]=add(f[i%(p-1)],f[i],mod);
f.a.resize(p-1);
}
while(f.a.size()>1 && !f.a.back()) f.a.pop_back();
base=base*p;
}
return f;
}
int main() {
cin>>n>>p;
cb.init();
ORD::find_G();
for(int i=0,pw=1;i<p-1;i++) ind[pw]=i, ori[i]=pw, pw=mul(pw,G,p);
k=((n%p)==p-1) ? n/p+1 : n/p;
poly f1=calc_fac(max(n-k*p,0ll));
poly f2=calc_pow(k);
poly f3(p-2,0);
for(int i=0;i<=f1.deg();i++) if(f1[i]) ++f3[ind[f1[i]]];
poly f=mul(f3,f2,mod);
if(f.deg()>=p-1) {
for(int i=p-1;i<=f.deg();i++)
f[i%(p-1)]=add(f[i%(p-1)],f[i],mod);
} f.a.resize(p-1);
int total=(n+1)%mod;
for(int i=0;i<p-1;i++) total=dec(total,f[i],mod);
cout<<total<<'\n';
for(int i=1;i<p;i++) cout<<f[ind[i]]<<'\n';
}