• BZOJ4305 数列的GCD


    也许更好的阅读体验

    (mathcal{Description})


    (mathcal{Solution})

    这里就不用(N,M),还是(n,m)写的习惯些
    直接计算一个方案是十分不方便的
    所以考虑容斥
    (gleft(d ight))表示(d | gcd)的方案数
    (a)中有(cnt)个数是(d)的倍数
    那么有
    (gleft( d ight) =left( dfrac {m}{d} ight) ^{n-cnt}egin{pmatrix} cnt \ n-k end{pmatrix}left( dfrac {m}{d}-1 ight) ^{k-left(n-cnt ight)})

    (left( dfrac {m}{d} ight) ^{n-cnt})表示有(n-cnt)个数是必须修改的,每个有(frac{m}{d})种数选择

    那么还剩(k-left(n-cnt ight))个数必须要修改,我们可将其写为(cnt-(n-k)),这样就是等价于要选(n-k)个数出来

    (egin{pmatrix} cnt \ n-k end{pmatrix}left( dfrac {m}{d}-1 ight) ^{k-left(n-cnt ight)})表示将这(cnt)个数修改(k-(n-cnt))个数,每个数因为自己本身是(d)的一个倍数,所以只有(frac{m}{d})种选择

    (fleft(d ight))表示(gcd=d)的方案数

    然后可以考虑莫比乌斯反演
    显然有
    (egin{aligned}gleft( n ight) =sum _{n | d}fleft( d ight)end{aligned})
    则根据莫比乌斯反演,有
    (egin{aligned}fleft( n ight) =sum _{n | d}mu left( dfrac {d}{n} ight) gleft( d ight)end{aligned})

    当然,莫比乌斯反演什么的是不可能莫比乌斯反演的
    直接容斥就可以啦
    (egin{aligned}fleft(n ight)=gleft(n ight)-sum_{n|d}f(d)end{aligned})
    从大到小枚举(d),直接计算即可

    (mathcal{Code})

    /*******************************
    Author:Morning_Glory
    LANG:C++
    Created Time:2019年08月23日 星期五 08时14分25秒
    *******************************/
    #include <cstdio>
    #include <fstream>
    #define ll long long
    using namespace std;
    const int maxn = 300005;
    const int mod = 1000000007;
    //{{{cin
    struct IO{
    	template<typename T>
    	IO & operator>>(T&res){
    		res=0;
    		bool flag=false;
    		char ch;
    		while((ch=getchar())>'9'||ch<'0')	 flag|=ch=='-';
    		while(ch>='0'&&ch<='9') res=(res<<1)+(res<<3)+(ch^'0'),ch=getchar();
    		if (flag)	 res=~res+1;
    		return *this;
    	}
    }cin;
    //}}}
    int n,m,k;
    int a[maxn],num[maxn];
    ll fac[maxn],inv[maxn],f[maxn];
    //{{{ksm
    int ksm (int a,int b)
    {
    	a%=mod;
    	int s=1;
    	for (;b;b>>=1,a=1ll*a*a%mod)
    		if (b&1)	s=1ll*s*a%mod;
    	return s;
    }
    //}}}
    inline ll C (int n,int m)
    {
    	return fac[n]*inv[m]%mod*inv[n-m]%mod;
    }
    int main()
    {
    	cin>>n>>m>>k;
    	fac[0]=fac[1]=inv[0]=inv[1]=1;
    	for (int i=1;i<=n;++i)	cin>>a[i],++num[a[i]];
    	for (int i=2;i<=n;++i)	fac[i]=1ll*fac[i-1]*i%mod;
    	for (int i=2;i<=n;++i)	inv[i]=(-mod/i*inv[mod%i]%mod+mod)%mod;
    	for (int i=2;i<=n;++i)	inv[i]=1ll*inv[i]*inv[i-1]%mod;
    	for (int i=m;i>=1;--i){
    		int cnt=0;
    		for (int j=i;j<=m;j+=i)	cnt+=num[j];
    		if (cnt-n+k<0)	f[i]=0;
    		else	f[i]=C(cnt,n-k)*ksm(m/i-1,cnt-n+k)%mod*ksm(m/i,n-cnt)%mod;
    		for (int j=i<<1;j<=m;j+=i)	f[i]=1ll*(f[i]-f[j]+mod)%mod;
    	}
    	for (int i=1;i<=m;++i)	printf("%lld ",f[i]);
    	return 0;
    }
    
    

    如有哪里讲得不是很明白或是有错误,欢迎指正
    如您喜欢的话不妨点个赞收藏一下吧

  • 相关阅读:
    从MSFT Project中同步数据到PSA
    TracingService. 码农debug的救星
    Business Workflow Flow 最后一个stage阶段触发Exist
    BPF form 中的readonly 只读字段自动unblock
    Windows Server 使用fiddler中抓取IIS的请求
    Dynamics 365 CE 的快捷键 Shortcut
    Dynamics 365 online 服务器保护机制 server limitation
    在后端C#中 call web api 关联lookup 和 GUID
    PSA 需要 Sales Order中 读取的最低权限 read minimal privilege
    python学习笔记10:分析程序性能cProfile
  • 原文地址:https://www.cnblogs.com/Morning-Glory/p/11398479.html
Copyright © 2020-2023  润新知