题目求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\n",sum); } return 0; }
原文:http://www.cnblogs.com/Ritchie/p/5372369.html