[Luogu3768]简单的数学题

最新推荐文章于 2022-12-22 11:30:41 发布

DSL_HN_2002

最新推荐文章于 2022-12-22 11:30:41 发布

阅读量326

点赞数

CC 4.0 BY-SA版权

分类专栏：数论——杜教筛文章类型——题解

本文链接：https://blog.youkuaiyun.com/DSL_HN_2002/article/details/80325853

文章类型——题解同时被 2 个专栏收录

57 篇文章

订阅专栏

数论——杜教筛

2 篇文章

订阅专栏

本文探讨了给定整数n时，如何求解特定形式的数论求和问题，并提供了一种有效的算法实现，利用杜教筛法来优化计算过程。

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

Decription

给定 $n$ ，求

\sum_{i = 1}^{n} \sum_{j = 1}^{n} i j g c d (i, j)

$\sum_{i=1}^n\sum_{j=1}^nijgcd(i,j)$

Solution

\sum i = 1 n \sum j = 1 n i j g c d (i, j)

$\sum_{i=1}^n\sum_{j=1}^nijgcd(i,j)$

= \sum d 3 \sum i = 1 ⌊ n d ⌋ \sum j = 1 ⌊ n d ⌋ i j [g c d (i, j = 1)]

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

设 $S(n)=\sum_{i=1}^n\sum_{j=1}^nij[gcd(i,j=1)]$

则原式

= \sum d 3 G (⌊ n d ⌋)

$=\sum d^3G(\lfloor \frac{n}{d}\rfloor)$

S (n) = \sum i = 1 n \sum j = 1 n i j [g c d (i, j = 1)]

$S(n)=\sum_{i=1}^n\sum_{j=1}^nij[gcd(i,j=1)]$

= \sum i = 1 n i \sum j = 1 n j [g c d (i, j = 1)]

$=\sum_{i=1}^ni\sum_{j=1}^nj[gcd(i,j=1)]$

= 2 \sum i = 1 n i \sum j = 1 i j [g c d (i, j = 1)] - 1

$=2\sum_{i=1}^ni\sum_{j=1}^ij[gcd(i,j=1)]-1$

= 2 \sum i = 1 n i i ϕ ( i ) + [ i = 1 ] 2 - 1

$=2\sum_{i=1}^ni\frac{i\phi(i)+[i=1]}{2}-1$

= \sum i = 1 n i 2 i ϕ (i)

$=\sum_{i=1}^ni^2i\phi(i)$

设 $f(i)=i^2\phi(i),g(i)=i^2$ ,则 $S(n)=\sum_{i=1}^nf(i)$

最后用杜教筛求出 $S(n)$

S (n) = \sum i = 1 n (f * g) (i) - \sum i = 2 n g (i) S (i) = \sum i = 1 n i 3 - \sum i = 2 n g (i) S (i) = (n * ( n + 1 ) 2) 2 - \sum i = 2 n g (i) S (i)

$S(n)=\sum_{i=1}^n(f*g)(i)-\sum_{i=2}^ng(i)S(i)=\sum_{i=1}^ni^3-\sum_{i=2}^ng(i)S(i)=(\frac{n*(n+1)}{2})^2-\sum_{i=2}^ng(i)S(i)$

#include <bits/stdc++.h>
using namespace std;

const int maxn = 5000005, Max = 5000000;
typedef long long lint;

int p, inv1, inv2;
lint n;

inline int pow(int x, int k)
{
    int ret = 1;
    while (k) {
        if (k & 1) ret = (lint)ret * x % p;
        x = (lint)x * x % p; k >>= 1;
    }
    return ret;
}

int s[maxn], phi[maxn], pr[maxn], pcnt;
bool isp[maxn];
void prepare(int n)
{
    memset(isp, -1, sizeof(isp));
    isp[1] = 0; phi[1] = 1;
    for (int i = 2; i <= n; ++i) {
        if (isp[i]) {phi[i] = i - 1; pr[++pcnt] = i;}
        for (int j = 1; i * pr[j] <= n && j <= pcnt; ++j) {
            isp[i * pr[j]] = 0;
            if (i % pr[j]) phi[i * pr[j]] = (lint)phi[i] * phi[pr[j]] % p;
            else phi[i * pr[j]] = (lint)phi[i] * pr[j] % p;
        }
    }
    for (int i = 1; i <= n; ++i) {
        s[i] = (lint)i * i % p * phi[i] % p;
        s[i] += s[i - 1];
        if (s[i] >= p) s[i] -= p;
    }
}

bool vis[maxn];
lint S[maxn];

inline lint calc1(lint x)
{
    x %= p;
    return (lint)x * (x + 1) % p * (2 * x + 1) % p * inv1 % p;
}

inline lint calc2(lint x)
{
    x %= p;
    x = x * (x + 1) % p * inv2 % p;
    return x * x % p;
}

lint dfs(lint x)
{
    if (x < Max) return s[x];
    int m = n / x;
    if (vis[m]) return S[m];
    vis[m] = true;
    lint sum = calc2(x);
    for (lint i = 2, j; i <= x; i = j + 1) {
        j = x / (x / i);
        (sum -= (calc1(j) - calc1(i - 1)) * dfs(x / j) % p) %= p;
    }
    if (sum < 0) sum += p;
    return S[m] = sum;
}

int main()
{
    scanf("%d%lld", &p, &n); inv1 = pow(6, p - 2); inv2 = pow(2, p - 2);

    prepare(Max);
    dfs(n);

    lint ans = 0;
    for (lint i = 1, j; i <= n; i = j + 1) {
        j = n / (n / i);
        if (n / i <= Max) ans += (lint)(calc2(j) - calc2(i - 1) + p) % p * s[n / i] % p;
        else ans += (lint)(calc2(j) - calc2(i - 1) + p) % p * S[i] % p;
        if (ans >= p) ans -= p;
    }

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

    return 0;
}