BZOJ 4407: 于神之怒加强版

简单反演题,都是套路,直接写式子了……

[sum_{i=1}^nsum_{j=1}^m gcd(i,j)^k\=sum_{g=1}^{min(n,m)} g^ksum_{i=1}^{lfloorfrac{n}{g} floor}sum_{j=1}^{lfloorfrac{m}{g} floor} [gcd(i,j)=1]\=sum_{g=1}^{min(n,m)} g^ksum_{i=1}^{lfloorfrac{n}{g} floor}sum_{j=1}^{lfloorfrac{m}{g} floor} sum_{d|gcd(i,j)} mu(d)\=sum_{g=1}^{min(n,m)} g^ksum_{d=1}^{min(lfloorfrac{n}{g} floor,lfloorfrac{m}{g} floor)} mu(d) imeslfloorfrac{n}{dg} floorlfloorfrac{m}{dg} floor)\=sum_{D=1}^{min(n,m)} lfloorfrac{n}{D} floorlfloorfrac{m}{D} floor sum_{g|D} mu(frac{D}{g}) imes g^k\ ]

推导这里我们发现只要设(f(D)=sum_{g|D} mu(frac{D}{g}) imes g^k),预处理出(f)就可以(O(sqrt n+sqrt m))除法分块得到答案了

那么(f)怎么搞呢,首先这是个狄利克雷卷积的形式,同时里面的(mu,g^k)都是积性函数,因此(f)也是积性函数

既然是积性函数就考虑当(D=p_i^{a_i})时的情况,显然只有两种取值方法(注意(mu)有平方质因子时值为(0)),我们容易算出:

[f(n)=egin{cases}1&n=1\n^k-1&nin Primeend{cases} ]

然后容易推出线性筛的时候:

[f(i imes P_j)=egin{cases}f(i) imes f(P_j)&P_j ot| i\f(i) imes P_j^k&P_j|iend{cases} ]

然后就做完了,预处理是(O(n))

#include<cstdio>
#include<iostream>
#define RI register int
#define CI const int&
using namespace std;
const int N=5000000,mod=1e9+7;
int n,m,k,t,prime[N+5],cnt,f[N+5],ans; bool vis[N+5];
inline int quick_pow(int x,int p,int mul=1)
{
    for (;p;p>>=1,x=1LL*x*x%mod) if (p&1) mul=1LL*mul*x%mod; return mul;
}
inline int sum(CI x,CI y)
{
    int t=x+y; return t>=mod?t-mod:t;
}
inline int sub(CI x,CI y)
{
    int t=x-y; return t<0?t+mod:t;
}
#define P(x) prime[x]
inline void init(CI n)
{
    RI i,j; for (vis[1]=f[1]=1,i=2;i<=n;++i)
    {
        if (!vis[i]) prime[++cnt]=i,f[i]=sub(quick_pow(i,k),1);
        int tp; for (j=1;j<=cnt&&(tp=i*P(j))<=n;++j)
        {
            vis[tp]=1; if (i%P(j)) f[tp]=1LL*f[i]*f[P(j)]%mod;
            else f[tp]=1LL*f[i]*quick_pow(P(j),k)%mod;
        }
    }
    for (i=2;i<=n;++i) f[i]=sum(f[i],f[i-1]);
}
#undef P
int main()
{
    for (scanf("%d%d",&t,&k),init(N);t;--t)
    {
        scanf("%d%d",&n,&m); ans=0; if (n>m) swap(n,m);
        for (RI l=1,r;l<=n;l=r+1) r=min(n/(n/l),m/(m/l)),
        ans=sum(ans,1LL*(n/l)*(m/l)%mod*sub(f[r],f[l-1])%mod);
        printf("%d
",ans);
    }
    return 0;
}
原文地址:https://www.cnblogs.com/cjjsb/p/12261749.html