hdu 4059 The Boss on Mars

题意

  求1到之间所有与n互质的数的四次幂之和。

分析

  分解n^5-(n-1)^5..可以求出1^4+2^4+...+n^4的公式S(n)=1/30*(n^5+15*n^2*(n+1)^2-10*n*(n+1)*(2*n+1)+15*n*(n+1)-6*n)

  先求出S(n),利用容斥原理去掉是n的质因子的倍数的数的四次幂之和。。。

  1 #include<iostream>
  2 #include<algorithm>
  3 #include<cstring>
  4 #include<vector>
  5 #define N 100100
  6 #include<cstdio>
  7 using namespace std;
  8 typedef long long ll;
  9 const ll mod=1000000007;
 10 bool p[N];
 11 int s[N];
 12 vector<int>factor;
 13 void prime(){
 14     memset(p,1,sizeof(p));
 15     for(int i=2;i*i<N;i++)
 16         for(int j=i;i*j<N;j++)
 17             p[i*j]=0;
 18     s[0]=0;
 19     for(int i=2;i<N;i++)
 20         if(p[i])s[++s[0]]=i;
 21 }
 22 
 23 ll exgcd(int a,int b,ll &x,ll &y){
 24     if(b==0){
 25         x=1;
 26         y=0;
 27         return a;
 28     }
 29     ll d=exgcd(b,a%b,x,y);
 30     ll t=x;
 31     x=y;
 32     y=t-(a/b)*y;
 33     return d;
 34 }
 35 void div(ll n){          //分解因式
 36     factor.clear();
 37     for(int i=1;s[i]*s[i]<=n;i++){
 38         if(n%s[i]==0){
 39             factor.push_back(s[i]);
 40             while(n%s[i]==0)n/=s[i];
 41         }
 42     }
 43     if(n>1)factor.push_back(n);
 44 }
 45 
 46 ll calsum(ll n){
 47     ll a1=6*(n*n%mod)*(n*n%mod)%mod*n%mod;
 48     ll a2=15*n%mod*n%mod*(n+1)%mod*(n+1)%mod;
 49     ll a3=10*n%mod*(n+1)%mod*(2*n+1)%mod;
 50     ll a4=15*n%mod*(n+1)%mod;
 51     ll a5=6*n%mod;
 52     ll ans=(a1+a2-a3+a4-a5)%mod;
 53     return ans;
 54 }
 55 ll sum(ll n){     //求和
 56     ll ans=calsum(n);
 57     ll x,y;
 58     ll d=exgcd(30,mod,x,y);
 59     x*=ans;
 60     x=(x%mod+mod)%mod;
 61     return x;
 62 }
 63 
 64 ll Mod(ll a,ll b){
 65     return a*a%mod*a%mod*a%mod*b%mod;
 66 }
 67 
 68 void solve(ll n){
 69     div(n);
 70     ll ans=sum(n);
 71     int c=factor.size();
 72     for(int i=1;i<(1<<c);i++){//容斥原理
 73         ll num=1;
 74         int a=0;
 75         int cnt=0;
 76         int ii=i;
 77         while(ii>0){
 78             if(ii&1){
 79                 num*=factor[cnt];
 80                 ++a;
 81             }
 82             ii>>=1;
 83             ++cnt;
 84         }
 85         ll cc=Mod(num,sum(n/num));
 86         if(a&1){
 87             ans=(ans-cc+mod)%mod;
 88         }else{
 89             ans=(ans+cc)%mod;
 90         }
 91     }
 92     cout<<ans<<endl;
 93 }
 94 
 95 int main(){
 96 //    freopen("test.out","w",stdout);
 97     ll n,m;
 98     int t;
 99     prime();
100     cin>>t;
101     while(t--){
102         cin>>n;
103         solve(n);
104     }
105     return 0;
106 }
posted @ 2012-10-11 23:15  silver__bullet  阅读(104)  评论(0编辑  收藏  举报