「洛谷 3768」简单的数学题

传送门

---

problem

给定 $n$ 和 $p$，求：

$$\sum _{i=1}^n\sum _{j=1}^{n}ij\gcd (i,j)$$

答案对 $p$ 取模。

数据范围：$n\le 10^{10}$，$5\times 10^8\le p\le 1.1\times 10^9$。

---

solution

枚举 $\gcd$：

$$\begin{array}{rl}ans& =\sum _{k=1}^{n}k\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{k}\rfloor }(ik)(jk)[\gcd (i,j)=1]\\ & =\sum _{k=1}^n{k}^3\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{k}\rfloor }ij[\gcd (i,j)=1]\end{array}$$

这里其实可以直接反演了，但我们有一种更巧妙的做法。

把后面那一坨直接拿出来，即：

$$\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{k}\rfloor }ij[\gcd (i,j)=1]$$

把它记为 $res$ 吧，那么：

$$\begin{array}{rl}res& =2\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{i}ij[\gcd (i,j)=1]-\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }{i}^2[\gcd (i,i)=1]\\ & =2\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }i\left(\sum _{j=1}^{i}j[\gcd (i,j)=1]\right)-1\end{array}$$

由于 ${\sum }_{i=1}^{n}i[\gcd (i,n)=1]=\frac{n\phi (n)+[n=1]}{2}$，证明可以看这里，那么有：

$$\begin{array}{rl}res& =2\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }i\frac{i\phi (i)+[i=1]}{2}-1\\ & =\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }{i}^2\phi (i)\end{array}$$

那么回到我们原来的式子：

$$ans=\sum _{k=1}^n{k}^3\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }{i}^2\phi (i)$$

其实这剩下的部分就和这道题的末尾部分很像了，我就不写了。

---

code

#include<bits/stdc++.h>
#define N 5000005
#define ll long long
using namespace std;
ll n,P,inv2,inv4,inv6,f[N];
int sum,prime[N],mark[N],phi[N];
ll power(ll a,ll b,ll ans=1){
	for(;b;b>>=1,a=a*a%P)
		if(b&1)  ans=ans*a%P;
	return ans;
}
void linear_sieves(){
	f[1]=1;
	for(int i=2;i<N;++i){
		if(!mark[i])  prime[++sum]=i,phi[i]=i-1;
		for(int j=1;j<=sum&&i*prime[j]<N;++j){
			mark[i*prime[j]]=1;
			if(i%prime[j]==0)  {phi[i*prime[j]]=phi[i]*prime[j];break;}
			phi[i*prime[j]]=phi[i]*(prime[j]-1);
		}
		f[i]=(ll)i*i%P*phi[i]%P,f[i]=(f[i]+f[i-1])%P;
	}
	inv2=power(2,P-2),inv4=power(4,P-2),inv6=power(6,P-2);
}
ll Squ(ll x)  {return x%P*x%P;}
ll sum2(ll n)  {return n%P*((n+1)%P)%P*((2*n%P+1)%P)%P*inv6%P;}
ll sum3(ll n)  {return Squ(n%P)*Squ((n+1)%P)%P*inv4%P;}
unordered_map<ll,ll>Sum;
ll F(ll n){
	if(n<N)  return f[n];
	if(Sum[n])  return Sum[n];
	ll ans=sum3(n);
	for(ll i=2,j;i<=n;i=j+1){
		j=n/(n/i);
		ans=(ans-((sum2(j)-sum2(i-1)+P)%P*F(n/i)%P)+P)%P;
	}
	return Sum[n]=ans;
}
int main(){
	int ans=0;
	scanf("%lld%lld",&P,&n),linear_sieves();
	for(ll i=1,j;i<=n;i=j+1){
		j=n/(n/i);
		ans=(ans+(sum3(j)%P-sum3(i-1)%P+P)%P*F(n/i)%P)%P;
	}
	printf("%d",ans);
	return 0;
}