Generating mesh with svmtk

mesh
,
1

I have trouble generating a subdomain mesh as outline in an svmtk tutorial. I have been following the svmtk tutorials here “Mathematical Modeling of the Human Brain From Magnetic Resonance Images to Finite Element Simulation“ (Mathematical Modeling of the Human Brain: From Magnetic Resonance Images to Finite Element Simulation | SpringerLink) and got the tutorial files from zenodo ( Software for Mathematical modeling of the human brain - from magnetic resonance images to finite element simulation ). When I run “two-domain-tagged.py“ from Chapter 4, I do not get the same mesh “ernie-gw.mesh“ as provided. Instead of having two domains (green and orange) it is a single domain (orange) (see attached screenshots). I installed svmtk using conda and as outlined here: GitHub - SVMTK/SVMTK .

# two-domain-tagged.py Chapter 4
# %%
import SVMTK as svmtk 

def create_gw_mesh(pial_stl, white_stl, output):
    # Load the surfaces into SVM-Tk and combine in list
    pial  = svmtk.Surface(pial_stl)
    white = svmtk.Surface(white_stl)
    surfaces = [pial, white]
    
    # Create a map for the subdomains with tags
    # 1 for inside the first and outside the second ("10")
    # 2 for inside the first and inside the second ("11")
    smap = svmtk.SubdomainMap()
    smap.add("10", 1)
    smap.add("11", 2)

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

create_gw_mesh("lh.pial.stl", "lh.white.stl", "ernie-gw.mesh")

# %%

code_rerun
code_rerun1920×1171 109 KB

provided_by_zenodo
provided_by_zenodo1920×1175 117 KB

2

I would suggest converting the .mesh file with msh into xdmf, and try to visualize it in Paraview.
As far as I am aware, GMSH doesn’t support all .mesh files.

1 Like
3

The CGAL/SVMTK mesh format is not compatible with GMSH. If you want to visualize it with gmsh,
you should convert it by using e.g. meshio (v5.3.5):

meshio convert ernie-gw.mesh ernie-gw.msh --output-format gmsh22

If you want to visualize it in Paraview, I recommend using the conversion in the book.

1 Like
4

Great that worked well! I am now trying to import into fenicsx but I noticed that the bounding surface aren’t properly associated when I have multidomain geometries. How do I get the geometry alongside all the surfaces ready for fenicsx use?

(base) chris@gondor:/mnt/shared/chris/Projects/2024_PhysicsAug/code/svmtk_meshing/stl_files$ meshio convert combined_final.mesh combined_final.msh --output-format gmsh22
Warning: Appending zeros to replace the missing physical tag data.
Warning: Appending zeros to replace the missing geometrical tag data.
5

I would have a look at:
mri2fem-ii-chapter-3-code/meshing at main · mhornkjol/mri2fem-ii-chapter-3-code · GitHub (Create and convert mesh) and mri2fem-ii-chapter-3-code/dolfinx_implementation/stokes_solver.py at main · mhornkjol/mri2fem-ii-chapter-3-code · GitHub
These scripts are from the upcoming book mri2fem -part 2

6

Thanks so much, looking forward to the new book!

Just FYI, I am getting couple syntax errors when trying to run create_mesh.py:

line 75
    aqueduct.clip(clp1, invert=True ))
                                     ^
SyntaxError: unmatched ')'

line 41, in <module>
    white   = svm.Surface(Z.white) 
NameError: name 'Z' is not defined
7

Thanks for the help. Here is the pipeline I used for others to follow

  1. Generate mesh from stl files using svmtk
  2. convert .mesh file to xdmf
  3. Use xdmf for fenics simulation

Step2: Conversion

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)


if __name__ == "__main__":
    parser = argparse.ArgumentParser(formatter_class=argparse.ArgumentDefaultsHelpFormatter)
    parser.add_argument(
        "-i",
        "--input",
        type=Path,
        dest="infile",
        help="Path to input mesh file",
        default="./stl_files/combined_final.mesh",
    )
    parser.add_argument(
        "-o",
        "--output",
        type=Path,
        dest="outfile",
        help="Path to output mesh file",
        default="./stl_files/combined_final.xdmf",
    )

    parsed_args = parser.parse_args()
    convert_mesh(parsed_args.infile, parsed_args.outfile)
8

Thanks for sharing, much appreciated ! Working on fixing the problems.

9

Thank you! I saw couple changes on the repo and tried running it, it all runs great (I think there was import argparse missing in one). The meshes come out well but I have trouble getting well formed boundaries, with boundaries not aligning with the volumes