【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=T p這個質因子時，你會有選與不選這兩個決策，不選的話貢獻為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 long
using namespace std;
const ll maxn=e5,mo=e9+,mo2=e8+;
int d[maxn],bz[maxn+],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=,k;
    if (y>b1[]){
        t=(a%x)?a/x+:a/x;
        k=b/x;
        ans+=(t+k)*(k-t+)*z;
        return;
    }
    dg(x,y+,z);
    for (k=;k<=c[y];k++)
        t*=b1[y],dg(x*t,y+,z*(-b1[y]));
}
int main(){
    scanf("%lld",&n);
    for (i=;i<=maxn;i++){
        if (!bz[i]) d[++d[]]=i;
        for (j=;j<=d[];j++){
            if (i*d[j]>maxn) break;
            bz[i*d[j]]=;
            if (i%d[j]==) break;
        }
    }
    for (i=;i<=n;i++){
        scanf("%lld%lld",&a,&b);
        z=b;b1[]=;
        for (j=;d[j]*d[j]<=z;j++){
            if (z%d[j]) continue;
            b1[++b1[]]=d[j];c[b1[]]=;
            while (z%d[j]==) z/=d[j],c[b1[]]++;
        }
        if (z>)b1[++b1[]]=z,c[b1[]]=;
        ans=;
        dg(,,);ans=(ans%mo+mo)%mo*mo2%mo;
        ans=b%mo*ans%mo;
        printf("%lld\n",ans);
    }
}