题解：P1998 DGF 等比求和

masterhuang · 2025-03-26 11:14:12 · 题解

前置知识：狄利克雷 k 次方

首先运用等比数列求和公式：g=H=\dfrac{f^{m+1}-f}{f-1}，其中乘法是狄利克雷卷积。

令 F=f^{m+1}-f,G=f-1，则 H=\dfrac{F}{G}，其中除法是狄利克雷求逆。

但是注意到 G(1)=0，于是无法直接求逆。

但是注意到 G(2)\neq 0，于是我们有一些手段：

我们只需把 f^{m+1} 保留 2n 项，然后枚举因数通过上述过程求得 H 即可。

std 的实现是枚举倍数，减小了常数。

code：

#include<bits/stdc++.h>
#define LL long long
#define fr(x) freopen(#x".in","r",stdin);freopen(#x".out","w",stdout);
using namespace std;
namespace IO
{
    const int _Pu=2e7+5,_d=32;
    char buf[_Pu],*p1=buf+_Pu,*p2=buf+_Pu;
    inline void fin()
    {
        memmove(buf,p1,p2-p1);
        int rlen=fread(buf+(p2-p1),1,p1-buf,stdin);
        if(p1-rlen>buf) buf[p2-p1+rlen]=EOF;p1=buf;
    }
    inline int rd()
    {
        if(p1+_d>p2) fin();int x=0;
        for(;!isdigit(*p1);++p1);x=(*p1++-'0');
        for(;isdigit(*p1);++p1) x=x*10+(*p1-'0');
        return x;
    }
}using IO::rd;//快读
const int N=2e6+5,mod=1e9+7;
int n,m,pr[N],w[N],a[N],f[N],g[N],h[N],I[35],ans;bool v[N];
inline int md(int x){return x>=mod?x-mod:x;}
inline int ksm(int x,int p){int s=1;for(;p;(p&1)&&(s=1ll*s*x%mod),x=1ll*x*x%mod,p>>=1);return s;}
inline void init(int M)
{
    for(int i=2;i<=M;i++)
    {
        if(!v[i]) pr[++pr[0]]=i,w[i]=1;
        for(int j=1;j<=pr[0]&&i*pr[j]<=M;j++)
        {
            v[i*pr[j]]=1;w[i*pr[j]]=w[i]+1;
            if(i%pr[j]==0) break;
        }
    }I[1]=1;
    for(int i=2;i<=30;i++) I[i]=mod-1ll*I[mod%i]*(mod/i)%mod;
}//预处理素因数个数与逆元
inline void pw(int *f,int *g,int n,int k){
    static int F[N];
    for(int i=1;i<=n;i++) F[i]=1ll*f[i]*w[i]%mod*k%mod,g[i]=0;
    for(int i=1;i<=n;i++)
    {
        int d=g[i];
        g[i]=(1ll*g[i]*I[w[i]]+(i==1))%mod;
        for(int k=(i<<1),j=2;k<=n;k+=i,j++)
            g[k]=(g[k]+1ll*g[i]*F[j]+1ll*(mod-f[j])*d)%mod;
    }
}//狄利克雷幂
int main()
{
    n=rd(),m=rd();init(n<<1);
    for(int i=1;i<=n;i++) a[i]=g[i]=rd();
    g[1]=0;pw(a,f,n<<1,m+1);
    for(int i=1;i<=n;i++) f[i]=md(f[i]+mod-a[i]);
    int _I=ksm(g[2],mod-2);
    for(int i=1;i<=n;i++)
    {
        h[i]=1ll*_I*(f[i<<1]+mod-h[i])%mod;ans^=h[i];
        for(int j=3,k=i*3;k<=2*n;j++,k+=i) if(k&1^1) h[k>>1]=(h[k>>1]+1ll*h[i]*g[j])%mod;
    }//枚举倍数贡献
    return cout<<ans,0;
}