1227 平均最小公倍数

Lcm(a,b)表示a和b的最小公倍数，A(n)表示Lcm(n,i)的平均数(1 <= i <= n)，
例如：A(4) = (Lcm(1,4) + Lcm(2,4) + Lcm(3,4) + Lcm(4,4)) / 4 = (4 + 4 + 12 + 4) / 4 = 6。
F(a, b) = A(a) + A(a + 1) + ...... A(b)。(F(a,b) = ∑A(k), a <= k <= b)
例如：F(2, 4) = A(2) + A(3) + A(4) = 2 + 4 + 6 = 12。
给出a,b，计算F(a, b)，由于结果可能很大，输出F(a, b) % 1000000007的结果即可。

输入

输入2个数a,b，中间用空格分隔(1 <= a <= b <= 10^9)。

输出

输出F(a, b) % 1000000007的结果。

输入样例

1 100

输出样例

122726

题解

\[ A(n)=\frac 1n\sum_{i=1}^n\textrm{lcm(i,n)}=\frac 1n\sum_{i=1}^n\frac{in}{\gcd(i,n)}\\ =\sum_{i=1}^n\frac i{\gcd(i,n)}\\ =\sum_{d|n}\sum_{i=1}^{\frac nd}i[\gcd(i,\frac nd)=1] \]

根据公式\(\sum_{i=1}^ni[\gcd(i,n)=1]=\frac{n\varphi(n)+e(n)}{2}\)，化简可得
\[ A(n)=\sum_{d|n}\frac{\frac nd\varphi(\frac nd)+e(\frac nd)}2\\ \frac 12(\sum_{d|n}d\varphi(d)+1)\\ F(n)=\sum_{i=1}^nA(n)=\sum_{i=1}^n\frac 12(\sum_{d|i}d\varphi(d)+1)\\ =\frac 12\left(\sum_{i=1}^n\sum_{d|i}d\varphi(d)+n\right) \]
\(d\varphi(d)\)就是\(id\cdot \varphi\)，看起来非常眼熟。所以使用杜教筛套路，把枚举约数改为枚举倍数
\[ F(n)=\frac 12\left(\sum_{i=1}^n\sum_{k=1}^{\lfloor\frac ni\rfloor}k\varphi(k)+n\right) \]
令\(f(n)=n\varphi(n)\)，它的前缀和\(S(n)=\sum_{i=1}^nf(i)\)。则
\[ F(n)=\frac 12\left(\sum_{i=1}^nS(\lfloor\frac ni\rfloor)+n\right) \]
只要能求出\(S\)就可以整除分块。于是构造\(g=id\)，则
\[ (f*g)(n)=\sum_{d|n}f(n)g(\frac nd)=\sum_{d|n}d\varphi(d)\frac nd\\ =n\sum_{d|n}\varphi(d)=n^2 \]
至此问题解决，上杜教筛即可。时间复杂度\(O(n^\frac 23)\)

#include<bits/stdc++.h>
#define il inline
#define co const
template<class T>T read(){
    T data=0,w=1;char ch=getchar();
    for(;!isdigit(ch);ch=getchar())if(ch=='-') w=-w;
    for(;isdigit(ch);ch=getchar()) data=data*10+ch-'0';
    return data*w;
}
template<class T>il T read(T&x) {return x=read<T>();}
typedef long long LL;
using namespace std;

co int N=1e6+1,mod=1e9+7,i2=500000004,i6=166666668;
il int add(int a,int b){
    return (a+=b)>=mod?a-mod:a;
}
il int mul(int a,int b){
    return (LL)a*b%mod;
}
int pri[N],tot,phi[N];
void init(){
    pri[1]=phi[1]=1;
    for(int i=2;i<N;++i){
        if(!pri[i]) pri[++tot]=i,phi[i]=i-1;
        for(int j=1;j<=tot&&i*pri[j]<N;++j){
            pri[i*pri[j]]=1;
            if(i%pri[j]==0){
                phi[i*pri[j]]=phi[i]*pri[j];
                break;
            }
            phi[i*pri[j]]=phi[i]*phi[pri[j]];
        }
    }
    for(int i=2;i<N;++i) phi[i]=add(phi[i-1],mul(phi[i],i));
}
map<int,int> sf;
int sum_f(int n){
    if(n<N) return phi[n];
    if(sf.count(n)) return sf[n];
    int ans=mul(i6,mul(n,mul(n+1,2*n+1)));
    for(int l=2,r;l<=n;l=r+1){
        r=n/(n/l);
        ans=add(ans,mod-mul(mul(i2,mul(l+r,r-l+1)),sum_f(n/l)));
    }
    return sf[n]=ans;
}
int solve(int n){
    int ans=n;
    for(int l=1,r;l<=n;l=r+1){
        r=n/(n/l);
        ans=add(ans,mul(r-l+1,sum_f(n/l)));
    }
    return mul(i2,ans);
}
int main(){
    init();
    int a=read<int>(),b=read<int>();
    printf("%d\n",add(solve(b),mod-solve(a-1))); // edit 1:mod
    return 0;
}