HDU 6706 huntian oy（杜教筛 + 一些定理）题解

题意:

已知\(f(n,a,b)=\sum_{i=1}^n\sum_{j=1}^igcd(i^a-j^a,i^b-j^b)[gcd(i,j)=1]\mod 1e9+7\)，\(n\leq1e9\)，且保证\(ab\)互质，求\(f(n,a,b)\)

思路：

由不知道什么得：当\(ab\)互质，则\(gcd(i^a-j^a,i^b-j^b)=i-j\)。
那么可转化为：

\[\begin{aligned} &\sum_{i=1}^n\sum_{j=1}^ii-j[gcd(i,j)=1]\mod 1e9+7\\ =&\sum_{i=1}^ni\varphi(i)-\frac{i*\varphi(i)}{2}\\ =&\frac{1}{2}(\sum_{i=1}^ni\varphi(i) \ \ \ - 1) \end{aligned} \]

那么问题转化为求：\(\sum_{i=1}^ni\varphi(i)\)，因为前缀和有点大，所以估计要用到杜教筛，所以先卷积一下：

\[\begin{aligned} (f*g)(n)=&\sum_{d|n}f(d)g(\frac{n}{d})\\ =&\sum_{d|n}i\varphi(i)g(\frac{n}{d}) \end{aligned} \]

为了让卷积完的\(h\)和\(g\)前缀更好求，那么这里不如令\(g=id\)，那么\((f*g)(n)=n*\sum_{d|n}\varphi(d)\)，又因为\(\varphi*I=id\)，故\((f*g)(n)=n^2\)。
所以可得终式：

\[S(n)=\sum_{i=1}^ni^2-\sum_{i=2}iS(\lfloor\frac{n}{i}\rfloor) \]

因为超过\(1e6\)只有\(10\)组，所以我这里直接手动哈希比\(map\)快一点。

代码：

#include<map>
#include<set>
#include<cmath>
#include<cstdio>
#include<stack>
#include<ctime>
#include<vector>
#include<queue>
#include<cstring>
#include<string>
#include<sstream>
#include<iostream>
#include<algorithm>
typedef long long ll;
typedef unsigned long long ull;
using namespace std;
const int maxn = 4e6 + 5;
const ll MOD = 1e9 + 7;
const ull seed = 131;
const int INF = 0x3f3f3f3f;

int mark[maxn], prime[maxn], tot;
ll phi[maxn];
void init(int n){
    tot = 0;
    int i, j;
    phi[1] = 1;
    for(i = 2; i <= n; i++){
        if(!mark[i]){
            prime[++tot] = i;
            phi[i] = i - 1;
        }
        for(j = 1; j <= tot; j++){
            if(i * prime[j] > n)  break;
            mark[i * prime[j]] = 1;
            if(i % prime[j] == 0){
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            }
            else  phi[i * prime[j]] = phi[i] * (prime[j] - 1);
        }
    }

    phi[0] = 0;
    for(int i = 1; i <= n; i++){
        phi[i] = (1LL * phi[i] * i + phi[i - 1]) % MOD;
    }
}

ll inv6 = 166666668, inv2 = 500000004;
ll mphi[100000], UP;
int vis[100000];
int N = 4e6;
ll getans(int x){
    return x <= N? phi[x] : mphi[UP / x];
}
ll ksc(ll a, ll b, ll p){
    ll t = a * b - (ll)((long double)a * b / p + 0.5) * p;
    return (t < 0)? t + p : t;
}
void solve(int n){
    int t = UP / n;
    if(n <= N) return;
    if(vis[t]) return;
    vis[t] = 1;
    mphi[t] = ksc(ksc(n, n + 1, MOD), 2 * n + 1, MOD) * inv6 % MOD;
    for(int l = 2, r; l <= n; l = r + 1){
        r = n / (n / l);
        solve(n / l);
        mphi[t] -= 1LL * (l + r) * (r - l + 1) % MOD * inv2 % MOD * getans(n / l) % MOD;
        mphi[t] %= MOD;
    }
}
int main(){
    int T;
    init(N);
//    init(N + 1);
//    printf("%lld\n", (phi[N + 1] - 1) * inv2 % MOD);
    scanf("%d", &T);
    while(T--){
        ll ans;
        int n, a, b;
        scanf("%d%d%d", &n, &a, &b);
        if(n <= N) ans = phi[n];
        else{
            UP = n;
            memset(vis, 0, sizeof(vis));
            solve(n), ans = mphi[1];
        }
        ans = (ans - 1) * inv2;
        ans = (ans % MOD + MOD) % MOD;
        printf("%lld\n", ans);
    }
    return 0;
}