Svmtk meshing, boundary tagging

mesh
,
1

I am able to generate meshes with svmtk, but I am not understanding how I can define boundary surface tags? I.e. if I have a mesh with subdomains, how do I set the boundary tag to 1 for subdomain A that boarders subdomain B, and then set boundary tag to 2 for subdomain A that boarders subdomain C? I would preferably be able to do this within svmtk or must I use other tools?

This is the code I am using. Thanks for the help!

# %%
import SVMTK as svmtk 
import meshio

# %%
def create_mesh(stl_files, output, resolution=16):
    # Load the surfaces into SVM-Tk and combine in list
    surfaces = []
    for stl_file in stl_files:
        surfaces.append(svmtk.Surface(stl_file))
    
    # Create a map for the subdomains with tags
    # Bitstring order: [body, bone, lung, tumor]
    # 0 = outside surface, 1 = inside surface
    # Order matters: tumor patterns come FIRST to override other tissues
    smap = svmtk.SubdomainMap()
    
    # Tag 4: Tumor tissue (overrides all other tissues where present)
    smap.add("1111", 4)  # Inside body, bone, lung, and tumor
    smap.add("1101", 4)  # Inside body, bone, and tumor (outside lung)
    smap.add("1011", 4)  # Inside body, lung, and tumor (outside bone)
    smap.add("1001", 4)  # Inside body and tumor (outside bone and lung)
    
    # Tag 2: Bone tissue (inside body and bone, outside tumor)
    smap.add("1100", 2)
    smap.add("1110", 2) # Inside body, bone, lung, outside tumor
    smap.add("0100", 2) # Outside body, lung, tumor
    
    # Tag 3: Lung tissue (inside body and lung, outside tumor)
    smap.add("1010", 3)
    smap.add("0010", 3) # Outside body, bone, tumor
    
    # Tag 1: Body tissue (inside body, outside bone, lung, tumor)
    smap.add("1000", 1)

    # Create a tagged domain from the list of surfaces
    # and the map
    domain = svmtk.Domain(surfaces, smap)
       
    # Create and save the volume mesh 
    domain.create_mesh(resolution)
    domain.save(output)

    mesh = meshio.read(output)      # or "input.msh" for Gmsh
    output_xdmf = output.replace(".mesh", ".xdmf")

    meshio.write(output_xdmf, mesh)

if __name__ == "__main__":
    stl_files = ["body.remesh.smooth.stl","bone.remesh.smooth.stl",
                    "lung.remesh.smooth.stl","tumor.remesh.smooth.stl"]
    stl_files = ["./stl_files/" + stl_file for stl_file in stl_files]
    create_mesh(stl_files, "./stl_files/combined_final.mesh", resolution=16)
2

The meshio conversion from .mesh to .xdmf does not handle different elements, like triangles and tetrahedrons. So you need to extract the subdomain and boundaries separately.

mesh = meshio.read(output) 

cell_mesh = meshio.Mesh( 
        points=mesh.points,
        cells={"tetra": mesh.cells_dict["tetra"]},
        cell_data={"subdomains": [mesh.cell_data_dict["medit:ref"]["tetra"]]},
    )
            
meshio.write("subdomains.xdmf", cell_mesh)
facet_mesh = meshio.Mesh(
        points=mesh.points,
        cells={"triangle": mesh.cells_dict["triangle"]},
        cell_data={"patches": [mesh.cell_data_dict["medit:ref"]["triangle"]]},
    )
    
meshio.write("boundaries.xdmf", facet_mesh)
3

Thanks for the help! Sorry I posted a prior version and yes I had changed that conversion to your code above. But even if I do it like that, how do I more finely control the boundary labels? I.e. if want to change the label of Subdomain A from boundary value 1 to boundary value 2?

4

I tried a simplified mock with spheres, I get nicely meshed subdomains but subdomain boundary tags are all over the place. This seems also to happen when I try running online examples (subdomains mesh well but not the respective boundaries): mri2fem-ii-chapter-3-code/meshing at main · mhornkjol/mri2fem-ii-chapter-3-code · GitHub . Is there anything I am missing, I am wondering if there is any svmtk version inconsistency (currently using latest svmtk version)? Thanks for any insight!

Weird subdomain boundaries:

image
image493×372 68 KB

Subdomains themselves look fine (one big sphere containing 3 smaller spheres):

image
image492×369 43.9 KB

Stl file generation:

# %%
import gmsh

# %%
spheres = [((0,0,0), 100.0), # center, radius, sphere 1
           ((60,0,0), 20.0), # center, radius, sphere 2
           ((0,60,0), 20.0), # center, radius, sphere 3
           ((-60,0,0), 20.0)] # center, radius, sphere 4

gmsh.initialize()
for i, (c, r) in enumerate(spheres, 1):
    gmsh.model.add(f"sphere_{i}")
    gmsh.model.occ.addSphere(*c, r)
    gmsh.model.occ.synchronize()
    # Surface mesh:
    gmsh.model.mesh.generate(2)
    gmsh.write(f"sphere_{i}.stl")
    # Or volume mesh (tetrahedra) instead:
    # gmsh.model.mesh.generate(3)
    # gmsh.write(f"sphere_{i}.msh")
    gmsh.model.remove()  # start fresh for next file
gmsh.finalize()

# %%

Mesh generation:

# %%
import SVMTK as svmtk 
import meshio
from convert_mesh import convert_mesh
from pathlib import Path

# %%
def create_mesh(stl_files, output, resolution=16):
    # Load the surfaces into SVM-Tk and combine in list
    surfaces = []
    for stl_file in stl_files:
        surface = svmtk.Surface(stl_file)
        L = 1.0
        n = 3
        do_not_move_boundary_edges = False

        surface.isotropic_remeshing(L, n,do_not_move_boundary_edges)
        surfaces.append(surface)
        #surfaces.append(svmtk.Surface(stl_file))
    
    # Create a map for the subdomains with tags
    # Bitstring order: [body, bone, lung, tumor]
    # 0 = outside surface, 1 = inside surface
    # Order matters: tumor patterns come FIRST to override other tissues
    smap = svmtk.SubdomainMap()
    
    # Tag 4: Tumor tissue (overrides all other tissues where present)
    smap.add("1111", 4)  # Inside body, bone, lung, and tumor
    smap.add("1101", 4)  # Inside body, bone, and tumor (outside lung)
    smap.add("1011", 4)  # Inside body, lung, and tumor (outside bone)
    smap.add("1001", 4)  # Inside body and tumor (outside bone and lung)
    
    # Tag 2: Bone tissue (inside body and bone, outside tumor)
    smap.add("1100", 2)
    smap.add("1110", 2) # Inside body, bone, lung, outside tumor
    smap.add("0100", 2) # Outside body, lung, tumor
    
    # Tag 3: Lung tissue (inside body and lung, outside tumor)
    smap.add("1010", 3)
    smap.add("0010", 3) # Outside body, bone, tumor
    
    # Tag 1: Body tissue (inside body, outside bone, lung, tumor)
    smap.add("1000", 1)

    # Create a tagged domain from the list of surfaces
    # and the map
    domain = svmtk.Domain(surfaces, smap)
       
    # Create and save the volume mesh 
    domain.create_mesh(resolution)
    domain.save(output)

    convert_mesh(Path(output), Path(output.replace(".mesh", ".xdmf")))

Mesh conversion to xdmf:

import meshio
import argparse
from pathlib import Path


def convert_mesh(infile: Path, outfile: Path):
    msh = meshio.read(infile)

    # try:
    tmp_dir = Path(".tmp")
    tmp_dir.mkdir(exist_ok=True)

    cell_mesh = meshio.Mesh(
        points=msh.points,
        cells={"tetra": msh.cells_dict["tetra"]},
        cell_data={"subdomains": [msh.cell_data_dict["medit:ref"]["tetra"]]},
    )

    meshio.write(outfile.with_suffix(".xdmf"), cell_mesh)

    facet_mesh = meshio.Mesh(
        points=msh.points,
        cells={"triangle": msh.cells_dict["triangle"]},
        cell_data={"patches": [msh.cell_data_dict["medit:ref"]["triangle"]]},
    )
    facet_file = (
        outfile.with_stem(outfile.stem + "_facets").with_suffix(".xdmf").absolute().as_posix()
    )
    meshio.write(facet_file, facet_mesh)

All files: