Solving $x^2 + y^2 = 1 + z^4$ with (x,y,z) = 1 and z < x < y

Question

I have a computer programming problem where I need to find n many sets of integers that meet the condition $x^2 + y^2 = 1 + z^4$ with (x,y,z) = 1 and z < x < y

I can do this relatively easily with a brute force algorithm that increments z, then finds an x which has a GCD of 1 with z, then finds a y that has a GCD of 1 with z & x. Which then checks to see if the co-prime tripple of x, y & z satisfy the above equation.

However the program is unbearably slow making millions of function calls once I've discovered the first 5 - 10 sets.

How can I restate the problem to quickly find a candidate x & y given some z? Is there some other approach I could consider over this slow brute force method?

Answer 1

You don't have to worry about $(x,y,z)=1$ because it's automatic: if $x\equiv y\equiv z\equiv 0 \pmod{p}$ then your condition $x^2+y^2=z^4+1$ becomes $0\equiv1\pmod{p}$.

Edit: Modified code to cater for missed constraint $z<x$ as per @alex.jordan's comment.

This brute force python fragment produces about 380 solutions in a few seconds on my laptop:

for z in range(1,300+1):
  target = z*z*z*z+1
  for x in range(z+1,int((target/2)**.5)+1):
    y = int((target-x*x)**.5)
    if x*x + y*y == target:
      print("Solution: z=",z,"x=",x,"y=",y)

(Here (blah)**.5 is python for $\sqrt{\text{blah}}$.) And int(blah) is floor(blah).

Basically, for a given $z$ and viable $x$ values, you test the single plausible $y$ value of $\left\lfloor\sqrt{z^4-1-x^2}\right\rfloor$.

If you need to find even more solutions efficiently, I think you'll have to factor $z^4+1$ into primes of the form $4k+1$ (plus an optional factor of 2), decompose the primes into sums of squares, then recombine those solutions. It'd be faster, but take a lot more code and you'd still be limited by the size of the numbers you can factor which probably won't go much beyond the $300^4+1$ in the code fragment above.

Incidentally, there will always be at least one solution to $x^2+y^2=z^4+1$ for each $z$ since $z^4+1$ is a product of primes of the form $4k+1$, plus an optional factor of 2. For a prime $p$ of the form $4k+3$, $z^4+1\equiv0 \pmod{p}$ has no solutions since $x^2\equiv-1\pmod{p}$ has no solutions. And $z^4+1\equiv0\pmod{4}$ has no solutions since $x^2\equiv3\pmod{4}$ has no solutions, so there is at most 1 factor of 2 in $z^4+1$. Any number can be written as the sum of two squares iff is it the product of powers of 2 and $4k+1$ type primes and even powers of $4k+3$ primes. The constraint $z<x$ doesn't always hold in the solutions, though (thanks to alex for remark in comment below).

Answer 2

Write for the equation.

$$x^2+y^2=z^4+1$$

One simple formula.

$$x=(\frac{a^4+1}{2}t+2a)ta^2+1$$

$$y=\frac{a^8-1}{4}t^2+a(a^4-1)t+a^2$$

$$z=\frac{a^4+1}{2}t+a$$

If you use the solutions of the equation Pell. Where $t -$ ask yourself.

$$p^2-2(t^2+1)s^2=1$$

Make the change. $$a=2t(p+ts)s$$

Then decisions can be recorded.

$$x=(2a+1)(t^2+1)-1$$

$$y=2a(a+1)(t^2+1)-2a-1$$

$$z=t(p^2+2tps+2(t^2+1)s^2)$$

Answer 3

Here is a quick and dirty method for generating infinitely many solutions. The idea is very crude in the sense that it only produces solutions, where $y=z^2-3$. Restricting ourselves in this way reduces the problem to a Pell equation, and we can solve those efficiently. My write-up does not assume any knowledge of Pell equations, but those familiar with their theory will be bored.

Given all this the algorithm below only finds a tiny fraction of all the solutions of your equation, but it does find this subset very quickly.

---

If plug in $y=z^2-3$ we get the following equation $$ x^2-6z^2=-8.\qquad(*) $$ So the question is how to find $(x,z)$ pairs satisfying this. To that end let me introduce the following functions:

• The function $N:\Bbb{Z}^2\to\Bbb{Z}$ that takes two integers as inputs and spews out a single integer given by $$N(a,b)=a^2-6b^2$$
• The function $M:(\Bbb{Z}^2)^2\to\Bbb{Z}^2$ that take two pairs of integers as inputs and outputs one pair of integers. The formula is $$M((a,b);(c,d))=(ac+6bd,ad+bc).$$

Remark: The function $N$ gives the norm of the element $a+b\sqrt6$ in the ring $\Bbb{Z}[\sqrt6]$, and the function $M$ calculates the product of the numbers $(a+b\sqrt6)(c+d\sqrt6)$ in that ring.

Key Lemma. We have for all integers $a,b,c,d$ the identity $$ N(M((a,b);(c,d)))=N(a,b) N(c,d). $$

Proof. You do this. Or you can just take my word for it.

The Key Lemma allows us to find solutions to $(*)$, i.e. $N(x,z)=-8$, if we can locate one solution to it, and if we can also locate a solution to the equation $N(a,b)=1$. These are easy to find with a bit of searching among small integers. Namely $$ N(4,2)=4^2-6\cdot2^2=-8 $$ and $$ N(5,2)=5^2-6\cdot2^2=1. $$

Here is then an algorithm for producing infinitely many solutions to $(*)$. We first produce a sequence of pairs $(a_n,b_n)$, $n=1,2,\ldots$, such that $N(a_n,b_n)=1$:

1. The starting point is $(a_1,b_1)=(5,2)$, i.e. the solution found above.
2. Given a solution $(a_k,b_k)$ we then find another solution $(a_{k+1},b_{k+1})$ with the recurrence formula $$ (a_{k+1},b_{k+1})=M((a_k,b_k),(5,2))=(5a_k+12b_k,5b_k+2a_k). $$

The point is that the Key Lemma implies that if $a^2-6b^2=1$ and $c^2-6d^2=1$, then also $(ac+6bd)^2-6(ad+bc)^2=1$.

On with producing solutions to $(*)$. We saw that $(x_0,z_0)=(4,2)$ is one solution. We can then use the pairs produced above to get others by setting $$ (x_k,z_k)=M((x_0,z_0);(a_k,b_k))=(4a_k+12b_k,2a_k+4b_k) $$ for $k=1,2,\ldots$. Observe that $(*)$ implies $x\approx z\sqrt6$, so the inequality $x_k>z_k$ is more or less given. Because $y_k=z_k^2-3$, the other requirement $y_k>x_k$ is also more or less clear, at least for all sufficiently large values of $k$. A quadratic polynomial grows faster than a linear one.

---

Here's how the list of solutions begins. I also present the calculations using $\Bbb{Z}[\sqrt6]$ arithmetic. The numbers will grow very quickly here. That is a price for it being deterministic and fast (no search is involved). You have to judge whether this is a problem.

1. The starting solution to $(*)$, namely $(x_0,z_0)=(4,2)$ must be discarded because here $y_0=z_0^2-3=1$ is too small.
2. Because $$ (4+2\sqrt6)(5+2\sqrt6)=44+18\sqrt{6} $$ we get as the next solution $(x_1,z_1)=(44,18)$. This gives $y_1=z_1^2-3=321$ and a valid solution $(x,y,z)=(44,321,18)$.
3. Because $$ (44+18\sqrt6)(5+2\sqrt6)=436+178\sqrt6 $$ we get $(x_2,z_2)=(436,178)$. Here $y_2=178^2-3=31681$, and we have the solution $(x,y,z)=(436,31681,178)$.

Essentially we keep multiplying by $5+2\sqrt6$ to get to the next solution of this form.