Puede que esto no sea lo más fácil (tal vez por mucho no la más fácil) para llegar al resultado correcto, pero a mí siempre me ha parecido la forma más intuitiva de derivar cualquier ecuación en cualquier sistema de coordenadas, así que he pensado que esto puede ayudarte tanto como me ayudó a mí.
En primer lugar, algunas aclaraciones:
$$ \vec{x} = x(r, \phi) \vec{e}_x + y(r, \phi)\vec{e}_y = r(x,y)\vec{e}_r+\phi(x,y)\vec{e}_{\phi}, $$ donde $$ \begin{align} x(r, \phi) &= r\cos\phi\\ y(r, \phi) &= r\sin\phi\\ r(x,y) &=\sqrt{x^2+y^2}\\ \phi(x,y) &=\arctan\left( \frac{y}{x} \right). \end{align} $$ Primero queremos encontrar los vectores base $\vec{e}_{\phi}$ y $\vec{e}_{r}$ como funciones de $r$ , $\phi$ , $\vec{e}_x$ y $\vec{e}_y$ . Los encontramos por:
$$ \begin{align} \vec{e}_r&:=\frac{ \frac{\partial \vec{x}}{\partial r} }{ \left| \left| \frac{\partial \vec{x}}{ \partial r} \right| \right| } = \cos\phi\vec{e}_x + \sin\phi\vec{e}_y \\ \vec{e}_{\phi}&:=\frac{ \frac{\partial \vec{x}}{\partial \phi} }{ \left| \left| \frac{\partial \vec{x}}{ \partial \phi} \right| \right| } = -\sin\phi\vec{e}_x + \cos\phi\vec{e}_y. \end{align} $$
Con esto podemos encontrar los vectores base $\vec{e}_x$ y $\vec{e}_y$ como funciones de $r$ , $\phi$ , $\vec{e}_r$ y $\vec{e}_{\phi}$ . (Por ejemplo, invirtiendo la matriz de coeficientes del sistema de ecuaciones lineales anterior, lo cual es muy fácil, ya que sólo consiste en una matriz de rotación). Esto nos da: $$ \begin{align} \vec{e}_x &= \cos\phi\vec{e}_r - \sin\phi\vec{e}_{\phi} \\ \vec{e}_y &= \sin\phi\vec{e}_r + \cos\phi\vec{e}_{\phi}. \end{align} $$
El siguiente paso es escribir las derivadas parciales $\frac{\partial}{\partial i} =: \partial_i$ en términos de las nuevas coordenadas. Esto se logra mediante: $$ \begin{align} \partial_x &= \frac{\partial r}{ \partial x} \partial_r + \frac{\partial \phi}{\partial x}\partial_{\phi} = \cos\phi\partial_r - \frac{\sin\phi}{r}\partial_{\phi} \\ \partial_y &= \frac{\partial r}{ \partial y} \partial_r + \frac{\partial \phi}{\partial y} \partial_{\phi} = \sin\phi\partial_r + \frac{\cos\phi}{r}\partial_{\phi}. \end{align} $$
Ahora debería ser fácil calcular el gradiente a ser: $$ \vec{\nabla} = \vec{e}_x\partial_x + \vec{e}_y\partial_y = \vec{e}_r\partial_r + \vec{e}_{\phi}\frac{1}{r}\partial_{\phi}. $$
Pero como tu pregunta incluye el laplaciano $\Delta= \mathrm{div}(\mathrm{grad}) $ también debemos saber qué significa tomar la divergencia de un campo vectorial $\vec{A}=A_r(r,\phi)\vec{e}_r+A_{\phi}(r,\phi)\vec{e}_{\phi} = A_x(x,y)\vec{e}_x + A_{y}(x,y)\vec{e}_{y}$ dado en coordenadas polares. Como $\mathrm{div}(\vec{A})=\partial_xA_x+\partial_yA_y$ sólo tenemos que encontrar $A_x$ y $A_y$ como funciones de $r$ , $\phi$ , $A_r$ y $A_{\phi}$ : $$ A_x(r, \phi)=\langle \vec{e}_x, \vec{A} \rangle = \cos\phi A_r - \sin\phi A_{\phi} \\ A_y(r, \phi)=\langle \vec{e}_y, \vec{A} \rangle = \sin\phi A_r + \cos\phi A_{\phi}. $$
Ahora podemos calcular $$ \mathrm{div}(\vec{A})=\partial_xA_x+\partial_yA_y=\frac{1}{r}\partial_{\phi}A_{\phi} + \frac{1}{r}A_r + \partial_r A_r. $$
Si ahora aplicamos este resultado al gradiente, que ya hemos calculado, llegamos a $$ \begin{align} \Delta=\mathrm{div}(\mathrm{grad})&=\partial_x \vec{\nabla}_x + \partial_y \vec{\nabla}_y\\ &=\frac{1}{r}\partial_{\phi}\vec{\nabla}_{\phi} + \frac{1}{r}\vec{\nabla}_r + \partial_r \vec{\nabla}_r\\ &=\frac{1}{r}\partial_{\phi}\left( \frac{1}{r} \partial_{\phi} \right) + \frac{1}{r}\partial_r+\partial_r(\partial_r)\\ &=\frac{1}{r^2}\partial_{\phi}^2+\frac{1}{r}\partial_r+\partial_r^2. \end{align} $$
Ahora podemos ver que la ecuación de onda $\Delta u(x,y,t)=\frac{1}{c^2}\partial_t^2u(x,y,t)$ se convierte de hecho en $$ \begin{align} \Delta u(r,\phi,t)&=\left[\frac{1}{r^2}\partial_{\phi}^2+\frac{1}{r}\partial_r+\partial_r^2\right] u(r,\phi,t)\\ &=\frac{1}{r^2}\partial_{\phi}^2u+\frac{1}{r}\partial_ru+\partial_r^2u\\ &=\frac{1}{r^2}\partial_{\phi}^2u+\frac{1}{r}\partial_r\left(r \partial_ru\right) = \frac{1}{c^2} \partial_t^2 u. \end{align} $$
