Is it possible to build this kind of subdomain?

Greetings everyone.

I would like to build a subdomain in dolfin legacy.
This subdomain cannot be described easily with an algebraic definition, but I can have it from a physical group in the .msh file.
I achieved to open an XDMF file as a MeshValueCollection with correct tag in it, and I don’t know how to use it properly.

And I would like this subdomain to be defined on which a specific function has its gradient pointing outward the domain.
So for a given function u (let’s say piece-wise linear), condition is grad(u).n>0 with n a normal pointing outward the domain.
And n cannot be written by hand.

Could this hypothetic subdomain exist ?
If yes, would someone know the syntax for it (in c++)?

Here’s a link with a cpp file to fill, accompanied by .ufl, .geo, .msh and .xdmf files: mycore.core-cloud.net/index.php/s/HYulKlbDe88oleX

Thank you for your consideration.
Best regards.

There are plenty of posts on how to read MeshValueCollection objects, see for instance:How to read a physical object from a xdmf file - #7 by dokken
It is not very different in C++.

I guess n is defined by the surface of the domain, i.e. on the facets where the interface between marked subdomain and other domains are? Then you can use FacetNormal(mesh) + the dS measure, see: Integrating over an interior surface - #4 by MiroK

Do you want to solve a PDE where u is the unknown?
If so, what kind of PDE are you trying to solve?

Hi !
Thank you for the quick reply.

I achieved to build the MeshValueCollection but I don’t how to use it in the “inside” method of a subdomain, that’s the point.

My goal is to solve a contact problem in mechanic with time dependency.
And I would like to set a Dirichlet Boundary condition, at instant t+1, on boundary on which my contact condition was violated in term of flux at instant t (ie a body is unphysically sticked to a wall).

The geometry of the mesh isn’t well defined because it’s made of hundred of tiny elements in various orientations.

I know how to compute the flux, but the normal flux is much harder for me.

Mainly because I cannot manually write the coordinates of the normal vector, the geometry is not describable with a simple algebraic formula like “x² + y² = 1”.

There’s a whole example in the zip folder, with mesh, tag, u function, only the “inside” method lack : mycore.core-cloud.net/index.php/s/HYulKlbDe88oleX