Extract marked boundary and save as xdmf file for another simulation

mesh
1

Hi,

I have a 3D mesh that consists of mesh.xdmf (volume, tetra) and mf.xdmf (surface, triangle) where mf.xdmf has boundary markers.

What I want to do is to extract the boundary which is marked in mf.xdmf and save as xdmf file for another simulation, but I’m a bit confused how to proceed.

The way I read the mesh is as follows and from there I would like to use mf to save the marked boundary

from fenics import *

mesh = Mesh()
mvc = MeshValueCollection("size_t", mesh, mesh.topology().dim())
with XDMFFile(MPI.comm_world, "mesh.xdmf") as infile:
   infile.read(mesh)
   infile.read(mvc, “name_to_read")

mvc = MeshValueCollection("size_t", mesh, mesh.topology().dim()-1)
with XDMFFile(MPI.comm_world, "mf.xdmf") as infile:
    infile.read(mvc, “name_to_read")

mf = cpp.mesh.MeshFunctionSizet(mesh, mvc)

I attached the screen shot of how it looks when you load mf.xdmf. The inlet boundary is marked as 1, the wall part is marked as 2, and the outlet part is marked as 3.

Mesh files can be accessed here mesh_coarse - Google Drive

Screen Shot 2021-09-25 at 13.57.22
Screen Shot 2021-09-25 at 13.57.222880×1800 483 KB

2

I am not sure what your issue with the current code is, can’t you simply do:

surface_mesh = Mesh()
with XDMFFile(MPI.comm_world, "mf.xdmf") as infile:
    infile.read(surface_mesh)
3

Hi @dokken thank you for the answer. Sorry it was not clear but what I want to do is to save only some parts of the boundary which is marked, NOT the whole surface mesh.
For example, in the case above, I want to save only the inlet part which is marked as 1.

I could do the same manipulation by, for example, as follows

bmesh = BoundaryMesh(mesh, "exterior")
facet_f = MeshFunction('size_t', bmesh, bmesh.topology().dim(), 0)
CompiledSubDomain('x[0] <= 0 + tol ) +tol’,tol=tol).mark(facet_f, 1) 
inlet = SubMesh(bmesh, facet_f, 1)

with XDMFFile("inlet.xdmf") as file:
    file.write(inlet)

However, it is not robust way of saving specific parts of boundary when the mesh becomes complicated, and thus I’m looking for the way to save specific boundary.

4

How about the following:

from dolfin import *

mesh = UnitSquareMesh(10, 10)

facet_f = MeshFunction("size_t", mesh, mesh.topology().dim() -1, 0)
CompiledSubDomain('x[0] <= 0 + tol',tol=1e-6).mark(facet_f, 1)
CompiledSubDomain('x[0] >= 1 - tol',tol=1e-6).mark(facet_f, 2)
with XDMFFile("mf.xdmf") as xdmf:
    xdmf.write(facet_f)


surface_mesh = Mesh()
mvc = MeshValueCollection("size_t", surface_mesh, surface_mesh.topology().dim())
with XDMFFile("mf.xdmf") as xdmf:
    xdmf.read(surface_mesh)
    xdmf.read(mvc)

mf = cpp.mesh.MeshFunctionSizet(surface_mesh, mvc)

mv = MeshView.create(mf, 1)
with XDMFFile("mv.xdmf") as xdmf:
    xdmf.write(mv)