bzoj 2154

最新推荐文章于 2020-07-16 21:03:57 发布

weixin_30820151

最新推荐文章于 2020-07-16 21:03:57 发布

阅读量82

点赞数

CC 4.0 BY-SA版权

原文链接：http://www.cnblogs.com/mrasd/p/9520543.html

bzoj 2154

题意

求 $ \sum_{x=1}^{n} \sum_{y=1}^{m} lcm(x,y) $

解法

详见笔记本。

代码如下：

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <cctype>
#define INF 2139062143
#define MAX 0x7ffffffffffffff
#define del(a,b) memset(a,b,sizeof(a))
using namespace std;
typedef long long ll;
template<typename T>
inline void read(T&x)
{
    x=0;T k=1;char c=getchar();
    while(!isdigit(c)){if(c=='-')k=-1;c=getchar();}
    while(isdigit(c)){x=x*10+c-'0';c=getchar();}x*=k;
}
const int maxn=1e7+5;
const int mod=20101009;
int add(int x,int y) {return (x+y)%mod;}
int mul(int x,int y) {return 1ll*x*y%mod;}
int pul(int x,int y) {return (x-y+mod)%mod;}
int poww(int a,int b){
    int ans=1;
    while(b){
        if(b&1) ans=mul(ans,a);
        a=mul(a,a);
        b>>=1;
    }
    return ans;
}
int inv(int x) {return poww(x,mod-2);}
int p[maxn];
int mo[maxn];
int sum[maxn];
bool nop[maxn];

void make_mo(int n){
    nop[1]=1;
    mo[1]=1;
    for(int i=2;i<=n;i++){
        if(!nop[i]) {
            mo[i]=-1;
            p[++p[0]]=i;
        }
        for(int j=1;j<=p[0]&&i*p[j]<=n;j++){
            nop[i*p[j]]=1;
            if(i%p[j]==0) {
                mo[i*p[j]]=0;
                break;
            }
            mo[i*p[j]]=-mo[i];
        }
    }
    for(int i=1;i<=n;i++) mo[i]=mul(add(mo[i],mod),mul(i,i));
    for(int i=1;i<=n;i++) sum[i]=add(sum[i-1],mo[i]);
}
int inv_2=inv(2);
int g_sum(int x,int y){
    x=mul(inv_2,mul(x,add(x,1)));
    y=mul(inv_2,mul(y,add(y,1)));
    return mul(x,y);
}
int get_g(int x,int y){
    int ans=0;
    for(int l=1,r;l<=x;l=r+1){
        r=min( x/(x/l) , y/(y/l) );
        ans=add( ans , mul(pul(sum[r],sum[l-1]), g_sum( x/l , y/l ) ) );
    }
    return ans;
}
int n,m;
int main()
{
    read(n),read(m);
    if(n>m) swap(n,m);
    make_mo(m);
    int ans=0;
    for(int l=1,r;l<=n;l=r+1){
        r=min( n/(n/l) , m/(m/l) );
        ans=add( ans , mul( 1LL*(l+r)*(r-l+1)/2%mod ,  get_g(n/l,m/l) ));
    }
    printf("%d",ans);
    return 0;
}

转载于:https://www.cnblogs.com/mrasd/p/9520543.html