Solution (source code)

= Solution

Set $\Delta=0$ and take $\xi_r=Cr^{\ell-1}$. The divergence equation then fixes
$$
\xi_h=\frac{Cr^\ell}{\ell}.
$$
Hydrostatic balance gives $\delta p=-\xi_rp'=\rho g\xi_r$. The horizontal momentum equation becomes
$$
\rho\omega^2\frac{Cr^\ell}{\ell}
=\rho\Omega^2r\,Cr^{\ell-1},
$$
so
$$
\boxed{\omega_\ell^2=\ell\Omega^2}.
$$
The radial equation gives the same result because $\omega^2\xi_r=g\xi_r'+g'\xi_r$. This is an <incompressible stellar surface mode>: it changes the shape of the free surface without compressing fluid elements. For $\ell=1$ it is a rigid displacement of the star in the fixed harmonic potential.