[BZOJ4407]于神之怒加强版(莫比乌斯反演,数论分块,线性筛)

发布于 2018-01-20  125 次阅读


本文章由SYCstudio或本站其它作者所有,请勿擅自转载

本文链接地址:[BZOJ4407]于神之怒加强版(莫比乌斯反演,数论分块,线性筛)

Description

给下N,M,K.求\(\sum_{i=1}^{n} \sum_{j=1}^{m} gcd(i,j)^k \mod (10^9+7) \)

Http

BZOJ

Tag

莫比乌斯反演,数论分块

题目大意

求\(\sum_{i=1}^{n} \sum_{j=1}^{m} gcd(i,j)^k\)

解决思路

先提取\(gcd(i,j)==d\)
\[\sum_{d=1}^{n} d^k \sum_{i=1}^{n} \sum_{j=1}^{m} [gcd(i,j)==d]\]
继续化为
\[\sum_{d=1}^{n} d^k \sum_{i=1}^{n/d} \sum_{j=1}^{m/d} [gcd(i,j)==1]\]
根据莫比乌斯反演,后面部分可以化简,所以得到
\[\sum_{d=1}^{n} d^k \sum_{i=1}^{n/d} \lfloor \frac{n}{id} \rfloor \lfloor \frac{m}{id} \rfloor \]
令\(T=id\),则有
\[\sum_{d=1}^{n} d^k \sum_{i=1}^{n/d} \lfloor \frac{n}{T} \rfloor \lfloor \frac{m}{T} \rfloor \]
把枚举\(d\)变成枚举\(T\),考虑每一个\(\mu(i)\)的贡献,可以得到
\[\sum_{T=1}^{n} \lfloor \frac{n}{T} \rfloor \lfloor \frac{n}{T} \rfloor \sum_{d|T} d^k \mu(\frac{T}{d}) \]
后面一块东西是积性函数,可以线性筛,具体来说,设\(g(x)=\sum_{d|x} d^k \mu(\frac{x}{d})\),则分情况讨论
设\(p\)为质数
\[g(p)=p^k-1\]
这个比较好理解,\(p\)是质数那么\(d\)就只有两种取值,一个得到\(p^k*\mu(1)=p^k\),另一个得到\(1^k*\mu(p)=-1\)。
\[g(p*i)=g(p)*g(i) \quad [gcd(i,p)==1]\]
这个直接由积性函数性质得到
\[g(p*i)=g(i)*p^k \quad [gcd(i,p)==p]\]
这个的话,考虑加进去一个原来已经存在的质因子\(p\),那么只会影响\(d^k\)变成\((d*p)^k\),所以补一个\(p^k\)进去就可以了
那么求出后面那一块之后,前缀和数论分块即可。

代码

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
using namespace std;

#define ll long long
#define mem(Arr,x) memset(Arr,x,sizeof(Arr))

const int maxNum=5000010;
const ll Mod=1e9+7;
const int inf=2147483647;

int pricnt=0,Prime[maxNum];
bool notprime[maxNum];
ll H[maxNum];
ll Pow[maxNum];

int K;

void Init();
int QPow(ll x,int cnt);

int main()
{
    int T;
    scanf("%d%d",&T,&K);
    Init();//线性筛初始化
    while (T--)
    {
        ll n,m;scanf("%lld%lld",&n,&m);if (n>m) swap(n,m);
        ll ans=0;
        for (ll i=1,last;i<=n;i=last+1)
        {//数论分块
            last=min(n/(n/i),m/(m/i));
            ans=(ans+(H[last]-H[i-1]+Mod)%Mod*(ll)((n/i))%Mod*(ll)((m/i))%Mod)%Mod;
        }
        printf("%lld\n",ans%Mod);
    }
    return 0;
}

void Init()
{
    notprime[1]=1;H[1]=1;
    for (int i=2;i<maxNum;i++)
    {
        if (notprime[i]==0) Prime[++pricnt]=i,Pow[pricnt]=QPow(i,K),H[i]=Pow[pricnt]-1;
        for (int j=1;(j<=pricnt)&&((ll)i*(ll)Prime[j]<maxNum);j++)
        {
            notprime[i*Prime[j]]=1;
            if (i%Prime[j]==0)
            {
                H[i*Prime[j]]=H[i]*Pow[j]%Mod;
                break;
            }
            H[i*Prime[j]]=H[i]*H[Prime[j]]%Mod;
        }
    }
    for (int i=1;i<maxNum;i++) H[i]=(H[i]+H[i-1])%Mod;
    return;
}

int QPow(ll x,int cnt)
{
    ll ret=1;
    while (cnt)
    {
        if (cnt&1) ret=ret*x%Mod;
        x=x*x%Mod;
        cnt=cnt>>1;
    }
    return ret;
}

本文章由SYCstudio或本站其它作者所有,请勿擅自转载

本文链接地址:[BZOJ4407]于神之怒加强版(莫比乌斯反演,数论分块,线性筛)


HNCJ OIer