We are given up to 50,000 green and red points, for a total of N <= 100,000 points. The algorithm proceeds by picking a pivot point, then dividing the points into three sets based on which triangle they fall in. This requires some geometry to determine if a point is inside a triangle. We must also carefully partition the list of points so that we don't test a point for membership in a triangle we've already ruled out. We could try to partition the points in place, but maintaining three partitions is tricky code so I find it easier to just create new lists at each stage. We have 32MB of memory, which should just fit O(N log N) points for N = 100,000.
Runtime Analysis: Just like quicksort, the random pivot point is key. It guarantees that, on average, we roughly divide the data into thirds with a linear time process. That's O(N log N), which is what we need for N = 100,000.
Below is a C++ solution using STL vectors.
#include <cstdio>
#include <cstdlib>
#include <cassert>
#include <vector>
using namespace std;
#define MAXN 50000
typedef long long i64;
enum {GREEN, RED};
typedef struct {
int x, y; /* coordinates */
int id; /* numerical id */
int color; /* GREEN | RED */
} point;
point all_pts[2*MAXN];
int num_all;
typedef vector vpoint;
FILE *in, *out;
/**********************************************************************/
vpoint *get_points(int color) {
int n;
fscanf(in, " %d", &n);
vpoint *vp = new vpoint();
for (int i = 0; i < n; i++) {
point *p = all_pts + num_all++;
fscanf(in, " %d %d", &(p->x), &(p->y));
p->id = i + 1;
p->color = color;
vp->push_back(p);
}
return vp;
}
/**********************************************************************/
int ccw(point *p0, point *p1, point *p2) {
i64 dx1, dx2, dy1, dy2;
dx1 = p1->x - p0->x; dy1 = p1->y - p0->y;
dx2 = p2->x - p0->x; dy2 = p2->y - p0->y;
if (dx1 * dy2 > dy1 * dx2) return +1;
if (dx1 * dy2 < dy1 * dx2) return -1;
return 0;
}
int in_triangle(point *p, point *a, point *b, point *c) {
return ccw(p, a, b) == ccw(p, b, c)
&& ccw(p, b, c) == ccw(p, c, a);
}
/* filter_by_tri: Take points in src which fall *strictly within* the
* triangle (a, b, c) and add them to dest.
*/
void filter_by_tri(vpoint *src, point *a, point *b, point *c, vpoint *dest) {
for (vpoint::iterator i = src->begin(); i != src->end(); i++) {
point *p = *i;
if (in_triangle(p, a, b, c))
dest->push_back(p);
}
}
/**********************************************************************/
void emit_edge(point *a, point *b) {
assert(a->color == b->color);
fprintf(out, "%d %d %c\n", a->id, b->id, "gr"[a->color]);
}
/* solve: Solve problems for points in vp, assumed to be within triangle
* (a, b, c). Two points are one color, third is another. Order doesn't
* matter.
*/
void solve(vpoint *vp, point *a, point *b, point *c) {
/* Make it so a->color == b->color */
if (a->color == b->color) {
assert (a->color != c->color);
} else if (a->color == c->color) {
solve(vp, a, c, b); return;
} else if (b->color == c->color) {
solve(vp, b, c, a); return;
}
/* No points in triangle? Done. */
if (vp->size() == 0) return;
int neach[2] = {0, 0};
/* Count number of green/red points. */
for (vpoint::iterator i = vp->begin(); i != vp->end(); i++)
neach[(*i)->color]++;
/* All the same? Connect to same color outer point. */
if (neach[0] == 0 || neach[1] == 0) {
for (vpoint::iterator i = vp->begin(); i != vp->end(); i++)
if ((*i)->color == a->color)
emit_edge(a, *i);
else
emit_edge(c, *i);
return;
}
/* General case: Pick a random point to partition the triangle into
* three subtriangles. This point is the same color as the odd color
* in the surrounding triangle (point c).
*/
point *pivot = NULL;
int n = 0;
for (vpoint::iterator i = vp->begin(); i != vp->end(); i++)
if ((*i)->color == c->color && (rand() % ++n) == 0)
pivot = *i;
assert(pivot != NULL);
/* Emit edge. */
emit_edge(pivot, c);
/* Partition */
vpoint *tri1, *tri2, *tri3;
filter_by_tri(vp, pivot, a, b, tri1 = new vpoint());
filter_by_tri(vp, pivot, b, c, tri2 = new vpoint());
filter_by_tri(vp, pivot, c, a, tri3 = new vpoint());
/* Recurse */
solve(tri1, pivot, a, b);
solve(tri2, pivot, b, c);
solve(tri3, pivot, c, a);
/* Clean up */
delete tri1; delete tri2; delete tri3;
}
/**********************************************************************/
int main() {
vpoint *green, *red;
vpoint *upper_left, *lower_right;
in = fopen("points.in", "r");
out = fopen("points.out", "w");
assert(in && out);
green = get_points(GREEN);
red = get_points(RED);
point *g0 = green->at(0), *g1 = green->at(1);
point *r0 = red->at(0), *r1 = red->at(1);
/* Partition original square into two triangles, solve each seperately. */
upper_left = new vpoint();
filter_by_tri(green, g0, g1, r1, upper_left);
filter_by_tri(red, g0, g1, r1, upper_left);
emit_edge(g0, g1);
solve(upper_left, g0, g1, r0);
delete upper_left;
lower_right = new vpoint();
filter_by_tri(green, r0, r1, g1, lower_right);
filter_by_tri(red, r0, r1, g1, lower_right);
emit_edge(r0, r1);
solve(lower_right, r0, r1, g1);
delete lower_right;
fclose(in); fclose(out);
return 0;
}