HDU5382 GCD?LCM!(数学公式推导)

题意：

$\ \ \ \ $定义$F(n) = \sum\limits_{i = 1}^{n}\sum\limits_{j = 1}^{n}[gcd(i,j)+lcm(i,j)\geqslant n]$，其中$[]$是艾弗森约定，现在要求出$ans_{n} = \sum\limits_{i = 1}^{n}F(i)$.

思路：

$\ \ \ \ $如果定义$T(n) = \sum\limits_{i = 1}^{n}\sum\limits_{j = 1}^{n}[gcd(i,j)+lcm(i,j)=n]$，那么$F_{n} = i * i-\sum\limits_{i = 1}^{n-1}T(i)$，因为不可能存在$i\geqslant n$或者$j\geqslant n$时有$gcd(i,j)+lcm(i,j)=n$，看$T(n)$怎么算:
$T(n) =\sum\limits_{d = 1}^{n}\sum\limits_{i = 1}^{\lfloor \frac{n}{d} \rfloor}\sum\limits_{j = 1}^{\lfloor \frac{n}{d} \rfloor}[gcd(i,j)=1][i*j*d+d=n]$
$\ \ \ \ \ \ \ \ \ $$=\sum\limits_{d = 1}^{n}\sum\limits_{i = 1}^{\lfloor \frac{n}{d} \rfloor}\sum\limits_{j = 1}^{\lfloor \frac{n}{d} \rfloor}[gcd(i,j)=1][(i*j+1)d=n]$
$\ \ \ \ \ \ \ \ \ $$=\sum\limits_{d|n}\sum\limits_{i = 1}^{\frac{n}{d}}\sum\limits_{j = 1}^{\frac{n}{d}}[gcd(i,j)=1][(i*j+1)d=n]$
$\ \ \ \ \ \ \ \ \ $$=\sum\limits_{d|n}\sum\limits_{i = 1}^{\frac{n}{d}}\sum\limits_{j = 1}^{\frac{n}{d}}[gcd(i,j)=1][i*j=\frac{n}{d}-1]$
$\ \ \ \ \ \ \ \ \ $$=\sum\limits_{d|n}\sum\limits_{i = 1}^{\frac{n}{d}}\sum\limits_{j = 1}^{\frac{n}{d}}[gcd(i,j)=1][j=\frac{\frac{n}{d}-1}{i}]$
$\ \ \ \ \ \ \ \ \ $$=\sum\limits_{d|n}\sum\limits_{i = 1}^{\frac{n}{d}}[gcd(i,\frac{\frac{n}{d}-1}{i})=1]$
$\ \ \ \ \ \ \ \ \ $$=\sum\limits_{d|n}\sum\limits_{i = 1}^{\frac{n}{d} - 1}[gcd(i,\frac{\frac{n}{d}-1}{i})=1]$
记$f(x) = \sum\limits_{i | x}[gcd(i,\frac{x}{i})=1]$，如果$i$和$x/i$互素，那么可以知道$i$要么完全拥有$x$的素因子，要么完全没有$x$的素因子，暴力搜索一下素因子产生$x$的因子即可，不用$gcd$去判断，因为那样时间就变成$O(n\lg^{2} n).....$，会超时的。

#include<bits/stdc++.h>
typedef long long ll;
const int maxn = 1e6 + 5;
const ll mod = 258280327;
using namespace std;

typedef pair<int, int> pa;
int a[maxn];

pa pri[1000]; int num;
void div(int x) {
    while(x > 1) {
        int now = a[x];
        if(!num || pri[num - 1].first != now) pri[num++] = pa(now, 1);
        else pri[num - 1].second++;
        x /= now;
    }
}

int vec[1000], cnt, T;
void dfs(ll &ans, int sol, int sum, int id, int flag) {
    if(id == num) {
        int y = sol / sum;
        vec[cnt++] = sum;
        if(!flag) ans++;
        return ;
    }
    int k = 1, g = flag;
    for(int i = 0; i <= pri[id].second; i++) {
        if(!i || i == pri[id].second) flag = 0;
        else flag = 1;
        dfs(ans, sol, sum * k, id + 1, g | flag);
        k *= pri[id].first;
    }
}

ll ans[maxn], rec[maxn];
ll seq[maxn], f[maxn], tot[maxn];

void solve() {
    int max_p = 0;
    for(ll i = 1; i <= 1000000; i++) {
        num = 0; div(i); cnt = 0;
        dfs(f[i], i, 1, 0, 0);
        rec[i] = i * i;
        rec[i] = rec[i] - tot[i - 1];
        for(int j = 0; j < cnt; j++) {
            int x = vec[j];
            seq[i] += f[x - 1];
        }
        tot[i] = tot[i - 1] + seq[i];
        ans[i] = ans[i - 1] + rec[i];
    }
}


int main() {
    for(ll i = 2; i < maxn; i++) {
        if(a[i]) continue;
        a[i] = (int) i;
        for(ll j = i * i; j < maxn; j += i) if(!a[j]) a[j] = i;
    }
    solve(); int n;
    scanf("%d", &T);
    while(T--) {
        scanf("%d", &n);
        printf("%lld\n", ans[n] % mod);
    }
    return 0;
}