51Nod 1238 - 最小公倍数之和 V3（毒瘤数学+杜教筛）

题目

戳这里

推导

∑i=1n∑j=1nlcm(i,j)~~~\sum_{i=1}^{n}\sum_{j=1}^{n}lcm(i,j)   ∑i=1n​∑j=1n​lcm(i,j)

=∑i=1n∑j=1nijgcd(i,j)=\sum_{i=1}^{n}\sum_{j=1}^{n}\frac{ij}{gcd(i,j)}=∑i=1n​∑j=1n​gcd(i,j)ij​

=∑i=1nd−1∑i=1n∑j=1nij[gcd(i,j)==d]=\sum_{i=1}^{n}d^{-1}\sum_{i=1}^{n}\sum_{j=1}^{n}ij[gcd(i,j)==d]=∑i=1n​d−1∑i=1n​∑j=1n​ij[gcd(i,j)==d]

=∑i=1nd∑i=1⌊nd⌋∑j=1⌊nd⌋ij[gcd(i,j)==1]=\sum_{i=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}ij[gcd(i,j)==1]=∑i=1n​d∑i=1⌊dn​⌋​∑j=1⌊dn​⌋​ij[gcd(i,j)==1]

=∑i=1nd∑i=1⌊nd⌋i∑j=1⌊nd⌋j[gcd(i,j)==1]=\sum_{i=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}i\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}j[gcd(i,j)==1]=∑i=1n​d∑i=1⌊dn​⌋​i∑j=1⌊dn​⌋​j[gcd(i,j)==1]

=∑i=1nd(2∑i=1⌊nd⌋i∑j=1ij[gcd(i,j)==1]−1)=\sum_{i=1}^{n}d(2\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}i\sum_{j=1}^{i}j[gcd(i,j)==1]-1)=∑i=1n​d(2∑i=1⌊dn​⌋​i∑j=1i​j[gcd(i,j)==1]−1)

=∑i=1nd(2∑i=1⌊nd⌋iiφ(i)+[i==1]2−1)=\sum_{i=1}^{n}d(2\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}i\frac{i\varphi(i)+[i==1]}{2}-1)=∑i=1n​d(2∑i=1⌊dn​⌋​i2iφ(i)+[i==1]​−1)

=∑i=1nd∑i=1⌊nd⌋i2φ(i)=\sum_{i=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}i^2\varphi(i)=∑i=1n​d∑i=1⌊dn​⌋​i2φ(i)

子问题：

求∑i=1ni2φ(i)\sum_{i=1}^{n}i^2\varphi(i)∑i=1n​i2φ(i)

令f(i)=i2φ(i)f(i)=i^2\varphi(i)f(i)=i2φ(i)

使用狄利克雷卷积，卷一个g(i)=i2g(i)=i^2g(i)=i2

那么：

∑i=1n(f∗g)(i)~~~~\sum_{i=1}^{n}(f*g)(i)    ∑i=1n​(f∗g)(i)

=∑i=1n∑d∣if(d)g(id)=\sum_{i=1}^{n}\sum_{d|i}^{}f(d)g(\frac{i}{d})=∑i=1n​∑d∣i​f(d)g(di​)

=∑i=1n∑d∣id2φ(d)(id)2=\sum_{i=1}^{n}\sum_{d|i}^{}d^2\varphi(d)(\frac{i}{d})^2=∑i=1n​∑d∣i​d2φ(d)(di​)2

=∑i=1ni2∑d∣iφ(d)=\sum_{i=1}^{n}i^2\sum_{d|i}^{}\varphi(d)=∑i=1n​i2∑d∣i​φ(d)

=∑i=1ni3=\sum_{i=1}^{n}i^3=∑i=1n​i3

=n2(n+1)24=\frac{n^2(n+1)^2}{4}=4n2(n+1)2​

又因为：

∑i=1n∑d∣id2φ(d)(id)2~~~~\sum_{i=1}^{n}\sum_{d|i}^{}d^2\varphi(d)(\frac{i}{d})^2    ∑i=1n​∑d∣i​d2φ(d)(di​)2

=∑i=1ni2∑d=1⌊ni⌋d2φ(d)=\sum_{i=1}^{n}i^2\sum_{d=1}^{\lfloor\frac{n}{i}\rfloor}d^2\varphi(d)=∑i=1n​i2∑d=1⌊in​⌋​d2φ(d)

=∑i=2ni2∑d=1⌊ni⌋d2φ(d)+∑i=1ni2φ(i)=\sum_{i=2}^{n}i^2\sum_{d=1}^{\lfloor\frac{n}{i}\rfloor}d^2\varphi(d)+\sum_{i=1}^{n}i^2\varphi(i)=∑i=2n​i2∑d=1⌊in​⌋​d2φ(d)+∑i=1n​i2φ(i)

=n2(n+1)24=\frac{n^2(n+1)^2}{4}=4n2(n+1)2​

所以：

∑i=1ni2φ(i)=n2(n+1)24−∑i=2ni2∑d=1⌊ni⌋d2φ(d)\sum_{i=1}^{n}i^2\varphi(i)=\frac{n^2(n+1)^2}{4}-\sum_{i=2}^{n}i^2\sum_{d=1}^{\lfloor\frac{n}{i}\rfloor}d^2\varphi(d)∑i=1n​i2φ(i)=4n2(n+1)2​−∑i=2n​i2∑d=1⌊in​⌋​d2φ(d)

使用杜教筛将时间复杂度降到O(n23)O(n^{\frac{2}{3}})O(n32​)

数学太难了QAQ

代码：

#include<cstdio>

#include<cmath>

#include<cstring>

#include<map>

#include<algorithm>

#define maxn 5000000

#define INF 0x3f3f3f3f

#define MOD 1000000007

#define two 500000004

#define six 166666668

using namespace std;

inline long long getint()

{

	long long num=0,flag=1;char c;

	while((c=getchar())<'0'||c>'9')if(c=='-')flag=-1;

	while(c>='0'&&c<='9')num=num*10+c-48,c=getchar();

	return num*flag;

}

long long n;

bool not_prime[maxn+5];

int prime[maxn+5],cnt;

long long phi[maxn+5];

map<long long,long long>M;

inline void init()

{

	phi[1]=1;

	for(int i=2;i<=maxn;i++)

	{

		if(!not_prime[i])prime[++cnt]=i,phi[i]=i-1;

		for(int j=1;j<=cnt&&i*prime[j]<=maxn;j++)

		{

			not_prime[i*prime[j]]=1;

			if(i%prime[j])phi[i*prime[j]]=phi[i]*phi[prime[j]];

			else{phi[i*prime[j]]=phi[i]*prime[j];break;}

		}

	}

	for(int i=1;i<=maxn;i++)(phi[i]*=1ll*i*i%MOD)%=MOD;

	for(int i=1;i<=maxn;i++)(phi[i]+=phi[i-1])%=MOD;

}

inline long long getsqr(long long x)

{return x%MOD*((x+1)%MOD)%MOD*((2*x+1)%MOD)%MOD*six%MOD;}

inline long long solve(long long x)

{

	if(x<=maxn)return phi[x];

	if(M.count(x))return M[x];

	long long sum=x%MOD*((x+1)%MOD)%MOD*two%MOD;

	(sum*=sum)%=MOD;

	for(long long i=2,j;i<=x;i=j+1)

	{

		j=x/(x/i);

		(sum-=(getsqr(j)-getsqr(i-1))%MOD*solve(x/i)%MOD)%=MOD;

		(sum+=MOD)%=MOD;

	}

	return M[x]=sum;

}

int main()

{

	init();

	n=getint();

	long long sum=0;

	for(long long i=1,j;i<=n;i=j+1)

	{

		j=n/(n/i);

		(sum+=1ll*(j+i)%MOD*(j-i+1)%MOD*two%MOD*solve(n/i)%MOD)%=MOD;

	}

	printf("%lld\n",sum);

}