分治法--最近点对问题

#include <iostream>
#include <math.h>
#include <iomanip>
using namespace std;
#define N 20

class Point {
	public:
	double x;
	double y;
	static double distance(const Point& p1, const Point& p2) {
		return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y));
	}
	void operator=(const Point& p) { 
       this->x=p.x;this->y=p.y; 
    }
};

int partition_x(Point* p, int low, int high) {
	double pivot = p[low].x;
	while(low<high) {
		//从右向左搜索 找第一个<=pivot的元素 
		while(low<high && p[high].x>=pivot)
			--high;
		p[low] = p[high];  //把比枢轴小的元素移到左子表 
		//从左向右搜索 找第一个>=pivot的元素
		while(low<high && p[low].x<=pivot)
			++low; 
		p[high] = p[low]; //把比枢轴大的元素移到右子表
	}
	p[low].x=pivot; //把枢轴移动到正确的位置 
	return low; //返回枢轴位置 
} 

int partition_y(Point* p, int low, int high) {
	double pivot = p[low].y;
	while(low<high) {
		//从右向左搜索 找第一个<=pivot的元素 
		while(low<high && p[high].y>=pivot)
			--high;
		p[low] = p[high];  //把比枢轴小的元素移到左子表 
		//从左向右搜索 找第一个>=pivot的元素
		while(low<high && p[low].y<=pivot)
			++low; 
		p[high] = p[low]; //把比枢轴大的元素移到右子表
	}
	p[low].y=pivot; //把枢轴移动到正确的位置 
	return low; //返回枢轴位置 
} 


void quicksort_x(Point* p, int low, int high) {
	int pivotloc;
	if(low<high) {
		pivotloc=partition_x(p,low,high);
		quicksort_x(p,low,pivotloc-1); //对左子表递归排序 
		quicksort_x(p,pivotloc+1,high); //对右子表递归排序 
	}
}

void quicksort_y(Point* p, int low, int high) {
	int pivotloc;
	if(low<high) {
		pivotloc=partition_y(p,low,high);
		quicksort_y(p,low,pivotloc-1); //对左子表递归排序 
		quicksort_y(p,pivotloc+1,high); //对右子表递归排序 
	}
}


//分治法求最近点对 
double Nearestpointpair(Point S[], int low, int high) {

	
	double d1, d2, d3, d;
	int mid, i, j, index;
	
	if(high-low == 1) {
		//cout << "1:\t" << Point::distance(S[low], S[high]);;
		return Point::distance(S[low], S[high]); //S[0] S[1]
	}
		
	if(high-low == 2) {
		d1 = Point::distance(S[low], S[low+1]); //S[0] S[1] 
		d2 = Point::distance(S[low+1], S[high]); //S[1] S[2]
		d3 = Point::distance(S[low], S[high]); //S[0] S[2]
	//	cout << "2:\t" << d1 << " " << d2 << " " << d3 << endl;
		 
		if (d1 < d2) d=d1;
		else d=d2;
		if (d < d3) return d;
		else return d3;
	}
	
	quicksort_x(S,low,high); // 按x升序 
	mid = (low+high) / 2;
	
	d1 =  Nearestpointpair(S, low, mid); //递归 在左半平面求解 
	d2 =  Nearestpointpair(S, mid+1, high); // 递归 在右半平面求解 
	
	if(d1<d2) d = d1;
	else d = d2;
	
	int index1=0;
	Point P1[N], P2[N];
	for (i = mid; (i >= low) && (S[mid].x - S[i].x < d); i--)			//建立点集合P1
		P1[index++] = S[i];
	quicksort_y(P1, 0, index-1); // 按y升序 
	
	int index2=0;
	for (i = mid + 1; (i <= high) && (S[i].x - S[mid].x < d); i++)		//建立点集合P2
		P2[index++] = S[i];
	quicksort_y(P2, 0, index-1); // 按y升序

	for(i=0; i<index1; i++) {
		for(j=0; j<index2; j++) {
			if(P2[j].y-P1[i].y >= d)
				break;
			else {
				d3 = Point::distance(P1[i], P2[j]);
				if(d3 < d)
					d = d3;
			}
		} 
	} 	
	return d;
}



int main()
{
	int n,i;
	cin >> n;
	Point p[N];
	for(i=0;i<n;i++)
		cin >> p[i].x >> p[i].y;
	double result = Nearestpointpair(p,0,n-1);
	cout.setf(ios::fixed);
    cout << setprecision(6) << result << endl;
	return 0;
}