# Integral on interface between two domains

Hello,

I have a domain that is split in two (using the code by Guosheng:
interface integrals - Kunena)

In fluid 1 I define a FacetFESpace Q1 and in fluid 2 I define a FacetFESpace Q2.

Let p1 and q1 be test and trial functions in Q1.
Let p2 and q2 be test and trial functions in Q2.

On the interface between the fluids I compute the integral

100p1q1 + 100p2q2.

I then print the matrix and I see that I get different matrices depending on how I define Q2. If I define

Q2 = FacetFESpace(mesh, order=1)

I will get, for example, that A(54,54) = 16.6667, but if I define

Q2 = FacetFESpace(mesh, order=1, definedon=mesh.Materials(“fluid2”))

I will get that A(54,54) = 0

I would prefer to use the second definition of Q2, but this gives the wrong matrix. Any suggestions how to fix this?
I have attached the code with Q2 defined in lines 52 and 53.

Thanks,
Sander

Attachment: question.py

Hi Sander,

using HDG plus additionally dgjumps looks like an expensive method - do you really want that ?
I can help with the HDG, but will not spend time on the dgjumps.

A note on the side: You can directly iterate over the mesh topology to mark edges etc.
There is a quite new section in the i-tutorials:
https://ngsolve.org/docu/nightly/i-tutorials/unit-1.8-meshtopology/meshtopology.html

``````for el in mesh.Elements(BND):
if el.mat == "gammaSD":
for edge in el.edges:
interface_indicator[edge.nr] = True``````

Joachim

Hi Joachim,

Thanks for the suggestion regarding the new i-tutorials. I will update my codes.

Regarding the dgjumps, the reason I’m using dgjumps is because of the Stokes-Darcy coupled code by Guosheng. I have attached a slightly modified version. Here `dgjumps=True` is added to the FESpace. Setting `dgjumps=False` leeds to errors of the form

catch in AssembleBilinearform 2: SparseMatrixTM::AddElementMatrix: illegal dnums
Traceback (most recent call last):
File “”, line 204, in
File “”, line 191, in SolveBVP
RuntimeError: SparseMatrixTM::AddElementMatrix: illegal dnumsin Assemble BilinearForm ‘biform_from_py’

Do you know why this is? Since this is HDG I would have expected dgjumps not to be necessary.

Regards,
Sander

Attachment: coupled_code.py

Hi Sander,

I think what you want is this:

``````interfaceTerm = 0.5*(1+4*pi*pi)*vbar.Trace()*ubar.Trace()
a += SymbolicBFI (interfaceTerm, definedon=mesh.Boundaries("gammaSD"))``````

The trace of the vector-facet variable can be evaluated on the boundary. It lives already in the tangential space, so you don’t need the projection. You have to provide the Trace() operator.

You don’t need the dgjumps, which would lead to unnecessarily many matrix entries.

Best, Joachim

Hi Joachim,

The .Trace() was indeed the missing piece. It also answered my original question. I attach the now working version with the updates you suggested.

Thanks!
Sander

Attachment: question_2018-08-31.py

overseeing the Trace() operator WAS a common pitfall.
I just added a test whether the trial and test-functions match with the form.
Joachim