I think this surface has a rational parameterization in terms of (a,b), given by:
$x = (1 + a + a^2)^2/(9 (3 + a (6 + (-1 + a)^2 a)) b^2) + ((-2 + (-2 + a) a) b)/(1 + a + a^2)$
$y = (1 + a + a^2)^2/(9 (3 + a (6 + (-1 + a)^2 a)) b^2) + (b + 2 a b)/(1 + a + a^2)$
$z = (1 + a + a^2)^2/(9 (3 + a (6 + (-1 + a)^2 a)) b^2) + (b - a^2 b)/(1 + a + a^2)$
For rational (a,b), this should give you a dense set of rational points...
Let me explain where this parameterization comes from, so it'll be clear that this indeed shows that rational solutions are dense.
Consider the (linear, rational) change of variables:
$x=p+r$,
$y=q+r$,
$z=r−p−q$.
The equation then simplifies to: $p^2+pq+q^2=1/(9r)$. Since the quadratic is psd, for real points we need $r>0$. Each slice (for fixed $r$) is just an ellipse. If $r$ is the square of a rational number (dense on $\mathbb{R}$), we can parameterize all the rational solutions for that slice. Taking the union over real slices (for all $r>0$), we're done.
