首页 > 其他 > 详细

BZOJ4407: 于神之怒加强版

时间:2019-02-22 21:33:17      阅读:179      评论:0      收藏:0      [点我收藏+]

# BZOJ4407: 于神之怒加强版

题目描述

传送门

题目分析

题目让求:
\[ Ans=\sum_{i=1}^{n}\sum_{j=1}^{m}gcd(i,j)^k \]

可以发现是一个比较正常的式子,我们直接开始化:

\[ \begin{aligned} Ans&=\sum_{i=1}^{n}\sum_{j=1}^{m}gcd(i,j)^k\&=\sum_{d=1}^{min(n,m)}\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)=d]d^k\&=\sum_{d=1}^{min(n,m)}d^k\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)=d]\&=\sum_{d=1}^{min(n,m)}d^k\sum_{d\mid T}\mu(\frac{T}{d})\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor\&=\sum_{T=1}^{min(n,m)}\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor\sum_{d\mid T}d^k\mu(\frac{T}{d}) \end{aligned} \]
这个式子明显前面可以整除分块。

可以发现后面的是个积性函数,又由于所有数据的\(k\)都是一样的,所以可以线筛出后面的部分。然后直接求解就可以了。

#include <bits/stdc++.h>
using namespace std;
const int MAXN=5e6+7;
const int mo=1e9+7;
#define ll long long
int prime[MAXN],n,m;
bool vis[MAXN];
int phi[MAXN],f[MAXN],g[MAXN],ans,k;
inline int power(int x,int k)
{
    int cnt=1;
    while(k){
        if(k&1) cnt=1ll*cnt*x%mo;
        x=(1ll*x*x)%mo;k>>=1;
    }
    return cnt%mo;
}
inline void get(int N)
{
    f[1]=1;
    for(int i=2;i<=N;i++){
        if(!vis[i]){
            prime[++prime[0]]=i;
            g[prime[0]]=power(i,k)%mo;
            f[i]=1ll*(1ll*g[prime[0]]-1+mo)%mo;
        }
        for(int j=1;j<=prime[0];j++){
            if(i*prime[j]>N) break;
            vis[i*prime[j]]=1;
            if(i%prime[j]==0){
                f[i*prime[j]]=1ll*f[i]*g[j]%mo;
                break;
            } else f[i*prime[j]]=1ll*f[i]*f[prime[j]]%mo;
        }
    }
    for(int i=1;i<=N;i++) f[i]=1ll*(1ll*f[i]+f[i-1])%mo;
}
inline int read()
{
    int x=0,c=1;
    char ch=' ';
    while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
    while(ch=='-')c*=-1,ch=getchar();
    while(ch>='0'&&ch<='9')x=x*10+ch-'0',ch=getchar();
    return x*c;
}
int main()
{
    int T=read();k=read();get(5000000);
    while(T--){
        n=read(),m=read();
        int mx=min(n,m);ans=0;
        for(int l=1,r;l<=mx;l=r+1){
            r=min(n/(n/l),m/(m/l));
            ans=(ans+1ll*(n/l)*(m/l)%mo*1ll*(1ll*f[r]-f[l-1]+mo)%mo)%mo;
        }
        printf("%d\n", ans);
    }
}

BZOJ4407: 于神之怒加强版

原文:https://www.cnblogs.com/victorique/p/10420629.html

(0)
(0)
   
举报
评论 一句话评论(0
关于我们 - 联系我们 - 留言反馈 - 联系我们:wmxa8@hotmail.com
© 2014 bubuko.com 版权所有
打开技术之扣,分享程序人生!