I'm confused about the following problem I encountered. Given a set $\Omega$, which has on its boundary $\partial\Omega$ a surface charge (per area) $q_0$, I want to show that $$\int_\Omega E^2=\int_\Omega \rho\Phi,$$ where $E$ is the electric field $E=-\nabla\Phi$, $\Phi$ the potential and $\nabla E=-\rho$, $\rho$ the charge density.
Here, $\rho = q_0\delta(\partial\Omega)$, where the $\delta$-distribution gives me the condition that the density is localized on the surface.
Now I need to be very clear that $\Omega\supseteq\partial\Omega$. Given this, lets start from the lhs and use $$\int_\Omega E^2=\int_\Omega (\nabla\Phi)^2=\int_{\partial\Omega}\vec{n}(\nabla \Phi)\Phi-\int_\Omega\Phi\nabla^2\Phi,$$ where in the second equality, I used Green's formula.
But now comes the problem. I can insert for the first term on the RHS $\nabla \Phi=-E$. Because of $\nabla E=-\rho$ we have $$\int_{\partial\Omega}\vec{n}E=-\int_\Omega\rho=-\int_{\partial\Omega}q_0.$$ Now suppose our problem is symmetric, such that all quantities are constant on $\partial\Omega$ (e.g. a sphere).
This implies then: $\vec{n}E=-q_0$.
On the other hand, the second term reads $-\int_\Omega\Phi\nabla^2\Phi=\int_\Omega\Phi\rho$. Again using the fact the charge is localized on the surface, our equation reads $$\int_\Omega E^2=-\int_{\partial\Omega}q_0\Phi+\int_\Omega\rho\Phi=0$$
What am I doing wrong here? Normally, the boundary terms drops out since we assume $\Phi=0$ on the boundary. But here, without that assumption, it seems to break this relation (the first one I wrote) completely. Is there some argument why the Dirichlet boundary condition is natural?