This problem pops up quite often in the search for primes. The number of integral solutions to an elliptical curve ax^2+by^2=n influences the chances that n is prime or not. This fact is exploited in the Sieve of Atkin – see my post about that. If both coefficients are positive, the curve is an ellipse. If one of the coefficients are negative, the curve is a hyperbola. Let’s look at both cases. We will consider only Quadrant I solutions.
Ellipses
My first idea was to find the maximum x and y values, then generate the Cartesian product of coordinate pairs. In effect, this overlays a grid of points atop the ellipse. But the vast majority of points don’t lie anywhere near the curve! These points would be wasteful to check. Rather than drawing a box around an ellipse, why not walk down the curve like stairs?
We start at the positive (x=1)-intercept of the ellipse, which is (1,\sqrt{\frac{n-a}{b}}). We then iterate x-values up to the floor of the x-intercept, which is (\sqrt{\frac{n}{a}},0). y-values are generated using the x-values, then rounded to the nearest integer. This set of coordinate pairs “steps” down the ellipse.
We then just test these coordinates by plugging into ax^2+by^2 to see if it spits out n. The ones that do are the integral solutions. Here’s my implementation in Python:
def find_quadrant1_lattice_points_of_ellipse(n, a, b):
"""Steps down an ellipse, generating rounded-off coordinates,
then checks them. Returns a list of lattice points if there
are any."""
lattice_points = []
x_intercept = math.sqrt(n/a)
for x in range(1, math.floor(x_intercept) + 1):
y = round(math.sqrt((n - a*x**2)/b))
if y != 0 and n == a*x**2 + b*y**2:
lattice_points.append((x, y))
return lattice_points
Hyperbolas
Hyperbolas are a bit different. For starts, they look a bit like a butterfly with wings stretching out into infinity. Therefore, finding lattice points need some kind of stopping condition. For the Sieve of Atkin, only integral solutions where x>y are considered. So, let’s do that. We will start at the (y=1)-intercept and climb the steps up by iterating y.
def find_quadrant1_lattice_points_of_hyperbola(n, a, b):
"""Steps up a hyperbola, generating rounded-off coordinates,
then checks them. Returns a list of lattice points satisfying
x > y, if there are any"""
lattice_points = []
stopping_y = math.sqrt(n/(a+b))
for y in range(1, math.floor(stopping_y) + 1):
x = round(math.sqrt((n - b*y**2)/a))
if x > y and n == a*x**2 + b*y**2:
lattice_points.append((x, y))
return lattice_points
If you want to use them, fetch them from my GitHub.