ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

于神之怒加强版 (莫反)

2019-03-16 18:40:24  阅读:301  来源: 互联网

标签:lfloor 神之怒 frac 加强版 int sum 莫反 mu rfloor


题目链接: 点击这里

首先哦,令n<=m,我们来化简式子:
\[ \sum_{i=1}^n\sum_{j=1}^mgcd(i,j)^k \]
设gcd(i,j)=d,则可以变成这样:
\[ \sum_{i=1}^n\sum_{j=1}^md^k\\ \]
把d提到前面去,来枚举d:
\[ \sum_{d=1}^nd^k\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor}[gcd(i,j)=1] \]
然后就可以开始套了
\[ \sum_{d=1}^nd^k\sum_{i=1}^{\lfloor{\frac{n}{d}}\rfloor}\sum_{j=1}^{\lfloor{\frac{m}{d}}\rfloor}\sum_{x|i,x|j}\mu(x) \]
然后找到满足条件的x的个数,就可以和前面两个sigma说再见了
\[ \sum_{d=1}^nd^k\sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor}\mu(x)\lfloor{\frac{n}{dx}}\rfloor\lfloor{\frac{m}{dx}}\rfloor\\ 设T=dx,再来枚举T\\ \sum_{d=1}^nd^k\sum_{x=1}^{\lfloor{\frac{n}{d}}\rfloor}\mu(\frac{T}{d})\lfloor{\frac{n}{T}}\rfloor\lfloor{\frac{m}{T}}\rfloor\\ \sum_{T=1}^n\lfloor{\frac{n}{T}}\rfloor\lfloor{\frac{m}{T}}\rfloor\sum_{d|T}d^k\mu(\frac{T}{d}) \]
一般的题到这一步就能直接搞了,然而会T,所以我们来继续推:

注:以下受@hby大佬的启发,与本人博客差不多


\[ f(x)=\sum_{d|T}d^k\mu(\frac{T}{d})\\ g(x)=x^k \]

看!\(f(x)\)是不是\(/mu\)和\(g\)的卷积形式!所以\(f\)是个积性函数!

我们来对\(f(x)\)进行讨论,首先我们将\(x\)分解质因数
\[ x=\Pi_{p_i^{a_i}},p_i\in\mathbb{P}\\ f(x)=\Pi_{f(p_i^{a_i})}\\ f(p_i^{a^i})=\sum_{d|p_i^{a_i}}d^k\mu(\frac{p_i^{a_i}}{d})\\ 因为p_i是个质数,所以\\ f(p_i^{a^i})=\sum_{i=0}^{a_i}{p_i}^k\mu(p_i^{a_i-i})\\ 再考虑一下\mu的性质,可以得出只有当i=a^i或i=a^i-1时,\mu不为0,于是\\ f(p_i^{a_i})=p_i^{a_ik}-p_i^{(a_i-1)k}\\ 于是对于f(xp),若p|x,则只是多了几个指数而已,那么f(xp)=f(x)*p^k\\ 否则,就是添加了几个质因数,则f(xp)=f(x)*(p^k-1) \]

Code:

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=5e6+1;
const int pps=1e9+7;
int n,m,k,t,tot;
int pri[N],mu[N],vis[N],p[N],f[N];
int quick(int a,int b){
    int sum=1;
    while(b){
        if(b&1) sum*=a,sum%=pps;
        a=(a*a)%pps,b>>=1;
    }return sum%pps;
}
void prepare(){
    f[1]=1;
    for(int i=2;i<=N;i++){
        if(!vis[i]) pri[++tot]=i,p[tot]=quick(i,k),f[i]=p[tot]-1;
        for(int j=1;j<=tot&&pri[j]*i<=N;j++){
            vis[i*pri[j]]=1;
            if(!(i%pri[j])){
                f[i*pri[j]]=f[i]*p[j]%pps;
                break;
            }f[i*pri[j]]=f[i]*(p[j]-1)%pps;
        }
    }for(int i=1;i<=N;i++) f[i]=(f[i]+f[i-1])%pps;
}
int calc(int n,int m){
    int d=1,sum=0;
    while(d<=n&&d<=m){
        int pre=d;d=min(n/(n/d),m/(m/d));
        sum=(sum+(n/d)*(m/d)%pps*(f[d]-f[pre-1])%pps)%pps;
        ++d;
    }return (sum+pps)%pps;
}
int read(){
    int x=0,f=1;char ch=getchar();
    while(!isdigit(ch)){if(ch=='-')f=-f;ch=getchar();}
    while(isdigit(ch)){x=x*10+ch-48;ch=getchar();}
    return x*f;
}
signed main(){
    t=read(),k=read();
    prepare();
    while(t--){
        n=read(),m=read();
        printf("%lld\n",calc(n,m));
    }
    return 0;
}

标签:lfloor,神之怒,frac,加强版,int,sum,莫反,mu,rfloor
来源: https://www.cnblogs.com/NLDQY/p/10543645.html

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有