原题链接:https://www.luogu.com.cn/problem/P3517
题意解读:有n个连续的点p1,p2,...,pn,将这n个点分成不超过m堆,每堆点连续,每一堆都缩成一个点qi,要使得原来的点p1~ps距离qi的最大值最小(最相似),求这个相似度,并计算一共分成几堆,以及每堆缩到的点qi的坐标。
解题思路:
要使得若干点缩到一个点,并且原来的点与新点距离的最大值最小,其实就是求以新点为圆心,能覆盖所有点的圆的最小半径。
1、最小圆覆盖
作为前置知识,先来了解一下最小圆覆盖问题:https://www.luogu.com.cn/problem/P1742
最小覆盖圆用到两个重要的性质:
性质一:最小覆盖圆唯一
性质二:若点 P不在点集 S的最小覆盖圆内部,则点 P必在P∪S的最小覆盖圆边上
由于过三点可以确定一个圆,可以枚举每三个点来确定圆是否包含其他点
但是直接枚举三个点会超时,为了规避掉最坏的情况,需要借助随机增量算法!
随机增量算法:
先将点的分布随机打散,可以借助random_shuffle函数,然后
取第一个点独立构成圆,圆心是第一个点,半径是0
for i:2->n
判断点i是否在当前圆内
如果不在圆内,以i点作为圆
for j:1->i-1
判断点j是否在当前圆内
如果不在圆内,以i,j两点构成圆
for k:1->j-1
判断点k是否在圆内
如果不在圆内,以i,j,k三点构成圆
两点构成圆:
点p1、p2构成圆的圆心是两点连线的中间点,半径是两点距离的一半
圆心坐标x = (p1.x + p2.x) / 2, y = (p1.y + p2.y) / 2
半径 = p1到圆心的距离
三点构成圆:
设三个点为(x1,y1)、(x2,y2)、(x3,y3),圆心为(x,y),半径为r,可得到方程组:
(x1-x)2 + (y1-y)2 = r2
(x2-x)2 + (y2-y)2 = r2
(x3-x)2 + (y3-y)2 = r2
展开后:
x12 - 2x1x + x2 + y12 - 2y1y + y2 = r2 (1)
x22 - 2x2x + x2 + y22 - 2y2y + y2 = r2 (2)
x32 - 2x3x + x2 + y32 - 2y3y + y2 = r2 (3)
用(1)-(2),(2)-(3),再移向后得到
2(x1-x2)x + 2(y1-y2)y = x12 - x22 + y12 - y22
2(x2-x3)x + 2(y2-y3)y = x22 - x32 + y22 - y32
设
a1 = (x1-x2),b1 = (y1-y2),c1 = x12 - x22 + y12 - y22
a2 = (x2-x3),b2 = (y2-y3),c2 = x22 - x32 + y22 - y32
解方程得到
x = (c1b2-c2b1)/(2a1b2-2a2b1)
y = (c1a2-c2a1)/(2a2b1-2a1b2)
半径是任意一个点到圆心的距离
100分代码(P1742):
#include <bits/stdc++.h>
using namespace std;
const int N = 100005;
const double eps = 1e-9;
struct Point
{
double x, y;
} a[N];
int n;
//计算两点之间距离
double distance(Point p1, Point p2)
{
return sqrt((p1.x-p2.x)*(p1.x-p2.x) + (p1.y-p2.y)*(p1.y-p2.y));
}
//两个点形成圆的圆心
Point circle(Point p1, Point p2)
{
return {(p1.x+p2.x) / 2, (p1.y+p2.y) / 2};
}
//三个点形成圆的圆心
Point circle(Point p1, Point p2, Point p3)
{
double x1 = p1.x, y1 = p1.y;
double x2 = p2.x, y2 = p2.y;
double x3 = p3.x, y3 = p3.y;
double a1 = x1 - x2, b1 = y1 - y2, c1 = x1*x1 - x2*x2 + y1*y1 - y2*y2;
double a2 = x2 - x3, b2 = y2 - y3, c2 = x2*x2 - x3*x3 + y2*y2 - y3*y3;
double x = (c1*b2 - c2*b1) / (2*a1*b2 - 2*a2*b1);
double y = (c1*a2 - c2*a1) / (2*a2*b1 - 2*a1*b2);
return {x, y};
}
int main()
{
cin >> n;
for(int i = 1; i <= n; i++) cin >> a[i].x >> a[i].y;
random_shuffle(a + 1, a + n + 1);
Point o = a[1]; //圆心
double r = 0; //半径
for(int i = 2; i <= n; i++)
{
if(distance(a[i], o) - r < eps) continue; //i点在圆内
o = a[i]; //更新圆为点i
for(int j = 1; j < i; j++)
{
if(distance(a[j], o) - r < eps) continue; //j点在圆内
o = circle(a[i], a[j]); //更新圆为i,j两点确定的圆
r = distance(a[i], o);
for(int k = 1; k < j; k++)
{
if(distance(a[k], o) - r < eps) continue; //k点在圆内
o = circle(a[i], a[j], a[k]); //更新圆为i,j,k三点确定的圆
r = distance(a[i], o);
}
}
}
cout << fixed << setprecision(10) << r << endl << o.x << " " << o.y;
}
2、二分答案求结果
回到本题,要将n个点分成不超过m堆,使得每堆的最小圆覆盖半径最小,可以借助二分。
二分答案mid,找尽可能多的点最小圆覆盖半径不超过mid,这些点即可分为一堆
如何找尽可能多的点分到一堆,可以一个一个枚举,直到最小圆覆盖半径超过mid,则分为一堆,最后统计堆数,不超过m则返回true,否则false
这样枚举+求最小圆覆盖总体复杂度就在O(n^2),再算上二分,复杂度一共是O(logn*n*n)
因此,还要继续优化!
3、倍增+二分找一堆的右端点
在找尽可能多的点分到一堆时,给定左端点l,可以通过倍增,先找到右端点r,使得最小圆覆盖半径超过mid
然后在l~r之间再二分准确的右端点R,这样l~R分为一堆
再更新l=R+1,继续找下一堆,最后记录堆数,看是否超过m
这样先倍增找一堆的范围,再在范围内二分找右端点,再计算最小圆覆盖半径,算法二分答案,总体复杂度在O(logn*logn*n)
100分代码:
#include <bits/stdc++.h>
using namespace std;
const int N = 1000005;
const double eps = 1e-9;
struct Point
{
double x, y;
} a[N], b[N], anso[N];
int cnt;
int n, m;
Point o; //圆心
double radius; //半径
//计算两点之间距离
double distance(Point p1, Point p2)
{
return sqrt((p1.x-p2.x)*(p1.x-p2.x) + (p1.y-p2.y)*(p1.y-p2.y));
}
//两个点形成圆的圆心
Point circle(Point p1, Point p2)
{
return {(p1.x+p2.x) / 2, (p1.y+p2.y) / 2};
}
//三个点形成圆的圆心
Point circle(Point p1, Point p2, Point p3)
{
double x1 = p1.x, y1 = p1.y;
double x2 = p2.x, y2 = p2.y;
double x3 = p3.x, y3 = p3.y;
double a1 = x1 - x2, b1 = y1 - y2, c1 = x1*x1 - x2*x2 + y1*y1 - y2*y2;
double a2 = x2 - x3, b2 = y2 - y3, c2 = x2*x2 - x3*x3 + y2*y2 - y3*y3;
double x = (c1*b2 - c2*b1) / (2*a1*b2 - 2*a2*b1);
double y = (c1*a2 - c2*a1) / (2*a2*b1 - 2*a1*b2);
return {x, y};
}
//找l~r之间点的最小覆盖圆圆心和半径
void find(int l, int r)
{
int idx = 0;
for(int i = l; i <= r; i++) b[++idx] = a[i];
random_shuffle(b + 1, b + idx + 1);
o = b[1]; //圆心
radius = 0; //半径
for(int i = 2; i <= idx; i++)
{
if(distance(b[i], o) - radius < eps) continue; //i点在圆内
o = b[i]; //更新圆为点i
for(int j = 1; j < i; j++)
{
if(distance(b[j], o) - radius < eps) continue; //j点在圆内
o = circle(b[i], b[j]); //更新圆为i,j两点确定的圆
radius = distance(b[i], o);
for(int k = 1; k < j; k++)
{
if(distance(b[k], o) - radius < eps) continue; //k点在圆内
o = circle(b[i], b[j], b[k]); //更新圆为i,j,k三点确定的圆
radius = distance(b[i], o);
}
}
}
}
bool check(double maxx)
{
cnt = 0;
for(int i = 1; i <= n; i++) //枚举左边界
{
int j;
for(j = 1; i + (1 << j) - 1 <= n; j++) //倍增找右边界
{
find(i, i + (1 << j) - 1);
if(radius - maxx > eps) break; //一段点的最小圆覆盖半径超过限定值
}
int l = i + (1 << j - 1) - 1, r = min(n, i + (1 << j) - 1), ans = i;
//二分找到准确的右端点,使得找到的一段最小圆覆盖半径都在maxx以内
while(l <= r)
{
int mid = (l + r) >> 1;
find(i, mid);
if(radius - maxx < eps) ans = mid, l = mid + 1;
else r = mid - 1;
}
find(i, ans); //更新圆心、半径
i = ans;
cnt++;
anso[cnt] = o;
if(cnt > m) return false;
}
return true;
}
int main()
{
scanf("%d%d", &n, &m);
for(int i = 1; i <= n; i++) scanf("%lf%lf", &a[i].x, &a[i].y);
find(1, n); //在所有点中计算最小圆覆盖的圆心、半径,得到的半径就是最大半径
double l = 0, r = radius;
//二分找答案
while(r - l > eps)
{
double mid = (l + r) / 2;
if(check(mid)) r = mid;
else l = mid;
}
check(r); //更新一次正确答案对应的圆心、半径等
printf("%.8lf\n", r);
printf("%d\n", cnt);
for(int i = 1; i <= cnt; i++) printf("%.8lf %.8lf\n", anso[i].x, anso[i].y);
return 0;
}
标签:p2,Plot,p1,Point,洛谷题,POI2011,double,x2,y2 From: https://www.cnblogs.com/jcwy/p/18429756