Skip to content

Instantly share code, notes, and snippets.

Created June 28, 2018 10:34
  • Star 4 You must be signed in to star a gist
  • Fork 3 You must be signed in to fork a gist
Star You must be signed in to star a gist
What would you like to do?
Terribly hacky MVC workaround for pygmsh-meshio-dolfin boundary tagging
from pygmsh.built_in.geometry import Geometry
from pygmsh import generate_mesh
import meshio
import dolfin
geom = Geometry()
lc = .1
p0 = geom.add_point([0, 0, 0], lcar=lc)
p1 = geom.add_point([1, 0, 0], lcar=lc)
p2 = geom.add_point([1, 1, 0], lcar=lc)
p3 = geom.add_point([0, 1, 0], lcar=lc)
l0 = geom.add_line(p0, p1)
l1 = geom.add_line(p1, p2)
l2 = geom.add_line(p2, p3)
l3 = geom.add_line(p3, p0)
ll = geom.add_line_loop(lines=[l0, l1, l2, l3])
ps = geom.add_plane_surface(ll)
# Tag line and surface
geom.add_physical_line(lines=l3, label=444)
geom.add_physical_surface(surfaces=ps, label=456)
points, cells, point_data, cell_data, field_data = generate_mesh(geom)
print("Writing 2d mesh for dolfin Mesh")
meshio.write("mesh_2d.xdmf", meshio.Mesh(
cells={"triangle": cells["triangle"]}))
print("Writing 1d line data for dolfin MVC")
meshio.write("mvc_1d.xdmf", meshio.Mesh(
cells={"line": cells["line"]},
cell_data={"line": {"name_to_read": cell_data["line"]["gmsh:physical"]}}
mesh_2d = dolfin.Mesh()
with dolfin.XDMFFile("mesh_2d.xdmf") as infile:
print("Reading 2d mesh into dolfin")
mvc = dolfin.MeshValueCollection("size_t", mesh_2d, 1)
with dolfin.XDMFFile("mvc_1d.xdmf") as infile:
print("Reading 1d line data into dolfin mvc"), "name_to_read")
print("Constructing MeshFunction from MeshValueCollection")
mf = dolfin.cpp.mesh.MeshFunctionSizet(mesh_2d, mvc)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment