【HAOI2011】problem b
数论好劲啊
原题：
对于给出的n个询问，每次求有多少个数对(x,y)，满足a≤x≤b，c≤y≤d，且gcd(x,y) = k，gcd(x,y)函数为x和y的最大公约数。
1≤n≤50000，1≤a≤b≤50000，1≤c≤d≤50000，1≤k≤50000
gcd出现，基本上可以确定是数论了
最开始脑补了一下感觉应该是容斥之类的东西，去网上搜题解，果然是容斥（反演算是建立在整除上的容斥？
然后补习po姐的ppt，发现这题居然是第一道例题
补习了一下反演，这次和第一次学的时候不一样，证明什么的都扔掉，直接记公式应用（应试套路暴力硬肛
然后实力提升了果然不一样，现在拿起笔推公式也有感觉了
题目中要求gcd(x,y)=k的xy对数，这个可以用一个函数表示：f(k)=cnt(gcd(x,y)=k)
然后gcd是建立在整除上的，所以可以往莫比乌斯反演的方向想，尝试构造函数F(i)=cnt(i|gcd(x,y))（我说的推导过程这么有代入感其实我靠自己是根本推不出来的hhh
于是
（我发现这个公式显示好像有点问题= =另存到本地即可查看）
然后就可以使用反演辣，直接从中反演出f(i)
F(i)比较容易计算
枚举i（实际上是k）的倍数即可算出答案
但是酱紫做是O(N)的，会T掉（还有n个询问）
因为可能会有很多
是相等的，所以可以直接搞一个miu的前缀和，把相等的F(i)乘上miu的前缀和即可
怎么求相等的F(i)呐，每次前缀和的右端点是min(n/(n/i),m/(m/i))，左端点是上一次的右端点即可
（我说不下去了具体还是看po姐的ppt理解吧
写代码的时候又出现傻逼错误了，不过这一次很快召唤了数学大神syq帮忙查看，所以很快就找到问题了
summiu[1]没设初值……妙啊……
代码：
1 #include<iostream> 2 #include<cstdio> 3 #include<algorithm> 4 #include<cstring> 5 #include<cmath> 6 using namespace std; 7 #define ll long long 8 int rd(){int z=0,mk=1; char ch=getchar(); 9 while(ch<'0'||ch>'9'){if(ch=='-')mk=-1; ch=getchar();} 10 while(ch>='0'&&ch<='9'){z=(z<<3)+(z<<1)+ch-'0'; ch=getchar();} 11 return z*mk; 12 } 13 int n; 14 int prm[51000],ptt=0; 15 bool flg[51000]; 16 int miu[51000],smmiu[51000]; 17 void slct(){ 18 memset(flg,0,sizeof(flg)); 19 miu[1]=1,smmiu[1]=1; 20 for(int i=2;i<=50000;++i){ 21 if(!flg[i]) prm[++ptt]=i,miu[i]=-1; 22 for(int j=1;prm[j]*i<=50000 && j<=ptt;++j){ 23 flg[prm[j]*i]=true; 24 if(!(i%prm[j])){ miu[i*prm[j]]=0; break;} 25 miu[prm[j]*i]=-miu[i]; 26 } 27 smmiu[i]=miu[i]+smmiu[i-1]; 28 } 29 } 30 ll cclt(int n,int m,int k){ 31 n/=k,m/=k; 32 if(n>m) swap(n,m); 33 ll bwl=0,tmp; 34 for(int i=1;i<=n;i=tmp+1){ 35 tmp=min(n/(n/i),m/(m/i)); 36 bwl+=(ll)(n/i)*(m/i)*(smmiu[tmp]-smmiu[i-1]); 37 } 38 return bwl; 39 } 40 int main(){//freopen("ddd.in","r",stdin); 41 slct(); 42 cin>>n; 43 int a,b,c,d,k; 44 while(n--){ 45 a=rd(),b=rd(),c=rd(),d=rd(),k=rd(); 46 printf("%I64d\n",cclt(b,d,k)+cclt(a-1,c-1,k)-cclt(b,c-1,k)-cclt(a-1,d,k)); 47 } 48 return 0; 49 }
