【题意】于给出的n个询问,每次求有多少个数对(x,y),满足a≤x≤b,c≤y≤d,且gcd(x,y) = k,gcd(x,y)函数为x和y的最大公约数。n,a,b,c,d,k<=50000。
【算法】数论(莫比乌斯反演)
【题解】差分转化为询问有多少数对(x,y)满足x,y互素,1<=x<=n/k,1<=y<=m/k。
令f[x]表示gcd(a,b)=x的数对个数,F[x]表示满足 x | gcd(a,b) 的数对个数,则F[x]=Σx|df(d)。
易得F[x]=(n/x)*(m/x),那么根据莫比乌斯反演定理,f(x)=Σx|dμ(d/n)*F(d)=Σx|dμ(d/n)*(n/d)*(m/d)。
当x=1时,f(1)=Σμ(d)*(n/d)*(m/d),d=1~min(n,m),单次询问复杂度O(n)。
继续优化,n/d至多只有2*√n个取值,只要枚举这些取值后运用μ的前缀和(预处理)快速计算。
具体方法是:当前取值为n/i时,最小为i,最大为pos=n/(n/i),这m/(m/i)取min即可。
复杂度O(n√n)。
#include<cstdio> #include<algorithm> #define ll long long using namespace std; const int maxn=50010; int miu[maxn],prime[maxn],tot,s[maxn],n; bool mark[maxn]; void pre(int n){ miu[1]=1; for(int i=2;i<=n;i++){ if(!mark[i])miu[i]=-1,prime[++tot]=i; for(int j=1;j<=tot&&i*prime[j]<=n;j++){ mark[i*prime[j]]=1; miu[i*prime[j]]=-miu[i]; if(i%prime[j]==0){miu[i*prime[j]]=0;break;} } } for(int i=1;i<=n;i++)s[i]=s[i-1]+miu[i]; } ll solve(int n,int m){ ll ans=0;int pos=0; for(int i=1;i<=min(n,m);i=pos+1){ pos=min(n/(n/i),m/(m/i)); ans+=1ll*(s[pos]-s[i-1])*(n/i)*(m/i); } return ans; } int main(){ scanf("%d",&n); pre(50000); for(int i=1;i<=n;i++){ int a,b,c,d,k; scanf("%d%d%d%d%d",&a,&b,&c,&d,&k); a--;c--;a/=k;b/=k;c/=k;d/=k; printf("%lld ",solve(b,d)-solve(b,c)-solve(a,d)+solve(a,c)); } return 0; }
尝试从套路的角度来推导ans=Σx|dμ(d/n)*(n/d)*(m/d)
★当x=1时,Σd|xμ(x)=1。所以gcd(a,b)=1等价于Σd|a&&d|bμ(d)。——①
$$ans=sum_{i=1}^{n}sum_{j=1}^{m}[gcd(i,j)=k]$$
由(i,j)=k等价于(i/k,j/k)=1可以得到:——②
$$ans=sum_{i=1}^{frac{n}{k}}sum_{j=1}^{frac{m}{k}}[gcd(i,j)=1]$$
下一步代入经典gcd转μ,得到:
$$ans=sum_{i=1}^{frac{n}{k}}sum_{j=1}^{frac{m}{k}}sum_{d|icap d|j}mu (d)$$
套路化地改为枚举gcd,得到:——③
$$ans=sum_{d=1}^{min(frac{n}{k},frac{m}{k})}mu (d)sum_{i=1}^{frac{n}{k}}sum_{j=1}^{frac{m}{k}}[d|icap d|j]$$
最后部分满足条件的数对都可以除以d,就可以压缩上标直接计算了,即:——④
$$ans=sum_{d=1}^{min(frac{n}{k},frac{m}{k})}mu (d)left lfloor frac{n}{kd} ight floorleft lfloor frac{m}{kd} ight floor$$