HDU 4059 容斥原理+快速幂+逆元

E - The Boss on Mars

Time Limit:1000MS     Memory Limit:32768KB     64bit IO Format:%I64d & %I64u

Submit Status Practice HDU 4059

Description

On Mars, there is a huge company called ACM (A huge Company on Mars), and it’s owned by a younger boss. 

Due to no moons around Mars, the employees can only get the salaries per-year. There are n employees in ACM, and it’s time for them to get salaries from their boss. All employees are numbered from 1 to n. With the unknown reasons, if the employee’s work number is k, he can get k^4 Mars dollars this year. So the employees working for the ACM are very rich. 

Because the number of employees is so large that the boss of ACM must distribute too much money, he wants to fire the people whose work number is co-prime with n next year. Now the boss wants to know how much he will save after the dismissal. 

 

Input

The first line contains an integer T indicating the number of test cases. (1 ≤ T ≤ 1000) Each test case, there is only one integer n, indicating the number of employees in ACM. (1 ≤ n ≤ 10^8) 

 

Output

For each test case, output an integer indicating the money the boss can save. Because the answer is so large, please module the answer with 1,000,000,007. 

 

Sample Input

2 4 5

 

Sample Output

82 354

Hint

 Case1: sum=1+3*3*3*3=82 Case2: sum=1+2*2*2*2+3*3*3*3+4*4*4*4=354 

题目求1-n中与n互质的数的4次方之和，即S=a1^4+a2^4+……; a1,a2……均小于等于n且与n互质。

先求出1^4+2^4+……n^4然后减去与n不互质的数的4次方。

必然要先要用到4次方的求和公式。接下来简单的证明一下，这里前提是你知道3次方的公式，如果不会照下面的模式可以利用2次公式推出3次公式

（x+1)^5=x^5+5*x^4+10*x^3+10*x^2+5*x+1;

则 1=1;

2^5=(1+1)^5=1^5+5*1^4+10*1^3+10*1^2+5*1^1+1;

3^5=(2+1)^5=2^5+5*2^4+10*2^3+10*2^2+5*2^1+1;

……

……

(n+1)^5=(n+1)^5=n^5+5*n^4+10*n^3+10*n^2+5*n^1+1;

全部叠加起来，则(n+1)^5=5*(1^4+2^4+……n^4)+10*(1^3+2^3+……+n^3)+10*(1^2+2^2+……+n^2)+5*(1+2+……+n)+n+1;

然后将(1^3+2^3+……n^4)=(n+1)^2*n^2/4; (1^2+2^2+……n^2)=(n*(n+1)*(2*n+1))/6; 代入。

化简后得到(1^4+2^4+……+n^4)=(n*(n+1)*(2n+1)*(3*n*n+3*n-1))/30;

公式证毕，这里用到除以30，还得算一下30对MOD的逆元，也就是30^(MOD-2),

接下来要减掉与n不互质的数4次方，将n质因子分解后运用容斥原理即可，就是减掉一个因子的倍数的4次方结果，加上两个因子乘积的倍

数的4次方结果，减去……以此类推。也可以通过状态压缩枚举,貌似最多9个质因子吧。

在运算的时候，注意各种相乘溢出就行了，类似计算几何的精度问题，数论的溢出也很纠结。

#include<iostream>  
#include<cstdio>  
#include<cstring>  
#include<vector>  
#define LL long long  
#define MOD 1000000007  
using namespace std;  
LL res;  //30对MOD的逆元  
int prime[10000],cnt=0,flag[10005]={0};  
vector<int>fact;  
LL PowMod(LL a,LL b){  
    LL ret=1;  
    while(b){  
        if(b&1)  
            ret=(ret*a)%MOD;  
        a=(a*a)%MOD;  
        b>>=1;  
    }  
    return ret;  
}  
LL Sum(LL n){ //求an=n^4,的前n项和  
    LL ans=n;  
    ans=(ans*(n+1))%MOD;  
    ans=(ans*((2*n+1)%MOD))%MOD;  
    ans=(ans*(((3*n*n)%MOD+(3*n)%MOD-1+MOD)%MOD))%MOD;  
    ans=(ans*res)%MOD;  
    return ans;  
}  
LL Pow(LL n){ //求n^4  
    LL ans=n;  
    ans=(((((ans*n)%MOD)*n)%MOD)*n)%MOD;  
    return ans;  
}  
int t;  
void Prime(){  //筛选素数，便于后面的分解  
    for(int i=2;i<=10000;i++){  
        if(flag[i]) continue;  
        prime[cnt++]=i;  
        for(int j=2;j*i<=10000;j++)  
            flag[i*j]=1;  
    }  
}  
void Init(){   
    res=PowMod(30,MOD-2);   //求30对MOD的逆元  
    Prime();  
    scanf("%d",&t);  
}  
LL dfs(int idx,LL n){    //容斥原理  
    LL ret=0,tmp;  
    for(int i=idx;i<fact.size();i++){  
        tmp=fact[i];  
        ret=(ret+(Sum(n/tmp)*Pow(tmp))%MOD)%MOD;  
        ret=((ret-dfs(i+1,n/tmp)*Pow(tmp))%MOD+MOD)%MOD;  
    }  
    return ret%MOD;  
}  
int main(){  
    LL n;  
    Init();  
    while(t--){  
        scanf("%I64d",&n);  
        fact.clear();  
        LL tmp=n;  
        for(int i=0;i<cnt&&prime[i]<=tmp;i++)  
            if(tmp%prime[i]==0){  
                fact.push_back(prime[i]);  
                while(tmp%prime[i]==0)  
                    tmp/=prime[i];  
            }  
        if(tmp!=1)  
            fact.push_back(tmp);          
        LL sum=((Sum(n)-dfs(0,n))%MOD+MOD)%MOD;  
        printf("%I64d
",sum);  
    }  
    return 0;  
}