51 Nod 1238 最小公倍数之和 V3 杜教筛

时间:2021-07-10 01:07:13

题目链接:http://www.51nod.com/Challenge/Problem.html#!#problemId=1238

题意:求$\sum_{i=1}^{n}\sum_{j=1}^{n}lcm(i,j)=\sum_{i=1}^{n}\sum_{j=1}^{n}{\frac{i*j}{gcd(i,j)}}$,$1\leq{n}\leq10^{10}$.

知识提要:小于等于n中与n互质的数总和为$\sum_{i=1}^{n}[(n,i)=1]i=\frac{\varphi(n)*n+[n=1]}{2}$

解析:

枚举最大公约数d,

$$Ans=\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}[(i,j)=1]i*j$$

我们先考虑 j<=i 的情况,

$$\quad\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{i}[(i,j)=1]i*j\\$$

$$=\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\frac{\varphi(i)*i+[i=1]}{2}*i$$

还有i<=j的情况没考虑,其实两者是对称的 ,上面的式子乘2就好了,然后(1,1)这一对多算了一次了,所以-1就好了,

$$Ans=\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\varphi(i)*i^2$$

令$F(n)=\sum_{i=1}^{n}\varphi(i)*i^2$

$$Ans=\sum_{d=1}^{n}d*F(\lfloor\frac{n}{d}\rfloor)$$

欧拉函数的前缀和$\phi(n)$之前博客里写过 按照类似的方法可以推出来

$$F(n)=\frac{n^2*(n+1)^2}{4}-\sum_{i=2}^{n}\sum_{j=1}^{\lfloor\frac{n}{i}\rfloor}\varphi(j)*i^2*j^2\\$$
$$=\frac{n^2*(n+1)^2}{4}-\sum_{i=2}^{n}i^2\sum_{j=1}^{\lfloor\frac{n}{i}\rfloor}\varphi(j)*j^2\\$$
$$=\frac{n^2*(n+1)^2}{4}-\sum_{i=2}^{n}i^2F(\lfloor\frac{n}{i}\rfloor)$$

到此为止可以$O(n^{\frac{2}{3}})$求出Ans

AC代码

#include <bits/stdc++.h>
#define pb push_back
#define mp make_pair
#define fi first
#define se second
#define all(a) (a).begin(), (a).end()
#define fillchar(a, x) memset(a, x, sizeof(a))
#define huan printf("\n");
#define debug(a,b) cout<<a<<" "<<b<<" ";
using namespace std;
const int maxn=1e6+,inf=0x3f3f3f3f;
typedef long long ll;
const ll mod = ;
typedef pair<int,int> pii;
int check[maxn],prime[maxn],phi[maxn],sum[maxn];
void Phi(int N)//线性筛
{
int pos=;sum[]=;
sum[]=phi[]=;
for(ll i = ; i <= N ; i++)
{
if (!check[i])
prime[pos++] = i,phi[i]=i-;
for (int j = ; j < pos && i*prime[j] <= N ; j++)
{
check[i*prime[j]] = ;
if (i % prime[j] == )
{
phi[i*prime[j]]=phi[i]*prime[j];
break;
}
else
phi[i*prime[j]]=phi[i]*(prime[j]-);
}
sum[i]=(sum[i-]+(phi[i]*i%mod)*i%mod)%mod;
}
}
unordered_map<ll,ll> ma;
ll inv2=;
ll inv4=;
ll inv6=;
ll solve(ll n)
{
if(n<=1e6)
return sum[n];
else if(ma.count(n))
return ma[n];
ll temp = ((n%mod)*((n+)%mod)%mod)*inv2%mod;
temp=temp*temp%mod;
for(ll i=,j;i<=n;i=j+)
{
j=n/(n/i);
ll r=(((j%mod)*((j+)%mod)%mod)*((*j+)%mod)%mod)*inv6%mod;
ll l=(((i%mod)*((i-)%mod)%mod)*((*i-)%mod)%mod)*inv6%mod;
r=(r-l+mod)%mod;
temp = (temp-solve(n/i)*r%mod+mod)%mod;
}
return ma[n]=temp;
}
int main()
{
ll n;
Phi(1e6);
scanf("%lld",&n);
ll ans=;
for(ll i=,j;i<=n;i=j+)
{
j=n/(n/i);
ll r=((j%mod)*((j+)%mod)%mod)*inv2%mod;
ll l=((i%mod)*((i-)%mod)%mod)*inv2%mod;
r=(r-l+mod)%mod;
ans=(ans+solve(n/i)*r%mod)%mod;
}
printf("%lld\n",ans);
}