74 lines
2.4 KiB
D
74 lines
2.4 KiB
D
import std.stdio, std.typecons, std.math, std.algorithm,
|
|
std.random, std.traits, std.range, std.complex;
|
|
|
|
auto bruteForceClosestPair(T)(in T[] points) pure nothrow @nogc {
|
|
// return pairwise(points.length.iota, points.length.iota)
|
|
// .reduce!(min!((i, j) => abs(points[i] - points[j])));
|
|
auto minD = Unqual!(typeof(T.re)).infinity;
|
|
T minI, minJ;
|
|
foreach (immutable i, const p1; points.dropBackOne)
|
|
foreach (const p2; points[i + 1 .. $]) {
|
|
immutable dist = abs(p1 - p2);
|
|
if (dist < minD) {
|
|
minD = dist;
|
|
minI = p1;
|
|
minJ = p2;
|
|
}
|
|
}
|
|
return tuple(minD, minI, minJ);
|
|
}
|
|
|
|
auto closestPair(T)(T[] points) pure nothrow {
|
|
static Tuple!(typeof(T.re), T, T) inner(in T[] xP, /*in*/ T[] yP)
|
|
pure nothrow {
|
|
if (xP.length <= 3)
|
|
return xP.bruteForceClosestPair;
|
|
const Pl = xP[0 .. $ / 2];
|
|
const Pr = xP[$ / 2 .. $];
|
|
immutable xDiv = Pl.back.re;
|
|
auto Yr = yP.partition!(p => p.re <= xDiv);
|
|
immutable dl_pairl = inner(Pl, yP[0 .. yP.length - Yr.length]);
|
|
immutable dr_pairr = inner(Pr, Yr);
|
|
immutable dm_pairm = dl_pairl[0]<dr_pairr[0] ? dl_pairl : dr_pairr;
|
|
immutable dm = dm_pairm[0];
|
|
const nextY = yP.filter!(p => abs(p.re - xDiv) < dm).array;
|
|
|
|
if (nextY.length > 1) {
|
|
auto minD = typeof(T.re).infinity;
|
|
size_t minI, minJ;
|
|
foreach (immutable i; 0 .. nextY.length - 1)
|
|
foreach (immutable j; i + 1 .. min(i + 8, nextY.length)) {
|
|
immutable double dist = abs(nextY[i] - nextY[j]);
|
|
if (dist < minD) {
|
|
minD = dist;
|
|
minI = i;
|
|
minJ = j;
|
|
}
|
|
}
|
|
return dm <= minD ? dm_pairm :
|
|
typeof(return)(minD, nextY[minI], nextY[minJ]);
|
|
} else
|
|
return dm_pairm;
|
|
}
|
|
|
|
points.sort!q{ a.re < b.re };
|
|
const xP = points.dup;
|
|
points.sort!q{ a.im < b.im };
|
|
return inner(xP, points);
|
|
}
|
|
|
|
void main() {
|
|
alias C = complex;
|
|
auto pts = [C(5,9), C(9,3), C(2), C(8,4), C(7,4), C(9,10), C(1,9),
|
|
C(8,2), C(0,10), C(9,6)];
|
|
pts.writeln;
|
|
writeln("bruteForceClosestPair: ", pts.bruteForceClosestPair);
|
|
writeln(" closestPair: ", pts.closestPair);
|
|
|
|
rndGen.seed = 1;
|
|
Complex!double[10_000] points;
|
|
foreach (ref p; points)
|
|
p = C(uniform(0.0, 1000.0) + uniform(0.0, 1000.0));
|
|
writeln("bruteForceClosestPair: ", points.bruteForceClosestPair);
|
|
writeln(" closestPair: ", points.closestPair);
|
|
}
|