摘要:
网上有很多关于分治方法求最近点对的讨论，但是没有完整的可运行代码，本文主要对于该问题介绍一完整的可运行代码，
供有兴趣者参考。
正文：
作为对比，我们也同时实现了最近点对的枚举解法，具体的主函数如下：
#include<stdio.h>
#include<stdlib.h>
#include "node.h"
void initList(node* p)
{
p[0].x= 2.0;
p[0].y= 1.0;
p[1].x= 1.0;
p[1].y= 2.0;
p[2].x= 1.2;
p[2].y= 3.0;
p[3].x= 3.0;
p[3].y= 4.0;
p[4].x= 5.0;
p[4].y= 5.0;
p[5].x= 1.5;
p[5].y= 5.5;
p[6].x= 2.5;
p[6].y= 7.0;
p[7].x= 3.5;
p[7].y= 8.0;
p[8].x= 4.0;
p[8].y= 9.0;
p[9].x = 3.9;
p[9].y = 8.8;
}
//测试分治和暴力;
int main()
{
node* p = (node*)malloc(sizeof(node)*10);
initList(p);
double ddf = force(p, 10); //9 is the number of elements;
printf("the output of force is %lf\n", ddf);
double dd = callMinDist(p, 10);
printf("the output of divide-conquer is %lf\n", dd);
getchar();
return 0;
}
上述中的force()是枚举的实现，callMinDist则是分治的实现，上述的initList主要对采用的测试案例进行初始化。
下面是node.h头文件的相关代码：
#ifndef __NODE__
#define __NODE__
#define SIZE 4
#define MAX 100000.0
typedef struct{
double x;
double y;
}node;
//排序部分;
void mergeX(node a[], node b[], int s, int m, int t);
void mergeY(node a[], node b[], int s, int m, int t);
void mergeSortX(node a[], node b[], int s, int t);
void mergeSortY(node a[], node b[], int s, int t);
//utility;
void show(node* a, int size);
void initList(node* list);
double dist(node* n1, node* n2);
//枚举部分;
double force(node* nodelist, int n);
//分治部分;
double combine(node* py, int n, double lx, double delta);
void copynode(node* dt, node* sr, int b, int n);
double minDist(node* px, node* py, int n);
double callMinDist(node*p, int n);
#endif
枚举部分的代码比较简单，具体看如下：
#include <math.h>
#include "node.h"
double dist(node* n1, node* n2)
{
double temp1 = n1->x - n2->x;
double temp2 = n1->y - n2->y;
double sum = temp1*temp1 + temp2*temp2; //pow(temp1, 2)+pow(temp2, 2);
return sqrt(sum);
}
double force(node* nodelist, int n) //n is number of elements
{
//这里d需要有一个初始的大值;
double d = MAX;
double t;
//函数的主体就是下面这个双层循环;
for(int i=0; i<n; i++){
for(int j=i+1; j<n; j++)
{
t = dist(&nodelist[i], &nodelist[j]);
if(t<d)
d = t;
}
}
//这里最后返回d;
return d;
}
下面是对于分治算法的调用部分，调用之前需要分别将其中的点按x轴和按y轴进行排序操作，并且
将排完序的点放置在新的存储空间中；
double callMinDist(node* p, int n){
node* px = (node*)malloc(n*sizeof(node)); //n主要是用于此处的空间申请;
node* py = (node*)malloc(n*sizeof(node));
//show(p, n);
mergeSortX(p, px, 0, n-1); //按点的x轴值排序；
copynode(px, p, 0, n-1);
//show(px, n);
//show(p, n);
mergeSortY(p, py, 0, n-1); //按点的y轴值排序；
copynode(py, p, 0, n-1);
//show(py, n);
double min = minDist(px, py, n);
free(px);
free(py);
return min;
}
下面就是分治算法的主体：
double minDist(node* px, node* py, int n) //这里n是number of elements呢？还是max of index？根据下面的空间申请，应该是number of elements.
{
//printf("n is %d\n", n);
if(n<=3){
//show(px, n); //n is number of elements;
double d = force(px, n); //n is number of elements;
//printf("n=%d is %lf\n",n, d);
return d;
}
int m=n/2;
double fx = px[m].x;
node* lx = (node*)malloc(m*sizeof(node));
node* ly = (node*)malloc(m*sizeof(node));
node* rx = (node*)malloc((n-m)*sizeof(node));
node* ry = (node*)malloc((n-m)*sizeof(node));
copynode(lx, px, 0, m-1); //对copy而言,这里的m应该是index;
//show(lx, m); //show中的n是number of elements;
//printf("lx :%x\n", lx);
copynode(rx, px, m, n-1); //copy这边是index;
//show(rx, n-m);
copynode(ly, py, 0, m-1);
//show(ly, m);
copynode(ry, py, m, n-1);
//show(ry, n-m);
double d1 = minDist(lx, ly, m); //m is number of elements;
double dr = minDist(rx, ry, n-m);
double delta = min(d1, dr);
double d = combine(py, n, fx, delta); //对combine而言,这里的n是number of elements;
//printf("lx :%x\n", lx);
free(lx);
free(ly);
free(rx);
free(ry);
return min(delta, d);
}
通过dl = minDist(lx, ly, m)来完成左半部分的计算；
dr = minDist(rx, ry, n-m)完成右半部分的计算；
最后通过combine(py, n, fx, delta)将两半部分的结果整合在一起；
这里的关键之处在于combine函数：
double combine(node* py, int n, double lx, double delta)
{
int num; double d=MAX;
double tempd;
node* temp = (node*)malloc(n*sizeof(node));
int j=0;
for(int i=0; i<n; i++)
{
if(fabs(py[i].x - lx)<= delta){ //求取在区间范围内的点；
temp[j].x = py[i].x;
temp[j].y = py[i].y;
j++;
}
}
num = j; //temp中的元素
for(int i=0; i<num; i++)
{
for(j=i+1; j<(i+8) && (j<num); j++)
{
tempd = dist(&temp[i], &temp[j]);
if(tempd < d)
d=tempd;
}
}
free(temp);
return d;
}
根据书本上的分析，在区间中求取时，只需要计算当前点后(按y轴的值排序)的6到7
个点即可，因此此处的语句表现为：
for(j=i+1; j<(i+8)&&(j<num); j+)....
最后，我们来看下上述代码中用到的一些周边函数：
void copynode(node* dt, node* sr, int b, int n) //n is max of index;
{
int k=0;
for(int i=b; i<=n; i++)
{
dt[k].x = sr[i].x;
dt[k].y = sr[i].y;
k++;
}
}
double min(double x, double y)
{
if(x<=y)
return x;
return y;
}
还有是通过归并排序对点集中的点进行排序的过程：
void mergeSortX(node a[], node b[], int s, int t)
{
if(s == t){
b[s].x = a[s].x;
b[s].y = a[s].y;
}
else{
int m = (s+t)/2;
mergeSortX(a, b, s, m);
mergeSortX(a, b, m+1, t);
mergeX(a, b, s, m, t);
}
}
void mergeSortY(node a[], node b[], int s, int t)
{
if(s == t){
b[s].x = a[s].x;
b[s].y = a[s].y;
}
else{
int m = (s+t)/2;
mergeSortY(a, b, s, m);
mergeSortY(a, b, m+1, t);
mergeY(a, b, s, m, t);
}
}
void mergeX(node* a, node* b, int s, int m, int t)
{
int i, j, n;
for(i=s, j=m+1, n=s; (i<=m)&&(j<=t); n++){
if(b[i].x<=b[j].x){
a[n].x = b[i].x;
a[n].y = b[i].y;
i++;
}else{
a[n].x = b[j].x;
a[n].y = b[j].y;
j++;
}
}
while(i<=m){
a[n].x = b[i].x;
a[n++].y = b[i++].y;
}
while(j<=t){
a[n].x = b[j].x;
a[n++].y = b[j++].y;
}
//这里需要将a中的数据重新拷贝到b中;
for(int i=s; i<=t; i++){
b[i].x = a[i].x;
b[i].y = a[i].y;
}
}
void mergeY(node* a, node* b, int s, int m, int t)
{
int i, j, n;
for(i=s, j=m+1, n=s; (i<=m)&&(j<=t); n++){
if(b[i].y<=b[j].y){
a[n].x = b[i].x;
a[n].y = b[i].y;
i++;
}else{
a[n].x = b[j].x;
a[n].y = b[j].y;
j++;
}
}
while(i<=m){
a[n].x = b[i].x;
a[n++].y = b[i++].y;
}
while(j<=t){
a[n].x = b[j].x;
a[n++].y = b[j++].y;
}
//这里需要将a中的数据重新拷贝到b中;
for(int i=s; i<=t; i++){
b[i].x = a[i].x;
b[i].y = a[i].y;
}
}
注意：上述的归并排序函数写得不好，没必要用到a这个参数，完全可以函数内部堆上分配局部变量进行
替换。
结论：
针对代码中的简单测试案例，分治案例结果正常；该算法的主要时间复杂度在于排序部分，复杂度为
O(nlogn)，而枚举版本的复杂度为O(n2)。
时间： 2024-10-11 05:35:14