Generating meshes with Wildmeshing

Generating meshes with Wildmeshing#

In this demo, we will show how to generate a mesh from STL/ply files based on MRI images. We use the files from the paper In-silico molecular enrichment and clearance of the human intracranial space

We start by importing all the packages required for this task.

from mpi4py import MPI
import wildmeshing
import dolfinx
from pathlib import Path
import json
import numpy as np
import pyvista
import ufl
import basix.ufl
import os

Next, we use the already downloaded stl files from this chapter. We will use WildMeshing to combine these. We create unique volume markers for the lateral ventricles (LV.ply) and the third and fourth ventricle (V34.ply). Next we embed the ventricular system within the brain parenchyma (parenchyma_incl_ventr.ply), before we add the dura interface (skull.ply), which is the interfacebetween the CSF and the skull.

We extract these files from the "WILDFENICS_DATA_PATH". If you have not set this variable, please change it to the appropriate location on your system.

folder = Path(os.environ["WILDFENICS_DATA_PATH"])
assert folder.exists(), "Could not find surface files"
tree = {
    "operation": "union",
    "left": str((folder / "skull.ply").absolute().as_posix()),
    "right": {
        "operation": "union",
        "left": str((folder / "parenchyma_incl_ventr.ply").absolute().as_posix()),
        "right": {
            "operation": "union",
            "left": str((folder / "LV.ply").absolute().as_posix()),
            "right": str((folder / "V34.ply").absolute().as_posix()),
        },
    },
}

Next, we can set up a 3D meshing instance and load in the instruction-set

tetra = wildmeshing.Tetrahedralizer(
    stop_quality=10, max_its=30, edge_length_r=0.03, epsilon=0.00225
)
tetra.load_csg_tree(json.dumps(tree))
TBB threads 4
o-[I/O         ] Loading file /root/data/mesh/standard/surfaces/skull.ply...
                 (FP64) nb_v:138154 nb_e:0 nb_f:276304 nb_b:0 tri:1 dim:3
                 Attributes on vertices: point[3]
                 Loading file /root/data/mesh/standard/surfaces/parenchyma_incl_ventr.ply...
                 (FP64) nb_v:317019 nb_e:0 nb_f:635594 nb_b:0 tri:1 dim:3
                 Attributes on vertices: point[3]
                 Loading file /root/data/mesh/standard/surfaces/LV.ply...
                 (FP64) nb_v:17572 nb_e:0 nb_f:35144 nb_b:0 tri:1 dim:3
                 Attributes on vertices: point[3]
                 Loading file /root/data/mesh/standard/surfaces/V34.ply...
                 (FP64) nb_v:4608 nb_e:0 nb_f:9212 nb_b:0 tri:1 dim:3
                 Attributes on vertices: point[3]
bbox_diag_length = 0.266221
ideal_edge_length = 0.00798662
stage = 2
eps_input = 0.000598997
eps = 0.000331599
eps_simplification = 0.000265279
eps_coplanar = 2.66221e-07
dd = 0.000399331
dd_simplification = 0.000319465

Then, we can generate the mesh

# Create mesh
tetra.tetrahedralize()
[2025-10-06 15:01:19.262] [float-tetwild] [info] remove duplicates: 
[2025-10-06 15:01:19.262] [float-tetwild] [info] #v: 477353 -> 477353
[2025-10-06 15:01:19.262] [float-tetwild] [info] #f: 956254 -> 956254
collapsing 21.2289
swapping 0.251291
[2025-10-06 15:01:40.979] [float-tetwild] [info] remove duplicates: 
[2025-10-06 15:01:40.980] [float-tetwild] [info] #v: 49442 -> 49442
[2025-10-06 15:01:40.980] [float-tetwild] [info] #f: 100369 -> 100369
[2025-10-06 15:01:40.981] [float-tetwild] [info] #v = 49442
[2025-10-06 15:01:40.981] [float-tetwild] [info] #f = 100256
#boundary_e1 = 309
#boundary_e2 = 120
[2025-10-06 15:01:41.338] [float-tetwild] [info] preprocessing 22.2886s
[2025-10-06 15:01:41.338] [float-tetwild] [info] 
[2025-10-06 15:01:41.569] [float-tetwild] [info] #v = 50230
[2025-10-06 15:01:41.570] [float-tetwild] [info] #t = 333917
[2025-10-06 15:01:41.570] [float-tetwild] [info] tetrahedralizing 0.232722s
[2025-10-06 15:01:41.570] [float-tetwild] [info] 
[2025-10-06 15:01:41.571] [float-tetwild] [info] triangle insertion start, #f = 100256, #v = 50230, #t = 333917
[2025-10-06 15:01:41.707] [float-tetwild] [info] matched #f = 74494, uninserted #f = 25762
[2025-10-06 15:01:50.708] [float-tetwild] [info] insert_one_triangle * n done, #v = 154009, #t = 869503
[2025-10-06 15:01:50.708] [float-tetwild] [info] uninserted #f = 0/25762
[2025-10-06 15:01:50.708] [float-tetwild] [info] total timing: 0s
[2025-10-06 15:01:51.157] [float-tetwild] [info] pair_track_surface_fs done
#boundary_e1 = 309
#boundary_e2 = 120
[2025-10-06 15:01:51.211] [float-tetwild] [info] find_boundary_edges done
[2025-10-06 15:01:51.266] [float-tetwild] [info] time1 = 0.055371
[2025-10-06 15:01:51.280] [float-tetwild] [info] uninsert boundary #e = 0/429
[2025-10-06 15:01:51.280] [float-tetwild] [info] time2 = 7.1e-05
[2025-10-06 15:01:51.280] [float-tetwild] [info] time3 = 0.012654
[2025-10-06 15:01:51.280] [float-tetwild] [info] time4 = 0.000317
[2025-10-06 15:01:51.280] [float-tetwild] [info] time5 = 9.5e-05
[2025-10-06 15:01:51.280] [float-tetwild] [info] time6 = 4e-06
[2025-10-06 15:01:51.288] [float-tetwild] [info] uninserted #f = 0/25762
known_surface_fs.size = 0
known_not_surface_fs.size = 0
[2025-10-06 15:01:53.999] [float-tetwild] [info] mark_surface_fs done
[2025-10-06 15:01:54.022] [float-tetwild] [info] #b_edge1 = 429, #b_edges2 = 0
[2025-10-06 15:01:54.056] [float-tetwild] [info] cutting 12.4843s
[2025-10-06 15:01:54.056] [float-tetwild] [info] 
initializing...
edge collapsing...
fixed 0 tangled element
success(env) = 85554
success = 100381(816258)
success(env) = 3278
success = 3524(534362)
success(env) = 288
success = 296(223781)
success(env) = 42
success = 43(29839)
success(env) = 11
success = 12(4894)
success(env) = 5
success = 5(1354)
success(env) = 3
success = 3(339)
success(env) = 0
success = 0(460)
edge collapsing done!
time = 14.5184s
#v = 44070
#t = 288346
max_energy = 4037.86
avg_energy = 7.75031
//////////////// pass 0 ////////////////
edge splitting...
fixed 0 tangled element
success = 43049(43049)
edge splitting done!
time = 0.327991s
#v = 87119
#t = 503700
max_energy = 4037.86
avg_energy = 7.66828
edge collapsing...
fixed 0 tangled element
success(e
nv) = 202
success = 29336(623550)
success(env) = 66
success = 1395(158154)
success(env) = 14
success = 202(46243)
success(env) = 4
success = 39(7659)
success(env) = 0
success = 11(1559)
success(env) = 0
success = 3(319)
success(env) = 0
success = 1(71)
success(env) = 0
success = 0(12)
edge collapsing done!
time = 5.25705s
#v = 56132
#t = 349178
max_energy = 3666.9
avg_energy = 6.66422
edge swapping...
fixed 0 tangled element
success3 = 33813
success4 = 39472
success5 = 3835
success = 77120(343547)
edge swapping done!
time = 2.49592s
#v = 56132
#t = 319200
max_energy = 1442.78
avg_energy = 5.5814
vertex smoothing...
success = 30140(52237)
vertex smoothing done!
time = 1.99946s
#v = 56132
#t = 319200
max_energy = 974.659
avg_energy = 5.22862
//////////////// pass 1 ////////////////
edge splitting...
fixed 0 tangled element
success = 5137(5137)
edge splitting done!
time = 0.193294s
#v = 61269
#t = 344334
max_energy = 974.659
avg_energy = 5.17587
edge collapsing...
fixed 0 tangled element
success(env) = 2732
success = 7002(519745)
success(env) = 166
success = 321(264748)
success(env) = 14
success = 26(25401)
success(env) = 1
success = 6(2252)
success(env) = 2
success = 4(229)
success(env) = 0
success = 0(200)
edge collapsing done!
time = 4.9521s
#v = 53910
#t = 306456
max_energy = 974.659
avg_energy = 5.12646
edge swapping...
fixed 0 tangled element
success3 = 4084
success4 = 14444
success5 = 1232
success = 19760(201287)
edge swapping done!
time = 1.38207s
#v = 53910
#t = 303604
max_energy = 974.659
avg_energy = 5.01462
vertex smoothing...
success = 25337(50014)
vertex smoothing done!
time = 1.7666s
#v = 53910
#t = 303604
max_energy = 974.658
avg_energy = 4.90411
updating sclaing field ...
filter_energy = 9.74658
is_hit_min_edge_length = 0
enlarge envelope, eps = 0.000368443
//////////////// pass 2 ////////////////
edge splitting...
fixed 0 tangled element
success = 58876(58876)
edge splitting done!
time = 0.434348s
#v = 112786
#t = 617593
max_energy = 974.658
avg_energy = 5.04069
edge collapsing...
fixed 0 tangled element
success(env) = 2359
success = 41996(544611)
success(env) = 314
success = 3085(301893)
success(env) = 31
success = 350(90203)
success(env) = 2
success = 79(12306)
success(env) = 0
success = 12(2437)
success(env) = 0
success = 0(377)
edge collapsing done!
time = 6.57755s
#v = 67264
#t = 384510
max_energy = 974.658
avg_energy = 4.49435
edge swapping...
fixed 0 tangled element
success3 = 5460
success4 = 22153
success5 = 2109
success = 29722(282602)
edge swapping done!
time = 2.53013s
#v = 67264
#t = 381159
max_energy = 974.658
avg_energy = 4.33056
vertex smoothing...
success = 41655(63358)
vertex smoothing done!
time = 2.04644s
#v = 67264
#t = 381159
max_energy = 974.658
avg_energy = 4.14123
//////////////// pass 3 ////////////////
edge splitting...
fixed 0 tangled element
success = 12718(12718)
edge splitting done!
time = 0.258573s
#v = 79982
#t = 442925
max_energy = 974.658
avg_energy = 4.166
edge collapsing...
fixed 0 tangled element
success(env) = 1277
success = 12617(440605)
success(env) = 86
success = 559(184020)
success(env) = 5
success = 37(17567)
success(env) = 2
success = 8(1318)
success(env) = 0
success = 2(293)
success(env) = 0
success = 1(31)
success(env) = 0
success = 0(23)
edge collapsing done!
time = 4.02853s
#v = 66758
#t = 378337
max_energy = 974.658
avg_energy = 4.13202
edge swapping...
fixed 0 tangled element
success3 = 1658
success4 = 11359
success5 = 776
success = 13793(247440)
edge swapping done!
time = 1.71702s
#v = 66758
#t = 377455
max_energy = 974.658
avg_energy = 4.09165
vertex smoothing...
success = 38829(62849)
vertex smoothing done!
time = 1.88821s
#v = 66758
#t = 377455
max_energy = 974.658
avg_energy = 4.01561
updating sclaing field ...
filter_energy = 9.74658
is_hit_min_edge_length = 0
//////////////// pass 4 ////////////////
edge splitting...
fixed 0 tangled element
success = 8191(8191)
edge splitting done!
time = 0.23631s
#v = 74949
#t = 419633
max_energy = 974.658
avg_energy = 4.07592
edge collapsing...
fixed 0 tangled element
success(env) = 1072
success = 6593(648827)
success(env) = 64
success = 434(136846)
success(env) = 5
success = 49(13041)
success(env) = 0
success = 4(1416)
success(env) = 0
success = 0(132)
edge collapsing done!
time = 4.42754s
#v = 67869
#t = 383884
max_energy = 974.658
avg_energy = 4.01525
edge swapping...
fixed 0 tangled element
success3 = 873
success4 = 6871
success5 = 434
success = 8178(241018)
edge swapping done!
time = 1.64698s
#v = 67869
#t = 383445
max_energy = 974.658
avg_energy = 3.99547
vertex smoothing...
success = 37683(63961)
vertex smoothing done!
time = 1.83436s
#v = 67869
#t = 383445
max_energy = 974.658
avg_energy = 3.95018
//////////////// pass 5 ////////////////
edge splitting...
fixed 0 tangled element
success = 3192(3192)
edge splitting done!
time = 0.222488s
#v = 71061
#t = 398858
max_energy = 974.658
avg_energy = 3.96589
edge collapsing...
fixed 0 tangled element
success(env) = 546
success = 3501(631450)
success(env) = 30
success = 142(76010)
success(env) = 1
success = 8(4538)
success(env) = 0
success = 1(433)
success(env) = 0
success = 0(43)
edge collapsing done!
time = 3.93022s
#v = 67409
#t = 380959
max_energy = 974.658
avg_energy = 3.95364
edge swapping...
fixed 0 tangled element
success3 = 372
success4 = 3706
success5 = 201
success = 4279(231392)
edge swapping done!
time = 1.42928s
#v = 67409
#t = 380788
max_energy = 974.658
avg_energy = 3.94678
vertex smoothing...
success = 35276(63494)
vertex smoothing done!
time = 1.73794s
#v = 67409
#t = 380788
max_energy = 974.658
avg_energy = 3.92505
updating sclaing field ...
filter_energy = 9.74658
is_hit_min_edge_length = 0
//////////////// pass 6 ////////////////
edge splitting...
fixed 0 tangled element
success = 717(717)
edge splitting done!
time = 0.212236s
#v = 68126
#t = 384145
max_energy = 974.658
avg_energy = 3.92864
edge collapsing...
fixed 0 tangled element
success(env) = 163
success = 1002(721458)
success(env) = 13
success = 42(40788)
success(env) = 2
success = 5(2050)
success(env) = 0
success = 1(275)
success(env) = 0
success = 0(105)
edge collapsing done!
time = 4.15419s
#v = 67076
#t = 378963
max_energy = 974.658
avg_energy = 3.9266
edge swapping...
fixed 0 tangled element
success3 = 156
success4 = 1983
success5 = 68
success = 2207(225858)
edge swapping done!
time = 1.35893s
#v = 67076
#t = 378875
max_energy = 974.658
avg_energy = 3.92382
vertex smoothing...
success = 33287(63169)
vertex smoothing done!
time = 1.67512s
#v = 67076
#t = 378875
max_energy = 974.658
avg_energy = 3.91376
//////////////// pass 7 ////////////////
edge splitting...
fixed 0 tangled element
success = 598(598)
edge splitting done!
time = 0.206554s
#v = 67674
#t = 381590
max_energy = 974.658
avg_energy = 3.91532
edge collapsing...
fixed 0 tangled element
success(env) = 69
success = 725(718328)
success(env) = 5
success = 16(20544)
success(env) = 0
success = 1(959)
success(env) = 0
success = 0(108)
edge collapsing done!
time = 3.73483s
#v = 66932
#t = 378047
max_energy = 974.658
avg_energy = 3.91486
edge swapping...
fixed 0 tangled element
success3 = 95
success4 = 1079
success5 = 54
success = 1228(223443)
edge swapping done!
time = 1.28825s
#v = 66932
#t = 378006
max_energy = 974.658
avg_energy = 3.91382
vertex smoothing...
success = 31883(63026)
vertex smoothing done!
time = 1.60436s
#v = 66932
#t = 378006
max_energy = 974.658
avg_energy = 3.90868
updating sclaing field ...
filter_energy = 9.74658
is_hit_min_edge_length = 1
//////////////// pass 8 ////////////////
edge splitting...
fixed 0 tangled element
success = 547(547)
edge splitting done!
time = 0.206473s
#v = 67479
#t = 380484
max_energy = 974.658
avg_energy = 3.91022
edge collapsing...
fixed 0 tangled element
success(env) = 36
success = 636(734638)
success(env) = 1
success = 19(15284)
success(env) = 0
success = 0(892)
edge collapsing done!
time = 3.70153s
#v = 66824
#t = 377387
max_energy = 974.658
avg_energy = 3.90965
edge swapping...
fixed 0 tangled element
success3 = 58
success4 = 646
success5 = 19
success = 723(221942)
edge swapping done!
time = 1.29271s
#v = 66824
#t = 377348
max_energy = 974.658
avg_energy = 3.90906
vertex smoothing...
success = 30746(62904)
vertex smoothing done!
time = 1.54847s
#v = 66824
#t = 377348
max_energy = 974.658
avg_energy = 3.90585
//////////////// pass 9 ////////////////
edge splitting...
fixed 0 tangled element
success = 523(523)
edge splitting done!
time = 0.205002s
#v = 67347
#t = 379694
max_energy = 974.658
avg_energy = 3.90709
edge collapsing...
fixed 0 tangled element
success(env) = 15
success = 560(733479)
success(env) = 0
success = 5(8561)
success(env) = 0
success = 0(294)
edge collapsing done!
time = 3.56293s
#v = 66782
#t = 377103
max_energy = 974.658
avg_energy = 3.90645
edge swapping...
fixed 0 tangled element
success3 = 23
success4 = 418
success5 = 18
success = 459(221210)
edge swapping done!
time = 1.29352s
#v = 66782
#t = 377098
max_energy = 974.658
avg_energy = 3.906
vertex smoothing...
success = 29980(62870)
vertex smoothing done!
time = 1.51964s
#v = 66782
#t = 377098
max_energy = 974.658
avg_energy = 3.90426
updating sclaing field ...
filter_energy = 9.74658
is_hit_min_edge_length = 1
//////////////// pass 10 ////////////////
edge splitting...
fixed 0 tangled element
success = 556(556)
edge splitting done!
time = 0.203775s
#v = 67338
#t = 379620
max_energy = 974.658
avg_energy = 3.90597
edge collapsing...
fixed 0 tangled element
success(env) = 14
success = 560(733941)
success(env) = 0
success = 4(6496)
success(env) = 0
success = 4(80)
success(env) = 0
success = 2(41)
success(env) = 0
success = 0(37)
edge collapsing done!
time = 3.90926s
#v = 66768
#t = 377017
max_energy = 974.658
avg_energy = 3.90454
edge swapping...
fixed 0 tangled element
success3 = 12
success4 = 263
success5 = 6
success = 281(220825)
edge swapping done!
time = 1.31279s
#v = 66768
#t = 377011
max_energy = 974.658
avg_energy = 3.90443
vertex smoothing...
success = 29471(62859)
vertex smoothing done!
time = 1.48971s
#v = 66768
#t = 377011
max_energy = 974.658
avg_energy = 3.90337
//////////////// pass 11 ////////////////
edge splitting...
fixed 0 tangled element
success = 541(541)
edge splitting done!
time = 0.206126s
#v = 67309
#t = 379448
max_energy = 974.658
avg_energy = 3.90476
edge collapsing...
fixed 0 tangled element
success(env) = 4
success = 540(733716)
success(env) = 0
success = 5(4384)
success(env) = 0
success = 0(86)
edge collapsing done!
time = 3.7363s
#v = 66764
#t = 376992
max_energy = 974.658
avg_energy = 3.90358
edge swapping...
fixed 0 tangled element
success3 = 8
success4 = 196
success5 = 3
success = 207(220680)
edge swapping done!
time = 1.25538s
#v = 66764
#t = 376987
max_energy = 974.658
avg_energy = 3.90337
vertex smoothing...
success = 29022(62859)
vertex smoothing done!
time = 1.46521s
#v = 66764
#t = 376987
max_energy = 974.658
avg_energy = 3.90274
updating sclaing field ...
filter_energy = 9.74658
is_hit_min_edge_length = 1
//////////////// pass 12 ////////////////
edge splitting...
fixed 0 tangled element
success = 574(574)
edge splitting done!
time = 0.206517s
#v = 67338
#t = 379593
max_energy = 974.658
avg_energy = 3.90474
edge collapsing...
fixed 0 tangled element
success(env) = 7
success = 569(733663)
success(env) = 0
success = 1(4943)
success(env) = 0
success = 0(9)
edge collapsing done!
time = 3.53614s
#v = 66768
#t = 376998
max_energy = 974.658
avg_energy = 3.90311
edge swapping...
fixed 0 tangled element
success3 = 4
success4 = 168
success5 = 7
success = 179(220660)
edge swapping done!
time = 1.26492s
#v = 66768
#t = 377001
max_energy = 974.658
avg_energy = 3.90296
vertex smoothing...
success = 28652(62855)
vertex smoothing done!
time = 1.43854s
#v = 66768
#t = 377001
max_energy = 974.658
avg_energy = 3.90234
//////////////// pass 13 ////////////////
edge splitting...
fixed 0 tangled element
success = 550(550)
edge splitting done!
time = 0.206317s
#v = 67318
#t = 379475
max_energy = 974.658
avg_energy = 3.90411
edge collapsing...
fixed 0 tangled element
success(env) = 3
success = 551(733583)
success(env) = 0
success = 3(4356)
success(env) = 0
success = 0(95)
edge collapsing done!
time = 3.36847s
#v = 66764
#t = 376982
max_energy = 974.658
avg_energy = 3.90269
edge swapping...
fixed 0 tangled element
success3 = 8
success4 = 148
success5 = 3
success = 159(220635)
edge swapping done!
time = 1.20379s
#v = 66764
#t = 376977
max_energy = 974.658
avg_energy = 3.90252
vertex smoothing...
success = 28454(62854)
vertex smoothing done!
time = 1.43652s
#v = 66764
#t = 376977
max_energy = 974.658
avg_energy = 3.90207
updating sclaing field ...
filter_energy = 9.74658
is_hit_min_edge_length = 1
//////////////// pass 14 ////////////////
edge splitting...
fixed 0 tangled element
success = 583(583)
edge splitting done!
time = 0.206728s
#v = 67347
#t = 379620
max_energy = 974.658
avg_energy = 3.90453
edge collapsing...
fixed 0 tangled element
success(env) = 3
success = 570(733643)
success(env) = 0
success = 12(3548)
success(env) = 0
success = 0(151)
edge collapsing done!
time = 3.37652s
#v = 66765
#t = 376983
max_energy = 974.658
avg_energy = 3.90244
edge swapping...
fixed 0 tangled element
success3 = 3
success4 = 132
success5 = 1
success = 136(220538)
edge swapping done!
time = 1.17794s
#v = 66765
#t = 376981
max_energy = 974.658
avg_energy = 3.90228
vertex smoothing...
success = 28200(62855)
vertex smoothing done!
time = 1.42366s
#v = 66765
#t = 376981
max_energy = 974.658
avg_energy = 3.90195
//////////////// pass 15 ////////////////
edge splitting...
fixed 0 tangled element
success = 595(595)
edge splitting done!
time = 0.204842s
#v = 67360
#t = 379681
max_energy = 974.658
avg_energy = 3.90423
edge collapsing...
fixed 0 tangled element
success(env) = 3
success = 577(733762)
success(env) = 0
success = 8(3923)
success(env) = 0
success = 1(510)
success(env) = 0
success = 0(32)
edge collapsing done!
time = 3.47169s
#v = 66774
#t = 377023
max_energy = 974.658
avg_energy = 3.90245
edge swapping...
fixed 0 tangled element
success3 = 0
success4 = 122
success5 = 3
success = 125(220541)
edge swapping done!
time = 1.2179s
#v = 66774
#t = 377026
max_energy = 974.658
avg_energy = 3.90232
vertex smoothing...
success = 28042(62865)
vertex smoothing done!
time = 1.43987s
#v = 66774
#t = 377026
max_energy = 974.658
avg_energy = 3.90189
updating sclaing field ...
filter_energy = 9.74658
is_hit_min_edge_length = 1
//////////////// pass 16 ////////////////
edge splitting...
fixed 0 tangled element
success = 595(595)
edge splitting done!
time = 0.20575s
#v = 67369
#t = 379721
max_energy = 974.658
avg_energy = 3.90467
edge collapsing...
fixed 0 tangled element
success(env) = 2
success = 582(733702)
success(env) = 0
success = 4(3559)
success(env) = 0
success = 2(132)
success(env) = 0
success = 1(27)
success(env) = 0
success = 1(26)
success(env) = 0
success = 0(12)
edge collapsing done!
time = 3.68698s
#v = 66779
#t = 377059
max_energy = 974.658
avg_energy = 3.90253
edge swapping...
fixed 0 tangled element
success3 = 7
success4 = 129
success5 = 3
success = 139(220627)
edge swapping done!
time = 1.21204s
#v = 66779
#t = 377055
max_energy = 974.658
avg_energy = 3.90234
vertex smoothing...
success = 27863(62872)
vertex smoothing done!
time = 1.43034s
#v = 66779
#t = 377055
max_energy = 974.658
avg_energy = 3.90194
//////////////// pass 17 ////////////////
edge splitting...
fixed 0 tangled element
success = 587(587)
edge splitting done!
time = 0.204838s
#v = 67366
#t = 379696
max_energy = 974.658
avg_energy = 3.90431
edge collapsing...
fixed 0 tangled element
success(env) = 3
success = 583(733622)
success(env) = 0
success = 1(3469)
success(env) = 0
success = 0(9)
edge collapsing done!
time = 3.49652s
#v = 66782
#t = 377067
max_energy = 974.658
avg_energy = 3.90231
edge swapping...
fixed 0 tangled element
success3 = 4
success4 = 113
success5 = 2
success = 119(220613)
edge swapping done!
time = 1.17513s
#v = 66782
#t = 377065
max_energy = 974.658
avg_energy = 3.9022
vertex smoothing...
success = 27697(62870)
vertex smoothing done!
time = 1.42799s
#v = 66782
#t = 377065
max_energy = 974.658
avg_energy = 3.90179
updating sclaing field ...
filter_energy = 9.74658
is_hit_min_edge_length = 1
//////////////// pass 18 ////////////////
edge splitting...
fixed 0 tangled element
success = 601(601)
edge splitting done!
time = 0.205337s
#v = 67383
#t = 379781
max_energy = 974.658
avg_energy = 3.90438
edge collapsing...
fixed 0 tangled element
success(env) = 2
success = 595(733722)
success(env) = 0
success = 9(3517)
success(env) = 0
success = 0(274)
edge collapsing done!
time = 3.4433s
#v = 66779
#t = 377050
max_energy = 974.658
avg_energy = 3.90209
edge swapping...
fixed 0 tangled element
success3 = 2
success4 = 105
success5 = 2
success = 109(220576)
edge swapping done!
time = 1.2573s
#v = 66779
#t = 377050
max_energy = 974.658
avg_energy = 3.90194
vertex smoothing...
success = 27489(62868)
vertex smoothing done!
time = 1.42578s
#v = 66779
#t = 377050
max_energy = 974.658
avg_energy = 3.9

We want to use this mesh within FEniCSx, so we extract the cells, nodes and volume markers from the mesh instance in wildmeshing.

point_array, cell_array, marker = tetra.get_tet_mesh()
[2025-10-06 15:06:02.802] [float-tetwild] [info] after winding number
[2025-10-06 15:06:02.803] [float-tetwild] [info] #v = 43289
[2025-10-06 15:06:02.805] [float-tetwild] [info] #t = 187028
[2025-10-06 15:06:02.805] [float-tetwild] [info] winding number 1.75976e+09s
[2025-10-06 15:06:02.805] [float-tetwild] [info] 

We can pass these into the DOLFINx mesh constructor called create_mesh.

mesh = dolfinx.mesh.create_mesh(
    MPI.COMM_WORLD,
    cells=cell_array.astype(np.int64),
    x=point_array,
    e=ufl.Mesh(basix.ufl.element("Lagrange", "tetrahedron", 1, shape=(3,))),
)

We add the cell markers as MeshTags with the following lines of code

local_entities, local_values = dolfinx.io.distribute_entity_data(
    mesh,
    mesh.topology.dim,
    cell_array.astype(np.int64),
    marker.flatten().astype(np.int32),
)
adj = dolfinx.graph.adjacencylist(local_entities)
ct = dolfinx.mesh.meshtags_from_entities(
    mesh,
    mesh.topology.dim,
    adj,
    local_values.astype(np.int32, copy=False),
)

We store the mesh and cell markers to file for usage in later tutorials

with dolfinx.io.XDMFFile(mesh.comm, folder / "brain.xdmf", "w") as xdmf:
    xdmf.write_mesh(mesh)
    xdmf.write_meshtags(ct, mesh.geometry)

We visualize the resulting mesh and cell tags with pyvista

pv_grid = pyvista.UnstructuredGrid(*dolfinx.plot.vtk_mesh(mesh))
pv_grid.cell_data["marker"] = ct.values
plotter = pyvista.Plotter()
plotter.add_mesh(pv_grid)
plotter.show()
2025-10-06 15:06:04.212 ( 286.442s) [    7F87F5599140]vtkXOpenGLRenderWindow.:1458  WARN| bad X server connection. DISPLAY=

We also create som slices to see the various markers within the mesh

slices = pv_grid.slice_along_axis(n=7, axis="y")
slices.plot()