# How to find a closest integer point to intersection of two lines?

Hello. Here's a question that originates from StackOverflow (and the SO crowd isn't really qualified to solve it).

We're given two lines on the plane, each of them has at least two integer points (i.e. points $(b,a)$ where $b$ and $a$ are in $\mathbb{Z}$). The lines aren't parallel, but their intersection is not necessarily in an integer point. We also are given a point $(x,y)$ (also integer) that doesn't belong to both lines.

The problem is to find a point $(x',y')$ in the same quarter as $(x,y)$ (as separated by the lines) closest to intersection of the lines. It may be $(x,y)$ itself, but most likely it's closer.

The $(x,y)$ point is seemingly needed just to make the problem NP-complete instead of NP-hard. The problem is claimed to have a polynomial solution, but I doubt that it really persists. Could anyone solve it or provide a proof that it's NP-complete?

-
Polynomial in what? – Qiaochu Yuan Apr 27 '10 at 20:51
@Quaochu, in the length of input - i.e. sum of logarithms of the integers that define each line (normal vectors, with GCD equal to one, for simplicity; and one integer point it crosses), and coordinates of (x,y). – user5674 Apr 27 '10 at 20:53
I'm very sorry that I didn;t learn how to tag and paste TeX properly, and I hope that the community would be kind to do it for me... :-/ – user5674 Apr 27 '10 at 20:58
"The (x,y) point is seemingly needed just to make the problem NP-complete instead of NP-hard." No, you can find such an (x,y) easily. Just go far enough away from the point of intersection so that the vertical or horizontal distance between the two lines is >= 1. – TonyK Apr 27 '10 at 21:11
Gaussian lattice reduction probably works. – François G. Dorais Apr 28 '10 at 2:38

If the point $(x,y)$ is in an obtuse quadrant between the lines the problem is easily solved by enumerating the lattice points in a closed sphere of radius $\sqrt{2}$ about the intersection point so I will only consider the case that $(x,y)$ lies in a acute quadrant.

The problem can be converted to one that is very similar to inhomogeneous Diophantine approximation. I know of at least one algorithm that will definitely find the correct answer in $O\left(\tfrac{1}{\theta}\right)$ operations where $2\theta$ is the acute angle between the lines. The algorithm is a minor modification of that given on page 19 of Vaughan Clarkson's thesis. I am very confident that Cassels' algorithm (or a minor modification of it) will solve the problem in $O\left(1 + \log\tfrac{1}{\theta}\right)$ operations (page 34 of 1), but I am not quite sure how to show it, or what modification needs to be made. Before I describe this, you need to know a little bit about inhomogeneous Diophantine approximation.

Let $\alpha$ and $\beta$ be real numbers. Define the function

$F(p,q) = |\alpha p - q - \beta|$.

The problem of inhomogeneous Diophantine approximation involves minimising $F(p,q)$ over integers $p$ and $q$ where $q$ is positive. The pair $(p, q)$ is called a best approximation if $F(p,q) < F(p',q')$ for all $q' < q$. The best approximations describe the minima found as $q$ and the magnitude of $p$ increase. There are various algorithms than can enumerate all of the best approximations. Two examples are the 'naive algorithm' (page 19 of 1) and Cassel's algorithm (page 34 of 1). The OP's problem is not exactly the same as this, but it is so similar that the algorithms (at least the naive algorithm) carry over.

Let our two lines be $\ell_1$ and $\ell_2$ and let $2\theta < \pi/2$ be the angle between them. It will be easier to describe the approach if we set the intersection of these lines to be at the origin and we look for the nearest point in the translated lattice $\mathbb{Z}^2 + t$ where $t$ is the appropriate translation. Define $\ell$ to be the unique line that passes through the origin and bisects the acute angle between $\ell_1$ and $\ell_2$. The angle between $\ell$ and $\ell_1$ and $\ell$ and $\ell_2$ is $\theta$. The problem can now be stated as:

Find the point $x \in \mathbb{Z}^2 + > t$ that is nearest to the origin such that the angle between $x$ and $\ell$ is less than or equal to $\theta$

Let $A(x)$ denote the angle between $\ell$ and $x$. Our motivation is now very similar to that in Diophantine approximation. That is, find all of the best approximations for $A(x)$, our problem is solved by the best approximation that first yields $A(x) \leq \theta$. It so happens that $A(x)$ is a very similar function to $F(p,q)$. To give this some context I will say that $F(p,q)$ is, in a sense, computing an inner product between two vectors, whereas $A(x)$ is computing an angle. In this context it is not surprising that the algorithms for Diophantine approximation can be used.

I will only consider the 'naive algorithm' and I'll just give some geometric insight as to its functionality, this should be enough to convince most people. Working through this rigorously is really beyond a typed answer on MO, but all the required machinery is in 1. The 'naive algorithm' enumerates every point in $\mathbb{Z}^2 + t$ that is a nearest lattice point to any point in the line $\ell$. In other words it consecutively locates (starting from the origin) every lattice point in $\mathbb{Z}^2 + t$ whose Voronoi cell (in this case squares) intersect $\ell$. A picture might be useful

It is not difficult to devise an algorithm which does this, just start at the origin and check where $\ell$ next crosses a boundary of a Voronoi cell. It is also easy to see that the points it locates are a super set of the best approximations for $A(x)$. The first point that the algorithm finds such that $A(x) < \theta$ is the solution to the problem (the blue circle).

This algorithm is called 'naive' because it checks a lot of lattice points that are not best approximations. Cassels' algorithm improves this substantially for the function $F(p,q)$. It's likely that a similar improvement is possible for $A(x)$ and someone might wish to work it out.

The OP (particularly on stack overflow, but also here) seems to have thrown quite a number of red herrings into the problem statement (rather annoying). For example, knowledge of the point (x,y) does nothing other than tell you which quadrant you are looking in. The statement about it converting the problem to NP-complete rather than NP-hard doesn't make any sense. Also, the fact that the lines pass through integer points appears to be irrelevant.

-
The naiive algorithm, as dscribed, is exponential in the initial data, and therefore, out of question. We need $O(\log\frac 1\theta)$, not $O(\frac 1\theta)$. On the other hand, the explanation of the meaning of the question and the picture are exemplary! – fedja Apr 29 '10 at 9:36
I don't disagree. The analogue of Cassels' algorithm needs to be found! It might be that Cassels' algorithm works 'as is'. I would be hugely surprised if there was not an analogue of Cassels' algorithm for this problem, but, at the moment there are other things I should be doing (like writing a thesis :)). If someone else can put it all together then all the brownie points to them. – Robby McKilliam Apr 29 '10 at 11:10
@Robby: You said "Also, the fact that the lines pass through integer points appears to be irrelevant." What's relevant is whether or not the intersection is at a lattice point. In diophantine approximation it's well known that there's a qualitative difference between the homogeneous and inhomogeneous problem. – Victor Miller Apr 29 '10 at 14:45
Thanks, Robby. I didn't realize you were typing your dissertation. – Will Jagy Apr 29 '10 at 16:21
@Victor: You are correct, if the intersection point is an integer then this becomes homogeneous Diophantine approximation. @Will: No problem. – Robby McKilliam Apr 29 '10 at 21:10

I don't believe that the problem is NP-complete, because you're working in a fixed dimension. Since most of us believe that the hardest case is when the angle between the two lines is very small, then you can take the line L half way in between the two lines, and orthogonally project onto it for the objective function. Then we have an integer programming problem in 2 dimensions -- the two inequalities specify the proper side of the two lines. Hendrik Lenstra in "Integer Programming with a fixed number of variables" in Mathematics of Operations Research, showed that when the dimension is fixed there is a polynomial time algorithm for IP (using a variant of the L^3 lattice basis reduction). There's also the paper http://www.math.uni-klu.ac.at/or/doctoralschool/deloera.pdf "Integer Polynomial Optimization in Fixed Dimension" which mentions that a convex polynomial objective function also has a polynomial time algorithm in fixed dimension, so that should do it for this problem.

[Added Comments] Using the two papers that I mentioned in my comments http://mpi-inf.mpg.de/~soeren/pubs/2ip.ps http://homepages.cwi.nl/~aardal/journal_rev.ps one can proceed as follows: We're going to have upper and lower supporting lines orthogonal to the midpoint line L. These will always have the property that we know that there is at least one lattice point in the quadrilateral bounded by the the original lines and the supporting line (though at the beginning, it has degenerated into a triangle). At the beginning the lower supporting line is the one passing through the intersection of the two original lines. The upper supporting line can be calculated by noticing that any disk of radius > sqrt(2)/2 must contain a lattice point, so you can find the smallest distance along the midline where you can place the center of such a disk -- it will be where the line from center orthogonal to each of the original lines has distance sqrt(2)/2. By simple trigonometry you can see that if the number of bits in the original number is N, then we need at most 2N bits to specify this point (i.e. if the denominators are around n for the originals, then the denominators for the above points are at most n^2). Now use the algorithm in first paper which tells you in time linear in the number of bits of the problem whether or not there are any lattice points in the quadrilateral. Do a binary search by looking at the a test line half way in between the supporting lines, and testing each of the two quadrilaterals. After only about N steps (remember that N is the log of the coefficients) you'll be down to a quadrilateral with a small area and width. At that point you can quickly enumerate all lattice points in it and test them for the minimum. This algorithm probably runs in time O(N^2) where N is the number of bits in the original coefficients.

The idea in solving the problem is to do the following:

1) Find a good enough approximation to the minimum distance, and the lattice point attaining that, so that you can enclose that region in a rectangle whose sides are parallel to the $x$ and $y$ axes, with a small enough area, $A$ and perimeter $P$. It's easy to see that you can enumerate all lattice points inside of such a rectangle in time $A+P$, and then check those to see which ones give you the minimum. [Changing notation] Let $n$ be the number of bits to specify the problem and $N=2^n$. Let the angle between the two lines be denoted by $2 \theta$, and the midline by $L$.

2) As I mentioned, any disk with radius $> \sqrt{2}/2$ must contain a lattice point. We use this by placing such a disk in between the two lines, and as close to the point of intersection as possible. We see that the distance from the intersection point to the center of the disk which just fits, is $\sqrt{2}/2 \csc \theta$. If $\theta \ge \pi/6$ (say) we can enclose the whole triangle bounded by the lines and the tangent line to the disk orthogonal to $L$ in a rectangle of constant size and perimeter, so we just try all of those points. Note that in any case $\theta \ge c/N$ for some absolute constant $c$ since when $\theta$ is small enough $\cos \theta \approx 1 - \theta^2/2$, and $\cos 2 \theta$ is given by a dot product between the coefficients of the lines and so has at most $2n$ bits.

3) Otherwise we call the algorithm described in http://mpi-inf.mpg.de/~soeren/pubs/2ip.ps just once, with the objective function the dot product between $(x,y)$ and a vector parallel to $L$, and the four constraints: between the two lines, and the dot product with $L$ is $\ge 0$ and $\le$ the bound we get by placing the disk. Because $\theta$ is small enough and bounded away from 0, we can now draw a rectangle enclosing the point produced by the algorithm whose area and perimeter are bounded, independent of the number bits, which must contain the answer, and we again enumerate all points in that rectangle.

-
Hi Will, In the meantime this paper is also quite relevant mpi-inf.mpg.de/~soeren/pubs/2ip.ps – Victor Miller Apr 29 '10 at 1:48
And here's another, perhaps even more relevant paper homepages.cwi.nl/~aardal/journal_rev.ps "Solving a system of linear diophantine equations with upper and lower bounds on the variables" – Victor Miller Apr 29 '10 at 2:21
That's good, Victor. I deleted my answer after somebody wrote "How is this an answer?" There is definitely a jump in expertise from my minimizing dot products with vectors through the intersection point. – Will Jagy Apr 29 '10 at 16:16

You can try to apply the algorithm from the paper

B. N. Delone, “An algorithm for the “divided cells” of a lattice”, Izv. Akad. Nauk SSSR Ser. Mat., 11:6 (1947), 505–538

-

Doesn't the following work? UPDATE: Sorry, no it doesn't. I missed the requirement that the lattice point be on the correct side of the two lines.

Find the exact coordinates, $(u,v)$ of the intersection. These are rational numbers, and can be found by solving two linear equations. Then the nearest lattice point is $(\mathrm{ROUND}(u), \mathrm{ROUND}(v))$ where $\mathrm{ROUND}$ rounds to the nearest integer.

To see that this is the closest point, translate $(\mathrm{ROUND}(u), \mathrm{ROUND}(v))$ to the origin. So we need to show that, if $|u|$ and $|v|$ are $< 1/2$, then $(u,v)$ is closer to $(0,0)$ than to any other lattice point $(a,b)$. The cases of $(a,b) = (0, \pm 1)$, $(\pm 1, 0)$ or $(\pm 1, \pm 1)$ can be checked by hand. For any other $a$, $b$, one of $|u-a|$ and $|v-b|$ is greater than $1$, so the distance from $(u,v)$ to $(a,b)$ is at least $1$; which is much greater than $\sqrt{u^2+v^2} \leq \sqrt{2}/2$.

For a formal complexity analysis, the size of your input is the log of the number of digits needed to specify the two lines. All the arithmetic operations need to solve linear equations can be done in time polynomial in this logarithm.

This question is much more interesting for lattices other than the square grid. See Voronoi decomposition for details.

-
It seems that you didn't notice the phrase "in the same quarter [as (x, y)]" in the problem statement... – Reid Barton Apr 27 '10 at 21:05
Sorry, it's most likey my fault that I didn't formulate the question well, but @Reid is correct: we need to find a point in a specific quarter; and of course we assume automatically (since it's MathOverflow) that the quarter is the hardest. – user5674 Apr 27 '10 at 21:08

This problem is a little bit obscure, but there are several efficient solutions which you may be interested in. I'd start by first trying to find a copy of Lee and Chang's paper, unless you're fluent in French, in which case Guigue's thesis contains an updated exposition.

Mehta, S., M. Mukherjee, and G. Nagy. "Constrained integer approximation to 2-D line intersections." Second Canadian Conf. on Computational Geometry. 1990.

Lee, H. S., and R. C. Chang. "Approximating vertices of a convex polygon with grid points in the polygon." Algorithms and Computation. Springer Berlin Heidelberg, 1992. 269-278.

Guigue, Philippe. Constructions géométriques à précision ﬁxée. Diss. Université de Nice Sophia-Antipolis, 2003.

-