51nod1222 最小公倍数计数

题目描述

定义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组不同的组合。
输入
输入数据包括2个数：a, b，中间用空格分隔（1 <= a <= b <= 10^11)。
输出
输出最小公倍数在这个区间的不同二元组的数量。
输入样例
4 6
输出样例
10

题解

由于是自己想的所以厚颜无耻打上原创

min25裸题（并不）
发现有lcm和10¹¹和函数前缀和，那么肯定是莫比乌斯反演 or min25
又发现反演不行，所以就是min25了
是不是很显然

---

发现答案的函数并不好求，所以设F(n)表示最小公倍数为n的不同二元组的数量（不考虑X和Y的大小关系）
可以发现若$n=\prod a{i}^{pi}$，则X和Y中一定有一个pi的指数为ai，每个pi就有2ai+1种分配方法（一个指数为0~ai，另一个为ai，两个都为ai的算重了一遍要减去）
那么$F(n)=\prod 2\ast ai+1$

min25三 要 素：
①$F(x)$是积性函数
显然
②$F(p)$能表示成多项式的形式（p是质数）
$F(p)=F(p^1)=3$
③$F(p^c)$能快速计算
$F(p^c)=2c+1$

最后为n的答案是(∑F(i) +n+1)/2（x=y=n只算了一次）
F(1)=1
没了

code

#include <algorithm>
#include <iostream>
#include <cstdlib>
#include <cstring>
#include <cstdio>
#include <cmath>
#define fo(a,b,c) for (a=b; a<=c; a++)
#define fd(a,b,c) for (a=b; a>=c; a--)
using namespace std;

long long b[316228];
int p[316228];
long long P[316228];
bool bz[316228];
long long g[316228];
long long G[316228];
int i,j,k,l,len;
long long A,B,N,s,S,I,I2;

void init(long long n,long long N)
{
	fo(i,1,N) g[i]=(i-1)*3,G[i]=(b[i]-1)*3;
	
	memset(bz,0,sizeof(bz));
	len=0;
	
	fo(I,2,N)
	{
		if (!bz[I])
		{
			p[++len]=I;
			P[len]=(long long)I*I;
			I2=I*I;
			
			fo(i,1,N)
			if (I2<=b[i])
			{
				s=b[i]/I;
				
				if (s>N)
				G[i]-=G[n/s]-g[I-1];
				else
				G[i]-=g[s]-g[I-1];
			}
			else
			break;
			fd(i,N,1)
			if (I2<=i)
			{
				s=i/I;
				g[i]-=g[s]-g[I-1];
			}
			else
			break;
		}
		
		fo(j,1,len)
		if ((long long)I*p[j]<=N)
		{
			bz[I*p[j]]=1;
			if (!(I%p[j])) break;
		}
		else
		break;
	}
}

long long F(long long p,int c)
{
	return c*2+1;
}

long long Work(long long n,long long _n,int j)
{
	long long ans,_n2;
	int i,e;
	
	if (_n>N)
	ans=G[n/_n];
	else
	ans=g[_n];
	ans-=(j-1)*3;
	
	fo(i,j,len)
	if (P[i]<=_n)
	{
		_n2=_n/p[i];
		e=1;
		
		while (_n2>=p[i])
		{
			ans+=(Work(n,_n2,i+1)+(e>1))*F(p[i],e);
			++e;
			_n2/=p[i];
		}
		ans+=(e>1)*F(p[i],e);
	}
	else
	break;
	
	return ans;
}

long long work(long long n)
{
	if (!n) return 0;
	
	N=floor(sqrt(n));
	fo(i,1,N) b[i]=n/i;
	init(n,N);
	
	return (Work(n,n,1)+1+n)/2;
}

int main()
{
//	freopen("51nod1222.in","r",stdin);
	
	scanf("%lld%lld",&A,&B);
	printf("%lld\n",work(B)-work(A-1));
}