常用距离算法

Liucf / 2024-11-17 / 原文

常用距离算法

对于两个点 \((x_1,y_1)\)\((x_2,y_2)\) 的距离大致有 \(3\) 种:

  1. 欧氏距离
  2. 曼哈顿距离
  3. 切比雪夫距离

三维情况下表示为 \((x_1,y_1,z_1)\)\((x_2,y_2,z_2)\)

多维情况下表示为 \((x_1,x_2,...,x_d)\)\((y_1,y_2,...,y_d)\),其中 \(d\) 表示维数(与二、三维表示有所不同)。

1、欧几里得距离

二维形式:\(E(A,B)=\sqrt{(x_1-x_2)^2+(y_1-y_2)^2}\) .

三维形式:\(E(A,B)=\sqrt{(x_1-x_2)^2+(y_1-y_2)^2+(z_1-z_2)^2}\).

多维形式:\(E(A,B)=\sqrt{\sum_{i=1}^{d}(x_i-y_i)^2}\).

P3958 [NOIP2017 提高组] 奶酪

更多的是并查集的应用,理论上 \((x_1-x_2)^2\) 可以达到 \(10^{20}\) 的大小,会爆 long long,但是这题数据过水?不会爆 。代码用了 __int128

#include <bits/stdc++.h>
using namespace std;
typedef __int128 int128;
const int N = 1e5 + 10;
int T, n, h, r, fa[N], down[N], up[N]; // 能否到下表面/上表面
struct node {
	int x, y, z;
} a[N];
int Find(int x) {
	return x == fa[x] ? x : fa[x] = Find(fa[x]);
}
void Union(int x, int y) {
	int fx = Find(x), fy = Find(y);
	up[fy] |= up[fx], down[fy] |= down[fx], fa[fx] = fy; // 合并状态
}
int128 dic(int i, int j) {
	int128 dx = a[i].x - a[j].x, dy = a[i].y - a[j].y, dz = a[i].z - a[j].z;
	return (int128) dx * dx + dy * dy + dz * dz; // 这里可能会爆 long long
}
int main() {
	ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
	cin >> T;
	while (T--) {
		cin >> n >> h >> r;
		int x, y, z;
		for (int i = 1; i <= n; i++) {
			up[i] = down[i] = 0, fa[i] = i;
			cin >> x >> y >> z;
			a[i] = {x, y, z};
			if (z <= r) down[i] = 1;
			if (h - z <= r) up[i] = 1;
		}
		for (int i = 1; i < n; i++)
			for (int j = i + 1; j <= n; j++)
				if (dic(i, j) <= (int128) 4 * r * r)
					Union(i, j);
		bool flag = false;
		for (int i = 1; i <= n; i++) {
			if (fa[i] == i && up[i] && down[i]) {
				flag = true;
				break;
			}
		}
		if (flag) cout << "Yes" << '\n';
		else cout << "No" << '\n';
	}
	return 0;
}

2、曼哈顿距离

二维形式:\(M(A,B)=|x_1-x_2|+|y_1-y_2|\) .

三维形式:\(M(A,B)=|x_1-x_2|+|y_1-y_2|+|z_1-z_2|\).

多维形式:\(M(A,B)=\sum_{i=1}^{d}|x_i-y_i|\).

代数含义

考虑二维形式,不妨设 \(x_1<x_2\)\(x_1>x_2\) 的情况是等价的:

  • \(y_1<y_2\) 时,\(M(A,B)=x_2-x_1+y_2-y_1=(x_2+y_2)-(x_1+y_1)\).
  • \(y_1>y_2\) 时,\(M(A,B)=x_2-x_1+y_1-y_2=(x_2-y_2)-(x_1-y_1)\).

发现只与 \(x\)\(y\)和与差有关,可以根据这一点做题。

P5098 [USACO04OPEN] Cave Cows 3 双倍经验 [ABC178E] Dist Max

给出 \(N\) 个点的坐标,求最大的曼哈顿距离。

还是分两种考虑,第一种情况时,\(d_{\max}=\max\{x+y\}-\min\{x+y\}\),第二种情况时,\(d_{\max}=\max\{x-y\}-\min\{x-y\}\),所以同时维护这 \(4\) 个最值,最后在两种情况下选更大的一个就行。

#include <bits/stdc++.h>
using namespace std;
const int INF = 1e9;
int n, x, y;
int main() {
	ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
	cin >> n;
	int a = -INF, b = -INF, c = INF, d = INF;
	for (int i = 1; i <= n; i++) {
		cin >> x >> y;
		a = max(a, x + y);
		b = max(b, x - y);
		c = min(c, x + y);
		d = min(d, x - y);
	}
	cout << max(a - c, b - d) << '\n';
	return 0;
}

几何含义

观察一下,不难发现点 \((2,3)\) 到这 \(2\) 条线上所有点的曼哈顿距离都为 \(2\)。不难证明:

  • 在平面直角坐标系中,与点 \(M(x_0,y_0)\) 的哈曼顿距离相等的点组成的集合为 \(k=\pm 1\) 的直线
曼哈顿距离1

回到这道题,可以发现答案一定由这两种情况构成:

  • 最左上的点和最右下的点的曼哈顿距离
  • 最左下的点和最右上的点的曼哈顿距离

其中最角落的点怎么定义?观察 \((5,1)\) 这个点到另外 \(2\) 个点的曼哈顿距离,发现点 \((2,5)\) 是最大的。

结合直线方程,发现最左上角的点所在直线就是 \(x-y=-3\),而另一个点的方程为 \(x-y=-2\),故记录 \(x-y\) 的最小值就可以表示出最左上的点 \((2,5)\) 所在直线。而对于其他点 \((5,1)\) 到这个点的哈曼顿距离,和到这条直线上的其他点的曼哈顿距离是一样的,所以只需要记录 \(\min\{x-y\}\) 就可以表示出这条直线,换句话说,最角落的点的坐标和它所在的 \(k=\pm 1\) 的直线是等价的,对于其他点,既可以用点的坐标计算曼哈顿距离,也可以用直线的方程计算曼哈顿距离。

回到这道题,只需要记录 \(4\) 条直线的方程:

  • 最左上角:\(x-y=k_1,k_1=\min\{x-y\}\),其中 \(k_1\) 对应上个代码的 \(d\)
  • 最左下角:\(x+y=k_2,k_2=\min\{x+y\}\),其中 \(k_2\) 对应上个代码的 \(c\)
  • 最右下角:\(x-y=k_3,k_3=\max\{x-y\}\),其中 \(k_3\) 对应上个代码的 \(b\)
  • 最右上角:\(x+y=k_4,k_4=\max\{x+y\}\),其中 \(k_4\) 对应上个代码的 \(a\)

这时考虑最左上的直线和最右下的直线的曼哈顿距离怎么计算:令 \(x_1=x_3\),则 \(y_1-y_3=k_3-k_1\) 即为其曼哈顿距离。同样的,最左下的直线和最右上的直线的曼哈顿距离为 \(k_4-k_2\),所以最后的答案为:

\[ans =\max(k_3-k_1,k_4-k_2)=\max(\max\{x-y\}-\min\{x-y\},\max\{x+y\}-\max\{x+y\}) \]

New York Hotel

和上一题一样,可以先找 \(4\) 个最角落的点所对的直线,然后再一个个比较打擂台就行。

#include <bits/stdc++.h>
#define int long long
using namespace std;
int n, m, a, b, c, d, x, y, pos, dic = 1e18;
signed main() {
	ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
	cin >> n >> m >> n;
	a = b = -1e18, c = d = 1e18;
	for (int i = 1; i <= n; i++) {
		cin >> x >> y;
		a = max(a, x + y), b = max(b, x - y), c = min(c, x + y), d = min(d, x - y);
	}
	cin >> m;
	for (int i = 1; i <= m; i++) {
		cin >> x >> y;
		int d1 = max(a - (x + y), (x + y) - c);
		int d2 = max(b - (x - y), (x - y) - d);
		if (max(d1, d2) < dic) {
			pos = i, dic = max(d1, d2);
		}
	}
	cout << dic << '\n' << pos << '\n';
	return 0;
}

3、切比雪夫距离

二维形式:\(C(A,B)=\max(|x_1-x_2|,|y_1-y_2|)\).

三维形式:\(C(A,B)=\max\{|x_1-x_2|,|y_1-y_2|,|z_1-z_2|\}\).

多维形式:\(C(A,B)=\max\{|x_i-y_i|\}\).

没找到模板题,但是其一个绝对值比曼哈顿距离的两个绝对值在排序、比较时显然更有优越性。但是在其他某些情况显然曼哈顿距离更好计算。所以这两者的相互转化对解题就极为关键了。

4、曼哈顿距离和切比雪夫距离的转换

转换指在将坐标系 \(P\) 中的各个点和坐标系 \(Q\) 中的各个点相互转换,以达到 \(P\) 中的各点之间的曼哈顿距离和 \(Q\) 中的各点之间的切比雪夫距离相等(或成倍数关系)。

代数证明

考虑化简 \(|a-b|\),其等于 \(\max(a-b,b-a)\)。因为这里面的 \(2\) 个数一个正,一个负,而其中的正数正好为绝对值的答案。

接下来拆开点 \(A(x_1,y_1)\)\(B(x_2,y_2)\) 两点的曼哈顿距离:

\[M(A,B)=|x_1-x_2|+|y_1-y_2|=\max(x_1-x_2,x_2-x_1)+\max(y_1-y_2,y_2-y_1) \]

\(x\)\(y\) 相互独立,可以将两个 \(\max\) 相加:

\[M(A,B)=\max\{x_1-x_2+y_1-y_2,x_2-x_1+y_1-y_2,x_1-x_2+y_2-y_1,x_2-x_1+y_2-y_1\}\tag1 \]

然后拆开点 \(C(x_3,y_3)\)\(D(x_4,y_4)\) 两点的切比雪夫距离(两个 \(\max\) 结合了):

\[C(C,D)=\max\{x_3-x_4,x_4-x_3,y_3-y_4,y_4-y_3\}\tag2 \]

观察到 \((1),(2)\) 式的形式很像,考虑能否进行代换,如果可以,就可以将两种距离相互转换了。首先观察到两式中的 \(4\) 个数两两为相反数,观察 \((1)\) 中的第 \(1,4\) 项和 \((2)\) 中的第 \(1,2\) 项, 有 \(C_4^2=6\) 种代换形式,但是为了找出一个固定的转换方式,即在两种坐标系中的一对点一一对应,采用这种转换方式:

\[\begin{cases} x_3=x_1+y_1\\ x_4=x_2+y_2\\ \end{cases} \]

这样的话 \(1\) 个坐标就仅由 \(1\) 个点的参数决定,同样的,可以这样代换剩下两项:

\[\begin{cases} y_3=x_1-y_1\\ y_4=x_2-y_2\\ \end{cases} \]

这样,就成功代换了两个坐标系中的点,总结一下结论:

在坐标系 \(P\) 中的 \(A(x_1,y_1)\)\(B(x_2,y_2)\) 两点的曼哈顿距离和在坐标系 \(Q\)\(C(x_1+y_1,x_1-y_1)\)\(D(x_2+y_2,x_2-y_2)\) 两点的切比雪夫距离相等。

同样的,\(Q\) 中的 \(C(x_3,y_3)\)\(D(x_4,y_4)\) 两点可以转换为 \(P\) 中的 \(A(\dfrac{x_3+y_3}{2},\dfrac{x_3-y_3}{2})\)\(B(\dfrac{x_4+y_4}{2},\dfrac{x_4-y_4}{2})\) .

观察到坐标除以 \(2\) 后可能会有小数,所以可以先不除以 \(2\),最后结果再除以 \(2\),注意可能会爆 long long

这种证明会比几何证明更加严谨(毕竟几何证明只是靠观察)。

几何证明

红色的正方形(坐标系 \(P\) 中)为与 \((0,0)\) 的曼哈顿距离为 \(1\),绿色的正方形(坐标系 \(Q\) 中)与 \((0,0)\) 的切比雪夫距离为 \(1\)

容易发现:将坐标系 \(P\) 旋转 \(45^{\circ}\),然后扩大到 \(\sqrt2\) 倍,就变成了 \(Q\) 坐标系。其中的点一一对应。

\(P\) (转换为极坐标系)中有点 \(A(r,\alpha)\),则(顺时针)旋转 \(\dfrac{\pi }{4}\) 后得到点 \(A'(r,\alpha-\dfrac{\pi }{4})\),扩大到 \(\sqrt2\) 倍后得到点 \(A''(\sqrt {2}r,\alpha-\dfrac{\pi }{4})\),观察 \(A\)\(A''\) 的直角坐标表示:

\[\begin{cases} A&(rcos\alpha,rsin\alpha)\\ A''&(\sqrt2rcos(\alpha-\dfrac{\pi }{4}),\sqrt2rsin(\alpha-\dfrac{\pi }{4}))\\ &(\sqrt2r(\dfrac{\sqrt2}{2}cos\alpha+\dfrac{\sqrt2}{2}sin\alpha),\sqrt2r(-\dfrac{\sqrt2}{2}cos\alpha+\dfrac{\sqrt2}{2}sin\alpha))\\ &(rcos\alpha+rsin\alpha,-rcos\alpha+rsin\alpha)\\ \end{cases} \]

可以发现就是上面的结论,只是 \(A''\) 关于 \(x\) 轴对称了一下。

P3964 [TJOI2013] 松鼠聚会

切比雪夫距离转为曼哈顿距离。可以直接转换,转换后就很简单了。按 \(x\)\(y\) 排序后,对于点 \(i\) 只要二分点 \(x_i\)\(y_i\) 同理)为第 \(k\)\(x\) 坐标,答案即为:

\[\sum_{j=1}^{k}(x_i-x_j)+\sum_{j=k}^{n}(x_j-x_i) \]

拆完式子之后用前缀和就可以 \(O(1)\) 出答案。总复杂度 \(O(n\log n)\)

#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N = 1e5 + 10;
int n, x[N], y[N], px[N], py[N];
struct node {
	int x, y;
} a[N];
signed main() {
	ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
	cin >> n;
	int r, t;
	for (int i = 1; i <= n; i++) {
		cin >> r >> t;
		x[i] = a[i].x = r + t, y[i] = a[i].y = r - t;
	}
	sort(x + 1, x + 1 + n);
	sort(y + 1, y + 1 + n);
	for (int i = 1; i <= n; i++) px[i] = px[i - 1] + x[i], py[i] = py[i - 1] + y[i]; // 处理前缀和
	int ans = 1e18;
	for (int i = 1; i <= n; i++) {
		int X = lower_bound(x + 1, x + 1 + n, a[i].x) - x;
		int Y = lower_bound(y + 1, y + 1 + n, a[i].y) - y;
		int wx = (2 * X - n) * a[i].x + px[n] - 2 * px[X];
		int wy = (2 * Y - n) * a[i].y + py[n] - 2 * py[Y];
		ans = min(ans, wx + wy);
	}
	cout << ans / 2 << '\n';
	return 0;
}

Four Coloring

一道很妙的曼哈顿转切比雪夫问题。转化后还要转化为 \(d\times d\) 的小正方形进行染色。这样边长的小正方形的八连通小正方形染成与之不同的颜色即可满足所有条件。

#include <bits/stdc++.h>
using namespace std;
int n, m, d;
char s[4] = {'R', 'Y', 'G', 'B'};
int main() {
	ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
	cin >> n >> m >> d;
	for (int i = 1; i <= n; i++) {
		for (int j = 1; j <= m; j++) {
			int x = (i + j) / d, y = (i - j + m) / d;
			cout << s[((y & 1) << 1) + (x & 1)];
		}
		cout << '\n';
	}
	return 0;
}

[JOISC 2021 Day2] 道路の建設案 (Road Construction)

考虑到曼哈顿距离的双绝对值很难进行排序等操作,转化为切比雪夫距离。整体二分距离 \(f\) 是否满足 \(k\) 对点距离小于 \(f\),按 \(x\) 排序后进行偏序问题,遍历 \(a_i\),双指针维护 \(x\) 方向的 \(\max\) 限制,然后 set 里二分 \(y\) 的范围,枚举 \([y-mid,y+mid]\) 范围内的点,每一个点都计入 \(cnt\),这样就有了 check 函数,最后 check(f - 1),不够的用 \(f\) 添上就行。

#include <bits/stdc++.h>
#define int long long
#define pii pair <int, int>
using namespace std;
const int N = 3e5 + 10;
int n, k, ans[N], cnt;
struct coor {
	int x, y;
	friend bool operator < (const coor &a, const coor &b) {return a.x == b.x ? a.y < b.y : a.x < b.x;}
} a[N];
struct node {
	int x, y;
	friend bool operator <(const node &a, const node &b) {return a.y < b.y;}
};
int check(int mid) {
	multiset <node> m;
	queue<int> q;
	cnt = 0;
	for (int r = 1; r <= n; r++) {
		while (!q.empty() && a[q.front()].x < a[r].x - mid) {
			auto it = m.find({a[q.front()].x, a[q.front()].y});
			m.erase(it), q.pop();
		}
		auto pos = m.lower_bound({0, a[r].y - mid});
		while (pos != m.end() && (*pos).y <= a[r].y + mid) {
			ans[++cnt] = max(abs(a[r].x - (*pos).x), abs((*pos).y - a[r].y));
			pos++;
		}
		m.insert({a[r].x, a[r].y}), q.push(r);
		if (cnt >= k) return cnt;
	}
	return cnt;
}
signed main() {
	ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
	cin >> n >> k;
	int x, y;
	for (int i = 1; i <= n; i++) {
		cin >> x >> y;
		a[i] = {x + y, x - y};
	}
	sort(a + 1, a + 1 + n);
	int l = -1e18, r = 1e18, mid = 0, f = 0;
	while (l <= r) {
		mid = (l + r) >> 1;
		if (check(mid) >= k) r = mid - 1, f = mid;
		else l = mid + 1;
	}
	check(f - 1);
	sort(ans + 1, ans + 1 + cnt);
	for (int i = 1; i <= cnt; i++) cout << ans[i] << '\n';
	for (int i = cnt + 1; i <= k; i++) cout << f << '\n';
	return 0;
}

P5193 [TJOI2012] 炸弹

这道题和上一道题一样的,转换后问题转变为求 \(r\times r\) 的正方形连通块个数,再偏序求个数就行,但是对于 \(y_i\) 只要找最接近的 \(y\) 即可,注意 set 的迭代器 pos 操作后是否为 s.begin()s.end(),以及 pos-- 的顺序。复杂度 \(O(\log n)\)

#include <bits/stdc++.h>
#define int long long
#define pii pair <int, int>
using namespace std;
const int N = 1e5 + 10;
int n, r, fa[N], j = 1, x, y, res;
struct noor {
	int x, y;
	friend bool operator < (const noor &a, const noor &b) {
		return a.x == b.x ? a.y < b.y : a.x < b.x;
	}
} a[N];
int Find(int x) {return x == fa[x] ? x : fa[x] = Find(fa[x]);}
void Union(int x, int y) {fa[Find(x)] = Find(fa[y]);}
signed main() {
	ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
	cin >> n >> r;
	for (int i = 1; i <= n; i++) {
		cin >> x >> y;
		a[i] = {x + y, x - y}, fa[i] = i;
	}
	sort(a + 1, a + 1 + n);
	multiset <pii> s;
	for (int i = 1; i <= n; i++) {
		while (j < i && a[j].x < a[i].x - r) s.erase(s.find({a[j].y, j})), j++;
		auto pos = s.lower_bound({a[i].y, i});
		if (pos != s.end() && (*pos).first <= a[i].y + r) Union((*pos).second, i);
		if (pos != s.begin()) {
			pos--;
			if ((*pos).first >= a[i].y - r) Union((*pos).second, i);
		}
		s.insert({a[i].y, i});
	}
	for (int i = 1; i <= n; i++) res += (Find(i) == i);
	cout << res << '\n';
	return 0;
}

P3439 [POI2006] MAG-Warehouse

这道题转换完后和这道 P10452 货仓选址 的思路是一模一样的,就是找中位数。只不过这道题找到的点在原坐标系中可能不是整数,要再判一下周围的点。

#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N = 1e5 + 10;
typedef unsigned long long ull;
int n, X[N], Y[N];
struct node {
	int x, y, t;
} u[N], a[N];
bool cmp1(node a, node b) {return a.x == b.x ? a.y < b.y : a.x < b.x;}
bool cmp2(node a, node b) {return a.y == b.y ? a.x < b.x : a.y < b.y;}
int init(int X[N], int cmp) {
	if (cmp == 1) sort(a + 1, a + 1 + n, cmp1);
	else sort(a + 1, a + 1 + n, cmp2);
	for (int i = 1; i <= n; i++) X[i] = X[i - 1] + a[i].t;
	int ax = upper_bound(X + 1, X + 1 + n, X[n] / 2) - X;
	if (cmp == 1) return a[ax].x;
	else return a[ax].y;
}
ull check(int x, int y) {
	ull ans = 0ull;
	for (int i = 1; i <= n; i++) ans += 1ull * max(abs(u[i].x - x), abs(u[i].y - y)) * u[i].t;
	return ans;
}
signed main() {
	ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
	cin >> n;
	int x, y, t;
	for (int i = 1; i <= n; i++) {
		cin >> x >> y >> t;
		u[i] = {x, y, t}, a[i] = {x + y, x - y, t};
	}
	int ax = init(X, 1), ay = init(Y, 2);
	x = (ax + ay) / 2, y = (ax - ay) / 2;
	ull s1 = check(x + 1, y + 1);
	ull s2 = check(x + 1, y);
	ull s3 = check(x, y + 1);
	ull s4 = check(x, y);
	ull m = min(min(s1, s2), min(s3, s4));
	if (s1 == m) cout << x + 1 << ' ' << y + 1 << '\n';
	else if (s2 == m) cout << x + 1 << ' ' << y << '\n';
	else if (s3 == m) cout << x << ' ' << y + 1 << '\n';
	else if (s4 == m) cout << x << ' ' << y << '\n';
	return 0;
}

5、以上问题的三维及更高维

坐标转换

三维形式

\(A(x_1,y_1,z_1),B(x_2,y_2,z_2)\) 为坐标系 \(P\) 中的两点,则:

\[\begin{align} M(A,B)=&|x_1-x_2|+|y_1-y_2|+|z_1-z_2| \\ =&\max(x_1-x_2,x_2-x_1)+\max(y_1-y_2,y_2-y_1)+\max(z_1-z_2,z_2-z_1) \\ \end{align} \]

又因为 \(x,y,z\) 相互独立,所以:

\[\begin{align} M(A,B)=\max\{&x_1-x_2+y_1-y_2+z_1-z_2,\\ -&x_1+x_2-y_1+y_2-z_1+z_2,\\ &x_1-x_2+y_1-y_2-z_1+z_2,\\ -&x_1+x_2-y_1+y_2+z_1-z_2,\\ &x_1-x_2-y_1+y_2+z_1-z_2,\\ -&x_1+x_2+y_1-y_2-z_1+z_2,\\ &x_1-x_2-y_1+y_2-z_1+z_2,\\ -&x_1+x_2+y_1-y_2-z_1+z_2\}\\ \end{align} \]

按照和二维一样的配凑法,找相反数,然后点对点配对,一共 \(4\) 对:

\[\begin{align} M(A,B)=\max\{&|x_1+y_1+z_1|-|x_2+y_2+z_2|,\\ &|x_1+y_1-z_1|-|x_2+y_2-z_2|,\\ &|x_1-y_1+z_1|-|x_2-y_2+z_2|,\\ &|x_1-y_1-z_1|-|x_2-y_2-z_2|\} \end{align} \]

而这就是 \(A’(x_1+y_1+z_1,x_1+y_1-z_1,x_1-y_1+z_1,x_1-y_1-z_1),B’(x_2+y_2+z_2,x_2+y_2-z_2,x_2-y_2+z_2,x_2-y_2-z_2)\) 两点的四维切比雪夫距离,所以三维形式的结论为:

\[P_A(x,y,z)\to Q_D(x+y+z,x+y-z,x-y+z,x-y-z) \]

容易发现为除了 \(x\)\((y,z)\)\(2^2\) 个排列组合。几何证明难以描述,不再提及。

高维形式

定义 \([i,j]\) 表示 \(i\) 的二进制表示 \((i)_2\) 的从右往左数第 \(j\) 位(不够的补 \(0\))。

\(d\) 维空间内有形如 \(X(t_1,t_2,...,t_d)\)\(A,B\) 两点,则其曼哈顿距离可以表示为:

\[\begin{align} M(A,B)&=\sum_{i=0}^d|a_i-b_i|\\ &=\max_{i=0}^{2^d-1}\{\sum_{j=1}^{d}{(-1)^{[i,j]}(a_j-b_j)}\}\\ &=\max_{i=0}^{2^d-1}\{ \sum_{j=1}^d{(-1)^{[i,j]}a_j} -\sum_{j=1}^d{(-1)^{[i,j]}b_j} \}\\ &=\max_{i=0}^{2^{d-1}-1}\{ \max(\sum_{j=1}^d{(-1)^{[i,j]}a_j} -\sum_{j=1}^d{(-1)^{[i,j]}b_j},\sum_{j=1}^d{(-1)^{[2^{d-1}-1-i,j]}a_j} -\sum_{j=1}^d{(-1)^{[2^{d-1}-1-i,j]}b_j})\}\\ &=\max_{i=0}^{2^{d-1}-1}\{ \max(\sum_{j=1}^d{(-1)^{[i,j]}a_j} -\sum_{j=1}^d{(-1)^{[i,j]}b_j},\sum_{j=1}^d{(-1)^{[i,j]+1}a_j} -\sum_{j=1}^d{(-1)^{[i,j]+1}b_j})\}\\ &=\max_{i=0}^{2^{d-1}-1}\{ \max(\sum_{j=1}^d{(-1)^{[i,j]}a_j} -\sum_{j=1}^d{(-1)^{[i,j]}b_j},-(\sum_{j=1}^d{(-1)^{[i,j]}a_j} -\sum_{j=1}^d{(-1)^{[i,j]}b_j}))\}\\ &=\max_{i=0}^{2^{d-1}-1}\{ |\sum_{j=1}^d{(-1)^{[i,j]}a_j} -\sum_{j=1}^d{(-1)^{[i,j]}b_j}|\}\\ &=\max_{i=0}^{2^{d-1}-1}\{ |\sum_{j=1}^d{(-1)^{[i,j]}{(a_j-b_j)}}|\}\\ &=C(A',B') \end{align} \]

其中,\(A',B'\) 为形如 \(X'(T_0,T_1,...,T_{2^{d-1}-1})\),其中 \(T_i=\sum_{j=1}^d{(-1)^{[i,j]}{t_j}}\)\(2^{d-1}\) 维点。

不难发现:由于 \(i\in[0,2^{d-1}-1]\),所以 \([i,d]=0\),所以 \(T_i\) 的第 \(d\)\(t_d\) 的系数一定为正。而其他的 \(t_j\) 的系数可正可负,有 \(2^{d-1}\) 种排列组合,也就对应共 \(2^{d-1}\) 个多项式 \(T_i\)

总结一下结论:\(d\) 维的曼哈顿距离可以转换为 \(2^{d-1}\) 维的切比雪夫距离,而转换方式就是上面的 \(X(t_1,t_2,...,t_d)\) 转为 \(X'(T_0,T_1,...,T_{2^{d-1}-1})\) 的方式。

统计切比雪夫距离小于 \(K\) 的点对数

可以发现上面处理二维问题的复杂度为 \(O(n\log n)\),思考一下为什么?首先是保持一维 \(x\) 的单调性,然后用数据结构处理 \([y_i-K,y_i+K]\) 范围内点的个数。不妨放弃一切,尝试仅仅维护以点 \(i\) 为中心的 \(2K\times 2K\) 的正方形内有共 \(cnt_i\) 个点,则总的答案为 \(\dfrac {\sum_{i=1}^{n}cnt_i-n}{2}\),而这要用前缀和维护,时间复杂度 \(O(A^2+n)\)(预处理复杂度和查询复杂度),空间复杂度 \(O(A^2)\) ,其中 \(A\) 为值域。

拓展到 \(d\) 维,则时间复杂度为 \(O(A^d+n)\),空间复杂度 \(O(A^d)\),都很容易爆炸。考虑优化,由于点数往往少于 \(10^6\),而值域很大,考虑用树状数组维护多维前缀和,query 时用容斥,这样时间复杂度可以降到 \(O(n\log^d A+n\log^dA)\)。再考虑优化空间复杂度,可以偏序一维,这样可以省去一维的空间,此时时间复杂度 \(O(n\log^{d-1}A+n\log^{d}A)\),空间复杂度 \(O(A^{d-1})\)

P4648 [IOI2007] pairs 动物对数

根据上面的结论对三个维度分别做就行了。

#include <bits/stdc++.h>
#define int long long
#define INF 0x3f3f3f3f
using namespace std;
const int N = 1e5 + 10, M = 75010, P = 80;
int B, n, d, m, a[N], c[M * 2], t[P * 3][P * 3][P * 3];
struct noor {
	int x, y;
	friend bool operator < (const noor &a, const noor &b) {
		return a.x == b.x ? a.y < b.y : a.x < b.x;
	}
} b[N];
struct node {
	int a, b, c, d;
	friend bool operator < (const node &a, const node &b) {
		return a.a < b.a;
	}
} u[N];
int lowbit(int x) {
	return x& -x;
}
void add(int i, int x) {
	for (; i < M * 2; i += lowbit(i)) c[i] += x;
}
int query(int i) {
	i = min(i, 2 * M - 1);
	int ans = 0;
	for (; i > 0; i -= lowbit(i)) ans += c[i];
	return ans;
}
inline void Add(int x, int y, int z, int v) {
	for (int i = x; i < P * 3; i += i & -i) {
		for (int j = y; j < P * 3; j += j & -j) {
			for (int k = z; k < P * 3; k += k & -k) {
				t[i][j][k] += v;
			}
		}
	}
}
inline int Query(int x, int y, int z) {
	x = min(x, P * 3 - 1), y = min(y, P * 3 - 1), z = min(z, P * 3 - 1);
	int res = 0;
	for (int i = x; i > 0; i -= i & -i) {
		for (int j = y; j > 0; j -= j & -j) {
			for (int k = z; k > 0; k -= k & -k) {
				res += t[i][j][k];
			}
		}
	}
	return res;
}
int ask(int x1, int x2, int y1, int y2, int z1, int z2) {
	int ans = Query(x2, y2, z2);
	ans -= Query(x1, y2, z2) + Query(x2, y1, z2) + Query(x2, y2, z1);
	ans += Query(x1, y1, z2) + Query(x1, y2, z1) + Query(x2, y1, z1);
	ans -= Query(x1, y1, z1);
	return ans;
}
signed main() {
	ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
	cin >> B >> n >> d >> m;
	if (B == 1) {
		for (int i = 1; i <= n; i++) cin >> a[i];
		sort(a + 1, a + 1 + n);
		int pos = 1, ans = 0;
		for (int i = 1; i <= n; i++) {
			while (pos <= n && a[pos] - a[i] <= d) pos++;
			ans += pos - i - 1;
		}
		cout << ans << '\n';
	} else if (B == 2) {
		int x, y, j = 1, cnt = 0;
		for (int i = 1; i <= n; i++) {
			cin >> x >> y, b[i] = {x + y, x - y};
		}
		sort(b + 1, b + 1 + n);
		for (int i = 1; i <= n; i++) {
			while (j < i && b[j].x < b[i].x - d) add(b[j].y + m, -1), j++;
			cnt += query(b[i].y + d + m) - query(b[i].y - d - 1 + m);
			add(b[i].y + m, 1);
		}
		cout << cnt << '\n';
	} else if (B == 3) {
		int x, y, z, j = 1, cnt = 0;
		for (int i = 1; i <= n; i++) {
			cin >> x >> y >> z;
			u[i] = {x + y + z, x + y - z, x - y + z, -x + y + z};
		}
		sort(u + 1, u + 1 + n);
		for (int i = 1; i <= n; i++) {
			while (j < i && u[j].a < u[i].a - d) {
				Add(u[j].b + m, u[j].c + m, u[j].d + m, -1);
				j++;
			}
			cnt += ask(u[i].b - d - 1 + m, u[i].b + d + m, u[i].c - d - 1 + m, u[i].c + d + m, u[i].d - d - 1 + m, u[i].d + d + m);
			Add(u[i].b + m, u[i].c + m, u[i].d + m, 1);
		}
		cout << cnt << '\n';
	}
	return 0;
}