bzoj3994 [SDOI2015]约数个数和

题目链接

solution

先(orz)一波大佬，然后扔上一个公式(233)

[d(ij)=sumlimits_{x | i}sumlimits_{y|j}[gcd(x,y)==1] ]

是不是看到这个公式瞬间就有思路了。

我们开始推柿子。

[sumlimits_{i=1}^nsumlimits_{j=1}^m d(ij)\ =sumlimits_{i=1}^nsumlimits_{j=1}^msumlimits_{x|i}sumlimits_{y|j}[gcd(x,y)==1]\ =sumlimits_{x=1}^nsumlimits_{y=1}^mlfloorfrac{n}{x} floorlfloorfrac{m}{y} floorsumlimits_{d|gcd(x,y)}mu(d)\ =sumlimits_{d=1}^nmu(d)sumlimits_{x=1}^{lfloorfrac{n}{d} floor}sumlimits_{y=1}^{lfloorfrac{m}{d} floor}lfloorfrac{n}{dx} floorlfloorfrac{m}{dy} floor\=sumlimits_{d=1}^nmu(d)sumlimits_{x=1}^{lfloorfrac{n}{d} floor}lfloorfrac{n}{dx} floorsumlimits_{y=1}^{lfloorfrac{m}{d} floor}lfloorfrac{m}{dy} floor ]

令(f(n)=sumlimits_{i=1}^nlfloorfrac{n}{i} floor)

我们可以利用整除分块在(O(nsqrt{n}))的复杂度内预处理所有的(f)

所以原式=

[sumlimits_{d=1}^nmu(d)f(lfloorfrac{n}{d} floor)f(lfloorfrac{m}{d} floor) ]

因为(f)都已经预处理出来，所以只要整除分块，就可以在(O(sqrt{n}))的复杂度内计算每次询问了。

code

/*
* @Author: wxyww
* @Date:   2020-04-23 19:36:37
* @Last Modified time: 2020-04-23 20:04:55
*/
#include<cstdio>
#include<iostream>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#include<queue>
#include<vector>
#include<ctime>
using namespace std;
typedef long long ll;
const int N = 50010;
ll read() {
	ll x = 0,f = 1;char c = getchar();
	while(c < '0' || c > '9') {
		if(c == '-') f = -1; c = getchar();
	}
	while(c >= '0' && c <= '9') {
		x = x * 10 + c - '0'; c = getchar();
	}
	return x * f;
}
int prime[N],tot,mu[N],vis[N];
void pre() {
	mu[1] = 1;
	for(int i = 2;i <= 50000;++i) {
		if(!vis[i]) {
			prime[++tot] = i;
			mu[i] = -1;
		}
		for(int j = 1;j <= tot && prime[j] * i <= 50000;++j) {
			vis[prime[j] * i] = 1;
			if(i % prime[j])
				mu[i * prime[j]] = -mu[i];
			else break;
		}
	}
	for(int i = 1;i <= 50000;++i) mu[i] += mu[i - 1];
}
ll f[N];
int main() {
	pre();
	
	for(int i = 1;i <= 50000;++i) {
		for(int l = 1,r;l <= i;l = r + 1) {
			r = i / (i / l);
			f[i] += (r - l + 1) * (i / l);
		}
	}


	int T = read();
	while(T--) {
		int n = read(),m = read();
		if(n > m) swap(n,m);
		ll ans = 0;
		for(int l = 1,r;l <= n;l = r + 1) {
			r = min(n / (n / l),m / (m / l));
			ans += 1ll * (mu[r] - mu[l - 1]) * f[n / l] * f[m / l];
		}
		printf("%lld
",ans);
	}
	return 0;
}