With
$f(x, y) = (x + y)(1 - xy) \tag 1$
and
$\nabla f = (f_x, f_y), \tag 2$
we see that
$f_x = 1 - xy + (x + y)(-y) = 1 - xy - xy - y^2 = 1 - 2xy - y^2, \tag 3$
and likewise
$f_y = 1 - xy + (x + y)(-x) = 1 - xy - x^2 - xy = 1 - 2xy - x^2; \tag 4$
at critical points of $f(x, y)$, we have
$\nabla f(x, y) = 0, \tag 5$
whence from (3) and (4), at the critical points,
$1 - 2xy - y^2 = 0 = 1 - 2xy - x^2; \tag 6$
we solve this system by observing that it implies
$y^2 = 1 - 2xy = x^2, \tag 7$
so that
$y = \pm x; \tag 8$
using (8) in (6) we may write an equation for $x$:
$0 = x^2 + 2xy - 1 = x^2 \pm 2x^2 - 1; \tag 9$
$x$ must thus obey
$3x^2 - 1 = 0 \tag{10}$
or
$x^2 + 1 = 0; \tag{11}$
we rule out (11) since $x$ is real; thus
$x = \pm \dfrac{1}{\sqrt 3} = \pm \dfrac{\sqrt 3}{3}; \tag{12}$
again from (6),
$y = \dfrac{1 - x^2}{2x} = \dfrac{\dfrac{2}{3}}{2x} = \dfrac{1}{3x}; \tag{13}$
therefore the critical points are
$(x, y) = \left ( \dfrac{\sqrt 3}{3}, \dfrac{\sqrt 3}{3} \right ), \; (x, y) = \left ( -\dfrac{\sqrt 3}{3}, -\dfrac{\sqrt 3}{3} \right ); \tag{14}$
the Hessian $H_f$ of $f(x, y)$ has been provided for us courtesy of our OP kronos:
$H_f = \begin{bmatrix} -2y & -2x-2y \\ -2x-2y & -2x \end{bmatrix}; \tag{15}$
we thus see that
$\det(H_f) = 4xy - 4(x + y)^2 = 4(xy - (x + y)^2) = -4(x^2 +xy + y^2); \tag{16}$
since at the critical points we have
$x = y = \pm \dfrac{\sqrt 3}{3}, \tag{17}$
it follows that at these points
$\det(H_f) = -4, \tag{18}$
which, as is well-known, implies that each critical point is a saddle.
Best Answer
To find the critical points, solve the simultaneous equations
$$f_x=3x^2+9y \tag{1}$$ $$f_y=-3y^2+9x \tag{2}$$
From (2), without introducing additional solutions or losing any, we get $x=\tfrac{1}{3}y^2$, so by substituting into (1):
$$\tfrac{1}{3}y^4+9y=0 \iff \tfrac{1}{3}y(y+3)(y^2-3y+9)=0$$
The quadratic term has no real roots, so the only solutions are $y=0,-3$, so solutions are:
$$(x,y)\in\{(0,0),(3,-3)\}$$
To find the nature of the critical points, you will need the partial derivatives of $f(x,y)$. These are:
$$\begin{align} f_{xx} &= \frac{\partial}{\partial x}f_x&=6x \\[2ex] f_{yy} &= \frac{\partial}{\partial y}f_y&=-6y \\[2ex] f_{xy} &= \frac{\partial}{\partial x}f_y&=9 \end{align}$$
Then conduct the second partial derivative test - as uniquesolution has commented, this involves the Hessian. In short, the sign of the determinant $D(x,y)$ and the sign of either $f_{xx}$ or $f_{yy}$ (your choice) determines the nature of the critical point. Wikipedia is a good source for this.
Now
$$D(x,y)=f_{xx}f_{yy}-(f_{xy})^2=-36xy-81$$
So for the two critical points:
Addendum: Other Cases
A. $D(x,y)=0$
Many of the ideas in this section (and the first two examples) are taken from https://www3.nd.edu/~eburkard/Teaching/20550%20F13/Second%20Derivatives%20Test.pdf.
When the determinant is zero, the Hessian matrix is said to be degenerate. I know of no higher-order (third derivatives or higher) test that can determine the nature of the critical points. Several examples will show that behaviour can vary markedly in cases where $f_x=f_y=D(x,y)=0$ at a critical point:
B. $D(x,y)>0\text{ and }f_{xx}=0$
When $f_{xx}=0$ we have $D(x,y)=f_{xx}f_{yy}-(f_{xy})^2=0-(f_{xy})^2\le 0$, so this case is not possible.