# Weak Formulation of Poisson's Equation for Electrostatics with Surface Charge

I am currently attempting to use FEnICS to solve an electrostatic problem with two materials of different permittivity $$\varepsilon_1$$ and $$\varepsilon_2$$ forming an interface:

Consider a domain $$\Omega_1 \cup \Omega_2 = \Omega$$ with $$\Gamma_T$$, $$\Gamma_B$$, $$\Gamma_R$$ and $$\Gamma_L$$ as the top, bottom, right and left exterior boundaries of a simple rectangular domain and $$\Gamma$$ as the interface between the two dielectrics, the boundary problem is:

$$-\nabla^2 \varphi = \frac{\rho+\sigma}{\varepsilon} \\ \varphi = U_0 \ \textrm{on} \ \Gamma_T \\ \varphi = U_1 \ \textrm{on} \ \Gamma_B \\ (\nabla \varphi)\cdot\hat{n} = 0 \ \textrm{on} \ \Gamma_R, \ \Gamma_L$$ and finally, the condition due to the surface charge accumulation at the interface: $$\hat{n}\cdot(D_1-D_2) = \sigma$$ With $$D = \varepsilon E$$.

Where $$\sigma$$ is the surface charge. Particularly, I am struggling with the weak formulation of the final boundary condition, I am able to derive, with $$w$$ being a test function and the general form of $$-\nabla^2\varphi = f$$ :$$\int_{\Omega} \nabla\varphi \cdot\nabla w \ d\Omega - \int_{\partial\Omega}w(\nabla\varphi)\cdot\hat{n} \ ds = \int_{\Omega}fw \ d\Omega$$

As I understand, the second integral on the left shall vanish on Dirichlet boundaries and, of course, where a zero Neumann boundary exists.

However, how can I consider the surface charge condition in the weak form on $$\Gamma$$?

First, I think that we should write the problem as

$$-\nabla(\varepsilon \nabla \phi) = \rho + \sigma\, .$$

Then, the weak form would be

$$\int\limits_\Omega \varepsilon \nabla \phi \nabla w d\Omega - \int\limits_{\partial\Omega} \varepsilon w\nabla\phi\cdot\hat{n}d\Gamma = \int\limits_\Omega f w d\Omega\, .$$

As you mentioned, the second term vanishes for the boundary conditions that you have and you end up with

$$\int\limits_\Omega \varepsilon \nabla \phi \nabla w d\Omega - \int\limits_{\Gamma} \varepsilon w\nabla\phi\cdot\hat{n}d\Gamma = \int\limits_\Omega f w d\Omega\, .$$

But, your interface condition can be rewritten as

$$(-\varepsilon_1 \nabla\phi + \varepsilon_2 \nabla \phi)\cdot\hat{n} = (\varepsilon_1 - \varepsilon_2)\nabla \phi\cdot\hat{n} = \sigma\, ,$$

or

$$\nabla \phi\cdot\hat{n} = \frac{\sigma}{\varepsilon_1 - \varepsilon_2}\, .$$

Plugging this into your weak form

$$\int\limits_\Omega \varepsilon \nabla \phi \nabla w d\Omega = \int\limits_\Omega f w d\Omega + \int\limits_{\Gamma}\frac{w \varepsilon \sigma}{\varepsilon_1 - \varepsilon_2} d\Gamma \, .$$

And I am not sure what value $$\varepsilon$$ should have for the interface.

### Update 20210609

If I understood @ConvexHull properly, you have the following

$$-\nabla^2 \phi = \frac{\rho + \sigma}{\varepsilon}\, .$$

Then, the weak form would be

$$\int\limits_\Omega \nabla \phi \nabla w d\Omega - \int\limits_{\partial\Omega} w\nabla\phi\cdot\hat{n}d\Gamma = \int\limits_\Omega f w d\Omega\, .$$

As you mentioned, the second term vanishes for the boundary conditions that you have and you end up with

$$\int\limits_\Omega \nabla \phi \nabla w d\Omega - \int\limits_{\Gamma} w\nabla\phi\cdot\hat{n}d\Gamma = \int\limits_\Omega f w d\Omega\, .$$

But, your interface condition can be rewritten as

$$(-\varepsilon_1 \nabla\phi + \varepsilon_2 \nabla \phi)\cdot\hat{n} = (\varepsilon_1 - \varepsilon_2)\nabla \phi\cdot\hat{n} = \sigma\, ,$$

or

$$\nabla \phi\cdot\hat{n} = \frac{\sigma}{\varepsilon_1 - \varepsilon_2}\, .$$

Plugging this into your weak form

$$\int\limits_\Omega \nabla \phi \nabla w d\Omega = \int\limits_\Omega f w d\Omega + \int\limits_{\Gamma}\frac{w \sigma}{\varepsilon_1 - \varepsilon_2} d\Gamma \, .$$

• I think $\epsilon$ is in fact not defined on the interface $\Gamma$. Therefore, the problem consists more of coupling two homogeneous regions via the interface condition. There are two possible solutions: (i) You either have a jump in your gradient and a kink in your potential (ii) or you have a kink in your gradient and a jump in your potential. I think the latter is what the questioner is interested in. Jun 9, 2021 at 10:49
• To be more presice: $\epsilon$ has in fact two values on $\Gamma$. It is not unique. Jun 9, 2021 at 11:33