SDOI2017遗忘的集合

题面链接

咕咕咕

题外话

为了这道题我敲了\(MTT\)、多项式求逆、多项式\(ln\)等模板,搞了将近一天。

sol

最近懒得写题解啊,随便搞搞吧。

看到这个就是生成函数套上去。

\[F(x)=\prod_{i=1}^{n}(\frac{1}{1-x^i})^{a_i}\]
\[-\ln F(x)=\sum_{i=1}^na_i\ln(1-x^i)\]
\[-\ln F(x)=-\sum_{i=1}^na_i\sum_{j=1}^{\infty}\frac{x^{ij}}{j}\]
常见莫比乌斯套路,令\(T=ij\),交换枚举顺序。
\[\ln F(x)=\sum_{T=1}^{\infty}x^T\sum_{i|T}a_i\times \frac iT\]
\(G(x)=lnF(x)\),于是我们得到了一个奇怪的生成函数,\(G(x)\)的每一项就是\(\sum_{i|T}a_i\times\frac iT\)
可以莫比乌斯反演,也可以调和级数暴力减。

感谢Cyhlnj巨sy

关于\[\ln(1-x^i)=-\sum_{j=1}^{\infty}\frac{x^{ij}}{j}\]洛谷上有人说是泰勒展开?害我证了好久没证出来
也许是我太弱了,但我不敢苟同,反正我没看出来泰勒展开怎么搞,如果有人证出来了欢迎留言。
\(zsy\)的博真舒服。
\[\ln F(x)=G(x)\\\frac{F'(x)}{F(x)}=G'(x)\\\frac{-ix^{i-1}}{1-x^i}=G'(x)\\-\sum_{j=0}^{\infty} ix^{i-1+ij}=G'(x)\\-\sum_{j=0}^{\infty}\frac{ix^{i+ij}}{i+ij}=G(x)\\-\sum_{j=1}^{\infty}\frac{x^{ij}}{j}=G(x)\]

#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define gt getchar()
#define ll long long
#define File(s) freopen(s".in","r",stdin),freopen(s".out","w",stdout)
inline int in()
{
    int k=0;char ch=gt;
    while(ch<'-')ch=gt;
    while(ch>'-')k=k*10+ch-'0',ch=gt;
    return k;
}

const int N=6e5+5;
const double PI=acos(-1);
int a[N],b[N],c[N],YL;
inline int MO(const int &a){return a>=YL?a-YL:a;}
inline int ksm(int a,int k){int r=1;while(k){if(k&1)r=1ll*r*a%YL;a=1ll*a*a%YL,k>>=1;}return r;}

int rev[N],Inv[N];
struct E
{
    double x,y;
    E(){}
    E(double a,double b):x(a),y(b){}
    E operator=(const int &a){x=a,y=0;return *this;}
    E conj(){return E(x,-y);}
}omg[N];
E operator+(const E &a,const E &b){return E(a.x+b.x,a.y+b.y);}
E operator-(const E &a,const E &b){return E(a.x-b.x,a.y-b.y);}
E operator*(const E &a,const E &b){return E(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}

inline void fft(E *a,int len)
{
    for(int i=0;i<len;++i)if(i<rev[i])std::swap(a[i],a[rev[i]]);
    for(int st=2,m=1;st<=len;st<<=1,m<<=1)
        for(E *p=a,x,y;p!=a+len;p+=st)
            for(int k=0;k<m;++k)
                x=p[k],y=omg[len/m*k]*p[k+m],p[k]=x+y,p[k+m]=x-y;
}

inline void mul(int *A,int *B,int *C,int len)
{
    int len2=len;len<<=1;int qwq=0;while((1<<qwq)<len)++qwq;--qwq;
    for(int i=0;i<len;++i)rev[i]=(rev[i>>1]>>1)|((i&1)<<qwq);
    for(int i=0;i<len;++i)omg[i]=E(cos(PI*i/len),sin(PI*i/len));
    for(int i=0;i<len;++i)A[i]=MO(A[i]+YL),B[i]=MO(B[i]+YL);
    static E mul_a[N],mul_b[N],dfta[N],dftb[N],dftc[N],dftd[N];
    for(int i=0;i<len2;++i)mul_a[i]=E(A[i]&32767,A[i]>>15);
    for(int i=0;i<len2;++i)mul_b[i]=E(B[i]&32767,B[i]>>15);
    for(int i=len2;i<len;++i)mul_a[i]=mul_b[i]=E(0,0);
    fft(mul_a,len),fft(mul_b,len);
    for(int i=0;i<len;++i)
    {
        int j=(len-i)&(len-1);E da,db,dc,dd;
        da=(mul_a[i]+mul_a[j].conj())*E( 0.5,0);
        db=(mul_a[i]-mul_a[j].conj())*E(0,-0.5);
        dc=(mul_b[i]+mul_b[j].conj())*E( 0.5,0);
        dd=(mul_b[i]-mul_b[j].conj())*E(0,-0.5);
        dfta[j]=da*dc,dftb[j]=da*dd,dftc[j]=db*dc,dftd[j]=db*dd;
    }
    for(int i=0;i<len;++i)mul_a[i]=dfta[i]+dftb[i]*E(0,1);
    for(int i=0;i<len;++i)mul_b[i]=dftc[i]+dftd[i]*E(0,1);
    fft(mul_a,len),fft(mul_b,len);
    for(int i=0;i<len;++i)
    {
        int da=(ll)(mul_a[i].x/len+0.5)%YL;
        int db=(ll)(mul_a[i].y/len+0.5)%YL;
        int dc=(ll)(mul_b[i].x/len+0.5)%YL;
        int dd=(ll)(mul_b[i].y/len+0.5)%YL;
        C[i]=(da+((ll)(db+dc)<<15)+((ll)dd<<30))%YL;
    }
}

void get_inv(int *a,int *b,int len)
{
    if(len==1){b[0]=ksm(a[0],YL-2);return;}
    static int tmp[N];get_inv(a,b,len>>1);mul(a,b,tmp,len);
    for(int i=0;i<len;++i)tmp[i]=YL-tmp[i];tmp[0]+=2;mul(tmp,b,b,len);
}

void get_ln(int *a,int *b,int len)
{
    static int d[N],inv[N];
    for(int i=1;i<len;++i)d[i-1]=1ll*i*a[i]%YL;
    get_inv(a,inv,len);mul(d,inv,b,len);
    for(int i=len-1;i;--i)b[i]=1ll*b[i-1]*ksm(i,YL-2)%YL;b[0]=0;
}

int main()
{
    int n=in();YL=in();int len=1,ans=0;while(len<=n)len<<=1;
    for(int i=1;i<=n;++i)a[i]=in();a[0]=1;get_ln(a,a,len);
    for(int i=1;i<=n;++i)a[i]=1ll*a[i]*i%YL;
    for(int i=1;i<=n;++i)
        for(int j=i<<1;j<=n;j+=i)
            a[j]=MO(a[j]-a[i]+YL);
    for(int i=1;i<=n;++i)if(a[i])++ans;printf("%d\n",ans);
    for(int i=1;i<=n;++i)if(a[i])printf("%d ",i);puts("");
    return 0;
}

转载于:https://www.cnblogs.com/cx233666/p/10030864.html

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值