【51nod 1190】最小公倍数之和 V2

Description

给出2个数a, b，求LCM(a,b) + LCM(a+1,b) + .. + LCM(b,b)。
例如：a = 1, b = 6，1,2,3,4,5,6 同6的最小公倍数分别为6,6,6,12,30,6，加在一起 = 66。
由于结果可能很大，输出Mod 10^9 + 7的结果。（测试数据为随机数据，没有构造特别坑人的Test）

Solution

ans=∑i=abi∗bgcd(i,b)=b∗∑d|b∑i=⌈a/d⌉⌊b/d⌋i∗∑d′|gcd(i,b/d)μ(d′)=b∗∑d|b∑d′|⌊b/d⌋μ(d′)∑i=⌈a/d⌉,d′|i⌊b/d⌋i

ans=b∗∑d|b∑d′|⌊b/d⌋μ(d′)(⌈add′⌉+⌊bdd′⌋)(⌊bdd′⌋−⌈add′⌉+1)2∗d′

令T=d*d’

ans=∑T|b(⌈aT⌉+⌊bT⌋)(⌊bT⌋−⌈aT⌉+1)2∑d′|Tμ(d′)∗d′

设f(T)=∑d′|Tμ(d′)∗d′。
由于μ(d)和d都是积性函数，所以当gcd（p,q）=1时，f(p*q)=f(p)*f(q)。
并且由μ(d)的性质可知，T里的每个质数指数最多为1，否则μ(d)=0。
我们可以对b分解质因数，枚举b的约数T，同时就能算出∑d′|Tμ(d′)∗d′，就如有x∗pk=Tp这个质因子时，你会有选与不选这两个决策，不选的话贡献为1，选的话贡献为-p,然后就f(T)=f(x)∗(1−p)。∑d′|Tμ(d′)∗d′的计算类似于(a+b)*(d+c)=ac+ad+bc+bd。时间复杂度O(28∗test).

Code

#include<iostream>#include<math.h>#include<string.h>#include<stdio.h>#include<algorithm>#define ll long longusing namespace std;const ll maxn=1e5,mo=1e9+7,mo2=5e8+4;int d[maxn],bz[maxn+5],b1[maxn],c[maxn];ll n,i,t,j,k,l,a,b,z,ans,p;void dg(int x,int y,ll z){    ll t=1,k;    if (y>b1[0]){        t=(a%x)?a/x+1:a/x;        k=b/x;        ans+=(t+k)*(k-t+1)*z;        return;    }    dg(x,y+1,z);    for (k=1;k<=c[y];k++)        t*=b1[y],dg(x*t,y+1,z*(1-b1[y]));}int main(){    scanf("%lld",&n);    for (i=2;i<=maxn;i++){        if (!bz[i]) d[++d[0]]=i;        for (j=1;j<=d[0];j++){            if (i*d[j]>maxn) break;            bz[i*d[j]]=1;            if (i%d[j]==0) break;        }    }    for (i=1;i<=n;i++){        scanf("%lld%lld",&a,&b);        z=b;b1[0]=0;        for (j=1;d[j]*d[j]<=z;j++){            if (z%d[j]) continue;            b1[++b1[0]]=d[j];c[b1[0]]=0;            while (z%d[j]==0) z/=d[j],c[b1[0]]++;        }        if (z>1)b1[++b1[0]]=z,c[b1[0]]=1;        ans=0;        dg(1,1,1);ans=(ans%mo+mo)%mo*mo2%mo;        ans=b%mo*ans%mo;        printf("%lld\n",ans);    }}