Consider a set of $n$ real-valued number pairs: $(x\_1,y\_1), (x\_2,y\_2), \dots, (x\_n,y\_n)$. I want to find a permutation $p$ of the indices which minimizes the sum of consecutive absolute differences: 


$$\sum\_{j=1}^{n-1} |x\_{p(j+1)} - x\_{p(j)}| + \sum\_{j=1}^{n-1} |y\_{p(j+1)} - y\_{p(j)}|.$$


I suspect this is reducible to a well known problem, so I'm looking for pointers to literature mainly, but would be happy to see a clever algorithm for doing this from scratch. 


Intuitively, I want to shuffle the observations so that the graph of $x$ elements against the index is smooth looking and the same graph of the $y$ elements is also smooth looking. If I cared only about one or the other, I could simply sort with respect to those elements. I want to shuffle in such a way that I compromise between the two coordinates.


My motivation is a statistical problem of estimating a smooth curve in the plane by assuming that the coordinate dimensions are each smooth functions of an unrecorded "time index". The above problem is maximizing the smoothness of the observed data under the assumption of evenly spaced observations in time.



sorting two paired lists of real numbers to minimize consecutive absolute differences
This is known in complexity circles as rectilinear TSP. Technically this is a path version TSP rather than the more common tour. In any of these specific settings, this is known to be NP-complete since the 70s (see [Papadimitriou's paper](http://www.sciencedirect.com/science/article/pii/0304397577900123)).


However, unlike some NP-hard problems, Euclidean and rectilinear TSP admit a polynomial-time approximation scheme: one can obtain a tour of cost at most $(1+\epsilon)$ times optimal in time $O(n (\log n)^{O(1/\epsilon)})$. See [Arora's web page](http://www.cs.princeton.edu/~arora/publist.html) for some of the key papers on this subject. You'll also find there a nice survey on approximation results ([postscript link](http://www.cs.princeton.edu/~arora/pubs/arorageo.ps)), outlining the key ideas.


The survey mentions a subexponential $2^{O(\sqrt{n})}$ exact algorithm for the Euclidean case by Smith (1988). While I would not be surprised if it adapts readily to $L^1$, I don't know enough about it to say for certain.


