Let a, b, c be a solution of the given equation. Putting
ba=X,cb=Y,
we have c/a=(1/XY). Note that X and Y are rational numbers. The equation takes the form
X2Y+XY2+1=3XY.(1)
This is a cubic curve. We observe that (1,1) is a point on the curve (1). If we know the rational points on (1), then we can get integer solutions of the given equation.
Suppose (x,y) is a rational point on the curve (1). Then the line joining (x,y) to (1,1) has rational slope. Conversely, any line through (1,1) having rational slope intersects the curve (1) in a point having rational coordinates.
Consider the line
x−1y−1=m,
where m=−1 is a non-zero rational number. Let its intersection with (1) be (X0,Y0)=(1,1). Then we have Y0=m(X0−1)+1. Putting this in (1), we get
X02(m(X0−1)+1)+X0(m(X0−1)+1)2+1=3X0(m(X0−1)+1).
This reduces to
(X0−1)2(m(1+m)X0+1)=0.
Since X0=1, we have
X0=−m(1+m)1,Y0=m(X0−1)+1=−(1+m)m2.
Thus
ba=−m(1+m)1,cb=−(1+m)m2,ca=(1+m)2m.
Suppose m=v/u, where v, u are integers and gcd(v,u)=1. We get
ca=v(v+u)2uv2.
This gives
uv2a=v(v+u)2c.
Similarly, the expression for a/b gives
v2a=−u(u+v)b.
We thus obtain
uv2a=−u2(u+v)b=v(v+u)2c.
Let each of these ratios be p/q for some integers p, q. We obtain
qa=puv2,qb=−pu2(u+v),qc=pv(u+v)2.
This shows that gcd(qb,qc)=±p(u+v). Thus
±q(b,c)=±p(u+v).
But qa=puv2. It follows that ±qgcd(a,b,c)=p. Hence
qp=±gcd(a,b,c)=k,
say. Finally, we have
a=kuv2,b=−ku2(u+v),c=kv(u+v)2.
It is easy to check that {a,b,c} satisfies the given relation.