Why does manually specifying the distance metric to `Nearest` slow down the calculation so much?

This is relevant to my previous [question](https://mathematica.stackexchange.com/questions/15580/compiling-map-over-expression-that-yields-a-ragged-array). I've been using `Select` to find sets of points close to a center point but I figured that an alternative (although uncompilable) solution might be to use `Nearest` with the `ChessboardDistance`.

However, I've found that if you specify the `DistanceFunction` to `Nearest` you tremendously slow down the calculation by orders of magnitude, with `ChessboardDistance` being the worst of the bunch by another order of magnitude. These results were observed in 1D, 2D and 3D with random data. I'm not sure what the effect, if any, of structured data would be.

Here is the example code for 3D points. I will give all the timings below. I've simply created a set of 1000 points and computed the `NearestFunction` and then computed closest point in the first set of 1000 points to each point in a second set of another 1000 points.

```
data1 = RandomReal[10, {10^3, 3}];
data2 = RandomReal[10, {10^3, 3}];

nf1 = Nearest[data1];
nf1 /@ data2; // AbsoluteTiming

nf2 = Nearest[data1, DistanceFunction -> EuclideanDistance];
nf2 /@ data2; // AbsoluteTiming

nf3 = Nearest[data1, DistanceFunction -> ManhattanDistance];
nf3 /@ data2; // AbsoluteTiming

nf4 = Nearest[data1, DistanceFunction -> ChessboardDistance];
nf4 /@ data2; // AbsoluteTiming

```

Here are the timings:

```
nf1 {0.0160009, Null}

nf2 {1.1730671, Null}

nf3 {1.1720671, Null}

nf4 {26.4165109, Null}

```

Very similar timings were found for 1D and 2D so I won't list all of them here.

The documentation says: If the elements are vectors or matrices of numbers, Nearest by default *in effect* uses the metric `Norm[#1-#2]&`. 

So I tried to manually use `Norm[#1-#2]&` as the `DistanceFunction`.

```
nf5 = Nearest[data1, DistanceFunction -> (Norm[#1 - #2] &)];
nf5 /@ data2; // AbsoluteTiming

```

This is, surprisingly, even slower.

```
nf5 {3.9752274, Null}

```

I figured that maybe the default distance transform, even though it is given in the documentation, is compiled. Lets try that.

```
euclideandistanceC = 
  Compile[{{pt1, _Real, 1}, {pt2, _Real, 1}}, Norm[pt1 - pt2], 
   CompilationTarget -> "C", RuntimeAttributes -> {Listable}, 
   Parallelization -> True, RuntimeOptions -> "Speed"];

nf6 = Nearest[data, DistanceFunction -> euclideandistanceC];
nf6 /@ data2; // AbsoluteTiming

nf6 {0.8990514, Null}

```

Better but still slower. What's going on? Why is the default method so much faster, and most importantly, how can I use other distance metrics without a huge performance penalty?