No sé si mi enfoque es correcto, pero encontré que:
$$f(x) = A \, x \pm \sqrt{1+A^2}, \quad A \in \mathbb{R},$$
es una solución de tu EDO. Usé el método de Lagrange-Charpit para resolver una "falsa" PDE de 1er orden para $f = f(x,y)$, asumiendo que no hay dependencia de $y$.
Si deseas más detalles, solo pídelos.
Editar
Las ecuaciones de Lagrange-Charpit establecen que, para una ecuación en derivadas parciales no lineal escrita en la forma:
$$G(f(x_1,x_2,\ldots,x_N),p_i,x_i) = 0,$
donde $f$ es una función desconocida (variable dependiente), $p_i = \frac{\partial f}{\partial x_i}$, o en el caso especial de $N = 2$, $x_i = (x,y)$, para simplificar la notación, se cumplen las siguientes relaciones:
$$ \frac{dx}{G_p} = \frac{dy}{G_q} = \frac{df}{pG_p + qG_q} = - \frac{dp}{G_x + pG_f} = -\frac{dq}{G_y + q G_f}, $$
Dado que:
$$\begin{align} G & = (f-xp)^2 -p^2-1 = 0, \\ G_x & = 2(f-xp)(-p), \\ G_y & = 0, \\ G_f & = 2 (f-xp), \\ G_p & = 2 (f-xp)(-x)-2p, \\ G_q & = 0, \end{align} $$ tienes:
$$\frac{dx}{2(f-xp)(-x)-2p} = \frac{dy}{0} = \frac{df}{2p(f-xp)(-x)-2p^2} =$$ $$= \frac{dp}{2(f-xp)(-p)+2p(f-xp) } = \frac{dp}{0} = \frac{dq}{2q(f-xp)}.$$
Esto te dice, entre otras cosas, que $dp = 0$, entonces $p = A$ es una constante. Dado que $p = f_x$ puedes integrar para obtener $f$, así que:
$$f(x,y) = A x + B(y),$$
para $B(y)$ alguna función arbitraria de $y$. Sustituye de nuevo en la ecuación original, $G = 0$ para obtener una relación entre $A$ y $B$, que resulta ser:
$$B^2(y) - A^2 - 1 = 0,$
así que $$B(y) = B = \pm \sqrt{1+A^2},$$ y $B$ no depende de $y, como se esperaba porque $f$ originalmente no dependía de $y$.
En resumen, he tratado tu EDO original como una PDE no lineal de 1er orden para la variable dependiente $f = f(x,y)$, pero, como ya sabíamos, $f = f(x)$, como se ha demostrado.
¡Saludos!