洛谷P4389 付公主的背包 [生成函数,NTT]

传送门


同样是回过头来发现不会做了,要加深一下记忆。


思路

只要听说过生成函数的人相信第一眼都可以想到生成函数。

所以我们要求

\[ans=\prod \sum_n x^{nV}=\prod \frac{1}{1-x^V}
\]

也就是\(\prod (1-x^V)\)。

但这玩意好像还是不会做,怎么办呢?

按照套路,可以先\(\ln\)一下,加起来,再\(\exp\)回去。

所以现在要求

\[\sum \ln(1-x^V)
\]

……

……

……

不会。

不会怎么办?

打表找规律!

经过打表,可以发现\(\ln(1-x^V)=\sum \frac{-1}{n}x^{nV}\)。

然后开个桶记录每种体积有多少,\(n\log n\)加一下系数,就做完了。


代码

#include<bits/stdc++.h>
namespace my_std{
using namespace std;
#define pii pair<int,int>
#define fir first
#define sec second
#define MP make_pair
#define rep(i,x,y) for (int i=(x);i<=(y);i++)
#define drep(i,x,y) for (int i=(x);i>=(y);i--)
#define go(x) for (int i=head[x];i;i=edge[i].nxt)
#define sz 404040
#define mod 998244353
typedef long long ll;
template<typename T>
inline void read(T& t)
{
t=0;char f=0,ch=getchar();
double d=0.1;
while(ch>'9'||ch<'0') f|=(ch=='-'),ch=getchar();
while(ch<='9'&&ch>='0') t=t*10+ch-48,ch=getchar();
if(ch=='.')
{
ch=getchar();
while(ch<='9'&&ch>='0') t+=d*(ch^48),d*=0.1,ch=getchar();
}
t=(f?-t:t);
}
template<typename T,typename... Args>
inline void read(T& t,Args&... args){read(t); read(args...);}
void file()
{
#ifndef ONLINE_JUDGE
freopen("a.txt","r",stdin);
#endif
}
// inline ll mul(ll a,ll b){ll d=(ll)(a*(double)b/mod+0.5);ll ret=a*b-d*mod;if (ret<0) ret+=mod;return ret;}
}
using namespace my_std; ll ksm(ll x,int y)
{
ll ret=1;
for (;y;y>>=1,x=x*x%mod) if (y&1) ret=ret*x%mod;
return ret;
}
ll inv(ll x){return ksm(x,mod-2);} int limit,r[sz];
void NTT_init(int n)
{
limit=1;int l=-1;
while (limit<=n+n) limit<<=1,++l;
rep(i,0,limit-1) r[i]=(r[i>>1]>>1)|((i&1)<<l);
}
void NTT(ll *a,int type)
{
rep(i,0,limit-1) if (i<r[i]) swap(a[i],a[r[i]]);
rep(i,0,limit-1) a[i]%=mod;
for (int mid=1;mid<limit;mid<<=1)
{
ll Wn=ksm(3,(mod-1)/mid>>1);if (type==-1) Wn=inv(Wn);
for (int j=0,len=mid<<1;j<limit;j+=len)
{
ll w=1;
for (int k=0;k<mid;k++,w=w*Wn%mod)
{
ll x=a[j+k],y=a[j+k+mid]*w;
a[j+k]=(x+y)%mod;a[j+k+mid]=(1ll*mod*mod-y+x)%mod;
}
}
}
if (type==1) return;
ll I=inv(limit);
rep(i,0,limit-1) a[i]=a[i]*I%mod;
}
ll tmp1[sz],tmp2[sz],tmp3[sz],tmp4[sz];
void PolyInv(ll *a,ll *f,int n) // f=a^{-1}
{
if (n==1) return (void)(f[0]=inv(a[0]));
int mid=(n+1)>>1;
PolyInv(a,f,mid);
NTT_init(n);
rep(i,0,mid-1) tmp1[i]=f[i];rep(i,mid,limit-1) tmp1[i]=0;
rep(i,0,n-1) tmp2[i]=a[i];rep(i,n,limit-1) tmp2[i]=0;
NTT(tmp1,1);NTT(tmp2,1);
rep(i,0,limit-1) tmp1[i]=tmp1[i]*(mod+2-tmp1[i]*tmp2[i]%mod)%mod;
NTT(tmp1,-1);
rep(i,0,n-1) f[i]=tmp1[i];
rep(i,0,limit-1) tmp1[i]=tmp2[i]=0;
}
void Derivative(ll *a,ll *b,int n){rep(i,0,n-2) b[i]=a[i+1]*(i+1)%mod;b[n-1]=0;}
void Integrate(ll *a,ll *b,int n){drep(i,n-1,1) b[i]=a[i-1]*inv(i)%mod;b[0]=0;}
void PolyLn(ll *a,ll *f,int n) // f=ln a
{
NTT_init(n);
PolyInv(a,tmp3,n);Derivative(a,tmp4,n);
NTT(tmp3,1);NTT(tmp4,1);
rep(i,0,limit-1) tmp1[i]=tmp3[i]*tmp4[i]%mod;
NTT(tmp1,-1);
Integrate(tmp1,f,n);
rep(i,0,limit-1) tmp1[i]=tmp3[i]=tmp4[i]=0;
}
void PolyExp(ll *a,ll *f,int n)
{
if (n==1) return (void)(f[0]=1);
int mid=(n+1)>>1;
PolyExp(a,f,mid);
rep(i,mid,n-1) f[i]=0;
PolyLn(f,tmp2,n);
rep(i,0,n-1) tmp1[i]=f[i];
rep(i,0,n-1) tmp2[i]=(a[i]-tmp2[i]+mod)%mod;
++tmp2[0];
NTT_init(n);
NTT(tmp1,1);NTT(tmp2,1);
rep(i,0,limit-1) tmp1[i]=tmp1[i]*tmp2[i]%mod;
NTT(tmp1,-1);
rep(i,0,n-1) f[i]=tmp1[i];
rep(i,0,limit-1) tmp1[i]=tmp2[i]=0;
} int n,m;
int cnt[sz];
ll f[sz],ans[sz]; int main()
{
file();
read(n,m);
int x;
rep(i,1,n) read(x),++cnt[x];
rep(i,1,m)
if (cnt[i])
for (int j=i;j<=m;j+=i)
(f[j]+=cnt[i]*inv(j/i)%mod)%=mod;
PolyExp(f,ans,m+1);
rep(i,1,m) printf("%lld\n",ans[i]);
}
上一篇:洛谷 P4389: 付公主的背包


下一篇:啊哈算法之巧用队列解密QQ号