hdu 4059 小学生容斥

啥都不说。。送上一个四次方求和公式:Sum(n)=n*(n+1)*(2n+1)(3n^2+3n-1)/30
再送上一个三次方求和公式: Sum(n) = (1+2+3+…+n)^2
抄的code

#include<iostream>
#include<algorithm>
#include<cstdio>
#include<cstring>
#define LL long long
#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("%lld\n",Solve());
        }
    }

    return 0;
}
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值