【 HDU 4059 The Boss on Mars】 差分序列+容斥

题目链接：http://acm.hdu.edu.cn/showproblem.php?pid=4059

题目大意：给你一个整数n，让你求前面n个数中与n互素的数E（x^4）的和。

解题思路：

    可以先求出Sum（n）：表示1^4+2^4+3^4+……+n^4;

    然后再求Sum（n）-E(x^4)  ,x表示与 n不互素的数。

   求与n不互素的数可以利用容斥定理，求Sum（n）可以利用差分。

简单介绍一下差分序列：

给你一列数 a[i][1]，a[i][2]，a[i][3]，a[i][4]，a[i][5]……

那么a[i][j]=a[i-1][j+1]-a[i-1][j]， 即后一行是上一行相邻两项的差（第一行除外）。

如果给你一个多项式， 比如 f（x）=(x+1)*(x+2)*……*（x+p）,即多项式最高项系数为p。

则这个差分表有如下性质：

1、将 0、1、2、……、p-1、p分别代入多项式中，则他们组成多项式的第一行，后面的差分表可以由上一行推出。第p+1行以后差分表的值都为0。

2、我们这里要用的差分序列是他的第0行对角线的数。 我们设他们为c0、c1、c2、……cp；

   ==>  第n项的值：f（n）=c0*C（n，0）+c1*C(n，1)+c2*C（n，2）+……+cp*（n，p）；

  ==>  前n项的值：Sum（n）=c0*C（n+1，1）+c1*C(n+1，2)+c2*C（n+1，3）+……+cp*（n+1，p+1）；

 

  这题藐视搞得我蛋疼，我是把求前n项和组合公式给化简出来（Sum（n）=(n^5)/5+(n^4)/2+(n^3)/3-n/30，最后发现求逆元取模的时候问题很多。最后又多写了几个函数，直接求组合数的值，分母和模求逆元的时候打表先求好。

 还有要注意的是：

1、这题不打表先求素数会超时。

2、运用容斥定理的时候取模遇减号要加mod，还有容斥定理求出的数可以用到Sum（n）：比如2^4+4^4+6^4+8^4 = 2^4*(1^4+2^4+3^4+4^4) 看懂了木有。

 

  1 #include<iostream>
  2 #include<cstdio>
  3 #include<cmath>
  4 #include<cstring>
  5 #include<vector>
  6 #define lld long long
  7 #define MOD 1000000007
  8 using namespace std;
  9 int prime[10000],cnt=0,flag[10005]= {0};
 10 vector<int>vt;
 11 
 12 lld ext_gcd(lld a, lld b, lld &x, lld &y)
 13 {
 14     if( !b ) return x=1,y=0,a;
 15     lld d= ext_gcd(b, a%b, y, x);
 16     y-= a/b*x;
 17     return d;
 18 }
 19 lld Inv(int a, int m)   ///求逆元
 20 {
 21     lld d,x,y;
 22     d= ext_gcd(a, m, x, y);
 23     return x<0 ? x+m : x;
 24 }
 25 
 26 lld base[8]= {0,0,1,14,36,24}; ///差分表第0行对角线的值
 27 lld dd[10]= {0, 1, 500000004, 166666668, 41666667, 808333339}; ///打表求出需要的逆元
 28 
 29 lld Comb(lld n, lld k)   ///组合
 30 {
 31     if(k==n) return 1;
 32     if(k>n) return 0;
 33     lld s= 1;
 34     for(lld i=n,j=1; j<=k; --i,++j)
 35         s= s*i%MOD;
 36     return s*dd[k]%MOD;
 37 }
 38 
 39 lld Sum(lld n)    ///求n^4的前n项和
 40 {
 41     if(n==1) return 1;
 42     lld sum= 0;
 43     for(int i=2; i<=5; ++i)    ///差分
 44     {
 45         sum+= ( Comb( n+1, i ) * base[i] )%MOD;
 46         if( sum>=MOD ) sum-= MOD;
 47     }
 48     return sum;
 49 }
 50 
 51 lld Pow(lld n)  ///求n^4
 52 {
 53     lld ans=n;
 54     ans=(((((ans*n)%MOD)*n)%MOD)*n)%MOD;
 55     return ans;
 56 }
 57 
 58 void Prime()   ///筛选素数，便于后面的分解
 59 {
 60     for(int i=2; i<=10000; i++)
 61     {
 62         if(flag[i]) continue;
 63         prime[cnt++]=i;
 64         for(int j=2; j*i<=10000; j++)
 65             flag[i*j]=1;
 66     }
 67 }
 68 
 69 lld dfs(int idx,lld n)     ///容斥原理
 70 {
 71     lld ret=0,tmp;
 72     for(int i=idx; i<vt.size(); i++)
 73     {
 74         tmp=vt[i];
 75         ret=(ret+(Sum(n/tmp)*Pow(tmp))%MOD)%MOD;
 76         ret=((ret-dfs(i+1,n/tmp)*Pow(tmp))%MOD+MOD)%MOD;
 77     }
 78     return ret%MOD;
 79 }
 80 int main()
 81 {
 82     int t;
 83     lld n;
 84     Prime();
 85     scanf("%d",&t);
 86     while(t--)
 87     {
 88         scanf("%I64d",&n);
 89         vt.clear();
 90         lld tmp=n;
 91         for(int i=0; i<cnt&&prime[i]<=tmp; i++)  ///求n的素因子
 92             if(tmp%prime[i]==0)
 93             {
 94                 vt.push_back(prime[i]);
 95                 while(tmp%prime[i]==0)
 96                     tmp/=prime[i];
 97             }
 98         if(tmp!=1)
 99             vt.push_back(tmp);
100 
101         lld sum=((Sum(n)-dfs(0,n))%MOD+MOD)%MOD;
102         printf("%I64d\n",sum);
103     }
104     return 0;
105 }

 

 

 

  

转载于:https://www.cnblogs.com/kane0526/archive/2012/12/25/2833248.html