Processing math: 0%

洛谷 3768 简单的数学题

发表于 | 分类于 |

zzq 的题啊……

直接开推 n∑i=1n∑j=1ijgcd

令 T = dk，更改枚举方式，有 \begin{aligned} & \sum\limits_{d = 1}^{n} d^3 \sum\limits_{k = 1}^{\lfloor \frac n d \rfloor} \mu(k) k^2 (\sum\limits_{idk \le n} i)^2 \\ = & \sum\limits_{T = 1}^{n} (\sum\limits_{iT \le n} i)^2 \sum\limits_{d | T} d^3 \mu(\dfrac T d) \dfrac {T^2} {d^2} \\ = & \sum\limits_{T = 1}^{n} (\sum\limits_{iT \le n} i)^2 T^2 \sum\limits_{d | T} d \mu(\dfrac T d) \\ = & \sum\limits_{T = 1}^{n} (\sum\limits_{iT \le n} i)^2 T^2 \varphi(T) \end{aligned}

杜教筛时，考虑设 f(i) = i^2 \varphi(i),g(i) = i^2,h(i) = (f * g)(i) = \sum\limits_{d | i} d^2 \varphi(d) \dfrac {i^2} {d^2} = d^3。
同时要知道 \sum\limits_{i = 1}^n i^2 = \dfrac {n (n + 1) (2n + 1)} 6,\sum\limits_{i = 1}^n i^3 = \dfrac {n^2 (n + 1)^2} 4。

代码： 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
#include <cstdio>
#include <tr1/unordered_map>
using namespace std;
using namespace tr1;
const long long N = 1e10;
const int CNT = 1e7;
long long n,mod,inv2,inv6;
int vis[CNT + 5],prime[CNT + 5],cnt;
long long phi[CNT + 5];
unordered_map<long long,long long> w;
inline long long fpow(long long a,long long b)
{
    long long ret = 1;
    for(;b;b >>= 1)
        (b & 1) && (ret = ret * a % mod),a = a * a % mod;
    return ret;
}
inline long long sum(long long n)
{
    n %= mod;
    return n * (n + 1) % mod * inv2 % mod;
}
inline long long sum1(long long n)
{
    n %= mod;
    return n * (n + 1) % mod * (2 * n + 1) % mod * inv6 % mod;
}
long long calc(long long n)
{
    if(n <= CNT)
        return phi[n];
    if(w.count(n))
        return w[n];
    long long ret = sum(n) * sum(n) % mod;
    for(long long l = 2,r;l <= n;l = r + 1)
    {
        r = n / (n / l);
        ret = (ret - (sum1(r) - sum1(l - 1)) * calc(n / l) % mod) % mod;
    }
    return w[n] = (ret + mod) % mod;
}
long long ans;
int main()
{
    scanf("%lld%lld",&mod,&n),inv2 = fpow(2,mod - 2),inv6 = fpow(6,mod - 2);
    phi[1] = 1;
    for(register int i = 2;i <= CNT;++i)
    {
        if(!vis[i])
            phi[prime[++cnt] = i] = i - 1;
        for(register int j = 1;j <= cnt && i * prime[j] <= CNT;++j)
        {
            vis[i * prime[j]] = 1;
            if(!(i % prime[j]))
            {
                phi[i * prime[j]] = phi[i] * prime[j] % mod;
                break;
            }
            else
                phi[i * prime[j]] = phi[i] * (prime[j] - 1) % mod;
        }
    }
    for(register int i = 1;i <= CNT;++i)
        phi[i] = (phi[i - 1] + phi[i] * i % mod * i % mod) % mod;
    for(register long long l = 1,r;l <= n;l = r + 1)
    {
        r = n / (n / l);
        ans = (ans + sum(n / l) * sum(n / l) % mod * (calc(r) - calc(l - 1)) % mod) % mod;
    }
    printf("%lld\n",(ans + mod) % mod);
}

Gitalking ...