Implementation of two-dimensional planar nearest point pair algorithm based on divide and conquer

Summary:

There is a lot of discussion on the method of divide and conquer on the Internet, but there is no complete running code, this article mainly introduces a complete and running code for this problem.

For those who are interested in the reference.

Body:

As a comparison, we also implement the nearest point pair enumeration solution, the main function is as follows:

#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;} Test division and violence; int main () {node* p = (node*) malloc (sizeof (node) *10), Initlist (p);d ouble DDF = Force (p, N);//9 is the number of elements;printf ("The output of force is%lf\n", DDF);d ouble dd = callmindist (P, ten);p rintf ("The output of Divide-conquer I S%lf\n ", DD); GetChar (); return 0;}

The above-mentioned force () is the implementation of the enumeration, Callmindist is the implementation of the Division, the above-mentioned initlist mainly for the use of the test case to initialize.

The following is the code for the node.h header file:

#ifndef __node__#define __node__#define SIZE 4#define MAX 100000.0typedef struct{double x;double y;}  node;//sort section; 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);d ouble Dist (node* N1, node* n2);//enumeration part; Double Force (node* nodelist, int n);//partition part; d Ouble Combine (node* py, int n, double lx, double delta), void Copynode (node* DT, node* sr, int b, int n);d ouble mindist (nod e* px, node* py, int n);d ouble callmindist (node*p, int n); #endif

The code for the enumeration section is relatively simple, as follows:

#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{//here D needs to have an initial large value; double d = The body of the max;double t;//function is the following double-layer loop; fo R (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;}} Here finally returns D;return D;}

The following is the invocation portion of the split algorithm, in which the points are sorted by the x-axis and the y-axis, respectively, and

Placing the ordered points in the new storage space;

Double Callmindist (node* p, int n) {node* px = (node*) malloc (n*sizeof (node)),//n is primarily used for space applications here; node* py = (node*) malloc (n*s Izeof (node));//show (P, N), Mergesortx (p, px, 0, n-1); Sorted by point x-axis values, Copynode (px, p, 0, n-1),//show (PX, n),//show (P, N), Mergesorty (p, py, 0, n-1); Sort by the y-axis value of the point; Copynode (py, p, 0, n-1),//show (PY, N);d ouble min = mindist (px, py, N), free (px), free (py), return min;}

Here is the main body of the divide-and-conquer algorithm:

Double mindist (node* px, node* py, int n)//where n is the number of elements? or Max of index? According to the space application below, it should be number of elements. {//printf ("N is%d\n", N), if (n<=3) {//show (px, n);//n are number of elements;double d = Force (px, n);//n is number of E lements;//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);  For copy, the m here should be index;//show (LX, m); The n in show is number of elements;//printf ("LX:%x\n", LX), Copynode (Rx, PX, M, n-1); Copy this is index;//show (Rx, n-m), Copynode (ly, py, 0, m-1);//show (ly, M); Copynode (ry, py, M, n-1);//show (ry, n-m);d ouble d1 = Mindist (lx, Ly, M); M is number of elements;double dr = Mindist (Rx, ry, n-m);d ouble delta = min (d1, DR);d ouble d = Combine (py, N, FX, Delta) ; For Combine, n is the number of elements;//printf ("LX:%x\n", LX), free (LX), free (ly), Free (RX), and Free (RY); rEturn min (delta, d);} 

The calculation of the left half is done by DL = Mindist (lx, Ly, M);

Dr = Mindist (Rx, Ry, n-m) completes the calculation of the right half;

Finally, the results of the two halves are integrated by combine (PY, N, FX, Delta);

The key point here is the Combine function:

Double Combine (node* py, int n, double lx, double delta) {int num; double d=max;double tempd;node* temp = (node*) malloc (n*s Izeof (node)); int j=0;for (int i=0; i<n; i++) {if (Fabs (PY[I].X-LX) <= delta) {//Find points within the interval range; temp[j].x = py[i].x;temp[ J].y = py[i].y;j++;}} num = j;  Elements in 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;}

According to the analysis in the book, when the interval is obtained, only the 6 to 7 of the current point (sorted by the value of the y-axis) is calculated.

Point, so the statement here shows:

for (j=i+1; j< (i+8) && (j<num); j+) ....

Finally, let's look at some of the surrounding functions used in the code above:

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;}

There is also the process of sorting points in a point set by merging:

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;} It is necessary to re-copy the data in A to 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;} It is necessary to re-copy the data in A to B; for (int i=s; i<=t; i++) {b[i].x = A[i].x;b[i].y = A[i].y;}}

Note: The above merge sort function is not well written, it is not necessary to use the A parameter, it can be assigned local variables on the internal heap of the function.

Replace.

Conclusion:

For the simple test cases in the code, the results of the divide-and-conquer case are normal; The main time complexity of the algorithm lies in the ordering part, the complexity is

O (Nlogn), and the complexity of the enumeration version is O (N2).

Implementation of two-dimensional planar nearest point pair algorithm based on divide and conquer