[ACM]【LCM/唯一分解定理】牛客算法周周练5 求幂

求幂

传送门
题意：寻找有多少对$a$、$b$、$c$、$d\in [1,n]$，满足$a^b=c^d$。

思路：

看似简单的题目往往都暗藏玄机。
本思路参考题解。
主要操作：
第一步，用唯一分解定理将原式简化为$(i^j{)}^x=(i^k{)}^y$。
第二步，按照范围枚举$i$、$j$、$k$。
第三步，找到$i^j$和$i^k$的最小公倍数，求得满足$(i^j{)}^x=(i^k{)}^y$且$x$、$y\in [1,n]$的$x$个数（ie.$y$的个数）。
具体原理：
很容易想到，一开始要将$A$和$C$用唯一分解定理表示。由于他们的某个幂相等，易得他们分解出来的质数是一一对应的。即$A=q_1^{a_1}+q_2^{a_2}+...+q_p^{a_p}$，$C=q_1^{b_1}+q_2^{b_2}+...+q_p^{b_p}$。并且有$a_i\times B=b_i\times D,i\in [1,p]$。我们设最简分数$\frac{m}{n}=\frac{d}{b}=\frac{a_i}{b_i}$，设$i=q_1^{\frac{a_1}{m}}+q_2^{\frac{a_2}{m}}+...+q_p^{\frac{a_p}{m}}=q_1^{\frac{b_1}{n}}+q_2^{\frac{b_2}{n}}+...+q_p^{\frac{b_p}{n}}$称为 “公共基底数”，原式就等价于$(i^m{)}^b=(i^n{)}^d$。
为了直观，把$m$换成$j$，$n$换成$k$，$b$换成$x$，$d$换成$y$。根据条件，要求$i^j$和$i^k$小于等于$n$，$x$、$y$小于等于$n$。
$(i^j{)}^x=(i^k{)}^y$意味着$j\times x=k\times y$，我们首先遍历$i$、$j$、$k$，然后要找到满足这个关系的$x$和$y$的对数。很好找，只需要求$j$与$k$的最小公倍数，在最小公倍数的基础上乘正整数，直到$x$和$y$任意一个大于$n$。设正整数为$q$，$q$的最大值就等于满足这个关系的$x$和$y$的对数。$q\times lcm(j,k)=j\times x=k\times y\le j\times n,k\times n$，于是，$q\le min(\frac{j\times n}{lcm},\frac{k\times n}{lcm})$。
实现：
有三个地方需要注意：
（1）$i=1$要特判，不然幂次运算陷入死循环。
（2）用vis数组判重，减少重复计算，寻找$i^j$和$i^k$的上限的时候，把经过的都标记一下就行了。
（3）如果直接用$j\times n$来算，会爆，又因为\符号自带floor功能，不能先除lcm再乘。于是选择直接拿gcd，$\frac{j\times n}{lcm}=\frac{gcd\times j\times n}{j\times k}=\frac{gcd\times n}{k}=\frac{n}{\frac{k}{gcd}}$。众所周知，$\frac{k}{gcd}$为一整数。

代码：

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=1000000007;
bool vis[1000004];
inline int read(){
	int k=0,j=1;char ch=getchar();
	while(ch<'0'||ch>'9'){if(ch=='-') j=-1;ch=getchar();}
	while(ch<='9'&&ch>='0'){k=(k<<1)+(k<<3)+ch-'0';ch=getchar();}
	return k*j;
}
ll gcd(ll a,ll b){
	if(b==0) return a;
	else return gcd(b,a%b);
}
int main(){
	ll n=(long long)read();
	ll ans=(n*n)%mod;//特判
	//开始遍历i
	for(ll i=2;i<=n;i++){
		if(vis[i]) continue;
		ll maxn=1;//不超过n的最大的幂次
		ll now=i;//幂次后的数
		while(now<=n/i){
			maxn++;
			now*=i;
			vis[now]=1;//确定可达，才标记
		}
		for(ll j=1;j<=maxn;j++){
			for(ll k=1;k<=maxn;k++){
				ll g=gcd(j,k);
				ll tmp=k/g,tmp2=j/g;
				ans=(ans+min(n/tmp,n/tmp2))%mod;
			}
		}
	}
	printf("%lld\n",ans);
}