BZOJ 3561 DZY Loves Math VI（莫比乌斯反演）

最新推荐文章于 2020-01-19 19:07:51 发布

v5zsq

最新推荐文章于 2020-01-19 19:07:51 发布

阅读量373

点赞数

分类专栏： BZOJ 莫比乌斯反演

本文链接：https://blog.csdn.net/V5ZSQ/article/details/78874149

版权

BZOJ 同时被 2 个专栏收录

50 篇文章 0 订阅

订阅专栏

莫比乌斯反演

28 篇文章 0 订阅

订阅专栏

Description

给定正整数 $n,m$ ，求 $\sum\limits_{i=1}^n\sum\limits_{j=1}^m lcm(i,j)^{gcd(i,j)}$

Input

一行两个整数 $n,m$ ， $1\le n,m\le 5e5$

Output

一个整数，为答案模 $10^9+7$ 后的值

Sample Input

5 4

Sample Output

424

Solution

不妨假设 $n\le m$

\sum i = 1 n \sum j = 1 m l c m (i, j) g c d (i, j) = = = = \sum i = 1 n \sum j = 1 m (i g c d ( i , j ) j g c d ( i , j ) g c d (i, j)) g c d (i, j) \sum d = 1 n d d \sum i = 1 n d \sum j = 1 m d i d j d [(i, j) = 1] \sum d = 1 n d d \sum i = 1 n d \sum j = 1 m d i d j d \sum k = 1 n d μ (k) \sum d = 1 n d d \sum k = 1 n d μ (k) k 2 d \sum i = 1 n d k i d \sum j = 1 m d k j d

$\begin{array}{rcl} \sum\limits_{i=1}^n\sum\limits_{j=1}^m lcm(i,j)^{gcd(i,j)}&=& \sum\limits_{i=1}^n\sum\limits_{j=1}^m(\frac{i}{gcd(i,j)}\frac{j}{gcd(i,j)}gcd(i,j))^{gcd(i,j)}\\ &=&\sum\limits_{d=1}^nd^d\sum\limits_{i=1}^{\frac{n}{d}}\sum\limits_{j=1}^{\frac{m}{d}}i^dj^d[(i,j)=1]\\ &=&\sum\limits_{d=1}^nd^d\sum\limits_{i=1}^{\frac{n}{d}}\sum\limits_{j=1}^{\frac{m}{d}}i^dj^d\sum\limits_{k=1}^{\frac{n}{d}}\mu(k)\\ &=&\sum\limits_{d=1}^nd^d\sum\limits_{k=1}^{\frac{n}{d}}\mu(k)k^{2d}\sum\limits_{i=1}^{\frac{n}{dk}}i^d\sum\limits_{j=1}^{\frac{m}{dk}}j^d \end{array}$
枚举

d $d$ ，预处理

Sum(i)=∑i=1ndid $Sum(i)=\sum\limits_{i=1}^{\frac{n}{d}}i^d$ ，总时间复杂度

O(nlogn) $O(nlogn)$

Code

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<vector>
#include<queue>
#include<map>
#include<set>
#include<ctime>
using namespace std;
typedef long long ll;
typedef pair<int,int>P;
const int INF=0x3f3f3f3f,maxn=500001;
#define mod 1000000007
int mu[maxn],p[maxn],mark[maxn],res=0;
void init(int n=5e5)
{
    mu[1]=1;
    for(int i=2;i<=n;i++)
    {
        if(!mark[i])p[res++]=i,mu[i]=-1;
        for(int j=0;j<res&&i*p[j]<=n;j++)
        {
            mark[i*p[j]]=1;
            if(i%p[j])mu[i*p[j]]=-mu[i];
            else
            {
                mu[i*p[j]]=0;
                break;
            }
        }
    }
}
int Pow(int a,int b)
{
    int ans=1;
    while(b)
    {
        if(b&1)ans=(ll)ans*a%mod;
        a=(ll)a*a%mod;
        b>>=1;
    }
    return ans;
}
void add(int &x,int y)
{
    x=x+y>=mod?x+y-mod:x+y;
}
int n,m,a[maxn],s[maxn];
int main()
{
    init();
    scanf("%d%d",&n,&m);
    if(n>m)swap(n,m);
    int ans=0;
    for(int i=1;i<=m;i++)a[i]=1;
    for(int d=1;d<=n;d++)
    {
        for(int i=1;i*d<=m;i++)a[i]=(ll)a[i]*i%mod,s[i]=(s[i-1]+a[i])%mod;
        int temp=0;
        for(int i=1;i*d<=n;i++)
        {
            if(mu[i]==1)add(temp,(ll)Pow(i,2*d)*s[n/d/i]%mod*s[m/d/i]%mod);
            else if(mu[i]==-1)add(temp,mod-(ll)Pow(i,2*d)*s[n/d/i]%mod*s[m/d/i]%mod);
        }
        add(ans,(ll)Pow(d,d)*temp%mod); 
    }
    printf("%d\n",ans);
    return 0;
}