[51NOD1237]最大公约数之和 V3

最新推荐文章于 2021-09-14 20:00:00 发布

a_crazy_czy

最新推荐文章于 2021-09-14 20:00:00 发布

阅读量1.1k

点赞数 1

CC 4.0 BY-SA版权

分类专栏：筛法分块 51NOD 文章标签： OI 杜教筛分块数论 51nod

本文链接：https://blog.youkuaiyun.com/a_crazy_czy/article/details/57149019

筛法同时被 3 个专栏收录

11 篇文章

订阅专栏

分块

10 篇文章

订阅专栏

51NOD

9 篇文章

订阅专栏

本文介绍了一种高效算法，用于解决特定的数论问题——在给定的大整数范围内求解两数最大公约数（GCD）之和，并通过数学转换简化问题，最终采用杜教筛算法进行前缀和的快速计算。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

题目大意

给定 $n$ ，试求

\sum i = 1 n \sum j = 1 n gcd (i, j)

$\sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)$
结果对

109+7 $10^9+7$ 取模。

$2\le n\le10^{10}$

题目分析

我们将题目改为求

\sum i = 1 n \sum j = 1 i gcd (i, j)

$\sum_{i=1}^n\sum_{j=1}^i\gcd(i,j)$
然后将答案乘二再减去

1 $1$ 到

n $n$ 的和即可。
那么上面这条式子是什么呢？

\sum i = 1 n \sum j = 1 i gcd (i, j) = \sum d = 1 n d \sum i = 1 ⌊ n d ⌋ φ (i) = \sum d = 1 n φ (d) ( 1 + ⌊ n d ⌋ ) ⌊ n d ⌋ 2

$\begin{align} \sum_{i=1}^n\sum_{j=1}^i\gcd(i,j)&=\sum_{d=1}^nd\sum_{i=1}^{\lfloor\frac nd\rfloor}\varphi(i)\\ &=\sum_{d=1}^n\varphi(d){\left(1+\lfloor\frac nd\rfloor\right)\lfloor\frac nd\rfloor\over2} \end{align}$
对

⌊nd⌋ $\lfloor\frac nd\rfloor$ 分块，对前面的欧拉函数使用杜教筛求前缀和。
时间复杂度？这个算起来有点麻烦，我来口胡一下。
假设我们预处理的部分大小为

a $a$

\sum i = 1 n \sqrt \sum j = 1 ⌊ i a ⌋ O ⎛ ⎝ i j - - \sqrt ⎞ ⎠ + \sum i = 1 n \sqrt \sum j = 1 ⌊ n i a ⌋ O (n i j - - \sqrt)

$\sum_{i=1}^\sqrt n\sum_{j=1}^{\lfloor\frac ia\rfloor}\mathrm O\left(\sqrt{\frac ij}\right)+\sum_{i=1}^\sqrt n\sum_{j=1}^{\lfloor\frac n{ia}\rfloor}\mathrm O\left(\sqrt\frac n{ij}\right)$
显然右边部分的复杂度要大于左边，使用积分来计算：

\sum i = 1 n \sqrt \sum j = 1 ⌊ n i a ⌋ O (n i j - - \sqrt) \approx O (\int n \sqrt 0 (\int n a x 0 n x y - - - \sqrt d y) d x) = O (\int n \sqrt 0 n x 1 a - - \sqrt) = O (1 a - - \sqrt n log n)

$\begin{align} \sum_{i=1}^\sqrt n\sum_{j=1}^{\lfloor\frac n{ia}\rfloor}\mathrm O\left(\sqrt\frac n{ij}\right)&\approx\mathrm O\left(\int_0^\sqrt n\left(\int_0^\frac n{ax}\sqrt\frac n{xy}\ \mathrm dy\right)\ \mathrm dx\right)\\ &=\mathrm O\left(\int_0^\sqrt n\frac nx\sqrt\frac1a\right)\\ &=\mathrm O\left(\sqrt\frac1an\log n\right) \end{align}$
然后

a $a$ 可以在空间限制内尽量取大一点，我取了

1.5×107 $1.5\times10^7$ 。

代码实现

我的代码莫名其妙常数巨大，卡着空间取的让我卡着时限过了这题。。。

#include <iostream>
#include <cstdio>
#include <cmath>

using namespace std;

typedef long long LL;

const int P=1000000007;
const int itwo=500000004;
const int L=15000000;
const int N=100050;

int pri[L+5],phi[L+5],f[L+5],sum[L+5];
int S[N],vis[N];
bool mark[L+5];
int ans,l,s,cnt;
LL n;

void pre()
{
    s=trunc(sqrt(n)),l=L;
    mark[1]=1,sum[1]=phi[1]=1,f[1]=1;
    for (int i=2;i<=l;++i)
    {
        if (!mark[i]) pri[++pri[0]]=i,phi[i]=i-1,f[i]=i;
        for (int j=1,k;j<=pri[0];++j)
        {
            if (1ll*pri[j]*i>l) break;
            mark[k=pri[j]*i]=1,f[k]=pri[j];
            phi[k]=phi[i]*(pri[j]-(f[i]!=f[k]));
            if (!(i%pri[j])) break;
        }
        sum[i]=(sum[i-1]+phi[i])%P;
    }
}

int id(LL x){return n/x;}

int sieve(LL n)
{
    if (n<=l) return sum[n];
    int idn=id(n);
    if (vis[idn]==cnt) return S[idn];
    int res=1ll*((n+1)%P)*(n%P)%P*itwo%P;
    for (LL st=2,en,x;st<n;st=en)
    {
        x=n/st,en=n/x+1;
        res=(res-1ll*sieve(x)*(en-st)%P+P)%P;
    }
    return vis[idn]=cnt,S[idn]=res;
}

int main()
{
    freopen("gcdsumIII.in","r",stdin),freopen("gcdsumIII.out","w",stdout);
    scanf("%lld",&n),pre();
    ans=0;
    for (LL st=1,en,x;st<n;st=en)
    {
        x=n/st,en=n/x+1;
        ans=(ans+1ll*((x+1)%P)%P*(x%P)%P*itwo%P*(((++cnt,sieve(en-1))-(++cnt,sieve(st-1))+P)%P)%P)%P;
    }
    ans=(ans*2%P-1ll*((n+1)%P)*(n%P)%P*itwo%P+P)%P;
    printf("%d\n",ans);
    fclose(stdin),fclose(stdout);
    return 0;
}