The only dubious thing I can see is the introduction of the equation $\cot x = \cot y$; this appears problematic since (i.) the cotangent functions have singularities in $(-\pi, \pi)$; and (ii.) I can see no way to infer the values of $x$ and $y$ given by $P_1$-$P_5$ from this equation; indeed, all it directly tells us is that $x = y$, $x = y + \pi$ or $x = y -\pi$ on $(-\pi, \pi) \setminus \{ 0 \}$. As far as issues with Wolfram Alpha are concerned, I can’t say much as I don’t know that much about it; but I did notice xsr’s input syntax appears different than Mary Star’s, so perhaps there is a clue to any mysteries.
These things being said, here’s my analysis of this problem:
As our O.P. xsr correctly showed,
$f_x =\dfrac{\partial f}{\partial x} = (\cos x)(\sin y), \tag{1}$
and similarly,
$f_y = (\sin x)(\cos y). \tag{2}$
Thus we seek those points
$(x, y) \in (-\pi, \pi) \times (-\pi, \pi) \tag{3}$
such that
$(\cos x)(\sin y) = 0, \tag{4}$
and
$(\sin x)(\cos y) = 0. \tag{5}$
Then from (4),
$\cos x = 0 \;\; \text{or} \;\; \sin y = 0, \tag{6}$
whilst from (5)
$\sin x = 0 \;\; \text{or} \;\; \cos y = 0. \tag{7}$
We recall that $\sin x$, $\cos x$ cannot simultaneously vanish, and likewise for $\sin y$, $\cos y$; thus if
$\cos x = 0, \tag{8}$
then $\sin x \ne 0$, so we must have
$\cos y = 0; \tag{9}$
it follows from (8) and (9) that
$x, y = \pm \dfrac{\pi}{2}. \tag{10}$
Next, if
$\cos x \ne 0, \tag{11}$
(4) yields
$\sin y = 0, \tag{12}$
hence
$\cos y \ne 0, \tag{13}$
hence, from (5),
$\sin x = 0. \tag{14}$
The only solutions to (12), (14) with $x, y \in (-\pi, \pi)$ are
$x = y = 0. \tag{15}$
We have thus covered all possible cases according to whether $\cos x = 0$ or not; we see the critical points are
$(x, y) = \left(\pm \dfrac{\pi}{2}, \pm \dfrac{\pi}{2}\right), (0, 0); \tag{16}$
we next proceed to classify them according to type. We do this by evaluating the Hessian of $f$, $H_f$:
$H_f(x, y) = \begin{bmatrix} f_{xx} & f_{xy} \\ f_{yx} & f_{yy} \end{bmatrix} = \begin{bmatrix} -(\sin x)(\sin y) & (\cos x)(\cos y) \\ (\cos x)(\cos y) & -(\sin x)(\sin y) \end{bmatrix}. \tag{17}$
From (17), we calculate that
$H_f\left(\dfrac{\pi}{2}, \dfrac{\pi}{2}\right) = H_f\left(-\dfrac{\pi}{2}, -\dfrac{\pi}{2}\right) = \begin{bmatrix} -1 & 0 \\ 0 & -1 \end{bmatrix}; \tag{18}$
$H_f\left(-\dfrac{\pi}{2}, \dfrac{\pi}{2}\right) = H_f\left(\dfrac{\pi}{2}, -\dfrac{\pi}{2}\right) = \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix}, \tag{19}$
and
$H_f(0, 0) = \begin{bmatrix} 0 & 1 \\ 1 & 0 \end{bmatrix}. \tag{20}$
We immediately see from (18), (19) that the points $P_1$, $P_4$ are local maxima, whilst $P_2$, $P_3$ are local minima; these assertions follow from the manifest fact that the eigenvalues of $H_f(P_1)$ and $H_f(P_4)$ are $-1, -1$, and those of $H_f(P_2)$, $H_f(P_3)$ are $1, 1$; as for $P_5 = (0, 0)$͵ the characteristic polynomial of $H_f(P_0)$ is
$\det(H_f(0, 0) - \lambda I) = \det\left(\begin{bmatrix} -\lambda & 1 \\ 1 & -\lambda \end{bmatrix}\right)$
$= \lambda^2 - 1; \tag{21}$
the eigenvalues are thus $\lambda = \pm 1$, being the roots of $\lambda^2 - 1 = 0$; $P_5$ is indeed a saddle point.