This is an ${\bf Edit}$ to my previous answers.
There are two questions here:
Is Fernando's reasoning correct?
Why are there problems in plotting the function $f$ near $(0,0)$?
Ad 1: As others have remarked the reasoning is correct and leads to the correct result that $f$ is differentiable at $(0,0)$ and has derivative $df(0,0)=0$ there. The cancelling of the factor $y^2$ (which could be $0$) is a little bit fishy, but this can be repaired by noticing that ${y^2\over y^2+r^3}\leq 1\qquad \bigl(r:=\sqrt{x^2+y^2}>0\bigr)\ .$
Ad 2: There is the possibility that $f$ is differentiable at $0$, but is not continuously differentiable there. Consider in this regard the one-variable example $f(t):=t^2\sin{1\over t}$, $\>f(0):=0$. In this case we have f'(0)=0, but f'(t)= 2t\sin{1\over t}-\cos{1\over t}\not\to 0 $\ (t\to 0)$.
In order to investigate the continuity of $df$ we are going to simplify our function $f$ considerably (this is already hinted at in Fernando's reasoning):
There are entire functions $z\mapsto g(z)$ and $z\mapsto h(z)$ with $g(0)=h(0)=1$ such that $\sin^2 z=z^2 \>g(z)\ ,\quad e^z-1=z \>h(z)\qquad\forall z\ .$
This means that up to a very smooth factor which is approximately $1$ in the neighborhood of $(0,0)$ we can replace $f$ by the following simpler function: $f(x,y):={x^2y^2\over y^2+r^3} \qquad \bigl(r:=\sqrt{x^2+y^2}>0\bigr)\ .$ We now look at the partial derivatives of this new $f$. The computation gives $f_x={r\sin^2\phi\cos\phi\bigl(r-(2+3r)\sin^2\phi\bigr)\over(r+\sin^2\phi)^2}\ ,$ where we have used the substitution $(x,y):=(r\cos\phi,r\sin\phi)$. As $r>0$ there is a $\tau\in{\mathbb R}$ with $\sin\phi=\tau\sqrt{r}$. It follows that $|f_x|\leq {r^3\tau^2(1+3 \tau^2)\over r^2(1+\tau^2)^2}\leq C r$ for some $C>0$ which does not depend on $r$ or $\tau$. Therefore $\lim_{(x,y)\to(0,0)}f_x(x,y)=0=f_x(0,0)$.
A similar computation gives $f_y={r^2\sin\phi\cos^2\phi (2-3\sin^2\phi)\over(r+\sin^2\phi)^2}\ .$ We again write $\sin\phi=\tau\sqrt{r}$ and then have $|f_y|\leq {2r^{5/2}|\tau|\over r^2(1+\tau^2)^2}\leq C\sqrt{r}$ for some $C>0$ which does not depend on $r$ or $\tau$. Therefore $\lim_{(x,y)\to(0,0)}f_y(x,y)=0=f_y(0,0)$ as well.
All in all we havce proven that the given function $f$ is in fact continuously differentiable in all of ${\mathbb R}^2$. So where are the problems in plotting $f$ coming from? Plotting $f_x$ and $f_y$ as a function of $\phi$ for very small $r>0$ shows sharp peaks (of small amplitude), which indicates that the second derivatives of $f$ are no longer bounded. Maybe the used surface plotting programs don't like that.