Literature reference for divergence theorem (Greens' first theorem) with surface discontinuity

General
1

In some applications (e.g. Poisson equations) multidimensional integration by parts, or Green’s first theorem (applying the divergence theorem), is helpful for deriving a weak formulation or handling functionals for variational calculus. Wikipedia helpfully provides a summary at Vector calculus identities - Wikipedia (I should send them some money)

\int_V \mathbf{A} \cdot \nabla \psi = \int_S \psi (\mathbf{A},n) ds - \int_V \psi \nabla \cdot \mathbf{A}

But we know from electrostatics that more attention to detail at the boundary is needed. The surface term is actually

\int_S \psi (\mathbf{A}_{out}-\mathbf{A}_{in},n) ds = \int_S \psi \sigma ds

involving a jump discontinuity in vector \mathbf{A} at the boundary, defining a scalar function on the boundary, \sigma=(\mathbf{A}_{out}-\mathbf{A}_{in},n) from the transverse components of vector \mathbf{A} on the inside and outside of the boundary. In electrostatics this is the surface charge, corresponding to a discontinuity in the electric displacement vector \mathbf{D}. (I may have the minus sign the wrong way around, which is the point of this request for reference)

The references I have at hand (Gradstein & Ryzhik §10.713, Riley, Hobson & Bence eq.11.19), that apply the divergence theorem to derive this vector calculus identity, all ignore the possibility of the jump condition, meaning they only allow for the case \sigma=0.

Does anyone know a good reference for the more general case involving the surface discontinuity?

Preferably something official with a doi that could be used for literature citation, but hey even a blog entry could be helpful for checking the minus signs are correct.

2

I guess what I would do is to divide the volume into two chunks, split at the discontinuity. Then you can apply Green’s theorem on either of them, where you get a common integral at the discontinuous boundary. Then you can choose whatever side the normal is defined at n_1 = -n_0 and get a jump condition.

This is typically how I arrive at DG formulations (before apply magics where jumps can become mixes of jumps and averages.

3

I guess I can try doing it manually if I have to. I’m missing some conceptual detail though. The two subdomains might be V_{in} and V_{out}, and V_{in} is the one I’m interested in. I can’t yet see how I should subtract out \int_{V_{out}}.

I think you meaning writing

\int_S \psi (\mathbf{A},n_{in}) ds = \int_{V_{in}} \mathbf{A} \cdot \nabla \psi + \int_{V_{in}} \psi \nabla \cdot \mathbf{A}
{} = \int_{V_{out}} \mathbf{A} \cdot \nabla \psi + \int_{V_{out}} \psi \nabla \cdot \mathbf{A} = \int_S \psi (\mathbf{A},n_{out}) ds

i.e. identifying \int_S \psi (\mathbf{A},n_{in}) ds =-\int_S \psi (\mathbf{A},n_{out}) ds (that is n_{out}=-n_{in}). But that’s not quite the same as a jump condition.

4

What I was thinking was split V=V_0\cup V_1, such that \Gamma = V_0\cap V_1

\int_V \mathbf{A}\cdot \nabla \psi ~\mathrm{d} x = \int_{V_0} \mathbf{A}\cdot \nabla \psi ~\mathrm{d} x + \int_{V_1} \mathbf{A}\cdot \nabla \psi ~\mathrm{d} x

Then apply the identity to each subdomain

\int_{V_0} \mathbf{A}\cdot \nabla \psi ~\mathrm{d} x = -\int_{V_0} \psi \nabla \cdot \mathbf{A}~\mathrm{d}x + \int_{S_0} \psi_0 \mathbf{A}_0\cdot n_0 ~\mathrm{d}s.
\int_{V_1} \mathbf{A}\cdot \nabla \psi ~\mathrm{d} x = -\int_{V_1} \psi \nabla \cdot \mathbf{A}~\mathrm{d}x + \int_{S_1} \psi_1 \mathbf{A}_1\cdot n_1 ~\mathrm{d}s.

Split S_0 into \Gamma \cup \Lambda_0 and S_1 into \Gamma \cup \Lambda_1.
This means \Lambda_0\cup\Lambda_1=\partial V

Pick n_0=-n_1
to get

\int_V \mathbf{A}\cdot \nabla \psi ~\mathrm{d} x = -\int_{V} \psi \nabla \cdot \mathbf{A}~\mathrm{d}x + \int_{\Lambda_0\cup\Lambda_1} \psi \mathbf{A}\cdot n ~\mathrm{d}s + \int_\Gamma (\psi_1 \mathbf{A_1} - \psi_0 \mathbf{A}_0)\cdot n_1\mathrm{d}s

If \psi is continuous you get your identity.
This means that the discontinuity cannot be at an outer boundary, as there it wouldn’t make sense (as you would always approach it from the inside).

5

I think I see. This way identifies an internal boundary \Gamma, distinct from the overall external boundary \Lambda=\partial V. So here \Lambda_0=\partial V_0/\Gamma is V_0's share of the external boundary, without the internal boundary. Likewise \Lambda_1=\partial V_1/\Gamma. And you’re defining the jump condition at the internal boundary \Gamma, not the external boundary.

It may be that the conceptual flaw in what I’m trying to do is to define my jump condition at an external boundary. Fair to say, it doesn’t really make sense to define an n_1 normal if V_1 does not exist. I will reflect on this a little longer.

In practice what I’m trying to do is perform my FEniCS calculation in one domain, which represents one piece of the universe. I’m trying to avoid having to do an explicit calculation of the rest of the universe outside the domain I’m interested in. So it’s not so much that V_1 does not exist, but that I’m trying to ignore it. I may need to invoke electroneutrality of my [sub]domain and its boundary, in which case the vector field outside might be zero (so the jump condition on the “external” boundary is not needed).