• 牛客小白月赛13-J小A的数学题 (莫比乌斯反演)


    链接:https://ac.nowcoder.com/acm/contest/549/J
    来源:牛客网
    题目描述
    小A最近开始研究数论题了,这一次他随手写出来一个式子,∑ni=1∑mj=1gcd(i,j)2∑i=1n∑j=1mgcd(i,j)2,但是他发现他并不太会计算这个式子,你可以告诉他这个结果吗,答案可能会比较大,请模上1000000007。
    输入描述:
    一行两个正整数n,m一行两个正整数n,m
    输出描述:
    一行一个整数表示输出结果一行一个整数表示输出结果
     
    输入:
    2 2
    输出:
    7
    1=<n,m<=1e6
    解题思路:这题应该算是一题莫比乌斯反演的套路题了,感觉莫比乌斯真的好难啊,学了好久感觉也没懂,这题算是它的一个简单应用。
    具体可以参考博客:https://blog.sengxian.com/algorithms/mobius-inversion-formula
    莫比乌斯反演经典套路:
    现在有个积性函数f(n),设n<m,则:

    于是原来的式子就变成了求fμ了,再用上整数分块就可以快速搞定了。

    自己推演了一遍:

    代码:
    #include<iostream>
    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #include<vector>
    #include<string>
    #include<set>
    #include<cmath>
    #include<list>
    #include<deque>
    #include<cstdlib>
    #include<bitset>
    #include<stack>
    #include<map>
    #include<cstdio>
    #include<queue>
    using namespace std;
    typedef long long ll;
    typedef unsigned long long ull;
    #define lson l,mid,rt<<1
    #define rson mid+1,r,rt<<1|1
    const int INF=0x3f3f3f3f;
    const double PI=acos(-1.0);
    const double eps=1e-6;
    const int maxn=1000005;
    ll gcd(ll a,ll b){return b?gcd(b,a%b):a;}
    ll lcm(ll a,ll b){return a/gcd(a,b)*b;}
    const int mod=1e9+7;
    const int dir[4][2]={{1,0},{-1,0},{0,1},{0,-1}};
    const int N=1e6+10;
    ll n,m,prime[N],mu[N],tot;
    void getMu(){
        for(int i=1;i<=1e6+5;i++) prime[i]=1;
        mu[1]=1;
        for(int i=2;i<=1e6+5;i++){
            if(prime[i]){
                prime[++tot]=i;
                mu[i]=-1;
            }
                for(int j=1;j<=tot&&prime[j]*i<=1e6+5;j++){
                    prime[i*prime[j]]=0;
                    if(i%prime[j]==0){
                        mu[i*prime[j]]=0;
                        break;
                    }else mu[i*prime[j]]=-mu[i];
                }
        }
    }
    int main(){
        cin>>n>>m;
        getMu();
        ll ans=0;
        for(ll i=1;i<=min(n,m);i++){
            ll tmp=0;
            for(ll j=i;j<=min(n,m);j+=i){
                tmp=(tmp+mu[j/i]*(n/j)*(m/j))%mod;
            }
            ans=(ans+tmp*i*i%mod)%mod;
        }
        cout<<ans<<endl;
        return 0;
    }

    整除分块优化:

    #include<iostream>
    #include<cstdio>
    using namespace std;
    typedef long long ll;
    const int maxn=1e6+7;
    const int mod=1e9+7;
    ll n,m,prime[maxn],mu[maxn],sum[maxn],tot,ans;
    void getMobius(int N){
        for(int i=1;i<=N;i++)prime[i]=1;
        mu[1]=1;
        for(int i=2;i<=N;i++){
            if(prime[i]){
                prime[tot++]=i;
                mu[i]=-1;
            }
            for(int j=0;j<tot&&i*prime[j]<=N;j++){
                prime[i*prime[j]]=0;
                if(i%prime[j]==0){
                    mu[i*prime[j]]=0;
                    break;
                }
                mu[i*prime[j]]=-mu[i];
            }
        }
    }
    ll solve(ll a,ll b){
        ll res=0;
        for(int l=1,r;l<=min(a,b);l=r+1){
            r=min(a/(a/l),b/(b/l));
            res=(res+(sum[r]-sum[l-1])%mod*(a/l)%mod*(b/l)%mod)%mod;
        }
        return res;
    }
    int main(){
        scanf("%lld%lld",&n,&m);
        if(n>m) swap(n,m);
        getMobius(1e6);
        sum[1]=1;
        for(int i=2;i<=1e6;i++) sum[i]=sum[i-1]+mu[i];
        for(ll l=1,r;l<=n;l=r+1){
            r=min(n/(n/l),m/(m/l));
            ll dd=(r*(r+1)*(2*r+1)/6-(l-1)*l*(2*l-1)/6)%mod;
            ans=(ans+dd*solve(n/l,m/l)%mod)%mod;
        }
        printf("%lld
    ",ans);
        return 0;
    }
  • 相关阅读:
    [LOJ537] DNA 序列
    yum安装jdk
    centos7添加tomcat自启服务
    ajax下载文件的方法 php下载图片的方法
    kali linux dpkg安装和删除的方法
    javascript 曲线流光动画
    javascript获取曲线路径每个像素的坐标
    抓取微信小程序源码的方法
    javascript canvas拖尾效果
    php错误execution timed out (118.236064 sec), terminating解决方法
  • 原文地址:https://www.cnblogs.com/zjl192628928/p/10758185.html
Copyright © 2020-2023  润新知