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互素,则-改变前的值。
感觉在处理容斥定理这一部分的代码写的很好。
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
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,实在。。。
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
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 }