importScript.py

import meshio
import numpy as np

msh = meshio.read("CTCSeparationExample2d03.msh")

print(msh.cells)
# cells={"tetra": msh.cells["tetra"]}))
meshio.write("mesh_CTCSep1.xdmf", meshio.Mesh(points=msh.points,
cells = [
        meshio.CellBlock(
            type="tetra",
            data=  np.concatenate([arr for (t, arr) in msh.cells if t == "tetra"]),
        )
    ])
)


meshio.write("mf_CTCSep1.xdmf", meshio.Mesh(points=msh.points, 
cells = [
        meshio.CellBlock(
            type="triangle",
            data=  np.concatenate([arr for (t, arr) in msh.cells if t == "triangle"]),
        )
    ],
cell_data = {
        "name_to_read": [
            np.concatenate(
                [
                    msh.cell_data["gmsh:physical"][i]
                    for i, cellBlock in enumerate(msh.cells)
                    if cellBlock.type == "triangle"
                ]
            )
        ]
    })
)

from dolfin import *
mesh = Mesh()
with XDMFFile("mesh_CTCSep1.xdmf") as infile:
    infile.read(mesh)
mvc = MeshValueCollection("size_t", mesh, 2)
with XDMFFile("mf_CTCSep1.xdmf") as infile: infile.read(mvc, "name_to_read")
mf = cpp.mesh.MeshFunctionSizet(mesh, mvc)
# Insert physical label number to test subdomain retrival
It_facet = SubsetIterator(mf,14)