hdu 4407 容斥定理

题意:有一个元素为 1~n 的数列{An},有2种操作(1000次):

1、求某段区间 [a,b] 中与 p 互质的数的和。

2、将数列中某个位置元素的值改变。

思路:acm数学方面实在很差~学习ing~看解题报告A的

容斥定理,看《组合数学》时看到过,但是从来没用过,求第x个数到第y个数中与p互素的数的和。

p=p1^a1 * p2^a2 * p3^a3 * ……其中p1,p2,p3……是n的素因子,ai是pi的指数。

求出 1~y的与p互素的数的和 减去 1~x的与p互素的数的和。

只要是pi的倍数的,则不可能与p互素。所以例如30=2*3*5;那么

1~y的与30互素的数的和=y*(y+1)-2的倍数和-3的倍数和-5的倍数和+2*3的倍数和+3*5的倍数和+2*5的倍数和-2*3*5的倍数和+(改变量)。

(改变量):如果改变的的大于y则无视,否则,如果 改变后的值 与p互素,则+改变后的值

                       如果 改变前的值 与p互素,则-改变前的值。

感觉在处理容斥定理这一部分的代码写的很好。

View Code
 1 int i,j,t,d;//j是第几个因子的下角标,t是临时
 2 int m=1<<n;//p有n个素因子
 3 int num;//记录第几次,决定加还是减 例如,30=2*3*5,-2-3-5 +2*3+3*5+2*5 -2*3*5
 4 long long ans=x*(x+1)/2;
 5 for(i=1;i<m;i++){
 6   t=i;j=num=0;d=1;
 7   while(t){
 8   if(t&1) {d*=factor[j];num++;}
 9   j++;t>>=1;
10   }
11   t=x/d;
12   if(num&1) ans-=(long long)d*t*(t+1)/2;//必须加强制类型转换,否则会WA
13   else ans+=(long long)d*t*(t+1)/2;//必须加强制类型转换,否则会WA
14 }

796msAC,实在。。。

View Code
 1 #include <map>
 2 #include <string.h>
 3 #include <stdio.h>
 4 #include <iostream>
 5 using namespace std;
 6 const int N=400005;
 7 int factor[30];//存因子的数组
 8 bool isPrime[N];//判断是否是素数表
 9 int prime[N],cnt;//存素数
10 map<int,int>M;//存第2种操作,改变的值
11 int gcd(int a,int b){//最大公约数
12     while(a=a%b)a^=b^=a^=b;
13     return b;
14 }
15 void IsPrime(){//筛选法打素数表
16     int i,j;
17     cnt=0;
18     memset(isPrime,1,sizeof(isPrime));
19     for(i=2;i<N;i++){
20         if(isPrime[i]){
21             for(j=i+i;j<N;j+=i)
22                 isPrime[j]=0;
23             prime[cnt++]=i;
24         }
25     }
26 }
27 int getPrime(int p){//求p的所有素数因子
28     if(isPrime[p]) {factor[0]=p;return 1;}
29     int i,k=0;
30     for(i=0;i<cnt;i++){
31         if(p%prime[i]==0) factor[k++]=prime[i];
32         while(p%prime[i]==0) p/=prime[i];
33         if(p!=1&&isPrime[p]) {factor[k++]=p;return k;}
34     }
35     return k;
36 }
37 long long find(int x,int n,int p){//find()函数求在1~x内与p互素的个数
38     int i,j,t,d;//j是第几个因子的下角标,t是临时
39     int m=1<<n;//p有n个素因子
40     int num;//记录第几次,决定加还是减 例如,30=2*3*5,-2-3-5 +2*3+3*5+2*5 -2*3*5
41     long long ans=x*(x+1)/2;
42     for(i=1;i<m;i++){
43         t=i;j=num=0;d=1;
44         while(t){
45             if(t&1) {d*=factor[j];num++;}
46             j++;t>>=1;
47         }
48         t=x/d;
49         if(num&1) ans-=(long long)d*t*(t+1)/2;//必须加强制类型转换,否则会WA
50         else ans+=(long long)d*t*(t+1)/2;//必须加强制类型转换,否则会WA
51     }
52     map<int,int>::iterator it;
53     for(it=M.begin();it!=M.end();it++){
54         if(it->first>x) continue;
55         if(gcd(it->first,p)==1) ans-=it->first;
56         if(gcd(it->second,p)==1) ans+=it->second;
57     }
58     return ans;
59 }
60 int main(){
61     int i;
62     int t,n,m;
63     int flag,x,y,p;
64     IsPrime();
65     scanf("%d",&t);
66     while(t--){
67         M.clear();
68         scanf("%d%d",&n,&m);
69         for(i=0;i<m;i++){
70             scanf("%d",&flag);
71             if(flag==1){
72                 scanf("%d%d%d",&x,&y,&p);
73                 int num=getPrime(p);//p素因子个数num
74                 printf("%I64d\n",find(y,num,p)-find(x-1,num,p));
75             }else{
76                 scanf("%d%d",&x,&p);
77                 M[x]=p;
78             }
79         }
80     }
81     return 0;
82 }

 

posted @ 2013-02-06 02:02  _sunshine  阅读(297)  评论(0编辑  收藏  举报