hdu4059 The Boss on Mars

求出(sum_{i = 1}^n [gcd(n, i) == 1] i ^ 4)
即求出([1,n])里面与(n)互质的数的四次方和

考虑互质,想到唯一分解求出n的所有质因子,然后利用这些质因子容斥求出所有不符合的数字的四次方和
然后求出所有的和,减去不符合的值就是答案。

[sum_{i=1}^{n} i^{4}=frac{n(n+1)(2 n+1)left(3 n^{2}+3 n-1 ight)}{30}=frac{1}{5} n^{5}+frac{1}{2} n^{4}+frac{1}{3} n^{3}-frac{1}{30} n ]

在求一个数字的倍数的四次方和时,比如2的倍数,把(2^4)提取出来就行了
( 2 4 6 8 \ 2^4 + (2*2)^4 + (2*3)^4 + (2*4)^4\ = 2^4 * (1^4 + 2^4 + 3^4 + 4^4) )

#include <iostream>
#include <cstdio>
#include <vector>
#include <cmath>
#define ll long long
using namespace std;
const int MOD = 1e9 + 7;
std::vector<int> divi(int n){
    std::vector<int> v;
    for(int i = 2; i * i <= n; i++) {
        if(n % i == 0) {
            v.push_back(i);
            while(n % i == 0) n /= i;
        }
    }
    if(n > 1) v.push_back(n);
    return v;
}
int gcd(int a, int b){
    return b == 0 ? a : gcd(b, a % b);
}
int lcm(int a, int b){
    return a / gcd(a, b) * b;
}
ll pow(ll a, ll b, ll p){
    ll ans = 1; a %= p;
    while(b){
        if(b & 1) ans = ans * a % p;
        a = a * a % p;
        b >>= 1;
    }
    return ans;
}
ll inv(ll a, ll p){
    return pow(a, p - 2, p);
}
ll cal(ll n){
    return n * (n + 1) % MOD * (2 * n + 1) % MOD * (3 * n * n % MOD + 3 * n - 1) % MOD * inv(30, MOD) % MOD;
}
namespace IO {
    template <typename T>
    inline void w(T x) { if (x > 9) w(x / 10); putchar(x % 10 + 48); }
    template <typename T>
    inline void write(T x, char c) { if(x < 0) putchar('-'), x = -x; w(x); putchar(c); }
    template <typename T>
    inline void read(T &x) {
        x = 0; T f = 1; char c = getchar();
        for (; !isdigit(c); c = getchar()) if (c == '-') f = -1;
        for (; isdigit(c); c = getchar()) x = (x << 1) + (x << 3) + (c ^ 48);
        x *= f;
    }
};
void solve(){
    int n; scanf("%d", &n);
    ll ans = 0, ans2 = cal(n);
    std::vector<int> v = divi(n);
    for(int i = 1; i < (1 << v.size()); i++) {
        int now = 1, k = 0;
        for(int j = 0; j < v.size(); j++) {
            if(i & (1 << j)) {
                now = now * v[j]; k++;
            }
        }
        now = pow(now, 4, MOD) * cal(n / now) % MOD;
        if(k % 2 == 0) {
            ans = ((ans - now) % MOD + MOD) % MOD;
        }else {
            ans = ((ans + now) % MOD + MOD) % MOD;
        }
    }
    IO::write(((ans2 - ans) % MOD + MOD) % MOD, '
');
}
int main(){
    int t; scanf("%d", &t);
    for(int i = 1; i <= t; i++) solve();
    return 0;
}
原文地址:https://www.cnblogs.com/Emcikem/p/14103238.html