solving diophantine problem $a^3+b^3=2019(1+ab)$ for coprime $a$ and $b$ We're interested in solving
$$\begin{cases} a^3+b^3=2019(1+ab) \\ \gcd(a,b)=1 \end{cases}$$
I'm stuck with the deduction

*

*Say why $a^3 \equiv -b^3 \pmod{2019}$ . (done)


*Using Fermat, prove that $a^{672} \equiv b^{672} \pmod{2019}$ . (done)


*Deduce that $a \equiv -b \pmod{2019}$ (note that $673 = 1 + 244*3$) (stuck)
thanks for any help.
 A: Here's a method that will work, though you'll need a computer like Wolfram Alpha to do the computations... unless you really want to do it by hand, which is also doable.
The idea is to use the factorization $$x^3+y^3+z^3-3xyz=(x+y+z)(x^2+y^2+z^2-xy-yz-zx).$$ Luckily, $$2019 = 3\cdot 673,$$ where the factor of $3$ is exceedingly important for the applicability of this technique. Then the equation $$a^3+b^3=2019(1+ab)$$ can be rewritten as $$a^3 + b^3+673^3 - 3ab\cdot 673 =3\cdot 673 +673^3$$ or equivalently $$(a+b+673)(a^2+b^2+673^2-ab-673a-673b)=2^2\cdot 673\cdot 113233.$$ The rest is easy, albeit computationally intensive. The number of positive divisors of the right side is $$(2+1)(1+1)(1+1)=12,$$ so there are $24$ ordered pairs of integers $(s,t)$ that can multiply to the right side, since we have to account for negative factors. So set the $a+b+673$ equal to each possible factor $s,$ isolate $b$ and substitute it into $$a^2+b^2+673^2-ab-673a-673b=t.$$ That will give you a quadratic in $a$ that you can solve using the quadratic formula. If you get an integer, find the corresponding $b$ using $b=s-a-673.$ Good luck with the details.
A: This is simple but not elementary method.
$$a^3+b^3=2019(1+ab)\tag{1}$$
Let $s=a+b$ and substitute $b=s-a$ to equation $(1)$, then we get
$(3s+2019)a^2+(-3s^2-2019s)a+s^3-2019=0\tag{2}$
In order to quadratic equation for $a$ has rational solution, the discriminant must be a rational number.
$v^2 = -3s^4+4038s^3+4076361s^2+24228s+16305444$
So, above quartic equation must have integral solutions.
Using online Magma calculator as follows.
IntegralQuarticPoints($[-3,4038,4076361,24228,16305444]$);
It says that all integral points are $[[-673, 0 ], [ 0, 4038 ], [ 2019, 8076 ]]$.
Hence all integral solutions are $(a,b)=(1,-1),(-1,1),(1010, 1009),(1009, 1010).$
