[BZOJ4916]神犇和蒟蒻

Posted 2020-10-25 xjr_01

[BZOJ4916]神犇和蒟蒻

试题描述

很久很久以前，有一只神犇叫yzy;

很久很久之后，有一只蒟蒻叫lty;

输入

请你读入一个整数 \(N\); \(1 \le N \le 10^9\), \(A\)、\(B\) 模 \(10^9+7\);

输出

请你输出一个整数 \(A=\sum_{i=1}^N{\mu (i^2)}\);

请你输出一个整数 \(B=\sum_{i=1}^N{\varphi (i^2)}\);

输入示例

1

输出示例

1
1

数据规模及约定

见“输入”

题解

知道莫比乌斯函数定义的人都知道 \(A\) 恒为 \(1\)。

然后我们把欧拉函数展开一下就会发现 \(\varphi(i^2) = i \cdot \varphi(i)\)。

那么就是杜教筛 \(i \cdot \varphi(i)\) 的前缀和。这个东西很神奇，你会发现将 \(i \cdot \varphi(i)\) 和 \(i\) 狄利克雷卷积一下就出来了。

\[ \sum_{i=1}^N \sum_{j|i} j \cdot \varphi(j) \cdot \frac{i}{j} \= \sum_{i=1}^N i \sum_{j|i} \varphi(j) \= \sum_{i=1}^N i^2 \]

目前还不知道这能干什么，令 \(F(N) = \sum_{i=1}^N i \cdot \varphi(i)\)，接着往下推就知道了

\[ \sum_{i=1}^N \sum_{j|i} j \cdot \varphi(j) \cdot \frac{i}{j} \= \sum_{j=1}^N { j \cdot \varphi(j) \sum_{i=1}^{\lfloor \frac{N}{j} \rfloor} i } \= \sum_{i=1}^N { i \sum_{j=1}^{\lfloor \frac{N}{i} \rfloor} j \cdot \varphi(j) } \= \sum_{i=1}^N i \cdot F \left( \lfloor \frac{N}{i} \rfloor \right) \]

于是有

\[ \sum_{i=1}^N i^2 = \sum_{i=1}^N i \cdot F \left( \lfloor \frac{N}{i} \rfloor \right) \F(n) = \sum_{i=1}^N i^2 - \sum_{i=2}^N i \cdot F \left( \lfloor \frac{N}{i} \rfloor \right) \]

解决！

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <algorithm>
#include <map>
using namespace std;
#define rep(i, s, t) for(int i = (s), mi = (t); i <= mi; i++)
#define dwn(i, s, t) for(int i = (s), mi = (t); i >= mi; i--)

int read() {
    int x = 0, f = 1; char c = getchar();
    while(!isdigit(c)){ if(c == ‘-‘) f = -1; c = getchar(); }
    while(isdigit(c)){ x = x * 10 + c - ‘0‘; c = getchar(); }
    return x * f;
}

#define maxn 1000010
#define MOD 1000000007
#define inv6 166666668
#define LL long long

bool vis[maxn];
int prime[maxn], cp, phi[maxn], sum[maxn];
void init() {
    int n = maxn - 1;
    phi[1] = sum[1] = 1;
    rep(i, 2, n) {
        if(!vis[i]) prime[++cp] = i, phi[i] = i - 1;
        for(int j = 1; j <= cp && i * prime[j] <= n; j++) {
            vis[i*prime[j]] = 1;
            if(i % prime[j] == 0) {
                phi[i*prime[j]] = phi[i] * prime[j];
                break;
            }
            phi[i*prime[j]] = phi[i] * (prime[j] - 1);
        }
        sum[i] = sum[i-1] + (LL)phi[i] * i % MOD;
        if(sum[i] >= MOD) sum[i] -= MOD;
    }
    return ;
}

map <int, int> hSum;
LL pre(int n) { return ((LL)n * (n + 1) >> 1) % MOD; }
int calc(int n) {
    if(n < maxn) return sum[n];
    if(hSum.count(n)) return hSum[n];
    int ans = (LL)n * (n + 1) % MOD * (n << 1 | 1) % MOD * inv6 % MOD;
    for(int i = 2; i <= n; ) {
        int r = min(n / (n / i), n);
        ans -= (LL)calc(n / i) * (pre(r) - pre(i - 1) + MOD) % MOD;
        if(ans < 0) ans += MOD;
        i = r + 1;
    }
    return hSum[n] = ans;
}

int main() {
    init();
    
    int n = read();

    printf("1\n%d\n", calc(n));
    
    return 0;
}