题面

DP学得好好的为什么突然想起写多项式呢？

这要从几天前说起

一开始是看到一道FWT的DP([SNOI2017]遗失的答案)

然后在洛谷上搜

没有搜到

但是搜到了这个题

仔细一看

MTT+莫比乌斯反演!

好像很有意思的样子

顺便写一写\(4\)次DFT的MTT

嘿嘿嘿……

于是LYC沉迷多项式无法自拔

说这个题

这个题是真的毒瘤

明明就有唯一解

偏要说什么输出字典序最小的

现在考虑计算方案数

一般来说会想到DP

然而总不可能DP回去吧

这时LYC耳边响起了语文老师的谆谆教诲

生命不能让阅(数)读(学)缺席

阅(数)读(学)让我们找到回家的路

所以只要找到方案数与集合的函数关系

反演将“让我们找到回家的路”

现在我们的问题是

已知集合求方案数的表达式

似乎除了生成函数也没什么好用的

设\(a_i\)表示\(i\)是否在集合中

设序列\(a\)的生成函数为\(A\)

方案数序列的生成函数为\(F\)
\[ F(x)=\prod_{i=1}^{n}{(\frac{1}{1-x^i})^{a_i}} \\ \ln(F(x))=\sum_{i=1}^{n}{(a_i*\ln(\frac{1}{1-x^i}))} \\-\ln(F(x))=\sum_{i=1}^{n}{(a_i*\ln({1-x^i}))} \]
考虑计算函数\(G_n(x)=\ln(1-x^n)\)

根据链式法则
\[ G_n'(x)=\frac{-nx^{n-1}}{1-x^n} \\\because \frac{1}{1-x^n}=\sum_{i=1}^{\infty}{x^{i*n}} \\\therefore G_n'(x)=\sum_{i=0}^{\infty}{-nx^{i*n+(n-1)}} \\\therefore G_n(x)=\sum_{i=0}^{\infty}\frac{-nx^{i*n+(n-1)+1}}{i*n+(n-1)+1} \\=-\sum_{i=0}^{\infty}{\frac{nx^{i*n+n}}{i*n+n}} \\=-\sum_{i=0}^{\infty}{\frac{x^{(i+1)*n}}{i+1}} \\=-\sum_{i=1}^{\infty}{\frac{x^{i*n}}{i}} \]
带入原式
\[ \\-\ln(F(x))=-\sum_{i=1}^{n}{(a_i*\sum_{j=0}^{\infty}{\frac{x^{ij}}{j}})} \\\ln(F(x))=\sum_{t=1}^{\infty}{x^t*\sum_{i|t}{a_i*\frac{i}{t}}} \]
令\(\ln(F(x))\)对应的数列为\(f\)
\[ n*f_n=\sum_{d|n}{a_d*d} \]
莫比乌斯反演即可求出\(a_i\)

注意只保证\(p\)为质数

要写MTT

话说洛谷上有个讨论问什么输出\(a_i\)是对的(本来要输出\(i\))

因为程序里反演求出的序列是\(\{a_i*i\}\)

而本来的\(a_i\)是零一串

当然求出的\(a_i\)要么是\(0\)要么是\(i\)啦

#include<bits/stdc++.h>using namespace std;#define gc c=getchar()
#define ll long long
#define db double template<typename T>
inline void read(T&x){x=0;T k=1;char gc;while(!isdigit(c)){if(c=='-')k=-1;gc;}while(isdigit(c)){x=x*10+c-'0';gc;}x*=k;
}struct comp{db r,i;comp(){}comp(const db &_r,const db &_i):r(_r),i(_i){}
};inline comp operator + (const comp &a,const comp &b){return comp(a.r+b.r,a.i+b.i);
}inline comp operator - (const comp &a,const comp &b){return comp(a.r-b.r,a.i-b.i);
}inline comp operator * (const comp &a,const comp &b){return comp(a.r*b.r-a.i*b.i,a.r*b.i+a.i*b.r);
}inline comp operator / (const comp &a,const db &b){return comp(a.r/b,a.i/b);
}inline comp conj(const comp &a){return comp(a.r,-a.i);
}const db PI=acos(-1.0);
const int N=1<<21|7;int rev[N];
comp Wn[N];
int Now_Len;
inline void revinit(int n){for(int i=0;i<n;++i)rev[i]=(rev[i>>1]>>1)|((i&1)*(n>>1));for(int i=0;i<n;++i)Wn[i]=comp(cos(PI*i/n),sin(PI*i/n));Now_Len=n;
}inline void DFT(comp *A,int n){if(Now_Len!=n)revinit(n);for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);for(int i=1,l=0;i<n;i<<=1,l++){for(int j=0;j<n;j+=i<<1){for(int k=0;k<i;++k){comp u=A[j+k],v=A[i+j+k]*Wn[(n>>l)*k];A[j+k]=u+v;A[i+j+k]=u-v;}}}
}inline void IDFT(comp *A,int n){if(Now_Len!=n)revinit(n);for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);for(int i=1,l=0;i<n;i<<=1,l++){for(int j=0;j<n;j+=i<<1){for(int k=0;k<i;++k){comp u=A[j+k],v=A[i+j+k]*conj(Wn[(n>>l)*k]);A[j+k]=u+v;A[i+j+k]=u-v;}}}for(int i=0;i<n;++i){A[i]=A[i]/n;}
}#define DFT(A) DFT(A,n)
#define IDFT(A) IDFT(A,n)int p;
comp a[N],b[N],c[N],d[N];
inline void MTT(int *A,int *B,int *C,int lenA,int lenB){int lenC=lenA+lenB-1,n;for(n=1;n<lenC;n<<=1);for(int i=0;i<lenA;++i)a[i]=comp(A[i]&32767,A[i]>>15);for(int i=0;i<lenB;++i)b[i]=comp(B[i]&32767,B[i]>>15);for(int i=lenA;i<n;++i)a[i].r=a[i].i=0;for(int i=lenB;i<n;++i)b[i].r=b[i].i=0;DFT(a);DFT(b);for(int i=0;i<n;++i){int j=(n-i)&(n-1);comp A0=(a[i]+conj(a[j]))*comp(0.5,0);comp A1=(a[i]-conj(a[j]))*comp(0,-0.5);comp B0=(b[i]+conj(b[j]))*comp(0.5,0);comp B1=(b[i]-conj(b[j]))*comp(0,-0.5);c[i]=A0*B0+A0*B1*comp(0,1);d[i]=A1*B0+A1*B1*comp(0,1);}IDFT(c);IDFT(d);for(int i=0;i<lenC;++i){ll s1=(ll)(c[i].r+0.5)%p;ll s2=(ll)(c[i].i+0.5)%p;ll s3=(ll)(d[i].r+0.5)%p;ll s4=(ll)(d[i].i+0.5)%p;C[i]=(s1+((s2+s3)<<15)+(s4<<30))%p;}
}inline ll qpow(ll a,ll b){ll ans=1;while(b){if(b&1)(ans*=a)%=p;(a*=a)%=p;b>>=1;}return ans;
}int Tmp_Inv[N];
inline void inverse(int *A,int *Inv,int lenA){Inv[0]=qpow(A[0],p-2);for(int i=2;i<=lenA;i<<=1){MTT(A,Inv,Tmp_Inv,i,i);MTT(Tmp_Inv,Inv,Tmp_Inv,i,i);for(int j=0;j<i;++j)Inv[j]=(2ll*Inv[j]-Tmp_Inv[j]+p)%p;}
}inline void differentiate(int *A,int *Dif,int len){for(int i=1;i<len;++i)Dif[i-1]=(ll)A[i]*i%p;Dif[len-1]=0;
}inline void integrate(int *A,int *Int,int len){for(int i=len-1;i>=1;--i)Int[i]=(ll)A[i-1]*qpow(i,p-2)%p;Int[0]=0;
}int Inv_ln[N];
inline void ln(int *A,int *Ln,int n){inverse(A,Inv_ln,n);differentiate(A,Ln,n); MTT(Ln,Inv_ln,Ln,n,n);integrate(Ln,Ln,n<<1);
}int tot;
int mu[N],pri[N];
bool mark[N];inline void init(int n){mu[1]=1;for(int i=2;i<=n;++i){if(!mark[i])pri[++tot]=i,mu[i]=-1;for(int j=1,tmp;j<=tot&&(tmp=i*pri[j])<=n;++j){mark[tmp]=1;if(i%pri[j]==0){mu[tmp]=0;break;}mu[tmp]=-mu[i];}}
}int F[N],G[N];
int main(){int n,m;read(n);read(p);for(int i=1;i<=n;++i)read(F[i]);F[0]=1;for(m=1;m<n;m<<=1);ln(F,F,m);for(int i=0;i<=n;++i)F[i]=(ll)F[i]*i%p;init(n);for(int i=1;i<=n;++i){for(int j=1;i*j<=n;++j){G[i*j]+=F[i]*mu[j];}}vector<int>ans;for(int i=1;i<=n;++i)if(G[i])ans.push_back(i);printf("%d\n",ans.size());for(int i=0;i<ans.size();++i){if(i)putchar(' ');printf("%d",ans[i]);}
}