Special Boundary-Condition

Dear Forum,
I am working with NGSolve on a nice application to plasma physics.

To make my model even better, I’m trying to add a non-trivial (and non-linear) boundary condition (BC) of the form E_t = grad(f(E)), where f(E) is computed externaly, to a mixed FEM Maxwell formulation.

I see several possibilities to impose such a condition but I’m always running into problems.

  1. I can compute grad(f(E)) by post-processing and impose E_t as a dirichlet BC. This would be the easiest solution.
    Problem 1: I did not find a clever way to only impose E_t (leaving E_n free) as the gfu.Set() function only works for the whole E field.
    Related problem: I did not find a way to select the vector En (only its amplitude by scalar product of E with the special normal coefficient function n). Is this possible?

  2. I can compute grad(f(E)) by post-processing and impose E_t through the surface term of my Maxwell equations of the form:
    (v \cross n) \vdot ( curl(E) ).
    However, as my BC only involves E_t, I end up with quite a tedious expression. My current attempt of implementing this method failed.

  3. I could add an additional test and trial function to express the weak form of my BC as a second equation and introduce this new E_t variable into the surface term of my weak Maxwell formulation. I did not test this yet but it probably suffers from the same difficulties found in (2).

Additional information:

  • As non-linear, all cases will need iterations.
  • I work with complex 2D geometries.

Is there a simpler technique that would make my life easier?
Can I avoid going through the Maxwell surface term appearing in the weak formulation?
Is there a better way to do this?

Looking for inputs and thanking you in advance for your kind help,
Vincent

Hi
I’d go with 1), if you discretize with hcurl then your natural trace is only the tangential component, so setting dirichlet bnd will automatically do what you want.
3) should also be possible in a second step to strongly couple the bc equation with the volume one.

You could also impose the equation as a penalty-robin type bc.

best Christopher

Hello Christopher,

For Hcurl, there is a catch (that I guess, as a FEM expert, will delight you):
The Dirichlet approach proposed in (1) is only working if the electric field is expressed using continuous elements and thus the normal component of the E fields behaves well. However, when using a basis function which defines a discontinuous En like Hcurl, the gradient is not well-defined and the whole technique fails.

For this reason of continuity, I implemented a Maxwell weak form on H1 (which can be regularized like in ERMES - CINME) which seems to behave quite well in my case (even without regularization).

In summary, if I want to implement (1), (2) or (3), I think that I need to access the normal vector En = (En \vdot \hat{n})\hat{n}. Is this possible in NGSolve? How would you do it?

I guess imposing a penalty-robin type BC would also do the work. Do you have an example of such a procedure with NGSolve/documentation?

Thanking you again for your help,
Vincent

Something like

n = specialcf.normal(mesh.dim)
a += penalty_factor * Cross(n,u) * Cross(n,v) * ds(boundary)
f += penalty_factor * Cross(n, u_dirichlet) * Cross(n, v) * ds(boundary)

should do the trick for penalty_factor large enough.