51nod 1222 最小公倍数计数

1222 最小公倍数计数
题目来源： Project Euler
基准时间限制：6 秒 空间限制：131072 KB 分值: 640 难度：8级算法题 收藏 关注
定义F(n)表示最小公倍数为n的二元组的数量。
即：如果存在两个数（二元组）X，Y(X <= Y)，它们的最小公倍数为N，则F(n)的计数加1。
例如：F(6) = 5，因为[2,3] [1,6] [2,6] [3,6] [6,6]的最小公倍数等于6。

给出一个区间[a,b]，求最小公倍数在这个区间的不同二元组的数量。
例如：a = 4，b = 6。符合条件的二元组包括：
[1,4] [2,4] [4,4] [1,5] [5,5] [2,3] [1,6] [2,6] [3,6] [6,6]，共10组不同的组合。
Input
输入数据包括2个数：a, b，中间用空格分隔（1 <= a <= b <= 10^11)。
Output
输出最小公倍数在这个区间的不同二元组的数量。
Input示例
4 6
Output示例
10

---

【分析】

出题人丧心病狂…奇葩反演层出不穷
$ans=\sum_{i=1}^{n}\sum_{j=1}^{n}[lcm(i,j)<=n]$
$=\sum_{d=1}^{n} \mu(d) \sum_{i=1}^{n}\sum_{j=1}^{n}\sum_{k=1}^{n}[ijk<=n/d^2]$

大概过程就是把 $lcm(i,j)$ 换成 $ij/gcd(i,j)$，然后枚举gcd搞一搞。
令 $f(x)=\sum_{i=1}^{n}\sum_{j=1}^{n}\sum_{k=1}^{n}[ijk<=x]$

观察到 $d$ 最多枚举至$\sqrt{n}，$暴算一波$f(x)$出解。

---

【代码】

//51nod 1222 最小公倍数计数 
#include<bits/stdc++.h>
#define N 1000000
#define ll long long
#define M(a) memset(a,0,sizeof a)
#define fo(i,j,k) for(i=j;i<=k;i++)
using namespace std;
const int mxn=1000005;
ll n,m,l,r;
bool vis[mxn];
int pri[mxn],miu[mxn];
inline void init()
{
    miu[1]=1;
    for(int i=2;i<=N;i++)
    {
        if(!vis[i]) pri[++pri[0]]=i,miu[i]=-1;
        for(int j=1;j<=pri[0]&&(ll)i*pri[j]<=N;j++)
        {
            vis[i*pri[j]]=1;
            if(i%pri[j]==0) break;
            miu[i*pri[j]]=-miu[i];
        }
    }
}
inline ll solve(ll n)
{
    ll ans=0;
    for(ll k=1;k*k<=n;k++) if(miu[k])
    {
        ll x=n/(k*k),tmp=0;
        for(ll i=1;i*i*i<=x;i++)
        {
            tmp+=(x/(i*i)-i)*3+1;   //前两项(三项)相同的情况 
            for(ll j=i+1;i*j*j<=x;j++)
              tmp+=(x/(i*j)-j)*6+3;  //三项不同(后两项相同)的情况 
        }
        ans+=tmp*miu[k];
    }
    return (ans+n)>>1;
}
int main()
{
    init();
    int i,j;
    scanf("%lld%lld",&l,&r);
    printf("%lld\n",solve(r)-solve(l-1));
    return 0;
}