51nod 1222 最小公倍数计数

题意

对于区间，我们求出\([1,b]\)的答案减去\([1,a-1]\)的即可。

求\([1,n]\)的答案：
\(ans=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^n[\frac{i*j}{\gcd(i,j)}\leqslant n]\)
\(=\sum\limits_{d=1}^n\sum\limits_{i=1}^{\frac{n}{d}}\sum\limits_{j=1}^{\frac{n}{d}}[i*j\leqslant \frac{n}{d}][\gcd(i,j)=1]\)
\(=\sum\limits_{d=1}^n\sum\limits_{k=1}^{\frac{n}{d}}\mu(k)\sum\limits_{i=1}^{\frac{n}{dk}}\sum\limits_{j=1}^{\frac{n}{dk}}[i*k*j*k\leqslant \frac{n}{d}]\)
\(=\sum\limits_{d=1}^n\sum\limits_{k=1}^{\frac{n}{d}}\mu(k)\sum\limits_{i=1}^{\frac{n}{dk}}\sum\limits_{j=1}^{\frac{n}{dk}}[i*j*d\leqslant \frac{n}{k^2}]\)
\(=\sum\limits_{k=1}^n\mu(k)\sum\limits_{d=1}^{\frac{n}{k}}\sum\limits_{i=1}^{\frac{n}{dk}}\sum\limits_{j=1}^{\frac{n}{dk}}[i*j*d\leqslant \frac{n}{k^2}]\)
\(=\sum\limits_{k=1}^{\sqrt{n}}\mu(k)\sum\limits_{d=1}^{\frac{n}{k^2}}\sum\limits_{i=1}^{\frac{n}{dk^2}}\sum\limits_{j=1}^{\frac{n}{dk^2}}[i*j*d\leqslant \frac{n}{k^2}]\)
于是我们枚举\(k\)，算出满足\([d*i*j\leqslant \frac{n}{k^2}]\)的三元组\((d,i,j)\)有多少。

我们发现\((d,i,j)\)的上界似乎没什么用，我们只要找到一个三元组\((a,b,c)\)满足\(a*b*c\leqslant \frac{n}{k^2}\)，就一定能满足上面的关系。

于是我们就是要求三元组\((a,b,c)\)满足\(a*b*c\leqslant \frac{n}{k^2}\)的个数。

我们分别统计三个都不相同，有两个相同和三个都相同的方案数。

我们先求各不相同的，枚举小的那个，设为\(x\)，\(x^3\leqslant \frac{n}{k^2}\)。之后枚举第二大的，设为\(y\)，\(x*y*y\leqslant \frac{n}{k^2}\)，最大的那个的范围就能确定。中途再算一下\((x,x,y),(x,y,y),(x,x,x)\)这种，统计一下就好了。

code:

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int maxn=1e6; 
int a,b;
int mu[maxn];
bool vis[maxn];
vector<int>prime;
inline int read()
{
	char c=getchar();int res=0,f=1;
	while(c<'0'||c>'9'){if(c=='-')f=-1;c=getchar();}
	while(c>='0'&&c<='9')res=res*10+c-'0',c=getchar();
	return res*f;
}
inline void prework(int n)
{
	vis[1]=mu[1]=1;
	for(int i=2;i<=n;i++)
	{
		if(!vis[i])prime.push_back(i),mu[i]=-1;
		for(unsigned int j=0;j<prime.size()&&i*prime[j]<=n;j++)
		{
			vis[i*prime[j]]=1;
			if(i%prime[j]==0)break;
			mu[i*prime[j]]=mu[i]*mu[prime[j]];
		}
	}
}
inline int calc(int n)
{
	int res=0;
	for(int k=1;k*k<=n;k++)
	{
		if(!mu[k])continue;
		int tmp=0,up=n/(k*k);
		for(int i=1;i*i*i<=up;i++)
		{
			for(int j=i+1;j*j*i<=up;j++)tmp+=(up/(i*j)-j)*6+3;//三个不同+(i,j,j)
			tmp+=(up/(i*i)-i)*3;//(i,i,j)
			tmp++;//(i,i,i)
		}
		res+=mu[k]*tmp;
	}
	return (res+n)/2;
}
signed main()
{
	prework(1e6);
	a=read(),b=read();
	printf("%lld",calc(b)-calc(a-1));
	return 0;
}