【LuoguP7439】[KrOI2021]Feux Follets 弱化版
【LuoguP7439】[KrOI2021]Feux Follets 弱化版
by AmanoKumiko
Description
设\(cyc_\pi\)表示将长度为\(n\)的排列\(\pi\)当成置换时能分解的循环置换个数。给定\(n,k\)和一个\(k-1\)次多项式\(F(x)\)
求：
\[\sum_{\pi}F(cyc_{\pi})
\]
其中\(\pi\)是长度为\(n\)的错排列
对\(998244353\)取模
Input
一行两个整数\(n,k\)
第二行\(k\)个整数读入\(F(x)\)
Output
一行一个整数表示答案
Sample Input
6 4
11 43 27 7
Sample Output
53070
Data Constraint
\(1\le n,k\le 10^5\)
Solution
一个长度\(\ge2\)的环\(EGF\)为\(G(x)=\sum_{i\ge2}\frac{x^i}{i}=-\ln(1-x)-x\)
那么我们要求的就是\(\sum_{i=0}^{+\infty}\frac{G^i(x)}{i}[x^n]n!F(i)\)
Sol1
考虑加入一元分离信息，直接求出每个次幂的答案，最后多点求值
令\(H(x)=\sum_{i=0}^{+\infty}u^iG^i(x)=\frac{1}{1-uG}\)
由于\(G(x)\)的最低次数为\(2\)，故复合逆不存在
考虑设\(\frac{g^2}{2}=-\ln(1-x)-x\)
即\(g=\sqrt{-2\ln(1-x)-2x}\)
显然复合逆存在，设为\(f\)
那么
\[\begin{align}
[x^n]H(x)&=\frac{1}{n}[x^{n-1}]\frac{ux}{(1-\frac{ux^2}{2})^2}(\frac{f}{x})^{-n}\\
&=\frac{1}{n}[x^{n-1}]ux\sum_{i=0}^{+\infty}(i+1)(\frac{ux^2}{2})^i(\frac{f}{x})^{-n}\\
[x^nu^k]&=\frac{k}{n2^{k-1}}[x^{n-2k}](\frac{f}{x})^{-n}
\end{align}
\]
Sol2
其实两种方法都要多点求值
由于次幂的存在，考虑把多项式转成牛顿级数，令其为\(F^*(x)\)
那么有
\[\begin{align}
\sum_{i=0}^{+\infty}\frac{G^i(x)}{i}[x^n]n!F(i)&=\sum_{i=0}^{+\infty}\frac{G^i(x)}{i}[x^n]n!\sum_{j=0}^{k-1}f_j^*\binom{i}{j}\\
&=\sum_{j=0}^{k-1}f_j^*[y^j]e^{(1+y)(-\ln(1-x)-x)}
\end{align}
\]
这里用组合意义，加入辅助因子变形
我们先处理\(e^{y(-\ln(1-x)-x)}\)，最后多项式平移弄回去
到这里就和上面差不多了，直接快进一下
\[[x^ny^k]e^{\frac{yg^2(x)}{2}}=\frac{1}{n2^{k-1}(k-1)!}[x^{n-2k}](\frac{f}{x})^{-n}
\]
由于一些项系数可能为\(0\)，在牛顿迭代里面有亿些偏移的细节
最后给出我的卡着\(1.50s\)过的垃圾代码
Code
#include<bits/stdc++.h>
using namespace std;
#define fo(i,a,b) for(int i=a;i<=b;i++)
#define fd(i,a,b) for(int i=a;i>=b;i--)
#define N 600010
#define mo 998244353
#define LL long long
#define ULL unsigned long long
namespace IO{
const int sz=1<<22;
char a[sz+5],b[sz+5],*p1=a,*p2=a,*t=b,p[105];
inline char gc(){
return p1==p2?(p2=(p1=a)+fread(a,1,sz,stdin),p1==p2?EOF:*p1++):*p1++;
}
template<class T>void read(T&x){
x=0;char c=gc();
for(;c<'0'||c>'9';c=gc());
for(;c>='0'&&c<='9';c=gc())x=x*10+(c-'0');
}
inline void flush(){fwrite(b,1,t-b,stdout),t=b;}
inline void pc(char x){*t++=x;if(t-b==sz)flush();}
template<class T>void write(T x,char c='\n'){
if(x==0)pc('0');int t=0;
for(;x;x/=10)p[++t]=x%10+'0';
for(;t;--t)pc(p[t]);pc(c);
}
struct F{~F(){flush();}}f;
}
using IO::read;
using IO::write;
int rev[N],G1[N],G2[N],fac[N],ifac[N],inv[N];
inline int getlen(int x){return 1<<32-__builtin_clz(x);}
inline int mod(int x){return x>=mo?x-mo:x;}
inline int mi(int x,int y){
if(!y)return 1;
if(y==1)return x;
return y%2?1ll*x*mi(1ll*x*x%mo,y/2)%mo:mi(1ll*x*x%mo,y/2);
}
void init(){
fac[0]=ifac[0]=1;
fo(i,1,N-10)fac[i]=1ll*fac[i-1]*i%mo,inv[i]=(i==1?1:1ll*mo/i*mod(mo-1ll*inv[mo%i]%mo)%mo);
ifac[N-10]=mi(fac[N-10],mo-2);
fd(i,N-11,1)ifac[i]=1ll*ifac[i+1]*(i+1)%mo;
for(int l=1;l<=N-10;l<<=1)G1[l]=mi(3,(mo-1)/(l*2)),G2[l]=mi(G1[l],mo-2);
}
inline void BRT(int x){fo(i,0,x-1)rev[i]=(rev[i>>1]>>1)|((i&1)?(x>>1):0);}
struct poly{
vector<int>val;
inline poly(int x=0){if(x)val.push_back(x);}
inline poly(const vector<int>&x){val=x;}
inline void Rev(){reverse(val.begin(),val.end());}
inline void ins(int x){val.push_back(x);}
inline void clear(){vector<int>().swap(val);}
inline int sz(){return val.size();}
inline void rsz(int x){val.resize(x);}
inline void shrink(){for(;sz()&&!val.back();val.pop_back());}
inline poly modxn(int x){
if(val.size()<=x)return poly(val);
else return poly(vector<int>(val.begin(),val.begin()+x));
}
inline int operator[](int x)const{
if(x<0||x>=val.size())return 0;
return val[x];
}
inline void NTT(int x){
static ULL f[N],w[N];
w[0]=1;
fo(i,0,sz()-1)f[i]=(((LL)mo<<5)+val[rev[i]])%mo;
for(int mid=1;mid<sz();mid<<=1){
int tmp=(x==1?G1[mid]:G2[mid]);
fo(i,1,mid-1)w[i]=w[i-1]*tmp%mo;
for(int i=0;i<sz();i+=(mid<<1)){
fo(j,0,mid-1){
int t=w[j]*f[i|j|mid]%mo;
f[i|j|mid]=f[i|j]+mo-t;f[i|j]+=t;
}
}
}
if(x==-1){int tmp=inv[sz()];fo(i,0,sz()-1)val[i]=f[i]%mo*tmp%mo;}
else{fo(i,0,sz()-1)val[i]=f[i]%mo;}
}
inline void DFT(){NTT(1);}
inline void IDFT(){NTT(-1);}
inline friend poly operator*(poly x,poly y){
if(x.sz()<30||y.sz()<30){
if(x.sz()>y.sz())swap(x,y);
poly ret;
ret.rsz(x.sz()+y.sz());
fo(i,0,ret.sz()-1){
for(int j=0;j<=i&&j<x.sz();j++)
ret.val[i]=mod(ret.val[i]+1ll*x[j]*y[i-j]%mo);
}
return ret;
}
int l=getlen(x.sz()+y.sz()-2);
x.rsz(l);y.rsz(l);BRT(l);
x.DFT();y.DFT();
fo(i,0,l-1)x.val[i]=1ll*x[i]*y[i]%mo;
x.IDFT();
return x;
}
inline friend poly operator+(poly x,poly y){
poly ret;
ret.rsz(max(x.sz(),y.sz()));
fo(i,0,ret.sz()-1)ret.val[i]=mod(x[i]+y[i]);
return ret;
}
inline friend poly operator-(poly x,poly y){
poly ret;
ret.rsz(max(x.sz(),y.sz()));
fo(i,0,ret.sz()-1)ret.val[i]=mod(x[i]-y[i]+mo);
return ret;
}
inline poly &operator*=(poly x){return (*this)=(*this)*x;}
inline poly &operator+=(poly x){return (*this)=(*this)+x;}
inline poly &operator-=(poly x){return (*this)=(*this)-x;}
inline poly deriv(){
poly f;
f.rsz(sz()-1);
fo(i,0,sz()-2)f.val[i]=1ll*(i+1)*val[i+1]%mo;
return f;
}
inline poly integ(){
poly f;
f.rsz(sz()+1);
fo(i,1,sz())f.val[i]=1ll*val[i-1]*inv[i]%mo;
return f;
}
inline poly inver(int Len){
poly f,g;
f.clear();g.clear();
g.ins(mi(val[0],mo-2));
for(int i=1;i<Len*2;i<<=1){
int Len=i<<1;
f.rsz(Len);
g.rsz(Len);
BRT(Len);
fo(j,0,i-1)f.val[j]=(j<val.size()?val[j]:0);
f.DFT();g.DFT();
fo(j,0,Len-1)g.val[j]=1ll*g[j]*mod(mo+2-1ll*f[j]*g[j]%mo)%mo;
g.IDFT();
g.rsz(i);
}
return g.modxn(Len);
}
inline poly Sqrt(int Len){
int i2=inv[2];
poly f,g,tmp;
f.clear();g.clear();
g.ins(1);
for(int i=1;i<Len*2;i<<=1){
int Len=i<<1;
f.rsz(Len);
tmp=g.inver(i);
tmp.rsz(Len);
BRT(Len);
fo(j,0,i-1)f.val[j]=(j<val.size()?val[j]:0);
f.DFT();tmp.DFT();
fo(j,0,Len-1)f.val[j]=1ll*f[j]*tmp[j]%mo;
f.IDFT();
g.rsz(i);
fo(j,0,i-1)g.val[j]=1ll*i2*mod(g[j]+f[j])%mo;
}
return g.modxn(Len);
}
inline poly Ln(int Len){
return (deriv()*inver(Len)).integ().modxn(Len);
}
inline poly Exp(int Len){
poly f;
f.clear();
f.ins(1);
for(int i=2;i<Len*2;i<<=1)f=(f*(1-f.Ln(i)+modxn(i))).modxn(i);
return f.modxn(Len);
}
inline poly Pow(int Len,int k){
poly f;
f.clear();
int tail=0;
while(val[tail]==0&&tail<sz())tail++;
if(tail>=sz())return f;
if(tail*k>=Len)return f;
f.rsz(Len);
int Mul=mi(val[tail],mo-2);
fo(i,0,min(Len-1,sz()-tail-1))f.val[i]=1ll*val[i+tail]*Mul%mo;
Mul=mi(val[tail],k);
f=f.Ln(Len);
fo(i,0,Len-1)f.val[i]=1ll*f[i]*(k%mo)%mo;
f=f.Exp(Len);
fd(i,Len-1,tail*k)f.val[i]=1ll*f[i-tail*k]*Mul%mo;
fo(i,0,tail*k-1)f.val[i]=0;
return f;
}
};
struct Evaluation{
#define ls x<<1
#define rs (x<<1)|1
poly f,g[N],h[N];
int b[N];
void solve1(int x,int l,int r){
if(l==r){g[x].rsz(2);g[x].val[0]=1;g[x].val[1]=mod(mo-l);return;}
int mid=l+r>>1,Len=1;
solve1(ls,l,mid);solve1(rs,mid+1,r);
while(Len<r-l+1+1)Len<<=1;BRT(Len);
poly A=g[ls],B=g[rs];
A.rsz(Len);B.rsz(Len);
A.DFT();B.DFT();
fo(i,0,Len-1)A.val[i]=1ll*A[i]*B[i]%mo;
A.IDFT();
g[x]=A.modxn(r-l+2);
}
void solve2(int x,int l,int r){
if(l==r){b[l]=h[x][0];return;}
int mid=l+r>>1,Len=1;
while(Len<(r-l+1<<1)+1)Len<<=1;BRT(Len);
g[ls].Rev();g[rs].Rev();
g[ls].rsz(Len);g[rs].rsz(Len);h[x].rsz(Len);
g[ls].DFT();g[rs].DFT();h[x].DFT();
fo(i,0,Len-1)g[ls].val[i]=1ll*g[ls][i]*h[x][i]%mo;
fo(i,0,Len-1)g[rs].val[i]=1ll*g[rs][i]*h[x][i]%mo;
g[ls].IDFT();g[rs].IDFT();
h[ls].rsz(mid-l+1+1);
fo(i,0,mid-l+1)h[ls].val[i]=g[rs][i+r-mid];
h[rs].rsz(r-mid+1);
fo(i,0,r-mid)h[rs].val[i]=g[ls][i+mid-l+1];
solve2(ls,l,mid);solve2(rs,mid+1,r);
}
void get(){
int n=f.sz()-1,m=n+1;
solve1(1,0,n);
g[1]=g[1].inver(n+1);
g[1].Rev();
int Len=1;
while(Len<(n<<1)+1)Len<<=1;BRT(Len);
g[1].rsz(Len);f.rsz(Len);
g[1].DFT();f.DFT();
fo(i,0,Len-1)g[1].val[i]=1ll*g[1][i]*f[i]%mo;
g[1].IDFT();
h[1].rsz(n+1);
fo(i,0,n)h[1].val[i]=g[1][i+n];
solve2(1,0,n);
}
}t;
inline poly Newton(int Len){
poly f,res;
f.ins(0);f.ins(1);
for(int i=2;i<Len;){
i<<=1;
poly tmp=1-f;tmp.rsz(i);
poly _Ln=(tmp.Ln(i+1)*(mo-2)-(2*f)).modxn(i+1);
poly _Sqrt;_Sqrt.clear();
fo(j,0,i-2)_Sqrt.ins(_Ln[j+2]);
_Sqrt=_Sqrt.Sqrt(i-1);
_Sqrt.rsz(i);
fd(j,i-1,1)_Sqrt.val[j]=_Sqrt[j-1];
_Sqrt.val[0]=0;
fo(j,0,i-1)_Ln.val[j]=_Ln[j+1];
_Ln.val.pop_back();
_Ln-=_Sqrt;
_Ln=(_Ln*tmp).modxn(i+1);
poly _inv;_inv.clear();
fo(j,0,i-1)_inv.ins(f[j+1]);
_Ln=(_inv.inver(i)*_Ln).modxn(i);
f.rsz(i);
f-=_Ln;
}
return f.modxn(Len);
}
poly f,g,h,fc;
int n,k,ans,st[N];
int main(){
init();
read(n);read(k);
fo(i,0,k-1){
int x;
read(x);
f.ins(x);
}
t.f=f;
t.get();
poly A,B;
fo(i,0,k-1)A.ins(1ll*t.b[i]*ifac[i]%mo);
fo(i,0,k-1)B.ins((i&1)?mo-ifac[i]:ifac[i]);
A*=B;
f=A.modxn(k);
fo(i,0,k-1)f.val[i]=1ll*f[i]*fac[i]%mo;
g=Newton(n+2);
fo(i,0,n)g.val[i]=g[i+1];
g.val.pop_back();
g=g.inver(n+1);
g=g.Pow(n+1,n);
h.rsz(k);
h.val[0]=0;
int tmp=1,i2=inv[2];
fo(i,1,k-1){
h.val[i]=1ll*tmp*i%mo*fac[n-1]%mo*g[n-2*i]%mo;
tmp=1ll*tmp*i2%mo;
}
h.Rev();
fo(i,0,k-1)fc.ins(ifac[i]);
fc=(fc*h).modxn(k);
fc.Rev();
fo(i,0,k-1)ans=mod(ans+1ll*f[i]*fc[i]%mo*ifac[i]%mo);
printf("%d",ans);
return 0;
}
