【51Nod 1222】最小公倍數計數

http://www.51nod.com/onlineJudge/questionCode.html#!problemId=1222
求\([a,b]\)中的個數轉化爲求\([1,b]\)中的個數減去\([1,a)\)中的個數。
\[ \begin{aligned} &\sum_{i=1}^n\sum_{j=1}^n\left[\frac{ij}{(i,j)}\leq n\right]\\ =&\sum_{d=1}^n\sum_{i=1}^{\left\lfloor\frac nd\right\rfloor}\sum_{j=1}^{\left\lfloor\frac nd\right\rfloor}[(i,j)=1][ijd\leq n]\\ =&\sum_{d=1}^n\sum_{i=1}^{\left\lfloor\frac nd\right\rfloor}\sum_{j=1}^{\left\lfloor\frac nd\right\rfloor}\sum_{d'|(i,j)}\mu(d')[ijd\leq n]\\ =&\sum_{d=1}^n\sum_{d'=1}^{\left\lfloor\frac nd\right\rfloor}\mu(d')\sum_{i=1}^{\left\lfloor\frac n{dd'}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac n{dd'}\right\rfloor}\left[ijdd'^2\leq n\right]\\ =&\sum_{d'=1}^{\left\lfloor\sqrt n\right\rfloor}\mu(d')\sum_{d=1}^{\left\lfloor\frac n{d'}\right\rfloor}\sum_{i=1}^{\left\lfloor\frac n{dd'}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac n{dd'}\right\rfloor}\left[ijd\leq\left\lfloor\frac n{d'^2}\right\rfloor\right]\\ \end{aligned} \]
枚舉\(\left\lfloor\sqrt n\right\rfloor\)個\(d'\)，對於每一個\(d'\)，先假定\(d\leq i\leq j\)，枚舉\(\left\lfloor\left(\frac n{d'}\right)^{\frac 13}\right\rfloor\)個\(d\)，\(i\)再從\(d\)枚舉到\(\left\lfloor\left(\frac n{dd'}\right)^{\frac 12}\right\rfloor\)，最後計算\(j\)的個數，乘上一個排列數便可。
\(f(x)\)表示對每一個\(d'\)計算的複雜度。
時間複雜度爲\(O\left(\sum\limits_{d'=1}^{\left\lfloor\sqrt n\right\rfloor}f(d')\right)\)。
\[ \begin{aligned} f(x)=&\int_{0}^{\left(\frac n{x^2}\right)^{\frac 13}}\left(\left(\frac n{x^2i}\right)^{\frac 12}-i\right)di\\ =&\left(\frac n{x^2}\right)^{\frac 12}\int_0^{\left(\frac n{x^2}\right)^{\frac 13}}i^{-\frac 12}di-\int_0^{\left(\frac n{x^2}\right)^{\frac 13}}idi\\ =&\left(\frac n{x^2}\right)^{\frac 12}\left(2\left({\left(\frac n{x^2}\right)^{\frac 13}}\right)^{\frac 12}-2\right)-\left(\frac 12\left({\left(\frac n{x^2}\right)^{\frac 13}}\right)^2\right)\\ =&\frac 34\left(\frac n{x^2}\right)^{\frac 23}-2\left(\frac n{x^2}\right)^{\frac 12} \end{aligned} \]
\[ \begin{aligned} &O\left(\sum_{d'=1}^{\left\lfloor\sqrt n\right\rfloor}f(d')\right)\\ =&O\left(\int_0^{\sqrt n}\left(\frac n{x^2}\right)^{\frac 23}dx\right)\\ =&O\left(n^{\frac 23}\int_0^{\sqrt n}x^{-\frac 43}dx\right)\\ =&O\left(n^{\frac 23}\left(\left(-3\sqrt{n}^{-\frac 13}\right)-\left(-3\times0^{-\frac 13}\right)\right)\right)\\ =&O\left(3n^{\frac 23}-3n^{\frac 12}\right)\\ =&O\left(n^{\frac 23}\right) \end{aligned} \]
因此總時間複雜度爲\(O\left(n^{\frac 23}\right)\)，好神奇啊qwq

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
typedef long long ll;

const int sq = 316228;

bool notp[sq + 1];
int mu[sq + 1], prime[sq], num = 0;

void Euler_shai() {
    mu[1] = 1;
    for (int i = 2; i <= sq; ++i) {
        if (!notp[i]) prime[++num] = i, mu[i] = -1;
        for (int j = 1; j <= num && prime[j] * i <= sq; ++j) {
            notp[prime[j] * i] = true;
            if (i % prime[j] == 0)
                break;
            mu[prime[j] * i] = -mu[i];
        }
    }
}

ll a, b;

ll cal(ll n) {
    ll ret = 0;
    for (int d = 1; 1ll * d * d <= n; ++d)
        if (mu[d]) {
            ll up = n / d / d, r = 0;
            for (int i = 1; 1ll * i * i * i <= up; ++i) {
                ll up2 = up / i;
                r += (up2 / i - i) * 3 + 1;
                for (int j = i + 1; 1ll * j * j <= up2; ++j)
                    r += (up2 / j - j) * 6 + 3;
            }
            mu[d] > 0 ? ret += r : ret -= r;
        }
    return (ret + n) >> 1;
}

int main() {
    scanf("%lld%lld", &a, &b);
    Euler_shai();
    printf("%lld\n", cal(b) - cal(a - 1));
    return 0;
}