Discontinuity at interface between two subdomains

1

Hi all,

I’m trying to simulate discontinuities at the interface between two subdomains but using CG elements everywhere else in the bulk.
I reckon I need some way to create a FunctionSpace based on a MixedElement with DG elements restricted on the subdomain of the interface.

I couldn’t find a way to do it in the docs though. Does that ring a bell ?

from fenics import *


class Left(SubDomain):
    def inside(self, x, on_boundary):
        return near(x[0], 0.0)


class Middle(SubDomain):
    def inside(self, x, on_boundary):
        return near(x[0], 0.5)


class Right(SubDomain):
    def inside(self, x, on_boundary):
        return near(x[0], 1.0)


mesh = UnitIntervalMesh(10)

boundaries = MeshFunction("size_t", mesh, 0)
boundaries.set_all(0)
left = Left()
left.mark(boundaries, 1)
middle = Middle()
middle.mark(boundaries, 2)
right = Right()
right.mark(boundaries, 3)

P1 = FiniteElement("CG", mesh.ufl_cell(), 1)
DG1 = FiniteElement("DG", mesh.ufl_cell(), 1)  # DG1 must be restricted to subdomain x=0.5
element = P1 + DG1
V = FunctionSpace(mesh, element)
u = Function(V)
v = TestFunction(V)

bcs = [DirichletBC(V, 5.0, boundaries, 1), DirichletBC(V, 0, boundaries, 3)]
dS = Measure('dS', domain=mesh, subdomain_data=boundaries)

F = dot(grad(u), grad(v))*dx
F += -2*jump(u)*dS(2)

solve(F == 0, u, bcs)

Thanks for the help!

2

You might want to have a look at the multiphenics library https://gitlab.com/multiphenics/multiphenics (disclaimer: I am the author of that), and at the related resources in section 5 of its readme, most notably the mixed dimensional branch by @cdaversin

2 Likes
3

Hi thanks for the reply! If I understood correctly by using @cdaversin 's branch there could be a way to achieve the restriction of an element to a subdomain ?

Thanks!

4

Hi, In my branch (cecile/mixed-dimensional) you can indeed create a lower-dimensional submesh for your interface. You can then define a FunctionSpace from this submesh, using the appropriate FiniteElement. And you can define the product of your two function spaces as :

V1 = FunctionSpace(mesh, P1)
V2 = FunctionSpace(submesh, DG1)
V = MixedFunctionSpace(V1,V2)

If you want to give it a try, I suggest using my Docker container including all the features and a 2D-1D demo:

docker pull ceciledc/fenics_mixed_dimensional
1 Like
5

Hi @RemDelaporteMathurin,

did you finally figure out how to solve your problem? I think I have trouble with a similiar task (Internal boundary condition with mixed_dimensional branch). I would love to see your solution!

6

Hi sorry for the late answer.

I did try to fiddle around with @cdaversin implementation but that wasn’t successful unfortunately. I managed to solve my issue with post-processing but applying a just on a boundary is definitely the proper way to do it.
I shall wait for an implementation in the master branch then.

7

Hi @cdaversin,
Thank you for sharing. Is the demo also available outside the Docker container? I would like to try it for simulating a discontinuous interface boundary.

8

Hi,
You can find some examples in the supplementary material we provide with our preprint presenting the mixed-dimensional framework.

1 Like
9

Hello,
I am having the same issue and I am trying to solve it with post-processing. Would you mind guiding me as to how you did it? Are you excluding some points of the boundary out by hand?
I’d really appreciate any advice.
Thanks!

10

Hi @kamilova-maths , I’m afraid I was never able to solve it :smiley: