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-16 07:41:40.960] [float-tetwild] [info] remove duplicates: 
[2025-10-16 07:41:40.960] [float-tetwild] [info] #v: 477353 -> 477353
[2025-10-16 07:41:40.960] [float-tetwild] [info] #f: 956254 -> 956254
collapsing 21.5822
swapping 0.265202
[2025-10-16 07:42:03.060] [float-tetwild] [info] remove duplicates: 
[2025-10-16 07:42:03.060] [float-tetwild] [info] #v: 49442 -> 49442
[2025-10-16 07:42:03.060] [float-tetwild] [info] #f: 100369 -> 100369
[2025-10-16 07:42:03.062] [float-tetwild] [info] #v = 49442
[2025-10-16 07:42:03.062] [float-tetwild] [info] #f = 100256
#boundary_e1 = 309
#boundary_e2 = 120
[2025-10-16 07:42:03.443] [float-tetwild] [info] preprocessing 22.6912s
[2025-10-16 07:42:03.443] [float-tetwild] [info] 
[2025-10-16 07:42:03.678] [float-tetwild] [info] #v = 50230
[2025-10-16 07:42:03.679] [float-tetwild] [info] #t = 333917
[2025-10-16 07:42:03.679] [float-tetwild] [info] tetrahedralizing 0.235683s
[2025-10-16 07:42:03.679] [float-tetwild] [info] 
[2025-10-16 07:42:03.680] [float-tetwild] [info] triangle insertion start, #f = 100256, #v = 50230, #t = 333917
[2025-10-16 07:42:03.816] [float-tetwild] [info] matched #f = 74494, uninserted #f = 25762
[2025-10-16 07:42:12.687] [float-tetwild] [info] insert_one_triangle * n done, #v = 154916, #t = 873987
[2025-10-16 07:42:12.687] [float-tetwild] [info] uninserted #f = 0/25762
[2025-10-16 07:42:12.687] [float-tetwild] [info] total timing: 0s
[2025-10-16 07:42:13.124] [float-tetwild] [info] pair_track_surface_fs done
#boundary_e1 = 309
#boundary_e2 = 120
[2025-10-16 07:42:13.174] [float-tetwild] [info] find_boundary_edges done
[2025-10-16 07:42:13.223] [float-tetwild] [info] time1 = 0.048667
[2025-10-16 07:42:13.235] [float-tetwild] [info] uninsert boundary #e = 0/429
[2025-10-16 07:42:13.235] [float-tetwild] [info] time2 = 4.2e-05
[2025-10-16 07:42:13.235] [float-tetwild] [info] time3 = 0.01137
[2025-10-16 07:42:13.235] [float-tetwild] [info] time4 = 0.00027
[2025-10-16 07:42:13.235] [float-tetwild] [info] time5 = 7.3e-05
[2025-10-16 07:42:13.235] [float-tetwild] [info] time6 = 2e-06
[2025-10-16 07:42:13.241] [float-tetwild] [info] uninserted #f = 0/25762
known_surface_fs.size = 0
known_not_surface_fs.size = 0
[2025-10-16 07:42:15.921] [float-tetwild] [info] mark_surface_fs done
[2025-10-16 07:42:15.942] [float-tetwild] [info] #b_edge1 = 429, #b_edges2 = 0
[2025-10-16 07:42:15.975] [float-tetwild] [info] cutting 12.295s
[2025-10-16 07:42:15.975] [float-tetwild] [info] 
initializing...
edge collapsing...
fixed 0 tangled element
success(env) = 86306
success = 101277(818074)
success(env) = 3177
success = 3411(534304)
success(env) = 308
success = 315(223412)
success(env) = 47
success = 48(32677)
success(env) = 7
success = 8(5227)
success(env) = 4
success = 4(1305)
success(env) = 0
success = 0(338)
edge collapsing done!
time = 13.6975s
#v = 44070
#t = 288270
max_energy = 7950.47
avg_energy = 7.68187
//////////////// pass 0 ////////////////
edge splitting...
fixed 0 tangled element
success = 42921(42921)
edge splitting done!
time = 0.306722s
#v = 86991
#t = 502900
max_energy = 7950.47
avg_energy = 7.61054
edge collapsing...
fixed 0 tangled element
success(env) = 180
success = 29300(621789)
success(env) = 80
success = 1392(156675)
success(env) = 20
success = 173(48145)
success(env) = 7
success = 33(7490)
success(env) = 0
success = 8(1728)
success(env) = 0
success = 6(185)
success(env) = 0
success = 0(217)
edge collapsing done!
time = 4.97643s
#v = 56079
#t = 348723
max_energy = 1146.91
avg_energy = 6.5897
edge swapping...
fixed 0 tangled element
success3 = 33933
success4 = 39591
success5 = 3876
success = 77400(344248)
edge swapping done!
time = 2.37762s
#v = 56079
#t = 318666
max_energy = 447.194
avg_energy = 5.54611
vertex smoothing...
success = 30000(52185)
vertex smoothing done!
time = 1.98857s
#v = 56079
#t = 318666
max_energy = 447.194
avg_energy = 5.21101
//////////////// pass 1 ////////////////
edge splitting...
fixed 0 tangled element
success = 5103(5103)
edge splitting done!
time = 0.190332s
#v = 61182
#t = 343553
max_energy = 447.194
avg_energy = 5.15867
edge collapsing...
fixed 0 tangled element
success(env) = 2653
success = 6930(518641)
success(env) = 168
success = 350(256741)
success(env) = 13
success = 33(25354)
success(env) = 1
success = 5(1923)
success(env) = 1
success = 1(272)
success(env) = 0
success = 0(177)
edge collapsing done!
time = 4.45017s
#v = 53863
#t = 305933
max_energy = 103.261
avg_energy = 5.11842
edge swapping...
fixed 0 tangled element
success3 = 3806
success4 = 14036
success5 = 1228
success = 19070(199173)
edge swapping done!
time = 1.26808s
#v = 53863
#t = 303355
max_energy = 103.261
avg_energy = 5.01371
vertex smoothing...
success = 25339(49967)
vertex smoothing done!
time = 1.76588s
#v = 53863
#t = 303355
max_energy = 103.261
avg_energy = 4.90499
//////////////// pass 2 ////////////////
edge splitting...
fixed 0 tangled element
success = 2800(2800)
edge splitting done!
time = 0.172871s
#v = 56663
#t = 316636
max_energy = 103.261
avg_energy = 4.88194
edge collapsing...
fixed 0 tangled element
success(env) = 913
success = 3458(493917)
success(env) = 40
success = 119(130718)
success(env) = 4
success = 13(7148)
success(env) = 1
success = 2(963)
success(env) = 0
success = 0(216)
edge collapsing done!
time = 3.37851s
#v = 53071
#t = 298821
max_energy = 103.261
avg_energy = 4.88656
edge swapping...
fixed 0 tangled element
success3 = 1602
success4 = 7098
success5 = 626
success = 9326(174011)
edge swapping done!
time = 1.03717s
#v = 53071
#t = 297845
max_energy = 103.261
avg_energy = 4.8501
vertex smoothing...
success = 21855(49177)
vertex smoothing done!
time = 1.61958s
#v = 53071
#t = 297845
max_energy = 40.5
avg_energy = 4.80148
//////////////// pass 3 ////////////////
edge splitting...
fixed 0 tangled element
success = 2072(2072)
edge splitting done!
time = 0.165004s
#v = 55143
#t = 307525
max_energy = 40.5
avg_energy = 4.78592
edge collapsing...
fixed 0 tangled element
success(env) = 398
success = 2342(485537)
success(env) = 26
success = 71(74734)
success(env) = 2
success = 4(4173)
success(env) = 0
success = 1(221)
success(env) = 0
success = 0(1)
edge collapsing done!
time = 3.0556s
#v = 52725
#t = 295869
max_energy = 40.5
avg_energy = 4.7958
edge swapping...
fixed 0 tangled element
success3 = 787
success4 = 3930
success5 = 312
success = 5029(163073)
edge swapping done!
time = 0.948127s
#v = 52725
#t = 295394
max_energy = 40.5
avg_energy = 4.7798
vertex smoothing...
success = 19510(48831)
vertex smoothing done!
time = 1.53287s
#v = 52725
#t = 295394
max_energy = 40.5
avg_energy = 4.75405
updating sclaing field ...
filter_energy = 8
is_hit_min_edge_length = 0
enlarge envelope, eps = 0.000368443
//////////////// pass 4 ////////////////
edge splitting...
fixed 0 tangled element
success = 66880(66880)
edge splitting done!
time = 0.44835s
#v = 119605
#t = 652526
max_energy = 39.1422
avg_energy = 4.8583
edge collapsing...
fixed 0 tangled element
success(env) = 2274
success = 44962(542250)
success(env) = 282
success = 3459(298833)
success(env) = 33
success = 413(90755)
success(env) = 8
success = 95(13280)
success(env) = 1
success = 22(3420)
success(env) = 0
success = 3(774)
success(env) = 0
success = 0(136)
edge collapsing done!
time = 6.01919s
#v = 70651
#t = 403102
max_energy = 39.1422
avg_energy = 4.37553
edge swapping...
fixed 0 tangled element
success3 = 4863
success4 = 21619
success5 = 1985
success = 28467(293781)
edge swapping done!
time = 2.55173s
#v = 70651
#t = 400224
max_energy = 39.1422
avg_energy = 4.2263
vertex smoothing...
success = 44098(66734)
vertex smoothing done!
time = 2.0484s
#v = 70651
#t = 400224
max_energy = 39.1422
avg_energy = 4.04574
//////////////// pass 5 ////////////////
edge splitting...
fixed 0 tangled element
success = 14688(14688)
edge splitting done!
time = 0.270489s
#v = 85339
#t = 471542
max_energy = 39.1422
avg_energy = 4.08277
edge collapsing...
fixed 0 tangled element
success(env) = 1338
success = 14360(430102)
success(env) = 65
success = 570(186992)
success(env) = 6
success = 58(15303)
success(env) = 0
success = 12(1417)
success(env) = 1
success = 2(237)
success(env) = 0
success = 2(36)
success(env) = 0
success = 0(50)
edge collapsing done!
time = 3.66097s
#v = 70335
#t = 398673
max_energy = 39.1422
avg_energy = 4.04039
edge swapping...
fixed 0 tangled element
success3 = 1517
success4 = 11095
success5 = 829
success = 13441(263088)
edge swapping done!
time = 1.68276s
#v = 70335
#t = 397985
max_energy = 39.1422
avg_energy = 4.00499
vertex smoothing...
success = 42078(66418)
vertex smoothing done!
time = 1.9282s
#v = 70335
#t = 397985
max_energy = 37.9377
avg_energy = 3.93227
//////////////// pass 6 ////////////////
edge splitting...
fixed 0 tangled element
success = 9562(9562)
edge splitting done!
time = 0.245499s
#v = 79897
#t = 443458
max_energy = 37.9377
avg_energy = 3.96403
edge collapsing...
fixed 0 tangled element
success(env) = 878
success = 9567(408567)
success(env) = 33
success = 264(127733)
success(env) = 2
success = 21(7290)
success(env) = 0
success = 4(613)
success(env) = 0
success = 1(76)
success(env) = 0
success = 0(20)
edge collapsing done!
time = 3.11253s
#v = 70040
#t = 396372
max_energy = 37.9377
avg_energy = 3.93619
edge swapping...
fixed 0 tangled element
success3 = 667
success4 = 6550
success5 = 388
success = 7605(249745)
edge swapping done!
time = 1.48773s
#v = 70040
#t = 396093
max_energy = 37.9377
avg_energy = 3.92187
vertex smoothing...
success = 39621(66128)
vertex smoothing done!
time = 1.83412s
#v = 70040
#t = 396093
max_energy = 37.9377
avg_energy = 3.88457
updating sclaing field ...
filter_energy = 8
is_hit_min_edge_length = 0
//////////////// pass 7 ////////////////
edge splitting...
fixed 0 tangled element
success = 24495(24495)
edge splitting done!
time = 0.308076s
#v = 94535
#t = 523538
max_energy = 39.4097
avg_energy = 4.03232
edge collapsing...
fixed 0 tangled element
success(env) = 2207
success = 16505(635500)
success(env) = 150
success = 1210(166677)
success(env) = 37
success = 158(24727)
success(env) = 7
success = 30(4253)
success(env) = 2
success = 7(996)
success(env) = 0
success = 1(232)
success(env) = 0
success = 0(40)
edge collapsing done!
time = 4.51474s
#v = 76624
#t = 433521
max_energy = 22.3441
avg_energy = 3.89868
edge swapping...
fixed 0 tangled element
success3 = 964
success4 = 7913
success5 = 539
success = 9416(282415)
edge swapping done!
time = 2.03168s
#v = 76624
#t = 433096
max_energy = 21.7772
avg_energy = 3.87404
vertex smoothing...
success = 45625(72708)
vertex smoothing done!
time = 1.9618s
#v = 76624
#t = 433096
max_energy = 17.9944
avg_energy = 3.81287
//////////////// pass 8 ////////////////
edge splitting...
fixed 0 tangled element
success = 9220(9220)
edge splitting done!
time = 0.266302s
#v = 85844
#t = 478087
max_energy = 17.9944
avg_energy = 3.86258
edge collapsing...
fixed 0 tangled element
success(env) = 1269
success = 9103(596118)
success(env) = 50
success = 341(112518)
success(env) = 6
success = 31(7571)
success(env) = 0
success = 3(730)
success(env) = 0
success = 0(28)
edge collapsing done!
time = 3.63024s
#v = 76366
#t = 432101
max_energy = 14.0016
avg_energy = 3.8173
edge swapping...
fixed 0 tangled element
success3 = 375
success4 = 4649
success5 = 215
success = 5239(273831)
edge swapping done!
time = 1.58442s
#v = 76366
#t = 431941
max_energy = 14.0016
avg_energy = 3.80976
vertex smoothing...
success = 43545(72458)
vertex smoothing done!
time = 1.88152s
#v = 76366
#t = 431941
max_energy = 14.0016
avg_energy = 3.78205
//////////////// pass 9 ////////////////
edge splitting...
fixed 0 tangled element
success = 7005(7005)
edge splitting done!
time = 0.253016s
#v = 83371
#t = 465608
max_energy = 14.0016
avg_energy = 3.82155
edge collapsing...
fixed 0 tangled element
success(env) = 1003
success = 6950(584188)
success(env) = 32
success = 187(88041)
success(env) = 2
success = 13(4004)
success(env) = 0
success = 0(371)
edge collapsing done!
time = 3.31702s
#v = 76221
#t = 431257
max_energy = 14.0016
avg_energy = 3.78715
edge swapping...
fixed 0 tangled element
success3 = 169
success4 = 2837
success5 = 122
success = 3128(269258)
edge swapping done!
time = 1.51073s
#v = 76221
#t = 431210
max_energy = 14.0016
avg_energy = 3.78393
vertex smoothing...
success = 41825(72311)
vertex smoothing done!
time = 1.81061s
#v = 76221
#t = 431210
max_energy = 14.0016
avg_energy = 3.76832
updating sclaing field ...
filter_energy = 8
is_hit_min_edge_length = 0
//////////////// pass 10 ////////////////
edge splitting...
fixed 0 tangled element
success = 2382(2382)
edge splitting done!
time = 0.23486s
#v = 78603
#t = 443371
max_energy = 18.7462
avg_energy = 3.78539
edge collapsing...
fixed 0 tangled element
success(env) = 327
success = 2121(812171)
success(env) = 17
success = 122(50072)
success(env) = 4
success = 20(3879)
success(env) = 0
success = 8(815)
success(env) = 1
success = 3(118)
success(env) = 0
success = 1(44)
success(env) = 0
success = 0(11)
edge collapsing done!
time = 4.22865s
#v = 76328
#t = 431900
max_energy = 11.9257
avg_energy = 3.77267
edge swapping...
fixed 0 tangled element
success3 = 99
success4 = 1740
success5 = 61
success = 1900(267681)
edge swapping done!
time = 1.47563s
#v = 76328
#t = 431862
max_energy = 11.9257
avg_energy = 3.77034
vertex smoothing...
success = 40789(72417)
vertex smoothing done!
time = 1.76812s
#v = 76328
#t = 431862
max_energy = 11.2283
avg_energy = 3.76061
//////////////// pass 11 ////////////////
edge splitting...
fixed 0 tangled element
success = 1299(1299)
edge splitting done!
time = 0.229431s
#v = 77627
#t = 438065
max_energy = 13.7945
avg_energy = 3.76789
edge collapsing...
fixed 0 tangled element
success(env) = 161
success = 1320(807260)
success(env) = 4
success = 42(19172)
success(env) = 1
success = 5(1068)
success(env) = 0
success = 1(222)
success(env) = 0
success = 0(96)
edge collapsing done!
time = 3.91088s
#v = 76259
#t = 431516
max_energy = 11.2283
avg_energy = 3.76155
edge swapping...
fixed 0 tangled element
success3 = 52
success4 = 993
success5 = 39
success = 1084(266044)
edge swapping done!
time = 1.42539s
#v = 76259
#t = 431503
max_energy = 11.2283
avg_energy = 3.76083
vertex smoothing...
success = 39742(72350)
vertex smoothing done!
time = 1.71217s
#v = 76259
#t = 431503
max_energy = 11.1994
avg_energy = 3.7561
updating sclaing field ...
filter_energy = 8
is_hit_min_edge_length = 1
//////////////// pass 12 ////////////////
edge splitting...
fixed 0 tangled element
success = 534(534)
edge splitting done!
time = 0.224954s
#v = 76793
#t = 433889
max_energy = 17.5679
avg_energy = 3.75868
edge collapsing...
fixed 0 tangled element
success(env) = 48
success = 649(864966)
success(env) = 2
success = 11(17328)
success(env) = 0
success = 1(820)
success(env) = 0
success = 1(154)
success(env) = 0
success = 3(101)
success(env) = 0
success = 0(232)
edge collapsing done!
time = 4.21127s
#v = 76128
#t = 430780
max_energy = 11.7479
avg_energy = 3.75699
edge swapping...
fixed 0 tangled element
success3 = 24
success4 = 575
success5 = 15
success = 614(264735)
edge swapping done!
time = 1.40537s
#v = 76128
#t = 430771
max_energy = 11.7479
avg_energy = 3.75667
vertex smoothing...
success = 38804(72209)
vertex smoothing done!
time = 1.66572s
#v = 76128
#t = 430771
max_energy = 11.6955
avg_energy = 3.75417
//////////////// pass 13 ////////////////
edge splitting...
fixed 0 tangled element
success = 494(494)
edge splitting done!
time = 0.2246s
#v = 76622
#t = 432935
max_energy = 11.6955
avg_energy = 3.75573
edge collapsing...
fixed 0 tangled element
success(env) = 16
success = 535(863877)
success(env) = 0
success = 2(7548)
success(env) = 0
success = 0(130)
edge collapsing done!
time = 3.8827s
#v = 76085
#t = 430533
max_energy = 11.6955
avg_energy = 3.75437
edge swapping...
fixed 0 tangled element
success3 = 10
success4 = 343
success5 = 10
success = 363(264095)
edge swapping done!
time = 1.40827s
#v = 76085
#t = 430533
max_energy = 11.6955
avg_energy = 3.75417
vertex smoothing...
success = 38202(72173)
vertex smoothing done!
time = 1.63285s
#v = 76085
#t = 430533
max_energy = 12.4609
avg_energy = 3.75278
updating sclaing field ...
filter_energy = 8
is_hit_min_edge_length = 1
//////////////// pass 14 ////////////////
edge splitting...
fixed 0 tangled element
success = 528(528)
edge splitting done!
time = 0.222697s
#v = 76613
#t = 432863
max_energy = 12.6897
avg_energy = 3.75479
edge collapsing...
fixed 0 tangled element
success(env) = 19
success = 565(867550)
success(env) = 0
success = 1(6677)
success(env) = 0
success = 0(17)
edge collapsing done!
time = 3.88329s
#v = 76047
#t = 430327
max_energy = 10.1335
avg_energy = 3.75313
edge swapping...
fixed 0 tangled element
success3 = 10
success4 = 249
success5 = 1
success = 260(263790)
edge swapping done!
time = 1.39888s
#v = 76047
#t = 430318
max_energy = 9.05774
avg_energy = 3.75292
vertex smoothing...
success = 37649(72127)
vertex smoothing done!
time = 1.60813s
#v = 76047
#t = 430318
max_energy = 9.05774
avg_energy = 3.75212
//////////////// postprocessing ////////////////
edge collapsing...
fixed 0 tangled element
success(env) = 7
success = 16(866391)
success(env) = 0
success = 0(1665)
edge collapsing done!
time = 3.72921s
#v = 76031
#t = 430229
max_energy = 9.05774
avg_energy = 3.75218
edge collapsing...
fixed 0 tangled element
success(env) = 6682
success = 29330(565803)
success(env) = 393
success = 871(414733)
success(env) = 24
success = 58(87459)
success(env) = 0
success = 5(6162)
success(env) = 0
success = 2(592)
success(env) = 0
success = 0(167)
edge collapsing done!
time = 7.40013s
#v = 45765
#t = 257076
max_energy = 9.99987
avg_energy = 4.80594
edge swapping...
fixed 0 tangled element
success3 = 5234
success4 = 19105
success5 = 2274
success = 26613(186943)
edge swapping done!
time = 1.20406s
#v = 45765
#t = 254116
max_energy = 9.99805
avg_energy = 4.55913
edge collapsing...
fixed 0 tangled element
success(env) = 775
success = 2244(403755)
success(env) = 69
success = 111(183906)
success(env) = 4
success = 6(14147)
success(env) = 0
success = 0(605)
edge collapsing done!
time = 3.59133s
#v = 43404
#t = 239899
max_energy = 9.99805
avg_energy = 4.71496
edge swapping...
fixed 0 tangled element
success3 = 1020
success4 = 4596
success5 = 594
success = 6210(134112)
edge swapping done!
time = 0.796577s
#v = 43404
#t = 239473
max_energy = 9.99805
avg_energy = 4.66066
edge collapsing...
fixed 0 tangled element
success(env) = 82
success = 321(382584)
success(env) = 17
success = 25(41395)
success(env) = 1
success = 2(3353)
success(env) = 0
success = 0(238)
edge collapsing done!
time = 2.606s
#v = 43056
#t = 237355
max_energy = 9.99805
avg_energy = 4.68971
edge swapping...
fixed 0 tangled element
success3 = 205
success4 = 873
success5 = 145
success = 1223(121831)
edge swapping done!
time = 0.680206s
#v = 43056
#t = 237295
max_energy = 9.99555
avg_energy = 4.67922
edge collapsing...
fixed 0 tangled element
success(env) = 14
success = 36(379098)
success(env) = 1
success = 1(4941)
success(env) = 0
success = 0(159)
edge collapsing done!
time = 2.34716s
#v = 43019
#t = 237075
max_energy = 9.99555
avg_energy = 4.68242
edge swapping...
fixed 0 tangled element
success3 = 27
success4 = 98
success5 = 27
success = 152(119388)
edge swapping done!
time = 0.659459s
#v = 43019
#t = 237075
max_energy = 9.99555
avg_energy = 4.6812
[2025-10-16 07:44:47.816] [float-tetwild] [info] mesh optimization 151.833s
[2025-10-16 07:44:47.816] [float-tetwild] [info] 
[2025-10-16 07:44:48.048] [float-tetwild] [info] correct_tracked_surface_orientation done

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-16 07:44:49.346] [float-tetwild] [info] after winding number

[2025-10-16 07:44:49.346] [float-tetwild] [info] #v = 43019
[2025-10-16 07:44:49.349] [float-tetwild] [info] #t = 185467
[2025-10-16 07:44:49.349] [float-tetwild] [info] winding number 1.7606e+09s
[2025-10-16 07:44:49.349] [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-16 07:44:50.669 ( 191.153s) [    7F664A5FE140]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()