description

求

\[\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=1][\frac ij在k进制下是纯循环小数] \]

data range

\(n,m\le 10^9,k\le 2000\)

solution

\(\frac ij\)在\(k\)进制下为纯循环小数当且仅当\(\gcd(j,k)=1\)

证明戳这里

那么原式

\[=\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=1][\gcd(j,k)=1]\\ =\sum_{j=1}^m[\gcd(j,k)=1]\sum_{i=1}^n[\gcd(i,j)=1]\\ =\sum_{j=1}^m[\gcd(j,k)=1]\sum_{i=1}^n\sum_{d\mid i,d\mid j}\mu(d)\\ =\sum_{d=1}^n\mu(d)\sum_{d\mid j}^m[\gcd(j,k)=1]\sum_{d\mid i}^n1\\ =\sum_{d=1}^n\mu(d)\sum_{j=1}^{\lfloor \frac md\rfloor}[\gcd(jd,k)=1]\lfloor \frac nd\rfloor\\ =\sum_{d=1}^n\lfloor \frac nd\rfloor\mu(d)[\gcd(d,k)=1]\sum_{j=1}^{\lfloor \frac md\rfloor}[\gcd(j,k)=1] \]

不妨设

\[F(n)=\sum_{i=1}^n[\gcd(i,k)=1]\\ S(n,k)=\sum_{i=1}^n\mu(i)[\gcd(i,k)=1]\\ \]

那么原式

\[=\sum_{d=1}^n\lfloor \frac nd\rfloor F(\lfloor \frac md\rfloor)(S(d,k)-S(d-1,k)) \]

显然如果可以快速求出\(S(n,k),F(n)\) ，那么通过整除分块就可以快速求得原式了

先来康康\(F(n)\),我们有

\[F(n)=\lfloor\frac nk\rfloor F(k)+F(n\bmod k) \]

只用预处理处\(n\le k\)的\(F(n)\)就可以做到快速查询了

再来康康\(S(n,k)\)

\[S(n,k)=\sum_{i=1}^n\mu(i)[\gcd(i,k)=1]\\ =\sum_{i=1}^n\mu(i)\sum_{d\mid i,d\mid k}\mu(d)\\ =\sum_{d\mid k}\mu(d)\sum_{d\mid i}^n\mu(i)\\ =\sum_{d\mid k}\mu(d)\sum_{i=1}^{\lfloor \frac nd\rfloor}\mu(id) \]

显然若\(\gcd(i,d)>1\) 则\(\mu(id)=0\)

因此

\[=\sum_{d\mid k}\mu(d)\sum_{i=1}^{\lfloor \frac nd\rfloor}\mu(id)[\gcd(i,d)=1]\\ =\sum_{d\mid k}\mu(d)\sum_{i=1}^{\lfloor \frac nd\rfloor}\mu(i)\mu(d)[\gcd(i,d)=1]\\ =\sum_{d\mid k}\mu(d)^2\sum_{i=1}^{\lfloor \frac nd\rfloor}\mu(i)[\gcd(i,d)=1]\\ =\sum_{d\mid k}\mu(d)^2S(\lfloor \frac nd\rfloor,d) \]

递归求解即可（记得记忆化）

注意递归边界

\(n=0\) 时,\(S(0,k)=0\)

\(k=1\) 时,\(S(n,1)=\sum_{i=1}^n\mu(i)[\gcd(1,k)=1]=\sum_{i=1}^n\mu(i)\)

杜教筛之即可

完结撒花。。。

time complexity

\(\mathcal O(能过)\)

code

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1e6+5,Bas=2001;
bool flag[N];int pr[N],pcnt,mu[N],smu[N],f[N];
unordered_map<ll,ll>mp;
inline ll qs(int a,int b){return 1ll*a*Bas+b;}
ll ss(int n,int k)
{
	if(!n)return 0;
	ll now=qs(n,k); 
	if(mp.count(now))return mp[now];
	if(k==1)
	{
		if(n<=1e6)return smu[n];
		ll ret=1;
		for(int l=2,r;l<=n;l=r+1)
		{
			r=n/(n/l);
			ret-=ss(n/l,1)*(r-l+1);
		}return mp[now]=ret;
	}ll ret=0;
	for(int i=1;i*i<=k;++i)
	{
		if(k%i)continue;
		if(mu[i])ret+=ss(n/i,i);
		if(i*i!=k&&mu[k/i])ret+=ss(n/(k/i),k/i);
	}return mp[now]=ret;
}
int k;
int gcd(int x,int y){return y?gcd(y,x%y):x;}
inline void pre(int n)
{
	mu[1]=1;
	for(int i=2;i<=n;++i)
	{
		if(!flag[i])pr[++pcnt]=i,mu[i]=-1;
		for(int j=1;j<=pcnt;++j)
		{
			int num=i*pr[j];if(num>n)break;
			flag[num]=1;
			if(i%pr[j])mu[num]=-mu[i];
			else{mu[num]=0;break;}
		}
	}
	for(int i=1;i<=n;++i)smu[i]=smu[i-1]+mu[i];
	for(int i=1;i<=k;++i)f[i]=f[i-1]+(gcd(i,k)==1);
}
inline int sf(int n){return n/k*f[k]+f[n%k];}
int n,m;
int main()
{
	scanf("%d%d%d",&n,&m,&k);pre(1e6);
	int mn=min(n,m);ll ans=0;
	for(int l=1,r;l<=mn;l=r+1)
	{
		r=min(n/(n/l),m/(m/l));
		ans+=(ss(r,k)-ss(l-1,k))*(n/l)*sf(m/l);
	}printf("%lld\n",ans);
	return 0;
}