An interesting equation: x/(y+z)+y/(z+x)+z/(x+y)=4 x,y,z∈Z⁺

The author have seen a problem on the Internet: Find positive integers \(x,y,z\) that satisfies equation \(\frac{x}{y+z}+\frac{y}{z+x}+\frac{z}{x+y}=4\).

(Thinking this equation is easy to solve? Try it before continue reading!)

One might want to solve it by brute-forcing with a computer (for example, using Python programming language):

from fractions import Fraction as Frac
n = 100
print('Solutions:')
for x in range(1, n + 1):
    for y in range(1, n + 1):
        for z in range(1, n + 1):
            if Frac(x, y + z) + Frac(y, z + x) + Frac(z, x + y) == 4:
                print((x, y, z))

Solutions:

However, it does not help.

Try with Wolfram Mathematica:

FindInstance[x/(y + z) + y/(z + x) + z/(x + y) == 4, {x, y, z}, Integers]

{{x -> 11, y -> 9, z -> -5}}

It gives an integer solution, but not a positive integer solution. Adding more constraints may solve this problem:

FindInstance[x/(y + z) + y/(z + x) + z/(x + y) == 4 && x > 0 && y > 0 && z > 0, {x, y, z}, Integers]

FindInstance::nsmet: The methods available to FindInstance are insufficient to find the requested instances or prove they do not exist.
FindInstance[z/(x + y) + y/(x + z) + x/(y + z) == 4 && x > 0 && y > 0 && z > 0, {x, y, z}, Integers]

Unfortunately, it does not work either!

This equation is much harder than expected. Only experts in the field of number theory will have a chance to solve it. The smallest solution to this equation is:

x = 36875131794129999827197811565225474825492979968971970996283137471637224634055579
y = 154476802108746166441951315019919837485664325669565431700026634898253202035277999
z = 4373612677928697257861252602371390152816537558161613618621437993378423467772036
from fractions import Fraction as Frac
print(Frac(x, y + z) + Frac(y, z + x) + Frac(z, x + y) == 4)

True

They are indeed enormous numbers since the most significant number has 81 digits. Surely brute-forcing would take an unacceptably long time.

(If you know some mathematics lovers, share them with this enjoyable problem!)

(Image above is from xkcd)

How to obtain this enormous solution? The answer will be given below.

This equation is simple on the outside, but complicated on the inside. Equations like this have puzzled the best mathematicians in this world for thousands of years. Diophantus of Alexandria, a Hellenistic mathematician of the 3rd century, made a study of several polynomial equations of this form:

\[P(x_{1},x_{2},\dots,x_{k})=\sum _{{0\leq i_{j}\leq n_{j}}}a_{{i_{1}i_{2}\dots i_{k}}}x_{1}^{{i_{1}}}x_{2}^{{i_{2}}}\dots x_{k}^{{i_{k}}}=0\]

Equations of this type are Diophantine equations. For example, \(x^2+y^2=z^2\) (\(x,y,z\) are the unknowns) is a Diophantine equation. System of Diophantine equations can be reduced to a Diophantine equation by means of \({\displaystyle\sum_i (RHS_i-LHS_i)^2=0}\). So this article will only discuss situations that contain one equation.

After applying reducing to common denominator technique to the equation in the beginning, an equivalent equation is obtained:

\[x^3 + y^3 + z^3 - 3x^2(y+z) - 3y^2(z+x) - 3z^2(x+y) - 5xyz=0\]

This equation is indeed a Diophantine equation! More precisely, a degree-3 homogeneous Diophantine equation (all items in this equation are of degree \(3\)).

Since the equation is homogeneous, if \((x_0,y_0,z_0)\) is a solution, \((kx_0,ky_0,kz_0),\,k\neq 0\) is a solution too. So with a rational solution, an integer solution can be generated with ease (by letting \(k\) equals to least common multiple of denominators of \(x_0,y_0,z_0\)). This property also implies that one can assume \(z=1\), without losing any wanted solution.

Then the equation is reduced to a simpler one (with fewer variables):

\(x^3-3 x^2 y-3 x^2-3 x y^2-5 x y-3 x+y^3-3 y^2-3 y+1=0,\quad x,y\in\mathbb{Q}^{+{}}\)

But one can make it even simpler, by doing some (rational) linear transformation:

\[(X,Y,Z)^T=\left(\begin{array}{ccc} -\frac{95}{2} & -\frac{95}{2} & \frac{277}{12} \\ -\frac{91}{2} & \frac{91}{2} & 0 \\ -6 & -6 & 1 \\\end{array}\right) (x,y,z)^T\]

It is easy to verify that:

\[X^3-\frac{11209}{48}X Z^2+\frac{1185157}{864}Z^3-Y^2 Z=8281(x^3 + y^3 + z^3 - 3x^2(y+z) - 3y^2(z+x) - 3z^2(x+y) - 5xyz)\]

The new equation is \(x^3-\frac{11209}{48}x z^2+\frac{1185157}{864}z^3-y^2 z=0\). Since the transformation matrix is invertible (determinant is not zero), one can convert the solution of the new equation to the solution of the previous equation, without losing any potential solutions of the original equation.

Assuming \(z=1\), the equation can be reduced to \(x^3-\frac{11209}{48}x+\frac{1185157}{864}-y^2=0\). This type of equations (\(y^{2}=x^{3}+px+q\)) are Weierstrass equations. Since the discriminant \(\Delta =-16(4p^{3}+27q^{2})\) is non-singular (which implies that the curve of this equation has no cusps, self-intersections, or isolated points), it is, in fact, an elliptic curve.

Open-source software Sagemath can convert a cubic curve to an elliptic curve automatically:

sage: R.<x,y,z> = QQ[];
sage: F = x^3 + y^3 + z^3 - 3*x^2*(y+z) - 3*y^2*(z+x) - 3*z^2*(x+y) - 5*x*y*z;
sage: WeierstrassForm(F)
(-11209/48, 1185157/864)

Each integer solution (with the greatest common divisor \(1\)) of the original equation one-to-one corresponds to a rational point on the elliptic curve, so finding integer solutions is equivalent to finding rational points on the curve.

But not all rational points correspond to a positive integer solution. Only rational points in the red region correspond to a positive integer solution. There are two parts of the red region: the left part has three intersections with the curve, while the right part has none.

It is easy to verify that point \((x,y)\) on the curve lies in red region if and only if \(x\in(\frac{-39 \sqrt{65}-109}{24},\frac{-84 \sqrt{3}-59}{12})\cup(\frac{84 \sqrt{3}-59}{12},\frac{95}{12})\). One can write code to judge whether a point corresponds to a positive solution.

from sympy import sqrt

O = ("Inf", "Inf")


def is_valid(P):
    if P == O:
        return False
    else:
        return ((-39 * sqrt(65) - 109) / 24 < P[0] <
                (-84 * sqrt(3) - 59) / 12) or (
                    (84 * sqrt(3) - 59) / 12 < P[0] < 95 / 12)

After converting ranges to decimal form, one will get \(x\in(-17.64,-17.04)\cup(7.21,7.92)\), which is a small range. It is intuitive that there is little probability that a rational point locates in this range, thus finding a positive integer solution may be much harder than finding an integer solution.

In mathematics, an elliptic curve is a plane algebraic curve defined by an equation of the form \(y^{2}=x^{3}+px+q\), which is non-singular.

Since the curve is symmetrical about the x-axis, given any point \(P\), one can take \(-P\) to be the point opposite it. Define operation \(+\) on points: assume \(P, Q\) are two different points on the curve, draw the line that intersects \(P, Q\), which will generally intersect the cubic curve (uniquely) at a third point \(R\). Take \(P+Q\) to be \(-R\), the point opposite \(R\).

• Situation 1 (generally):

A curve will have an additional point at infinity, \(O\). Taking \(-O\) to be just \(O\), \(P+(-P)=O,\,O+P=P=P+O\) will hold for all points \(P\). (\(O\) on elliptic curve serves like zero in integer addition arithmetic)

There are some other special situations:

• Situation 2 (\(P=Q\)):

• Situation 3 (the line is a tangent line at point \(P\) or \(Q\)):

(Images of various situations above come from Andrea Corbellini's blog)

The reader can change values of \(p,q\), and move points \(P,Q\), to visualize addition operation on an elliptic curve. (you can also zoom, full screen, move, reset canvas)

Addition of points on elliptic curve can be defined not only geometrically, but also algebraically, which is more precise. Given the curve \(y^2 = x^3 + px + q\) over the field \(K\) whose characteristic is assumed to be neither \(2\) nor \(3\), and points \(P = (x_P, y_P)\) and \(Q = (x_Q, y_Q)\) on the curve. Let \(s\) be the slope of the line containing \(P\) and \(Q\); i.e., \(s={\frac {y_{P}-y_{Q}}{x_{P}-x_{Q}}}\).

• \({\displaystyle x_P\neq x_Q}\): \({\displaystyle x_R=s^2-x_P-x_Q,\,y_R=-y_P+s(x_P-x_R)}\)
• \({\displaystyle x_P=x_Q}\):
  • \({\displaystyle y_P=-y_Q}\): \({\displaystyle P+Q=O}\)
  • \({\displaystyle y_P=y_Q}\): \({\displaystyle R=2P,\,s=\frac{3x_P^2+p}{2y_P},\,x_R=s^2-2x_P,\,y_R=-y_P+s(x_P-x_R)}\)

Code for doing points addition:

from fractions import Fraction as Frac

O = ("Inf", "Inf")


def ec_plus(P, Q):
    if P == O:
        return Q
    if Q == O:
        return P
    R = [None] * 2
    if P[0] == Q[0]:
        if P[1] == -Q[1]:
            R = O
        else:
            s = (3 * P[0]**2 + -Frac(11209, 48)) / (2 * P[1])
            R[0] = s**2 - 2 * P[0]
            R[1] = -P[1] + s * (P[0] - R[0])
    else:
        s = (P[1] - Q[1]) / (P[0] - Q[0])
        R[0] = s**2 - P[0] - Q[0]
        R[1] = -P[1] + s * (P[0] - R[0])
    return tuple(R)

It's easy to prove \(P,Q,-(P+Q)\) are collinear with algebraic definition. Points set \(G\) on elliptic curve (include \(O\)) satisfies following properties:

• Closure: For all \(P,Q\) in \(G\), the result of the operation \(P+Q\) is also in \(G\).
• Associativity: For all \(P,Q,R\) in \(G\), the equation \((P+Q)+R=P+(Q+R)\) holds.
• Identity element: There exists an element \(O\) in \(G\), such that for all elements \(P\) in \(G\), the equation \(O+P=P=P+O\) holds.
• Inverse element: For each \(P\) in \(G\), there exists an element \(Q\) in \(G\) such that \(P+Q=O=Q+P\), where \(O\) is the identity element.
• Commutativity: For all \(P,Q\) in \(G\), the equation \(P+Q=Q+P\) holds.

These properties show that points and additions on elliptic curve \((G,{}+)\) forms an abelian group (also called commutative group), and \(O\) is identity element. It is easy to verify that when taking all rational points (and \(O\)) from \(G\), the above five properties will still hold, which means that after finding a rational point \(P\), one can generate infinitely many rational points by calculating \(2P=P+P,\,3P=2P+P,\,\cdots\). It is far more efficient at finding rational points than brute-forcing approach.

Elliptic curves are familiar in cryptography, due to its abelian group property, just like RSA. However, a quantum computer can crack cryptosystems based on an abelian group efficiently.

It is known that finding an integer solution (not necessary positive) is easy, therefore one can get several rational points on the curve by brute forcing integer solutions of original equation.

from fractions import Fraction as Frac
from math import gcd
n = 100
for x in range(-n, n + 1):
    for y in range(-n, n + 1):
        for z in range(0, n + 1):
            if gcd(x, gcd(y, z)) == 1:
                if x**3 + y**3 + z**3 - 3 * x**2 * (y + z) - 3 * y**2 * (
                        z + x) - 3 * z**2 * (x + y) - 5 * x * y * z == 0:
                    print((x, y, z))

(-11, -9, 5)
(-11, -4, 1)
(-9, -11, 5)
(-5, 9, 11)
(-5, 11, 9)
(-4, -11, 1)
(-1, -1, 1)
(-1, 0, 1)
(-1, 1, 0)
(-1, 1, 1)
(-1, 4, 11)
(-1, 11, 4)
(0, -1, 1)
(1, -1, 0)
(1, -1, 1)
(4, -1, 11)
(9, -5, 11)
(11, -5, 9)
(11, -1, 4)

Pick (-1, 4, 11) as an example, transform it with the matrix:

\[\left(\begin{array}{ccc} -\frac{95}{2} & -\frac{95}{2} & \frac{277}{12} \\ -\frac{91}{2} & \frac{91}{2} & 0 \\ -6 & -6 & 1 \\\end{array}\right) (-1,4,11)^T = (\frac{1337}{12},-\frac{455}{2},-7)^T = -7\cdot(-\frac{191}{12},\frac{65}{2},1)^T\]

This solution corresponds to point \(P=(-\frac{191}{12},\frac{65}{2})\).

This point lies outside the red region. However, do not worry, one can use point addition to get a lot of rational points (with code below) until a valid point occurs.

<<is_valid>>
<<ec_plus>>
P = S = (-Frac(191, 12), Frac(65, 2))
cnt = 1
while not is_valid(S):
    S = ec_plus(S, P)
    cnt += 1
    print("{}P = ({}, {})".format(cnt, *S))

2P = (29233/300, -237679/250)
3P = (-449023/71286, 263764935/5180116)
4P = (2595564718153/100280426700, -171959440568022731/1527855867796750)
5P = (378326742585987529/182304915530542572, 112103221774745790541183075/3745033395604274176822442)
6P = (733306806381056483585521/52901390736340463822400, -343764653760831645784970282294394569/12167479807190947672003083538368000)
7P = (12186898088049137933023410813097849/1960207391964596384635579845637932, 52847537824372995383455949500060960584963350540405/4175525694290461543876101536662682824637424706682)
8P = (184291903459371594568624235201241521313551833/17546136990312774398727858490642680612252300, -31167188094052845691703198380025832449772385839996621463599210371/3536147729595630697250285084695028312237650433064461603785465250)
9P = (637231263346693346010600868176842100660564986054619857/81086401756297731072425447670017660116155468002441526, 29400417578654041653688375863673590665042836425148365175935874356653994350305605/6284274674391913738225458829167585075352347781720125288787017561681121629385876)

The results above show that the denominators of coordinates of rational points are overgrowing.

Use inverse of original transformation matrix to transform \(9P\) back to correspond rational solution:

\[\left(\begin{array}{ccc} \frac{1}{182} & -\frac{1}{91} & -\frac{277}{2184} \\ \frac{1}{182} & \frac{1}{91} & -\frac{277}{2184} \\ \frac{6}{91} & 0 & -\frac{95}{182} \\\end{array}\right) \left(\begin{array}{c} x_{9P} \\ y_{9P} \\ 1\\\end{array}\right) = \left(\begin{array}{c} -\frac{36875131794129999827197811565225474825492979968971970996283137471637224634055579}{1143737990739328300357033506908500483714127296273062802559237196225964136548229432} \\ -\frac{154476802108746166441951315019919837485664325669565431700026634898253202035277999}{1143737990739328300357033506908500483714127296273062802559237196225964136548229432} \\ -\frac{9405501114660716625534518421595417184664937929417781}{2459620853274364509196905245990535690190049196074059622}\\\end{array}\right)\]

Then, by multiplying the lowest common multiple of denominators and taking the opposite number, one will finally get the positive integer solution!

x = 36875131794129999827197811565225474825492979968971970996283137471637224634055579
y = 154476802108746166441951315019919837485664325669565431700026634898253202035277999
z = 4373612677928697257861252602371390152816537558161613618621437993378423467772036

You may be wondering how to prove that this positive integer solution is the smallest one. If so, you will be interested in following facts.

What an awesome equation! However, one may be wondering a harder equation. So it comes to degree-4 Diophantine equations.

It seems that it is an endless effort, for one can always increase the degree of an equation to make it harder. However, degree-4 is hard enough, since Mathematicians have found a degree-4 Diophantine equation (with \(58\) variables, and some of them are parameters) that is universal.

That means, one can reduce any Diophantine equation to a degree-4 Diophantine equation, and there does not exist an algorithm that can determine whether this equation has a solution (given values of parameters). As a result, mathematics cannot do that too. We finally reached the limit of computation and mathematics!

This step is not an ordinary step, but a step towards the abyss of unknowns. It is related to Hilbert's tenth problem. The author will post an article on this later. (On Hilbert's tenth problem)

1. An unusual cubic representation problem Andrew Bremner, Allan Macleod
2. https://en.wikipedia.org/wiki/Elliptic_curve Elliptic curve - Wikipedia