La ecuación es soluble en enteros. Toma, por ejemplo, $$ x = -19578556686240310295378317903565, \\ y = -101658411567714319887, \\ z = 418962851513108789978912616277180591709694. $$
La verificación puede hacerse por sustitución.
Encontré esta solución utilizando las transformaciones descritas por Denis Shatrov en el comentario. Consideramos la ecuación $$ 1 + x^2 + x^3 + y^2 + x y z = 0 \quad\quad (1) $$ y buscamos una solución tal que $z$ sea divisible por $9$. Si partimos de una solución cualquiera $(x_0,y_0,z_0)$, entonces, como observó Denis, $$ (x_1,y_1,z_1)=\left(\frac{y_0^2+1}{x_0}, y_0, -\frac{1+x_0+x_0^3+y_0^2}{x_0y_0}\right) $$ resuelve la ecuación $$ 1 + x + x^3 + y^2 + x y z = 0 \quad\quad (2). $$ Luego $$ (x_2,y_2,z_2)=\left(x_1, \frac{x_1^3+x_1+1}{y_1}, -\frac{1+x_1+x_1^3+y_1^2}{x_1y_1}\right) $$ también es una solución de (2), mientras que $$ (x_3,y_3,z_3)=\left(\frac{y_2^2+1}{x_2}, y_2, -\frac{1+x_2^2+x_2^3+y_2^2}{x_2y_2}\right) $$ es nuevamente una solución de (1).
Haciendo un análisis módulo 9, observé que si $(x_0,y_0,z_0)$ es $(4,0,3)$ módulo $9$, entonces $(x_3,y_3,z_3)$ es $(4,6,0)$ módulo $9$. Una búsqueda sencilla en computadora devolvió la solución $(x_0,y_0,z_0) = (-3965, 1446687, 354)$ a (1) que es $(4,0,3)$ módulo $9$. Entonces el correspondiente $(x_3,y_3,z_3)$ es una solución de (1) con $z$ divisible por (9), por lo tanto $(x_3,y_3,z_3/9)$ es una solución entera de la ecuación original $1 + x^2 + x^3 + y^2 + 9x y z = 0$. Como se mencionó anteriormente, su corrección puede verificarse fácilmente mediante la sustitución directa.
