HDU4059 The Boss on Mars【容斥原理】【乘法逆元】【高次求和】

题目链接：

http://acm.hdu.edu.cn/showproblem.php?pid=4059

题目大意：

给一个整数N，求1~N中与N互质的数的4次方的和。

解题思路：

题目简单，过程有点复杂。理清思路就简单了。

利用公式1^4 + 2^4 + … + n^4 = n*(n+1)*(2*n+1)*(3*n^2+3*n-1)/30，可以求出

1^4 + 2^4 + … + n^4，除以30可以先求出30模M的逆元，然后将上式中除以30改为

乘以30的逆元。

再来求与n不互质的数的4次方和。将n质因数分解为 n = p1^a1* p2^a2 * … * pn^an

(其中p1、p2为质数)。

这样不互质的数就是n的质因数的倍数，用容斥原理求出不互质的数的4次方和，最终

(1^4 + 2^4 + … + n^4 - 不互质的数的4次方和)即为所求，具体参考代码。

AC代码：

#include<iostream>
#include<algorithm>
#include<cstdio>
#include<cstring>
#define LL __int64
#define M 1000000007
using namespace std;

LL d,x,y,Inv,N;

void ExGcd(LL a,LL b,LL &d,LL &x,LL &y)
{
    if(b == 0)
    {
        x = 1;
        y = 0;
        d = a;
    }
    else
    {
        ExGcd(b,a%b,d,y,x);
        y -= x*(a/b);
    }
}

LL ModInverse(LL a,LL n,LL &d,LL &x,LL &y)    //a*a-1 = 1(mod n) 用来计算30模M的逆元
{
    ExGcd(a,n,d,x,y);
    if(1 % d != 0)
        return -1;
    LL ans = x/d < 0 ? x/d + n : x/d;
    return ans;
}

LL Four(LL k)   //计算1^4 + 2^4 + … + k^4( 公式为 k*(k+1)*(2*k+1)*(3*k*k+3*k-1)/30 )
{
    LL ans = k;
    ans = ans * (k+1) % M;
    ans = ans * (2*k+1) % M;
    ans = ans * ((3*k*k%M+3*k-1)%M) % M;
    ans = (ans * Inv) % M;  // 除以30 等于 乘以30模M的逆元
    return ans;
//    return (k%M) * ((k+1)%M) % M * ((2*k+1)%M) % M * ((3*k*k%M + (3*k-1)%M)%M)%M * Inv;
}

LL Factor[20],ct;

void Divide()   //将N分解素因子
{
    ct = 0;
    LL n = N;
    for(LL i = 2; i <= sqrt(n*1.0); ++i)
    {
        if(n % i == 0)
        {
            Factor[ct++] = i;
            while(n % i == 0)
                n /= i;
        }
    }
    if(n != 1)
        Factor[ct++] = n;
}

LL Cal(LL a)
//计算 质因数乘积为a的所有小于n的倍数的前四次方和，
//即a^4 + (2a)^4 + (3a)^4 + … + (k*a)^4 = a^4 * (1^4 + 2^4 + … + k^4)
{
    LL k = N / a;
    LL ans = 1;
    for(int i = 1; i <= 4; ++i)
        ans = ans * a % M;
    ans = ans * Four(k) % M;
    return ans;
}

LL Solve()  //容斥原理求解与n不互质数的4次方和，最终结果为(1^4 + 2^4 + … + n^4 - 与n不互质数的4次方和)
{
    Divide();
    LL ans = 0;
    for(LL i = 1; i < (1 << ct); ++i)
    {
        LL odd = 0;
        LL tmp = 1;
        for(LL j = 0; j < ct; ++j)
        {
            if((1<<j) & i)
            {
                tmp *= Factor[j];
                odd++;
            }
        }
        if(odd & 1)
            ans = (ans%M + Cal(tmp)%M) % M;
        else
            ans = (ans%M - Cal(tmp)%M + M) % M;
    }
    return (Four(N) - ans + M)%M;
}

int main()
{
    int T;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%I64d",&N);
        if(N == 1)
            printf("0\n");
        else
        {
            Inv = ModInverse(30,M,d,x,y);   //求30模M的逆元
            printf("%I64d\n",Solve());
        }
    }

    return 0;
}

版权声明：本文为博主原创文章，未经博主允许不得转载。