P4357 [CQOI2016]K 远点对

题目链接

P4357 [CQOI2016]K 远点对

题目描述

已知平面内 \(N\) 个点的坐标，求欧氏距离下的第 \(K\) 远点对。

两个点 \(P(x_1,y_1)\) 和 \(Q(x_2,y_2)\) 的欧氏距离定义为 \(\sqrt{(x_1-x_2)^2+(y_1-y_2)^2}\)

输入格式

输入文件第一行为用空格隔开的两个整数 \(N,K\) 。

接下来 \(N\) 行，每行两个整数 \(X,Y\) ，表示一个点的坐标。

输出格式

输出文件第一行为一个整数，表示第 \(K\) 远点对的距离的平方（一定是个整数）。

样例 #1

样例输入 #1

10 5
0 0
0 1
1 0
1 1
2 0
2 1
1 2
0 2
3 0
3 1

样例输出 #1

9

提示

对于 \(100\%\) 的测试点，\(N \le 100000,1 \le K \le 100,K \le \dfrac {N(N-1)}{2},0 \le X,Y < 2^{31}\)

解题思路

kdtree

同 P7883 平面最近点对（加强加强版），枚举点时，只需在 kdtree 上维护点到矩形的最远距离，同时用小根堆维护前 \(2k\) 大的点（由于求得点对会重复计算），最后小根堆栈顶元素即为答案

• 时间复杂度：\(O(n\times \sqrt{n})\)

代码

// Problem: P4357 [CQOI2016]K 远点对
// Contest: Luogu
// URL: https://www.luogu.com.cn/problem/P4357
// Memory Limit: 500 MB
// Time Limit: 3000 ms
// 
// Powered by CP Editor (https://cpeditor.org)

// %%%Skyqwq
#include <bits/stdc++.h>
 
// #define int long long
#define help {cin.tie(NULL); cout.tie(NULL);}
#define pb push_back
#define fi first
#define se second
#define mkp make_pair
using namespace std;
 
typedef long long LL;
typedef pair<int, int> PII;
typedef pair<LL, LL> PLL;
 
template <typename T> bool chkMax(T &x, T y) { return (y > x) ? x = y, 1 : 0; }
template <typename T> bool chkMin(T &x, T y) { return (y < x) ? x = y, 1 : 0; }
 
template <typename T> void inline read(T &x) {
    int f = 1; x = 0; char s = getchar();
    while (s < '0' || s > '9') { if (s == '-') f = -1; s = getchar(); }
    while (s <= '9' && s >= '0') x = x * 10 + (s ^ 48), s = getchar();
    x *= f;
}

const int N=1e5+5;
int n,k;
priority_queue<LL,vector<LL>,greater<LL>> q;
struct Point
{
	int x,y;
}point[N];
struct Tr
{
	int s[2];
	int U,D,L,R;
}tr[N];
bool cmp1(const Point &x,const Point &y)
{
	return x.x<y.x;
}
bool cmp2(const Point &x,const Point &y)
{
	return x.y<y.y;
}
void pushup(int x)
{
	tr[x].L=tr[x].R=point[x].x;
	tr[x].U=tr[x].D=point[x].y;
	if(tr[x].s[0])
	{
		tr[x].L=min(tr[x].L,tr[tr[x].s[0]].L);
		tr[x].D=min(tr[x].D,tr[tr[x].s[0]].D);
		tr[x].R=max(tr[x].R,tr[tr[x].s[0]].R);
		tr[x].U=max(tr[x].U,tr[tr[x].s[0]].U);
	}
	if(tr[x].s[1])
	{
		tr[x].L=min(tr[x].L,tr[tr[x].s[1]].L);
		tr[x].D=min(tr[x].D,tr[tr[x].s[1]].D);
		tr[x].R=max(tr[x].R,tr[tr[x].s[1]].R);
		tr[x].U=max(tr[x].U,tr[tr[x].s[1]].U);
	}
}
int build(int l,int r)
{
	if(l>r)return 0;
	if(l==r)
	{
		pushup(l);
		return l;
	}
	int mid=l+r>>1;
	double avg[2]={0},va[2]={0};
	for(int i=l;i<=r;i++)avg[0]+=point[i].x,avg[1]+=point[i].y;
	avg[0]/=(double)(r-l+1),avg[1]/=(double)(r-l+1);
	for(int i=l;i<=r;i++)
	{
		va[0]+=(point[i].x-avg[0])*(point[i].x-avg[0]);
		va[1]+=(point[i].y-avg[1])*(point[i].y-avg[1]);
	}
	if(va[0]>=va[1])nth_element(point+l,point+mid,point+r+1,cmp1);
	else
		nth_element(point+l,point+mid,point+r+1,cmp2);
	tr[mid].s[0]=build(l,mid-1),tr[mid].s[1]=build(mid+1,r);
	pushup(mid);
	return mid;
}
LL dist(Point x,Point y)
{
	return (LL)(x.x-y.x)*(x.x-y.x)+(LL)(x.y-y.y)*(x.y-y.y);
}
LL get(int x,int y)
{
	int dx=max(abs(point[y].x-tr[x].L),abs(point[y].x-tr[x].R));
	int dy=max(abs(point[y].y-tr[x].U),abs(point[y].y-tr[x].D));
	return (LL)dx*dx+(LL)dy*dy;
}
void ask(int l,int r,int x)
{
	if(l>r)return ;
	int mid=l+r>>1;
	if(x!=mid)
		q.push(dist(point[x],point[mid]));
	if(q.size()>k)q.pop();
	LL distL=get(tr[mid].s[0],x),distR=get(tr[mid].s[1],x);
	if(distL>q.top()&&distR>q.top())
	{
		if(distL>distR)
		{
			ask(l,mid-1,x);
			if(q.size()<k||distR>q.top())ask(mid+1,r,x);
		}
		else
		{
			ask(mid+1,r,x);
			if(q.size()<k||distL>q.top())ask(l,mid-1,x);
		}
	}
	else
	{
		if(q.size()<k||distL>q.top())ask(l,mid-1,x);
		if(q.size()<k||distR>q.top())ask(mid+1,r,x);
	}
}
int main()
{
    scanf("%d%d",&n,&k);
    k<<=1;
    for(int i=1;i<=n;i++)scanf("%d%d",&point[i].x,&point[i].y);
    build(1,n);
    for(int i=1;i<=n;i++)ask(1,n,i);
    printf("%lld",q.top());
    return 0;
}