51NOD 1222 最小公倍数计数 [莫比乌斯反演 杜教筛]

1222 最小公倍数计数

题意：求有多少数对\((a,b):a<b\)满足\(lcm(a,b) \in [1, n]\)

\(n \le 10^{11}\)

---

卡内存！

枚举\(gcd, \frac{a}{gcd}, \frac{b}{gcd}\)，然后\(\mu\)代入，就是
\[ \sum_{d=1}^{\sqrt{n}}\mu(d) \sum_i \sum_j \sum_k [ijk \le \frac{n}{d^2}] \]
问题就是怎么求后面的式子了

一开始我是
\[ f(n) = \sum_i \sum_j \sum_k [ijk \le n] = \sum_{i=1}^n g(i) \\ g(n) = \sum_{i=1}^n \sigma(i) \]
杜教筛\(g\)，方法是卷上\(\mu\)，当然还要杜教筛\(\sum\mu\)，但不影响复杂度，还是\(O(n^{\frac{2}{3}})\)

本机6.5s，改小预处理的大小后当然T了

然后又用分块的方法算\(g\)，预处理前\(O(n^{2/3})\)的\(\sigma\)剩下的分块\(O({\sqrt{n}})\)计算，复杂度也是\(O(n^{\frac{2}{3}})\)

本机4.6s，改小预处理大小又T了...

最后还是用了tangjz的方法，统计\(abc \le n\)的数对个数，规定\(a\le b \le c\)，然后a枚举\(n^{\frac{1}{3}}\)，b枚举\(\sqrt{\frac{n}{a}}\)。最后乘上\(3!\)再减去\(a= b <c,\ a<b = c,\ a=b=c\)

这一部分的复杂度\(T'(n)=O(\sum_{i=1}^{n^\frac{1}{3}}{\sqrt{\frac{n}{i}}})=O(n^\frac{2}{3})\)

总体复杂度\(T(n)=O(\sum_{d=1}^{\sqrt{n}}{T'(\frac{n}{d^2})})=?O(n^\frac{2}{3})​\)

求积分好像是\(\sqrt{n}\)啊，我也不知道这是怎么算的

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#include <ctime>
using namespace std;
typedef long long ll;
const int N=320000;
int U=316230;
inline ll read(){
    char c=getchar(); ll x=0,f=1;
    while(c<'0' || c>'9') {if(c=='-')f=-1; c=getchar();}
    while(c>='0' && c<='9') {x=x*10+c-'0'; c=getchar();}
    return x*f;
}

bool notp[N]; int p[N], mu[N];
void sieve(int n) {
    mu[1]=1; 
    for(int i=2; i<=n; i++) {
        if(!notp[i]) p[++p[0]] = i, mu[i] = -1;
        for(int j=1; j <= p[0] && i*p[j] <= n; j++) {
            notp[ i*p[j] ] = 1;
            if(i % p[j] == 0) {mu[ i*p[j] ] = 0; break;}
            mu[ i*p[j] ] = -mu[i];
        }
    }
}

inline ll cal(ll n) {
    ll ans=0;
    for(ll i=1; i*i*i <= n; i++) {
        for(ll j=i; i*j*j <= n; j++) {
            ans += (n/i/j - j + 1) * 6;
            if(i == j) ans -= (n/i/j - j) * 3;
            if(i != j) ans -= 3;
        }
        ans -= 5;
    }
    return ans;
}

ll solve(ll n) {
    ll ans=0;
    int m = sqrt(n);
    for(int i=1; i<=m; i++) if(mu[i]) ans += mu[i]>0 ? cal(n/i/i) : -cal(n/i/i);
    return (ans + n) / 2;
}
int main() {
    freopen("in", "r", stdin);
    sieve(U);
    ll l=read(), r=read();
    printf("%lld", solve(r) - solve(l-1));
    printf("\n\n%lf\n",(double)clock()/CLOCKS_PER_SEC);
}

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#include <map>
#include <ctime>
using namespace std;
typedef long long ll;
const int N=22000005;
int U=22000000;
inline ll read(){
    char c=getchar(); ll x=0,f=1;
    while(c<'0' || c>'9') {if(c=='-')f=-1; c=getchar();}
    while(c>='0' && c<='9') {x=x*10+c-'0'; c=getchar();}
    return x*f;
}

bool notp[N]; int p[N/10], mu[N], lp[N]; ll si[N];
void sieve(int n) {
    mu[1]=1; si[1]=1;
    for(int i=2; i<=n; i++) {
        if(!notp[i]) p[++p[0]] = i, mu[i] = -1, si[i] = lp[i] = 2;
        for(int j=1; j <= p[0] && i*p[j] <= n; j++) {
            int t = i*p[j];
            notp[t] = 1;
            if(i%p[j] == 0) {
                mu[t] = 0;
                lp[t] = lp[i] + 1;
                si[t] = si[i] / lp[i] * lp[t];
                break;
            }
            mu[t] = -mu[i];
            lp[t] = 2;
            si[t] = si[i] * 2;
        }
        mu[i] += mu[i-1];
        si[i] += si[i-1];
    }
}

namespace ha {
    const int p = 1001001;
    struct ha {
        struct meow{int ne; ll val, r;} e[10000];
        int cnt, h[p];
        ha() {cnt=0; memset(h, 0, sizeof(h));}
        inline void insert(ll x, ll val) {
            int u = x%p;
            for(int i=h[u];i;i=e[i].ne) if(e[i].r == x) return;
            e[++cnt] = (meow){h[u], val, x}; h[u] = cnt;
        }
        inline ll quer(ll x) {
            int u = x%p;
            for(int i=h[u];i;i=e[i].ne) if(e[i].r == x) return e[i].val;
            return -1;
        }
    } hs, hu, me; 
} using ha::hs; using ha::hu; using ha::me;

ll dj_u(ll n) { //return 1;
    if(n <= U) return mu[n];
    if(hu.quer(n) != -1) return hu.quer(n);
    ll ans = 1, r;
    for(ll i=2; i<=n; i=r+1) {
        r = n/(n/i);
        ans -= (r-i+1) * dj_u(n/i);
    }
    hu.insert(n, ans);
    return ans;
}

ll dj_s(ll n) {
    if(n <= U) return si[n];
    if(hs.quer(n) != -1) return hs.quer(n);
    dj_u(n);
    ll ans = n, r, now, last = dj_u(1);
    for(ll i=2; i<=n; i=r+1, last=now) {
        r = n/(n/i); now = dj_u(r);
        ans -= (now - last) * dj_s(n/i);
    }
    hs.insert(n, ans);
    return ans;
}
ll cal(ll n) { 
    ll ans=0, r;
    for(ll i=1; i<=n; i=r+1) {
        r = n/(n/i);
        ans += (r-i+1) * dj_s(n/i);  
    }
    return ans;
}
ll solve(ll n) { 
    ll ans=0;
    int m = sqrt(n);
    for(ll i=1; i<=m; i++) if(mu[i] - mu[i-1]) ans += (mu[i] - mu[i-1]) * cal(n / (i*i));
    return (ans + n) / 2;
}
ll l, r;
int main() {
    freopen("in", "r", stdin);
    sieve(U);
    l=read(); r=read();
    //printf("%lld", dj_u(r));
    printf("%lld", solve(r) - solve(l-1));
    printf("\n\n%lf\n",(double)clock()/CLOCKS_PER_SEC);
}

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#include <ctime>
using namespace std;
typedef long long ll;
const int N=22000005;
int U=22000000;
inline ll read(){
    char c=getchar(); ll x=0,f=1;
    while(c<'0' || c>'9') {if(c=='-')f=-1; c=getchar();}
    while(c>='0' && c<='9') {x=x*10+c-'0'; c=getchar();}
    return x*f;
}

bool notp[N]; int p[N/10], mu[N], lp[N]; ll si[N];
void sieve(int n) {
    mu[1]=1; si[1]=1;
    for(int i=2; i<=n; i++) {
        if(!notp[i]) p[++p[0]] = i, mu[i] = -1, si[i] = lp[i] = 2;
        for(int j=1; j <= p[0] && i*p[j] <= n; j++) {
            int t = i*p[j];
            notp[t] = 1;
            if(i%p[j] == 0) {
                mu[t] = 0;
                lp[t] = lp[i] + 1;
                si[t] = si[i] / lp[i] * lp[t];
                break;
            }
            mu[t] = -mu[i];
            lp[t] = 2;
            si[t] = si[i] * 2;
        }
        //mu[i] += mu[i-1];
        si[i] += si[i-1];
    }
}
inline ll g(ll n) { if(n <= U) return si[n];
    ll ans=0, r;
    for(ll i=1; i<=n; i=r+1) {
        r = n/(n/i);
        ans += (r-i+1) * (n/i);
    }
    return ans;
}
ll cal(ll n) {
    ll ans=0, r;
    for(ll i=1; i<=n; i=r+1) {
        r = n/(n/i);
        ans += (r-i+1) * g(n/i); 
    }
    return ans;
}

ll solve(ll n) { 
    ll ans = 0;
    int m = sqrt(n);
    for(ll i=1; i<=m; i++) 
        if(mu[i]) ans += mu[i] * cal(n / (i*i));
    return (ans + n) / 2;
}

ll l, r;
int main() {
    freopen("in", "r", stdin);
    sieve(U);
    l=read(); r=read();
    printf("%lld", solve(r) - solve(l-1));
    printf("\n\n%lf\n",(double)clock()/CLOCKS_PER_SEC);
}