29.1 最近点对
“最近”是指通常意义上的欧几里德距离。《算法导论》中对这个问题进行了介绍。[i]问对于经典的分治法求最近点对有深入的介绍和详尽的伪代码。
29.1.1 实例
PKU JudgeOnline, 3714, Raid.
29.1.2 问题描述
给定n个岗位和n个战士的坐标(0 ≤ X ≤1000000000, 0 ≤ Y ≤ 1000000000)。问这些战士当中，离某个岗位最近的距离是多少？
29.1.3 输入
2
4
00
01
10
11
22
23
32
33
4
00
00
00
00
00
00
00
0 0
29.1.4 输出
1.414
0.000
29.1.5 分析
这也是一个最近距离的问题，只不过是两个集合中的元素的最近距离的问题。采用的方法不同于经典的最近距离点对的分治算法，但是思想上也有类似的地方。
对于每个战士(x0, y0)，首先找到一个点，这个点是x最接近x0的所有点中y最接近y0的点，求出(x, y)、(x0, y0)两点距离dis。若minDis > dis，令minDis = dis。在[x - minDis, x + minDis]的横坐标范围内移动x，并对每个x，寻找y最接近y0的点(x, y)，若minDis > dis，令minDis = dis。
上述描述中，“最接近x0的x”或者“最接近y0的y”都通过对有序数组的二分查找法来找出。
29.1.6 程序
#include<stdio.h>
#include<string.h>
#include <stdlib.h>
#include <math.h>
#define ONLINE_JUDGE 0
#define maxNum 100001
#define MAX 0x7FFFFFFF
#define __int64MAX 0x7FFFFFFFFFFFFFFF
struct position
{
int x;
int y;
};
int cmp(const void*p1, const void*p2)
{
int y1;
int y2;
int x1;
int x2;
y1 = (*(position *)p1).y;
y2 = (*(position *)p2).y;
x1 = (*(position *)p1).x;
x2 = (*(position *)p2).x;
if(x1 >x2)
{
return1;
}else if(x1 < x2){
return-1;
}else{
if(y2< y1){
return1;
}else{
return-1;
}
}
}
positionstation[maxNum];
positionagent[maxNum];
__int64 minDis;
struct coordinate{
int pos;
int x;
};
int findMinY(int begin, int end, int t)
{
int y;
int low,high, mid;
low = begin;
high = end - 1;
while(low< high){
mid = (low + high) / 2;
if(t> station[mid].y){
low = mid + 1;
}else{
high = mid;
}
}
//现在有关系: t< station[high].y && t > station[high - 1].y
y = station[high].y - t;
//如果直接check[high+ 1].pos可能会越界，因为high有可能等于top
if(high> begin && y > t - station[high - 1].y)
{
y = t - station[high - 1].y;
}
return y;
}
coordinatecheck[maxNum];
int top;
int N;
void cal(position now)
{
int low,high, mid;
int t;
int x;
int y;
int temp;
int begin;
int end;
int pos;
__int64bigTemp;
t = now.x;
low = 0;
high = top;
while(low< high){
mid = (low + high) / 2;
if(t> check[mid].x){
low = mid + 1;
}else{
high = mid;
}
}
//现在有关系: t< check[high].x && t > check[high - 1].x
x = check[high].x - t;
begin = check[high].pos;
end = high + 1;
pos = high;
//如果直接check[high+ 1].pos可能会越界，因为high有可能等于top
if(high> 0 && x > t - check[high - 1].x)
{
pos = high - 1;
x = t - check[pos].x;
begin = check[pos].pos;
end = high;
}
end = check[end].pos;
//这下即使士兵的横坐标大于所有电站的横坐标，check也不可能越界了，
//因为只要top不为（即电站数目不为，这是可以保证的），
//就存在一个电站，使得电站与士兵的横坐标之差小于MAX - t的距离。
y = findMinY(begin, end, now.y);
bigTemp = ((__int64)x)* ((__int64)x)+ ((__int64)y)* ((__int64)y);
if(minDis> bigTemp){
minDis = bigTemp;
}
temp = pos - 1;
while(temp>= 0){
x = now.x - check[temp].x;
bigTemp = ((__int64)x)* ((__int64)x);
if(bigTemp> minDis)
{
break;
}
begin = check[temp].pos;
end = check[temp + 1].pos;
y = findMinY(begin, end, now.y);
bigTemp += ((__int64)y)* ((__int64)y);
if(minDis> bigTemp){
minDis = bigTemp;
}
temp--;
}
temp = pos + 1;
while(temp< top){
x = check[temp].x - now.x;
bigTemp = ((__int64)x)* ((__int64)x);
if(bigTemp> minDis)
{
break;
}
begin = check[temp].pos;
end = check[temp + 1].pos;
y = findMinY(begin, end, now.y);
bigTemp += ((__int64)y)* ((__int64)y);
if(minDis> bigTemp){
minDis = bigTemp;
}
temp++;
}
}
int main()
{
#ifndef ONLINE_JUDGE
FILE *fin;
fin = freopen("test.txt", "r", stdin);
if( !fin )
{
printf("reopen in file failed...\n");
while(1){}
return 0;
}
freopen("out.txt", "w", stdout);
#endif
int cases;
int i;
int temp;
scanf("%d",&cases);
for(; cases> 0; cases--){
scanf("%d",&N);
for(i =0; i < N; i++){
scanf("%d%d",&station[i].x, &station[i].y);
}
for(i =0; i < N; i++){
scanf("%d%d",&agent[i].x, &agent[i].y);
}
qsort(station, N, sizeof(station[0]), cmp);
top = 0;
temp = -MAX;
for(i =0; i < N; i++){
if(temp!= station[i].x)
{
temp = station[i].x;
check[top].x = temp;
check[top].pos = i;
top++;
}
}
check[top].x = MAX;
check[top].pos = N;
//check设立的最后这个元素，有诸多好处
minDis = __int64MAX;
for(i =0; i < N; i++){
cal(agent[i]);
}
//printf("%d\n",minDis);
printf("%.3f\n",sqrt((double)minDis));
}
#ifndef ONLINE_JUDGE
// fclose(stdout);
fclose(stdin);
#endif
return 1;
}关于线段相交，《算法艺术与信息学竞赛》介绍得比较细致。
29.2.1 实例
PKU JudgeOnline, 1039, Pipe.
29.2.2 问题描述
如上图所示，有一条管子，拐角处的两点坐标为(x, y)和(x, y + 1)。有一束光从入口射入，遇到管壁则被阻挡。求光所能达到的最大的x。
29.2.3 输入
4
01
22
41
64
6
01
2-0.6
5-4.45
7-5.57
12-10.8
17-16.55
0
29.2.4 输出
4.67
Through all the pipe.
29.2.5 分析
《算法艺术与信息学竞赛》对这个题有结题报告。
29.2.6 程序
#include <stdio.h>
#include <math.h>
typedef struct {
double x,y;
}point;
int dblcmp (double d) //判断d的符号
{
if (fabs(d) < 1e-6) return 0;
return (d> 0) ? 1 : -1;
}
void getline (point a, point b, double&k, double &t, double&c)
{ //根据点a，b求直线kx+ty+c=0
k = a.y - b.y;
t = b.x - a.x;
c = a.x*b.y - b.x*a.y;
return ;
}
int cal (double k, doublet, double c, point a)
{//计算点a在直线kx+ty+c=0的上方还是下方还是直线上
returndblcmp (k*a.x + t*a.y + c);
}
int n;
const int size = 21;
pointup[size], down[size];
double max;
void updata (double k, double t, double c,point a, point b)
{//更新最大值
double k2,t2, c2;
getline (a, b, k2, t2, c2);
double xp =(t*c2 - t2*c) / (k*t2 - t*k2);
if (xp >max)
{
max = xp;
}
return ;
}
bool judge (point a, point b)
{ // 判断点a，b构成的直线是否能通过整条管道
double k,t, c;
getline (a, b, k, t, c);
//判断光线能够从第一个通道口射入
if (cal(k,t, c, up[0]) * cal(k, t, c, down[0]) > 0) returnfalse;
int i;
//寻找第一个拐角，在这个拐角前，光线已穿越管壁
for (i = 1;i < n; i++)
if(cal(k, t, c, up[i]) * cal(k, t, c, down[i]) > 0) break;
if (i == n)return true;
if (cal(k,t, c, up[i]) * cal(k, t, c, up[i-1]) <= 0) updata (k, t, c, up[i], up[i-1]);
if (cal(k,t, c, down[i]) * cal(k, t, c, down[i-1]) <= 0) updata (k, t, c, down[i],down[i-1]);
return false;
}
bool solve ()
{
for (int i = 0; i < n; i++)
for (int j = 0; j < n; j++){
if(i == j) continue;
if(judge (up[i], down[j])) return true;
}
return false;
}
int main ()
{
// freopen("in.txt", "r", stdin);
while(scanf ("%d", &n), n){
for (int i = 0; i < n; i++){
scanf ("%lf%lf",&up[i].x, &up[i].y);
down[i].x = up[i].x; down[i].y =up[i].y - 1;
}
max = up[0].x;
if(solve ()) printf ("Through all thepipe.\n");
elseprintf ("%.2f\n", max);
}
}29.3 Melkman在线凸包算法
《算法艺术与信息学竞赛》比较全面、生动地介绍了凸包算法及其应用。
Melkman凸包算法继承Graham-Scan算法的主要思想，并更近一步地采用双端队列，动态地在队列两头进行增删操作，维护“凸性”。
一个Melkman凸包算法的实例如下：
[ii]中描述的Melkman凸包算法如下所示：
1. t ← 1; b ← 0;//这里似乎Melkman弄反了，应该改成：t ← 0;b ←1;
vl← input; v2 ← input; v3 ← input;
if (Vl, v2, v3) > 0
then begin push vl; push v2; end;
else begin push v2; push v1; end;
pushv3 ; insert v3;
2. v ← input;
until (v, db, db+l) < 0 or (dt-1, dt, v) < 0
do v ← inputend;
3. until (dt-l,dt, v) > 0 do pop dt end;
pushv;
4. until(v, db, db+l) > 0 do remove db end;
insertv;
goto 2.
其中定义：
(x, y, z)为根据z在x到y的向量右边、同一条直线上、或左边，而返回1, 0, -1的函数。
push操作就是使得t = t + 1; dt = v。
pop操作就是使得 t = t - 1。
insert操作就是使得b = b - 1; db = v。
remove操作就是使得b = b + 1。
29.3.1 实例
PKU JudgeOnline, 1113, Wall.
29.3.2 问题描述
已知N个点的坐标，要将这些点用围墙围起来，并使得围墙距离点的距离不超过L。求围墙的周长。输出四舍五入。
先输入N和L，然后是N个点的坐标。
1.3.3 输入
9100
200400
300400
300300
400300
400400
500400
500200
350200
200 200
29.3.4 输出
1628
29.3.5 分析
先求这些点的凸包。
可以很容易地证明：周长最短的围墙是由由于凸包的每条边，垂直向外平移L的距离，再将这些线段用半径为L的圆弧连接起来而构成的。
同时又很容易证明：这些圆弧的角度之和为360度。
故此这个问题直接转化为求凸包的周长问题。
凸包的求法中三点共线是十分需要注意的。Melkman算法描述中略去了对三点共线的许多细节的处理。
29.3.6 程序
#include<stdio.h>
#include<string.h>
#include<math.h>
#include<iostream>
using namespace std;
#define maxNum 2010
#define BEGIN 1005
#define PI 3.1415926
struct Vector{
int x;
int y;
};
Vectord[maxNum];
inline int crossProduct(Vector p1,Vector p2)
{
return p1.x* p2.y - p1.y * p2.x;
}
inline int fun(Vector v1,Vectorv2,Vector v3)
{
Vector tempv1;
Vector tempv2;
tempv1.x = v2.x - v1.x;
tempv1.y = v2.y - v1.y;
tempv2.x = v3.x - v1.x;
tempv2.y = v3.y - v1.y;
returncrossProduct(tempv1, tempv2);
}
inline double dis(Vector p1, Vectorp2)
{
double x;
double y;
x = p2.x - p1.x;
y = p2.y - p1.y;
x = x * x;
y = y * y;
returnsqrt(x + y);
}
inline int between(Vector p0,Vector p1, Vector p2)
{
int x1;
int x2;
if(p1.x> p2.x){
x1 = p1.x;
x2 = p2.x;
}else{
x1 = p2.x;
x2 = p1.x;
}
if(p0.x< x1 && p0.x > x2){
return1;
}else{
return0;
}
}
#define push(v) t++;d[t]= v;
#define pop() t--;
#define insert(v) b--;d[b] =v;
#define remove() b++;
int N;
void Melkman()
{
int t;
int b;
int i;
double sum;
int temp;
Vector v, v1, v2, v3;
int L;
cin >> L;
cin >> v1.x >> v1.y;
cin >> v2.x >> v2.y;
cin >> v3.x >> v3.y;
t = 0;
b = 1;
t += BEGIN;
b += BEGIN;
if(fun(v1,v2, v3) > 0){
push(v1);
push(v2);
}else{
push(v2);
push(v1);
}
push(v3);
insert(v3);
for(i = 3;i < N; i++){
cin >> v.x >> v.y;
/*
//如果有坐标相同的点，需要加入如下的判断
for(j = b; j< t; j++){
if(v.x ==d[j].x && v.y == d[j].y)
{
//int*a = (int *)0;
//a =0;
break;
}
}
if(j < t)
{
continue;
}
*/
//if(!(fun(v,d[b], d[b + 1]) < 0 || fun(d[t - 1], d[t], v) < 0))
//把共点的连同共线的都过滤了，过滤点的条件太松
if((fun(v,d[b], d[b + 1]) > 0 && fun(d[t - 1], d[t], v) > 0))
//把共点的连同共线的都没过滤，过滤点的条件太紧
{
//这里必须如此，不能是上面的，否则共线的出错
continue;
}
if(fun(v,d[b], d[b + 1]) == 0 && fun(d[t - 1], d[t], v) == 0)
{
if(between(v,d[b], d[b + 1]) || between(v, d[t - 1], d[t]))
{
//这里进行的是将共线的，但是处在已得凸包边上的点滤除
continue;
}
}
while(t- b >= 2 && fun(d[t - 1], d[t], v) <= 0){
//只要共线的都删除重新加了，直到只剩下两个点为止
pop();
// cout<< "poped"<< endl;
}
push(v);
while(t- b >= 2 && fun(v, d[b], d[b + 1]) <= 0){
//只要共线的都删除重新加了，直到只剩下两个点为止
remove();
// cout<< "removed"<<endl;
}
insert(v);
}
sum = 0;
for(i = b;i < t; i++){
sum += dis(d[i], d[i + 1]);
}
sum += 2 * PI * L;
temp = (int)sum;
sum -= temp;
if(sum>= 0.5)
{
temp++;
}
cout << temp << endl;
}
int main()
{
while(cin>> N){
Melkman();
}
}29.4 实例
29.4.1 线段相交的实例
PKU JudgeOnline, 1039, Pipe.
29.4.2 凸包的实例
PKU JudgeOnline, 1113, Wall.
本文章欢迎转载，请保留原始博客链接http://blog.csdn.net/fsdev/article
154
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
