Define Subdomains created from .geo into XDMF file into dolfin legacy

General
4

This not current, as the last argument should be 1.
Your new mesh runs fine on my system with:

import numpy as np
from IPython import embed
import dolfin
import numpy
import meshio


def create_mesh(mesh, cell_type, prune_z=False):
    cells = mesh.get_cells_type(cell_type)
    cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
    points = mesh.points[:, :2] if prune_z else mesh.points
    out_mesh = meshio.Mesh(points=points, cells={cell_type: cells}, cell_data={
                           "name_to_read": [cell_data]})
    return out_mesh


mesh_from_file = meshio.read("test_mesh.msh")

line_mesh = create_mesh(mesh_from_file, "line", prune_z=True)
meshio.write("mf.xdmf", line_mesh)
triangle_mesh = create_mesh(mesh_from_file, "triangle", prune_z=True)
meshio.write("mesh.xdmf", triangle_mesh)

assert np.allclose(line_mesh.points, triangle_mesh.points)
vertices_in_triangles = np.unique(
    triangle_mesh.cells_dict["triangle"].reshape(-1))

mesh = dolfin.Mesh()
with dolfin.XDMFFile("mesh.xdmf") as infile:
    infile.read(mesh)
mvc = dolfin.MeshValueCollection("size_t", mesh, 1)
with dolfin.XDMFFile("mf.xdmf") as infile:
    infile.read(mvc, "name_to_read")
mf = dolfin.cpp.mesh.MeshFunctionSizet(mesh, mvc)