1 s2.0 S0169260722003911 Main

Computer Methods and Programs in Biomedicine 224 (2022) 107009
Contents lists available at ScienceDirect
Computer Methods and Programs in Biomedicine

journal homepage: www.elsevier.com/locate/cmpb
Open-Full-Jaw: An open-access dataset and pipeline for finite element

models of human jaw
Torkan Gholamalizadeh a,b,∗, Faezeh Moshfeghifar a, Zachary Ferguson c, Teseo Schneider d,
Daniele Panozzo c, Sune Darkner a, Masrour Makaremi e,f, François Chan f,
Peter Lampel Søndergaard b, Kenny Erleben a
a
Department of Computer Science, University of Copenhagen, Copenhagen 2100, Denmark
b
3Shape A/S, Copenhagen 1060, Denmark
c
Courant Institute of Mathematical Sciences, New York University, 60 5th Ave, New York NY 10011, USA
d
Department of Computer Science, University of Victoria, Victoria BC V8P 5C2, Canada
e
Dentofacial Orthopedics Department, University of Bordeaux, Bordeaux, France
f
Orthodontie clinic, 2 Rue des 2 Conils, Bergerac 24100, France
a r t i c l e i n f o a b s t r a c t
Article history: Background: State-of-the-art finite element studies on human jaws are mostly limited to the geometry of
Received 23 March 2022 a single patient. In general, developing accurate patient-specific computational models of the human jaw
Revised 2 July 2022
acquired from cone-beam computed tomography (CBCT) scans is labor-intensive and non-trivial, which
Accepted 5 July 2022
involves time-consuming human-in-the-loop procedures, such as segmentation, geometry reconstruction,
and re-meshing tasks. Therefore, with the current practice, researchers need to spend considerable time
Keywords: and effort to produce finite element models (FEMs) to get to the point where they can use the models
Finite element to answer clinically-interesting questions. Besides, any manual task involved in the process makes it dif-
Open-access dataset ficult for the researchers to reproduce identical models generated in the literature. Hence, a quantitative
CBCT scan
comparison is not attainable due to the lack of surface/volumetric meshes and FEMs. Methods: We share
Human jaw
an open-access repository composed of 17 patient-specific computational models of human jaws and the
Geometry reconstruction
Conformal mesh utilized pipeline for generating them for reproducibility of our work. The used pipeline minimizes the
required time for processing and any potential biases in the model generation process caused by human
intervention. It gets the segmented geometries with irregular and dense surface meshes and provides
reduced, adaptive, watertight, and conformal surface/volumetric meshes, which can directly be used in
finite element (FE) analysis. Results: We have quantified the variability of our 17 models and assessed the
accuracy of the developed models from three different aspects; (1) the maximum deviations from the
input meshes using the Hausdorff distance as an error measurement, (2) the quality of the developed
volumetric meshes, and (3) the stability of the FE models under two different scenarios of tipping and
biting. Conclusions: The obtained results indicate that the developed computational models are precise,
and they consist of quality meshes suitable for various FE scenarios. We believe the provided dataset of
models including a high geometrical variation obtained from 17 different models will pave the way for
population studies focusing on the biomechanical behavior of human jaws.
© 2022 The Authors. Published by Elsevier B.V.
This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/)
1. Introduction and cardiovascular surgeries [2]. More specifically, in the orthodon-

tic and dental fields of studies, FEM is utilized to predict teeth
Finite element modeling (FEM) is a numerical approach for pre- movements, stress/strain distribution in different tissues (e.g., peri-
dicting responses of different tissues under physical loads, which odontal ligament, gingiva, and alveolar bone) or orthodontic appli-
can be difficult or impossible to measure directly in vivo [1]. It is a ances [3]. Except for a few recent studies [4–6], almost all of the
widely used tool as a pre-operative protocol in different medical previous studies in the field are limited to single-patient analysis
applications such as orthopedic surgery, orthodontic treatments, [2,7–11], in which the results might not be generalized to a larger
population with high geometrical variations in the teeth, periodon-
∗
tal ligament (PDL), and bone anatomies [4,6].
Corresponding author.
E-mail address: torkan@di.ku.dk (T. Gholamalizadeh).
https://doi.org/10.1016/j.cmpb.2022.107009
0169-2607/© 2022 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/)
T. Gholamalizadeh, F. Moshfeghifar, Z. Ferguson et al. Computer Methods and Programs in Biomedicine 224 (2022) 107009
The main reason for using a single model in the literature is

that developing accurate computational models of the human jaw
FE model file
Related studies in the literature and their details on the data availability and mesh conformity. The few studies with public models are listed below the dashed line. Note that there are empty cells in some studies, as
See Table 2
is challenging and involves time-consuming and labor-intensive
processes such as segmentation, geometry reconstruction, geom-
etry processing, re-meshing, and mesh simplification tasks. For in-
Yes
Yes
No
No
No
No
No
No
No
No
No
No
No
No
stance, generating a complex and highly detailed finite element
∗
CG: Congruent interface (no gap/penetration in the contacting interfaces); Non-CG: Non-congruent interface; CM: Conformal Mesh; Non-CM: Non-Cnformal Mesh. Obtained from cadaver or dried skull.
Volumetric Meshes
(FE) model of the entire human jaw can take up to several months
per scan [10]. Therefore, developing several patient-specific FE
models may not be feasible for many researchers. In addition,
currently, there are no publicly available datasets of full denti-
tion human jaw to be used by researchers, except for two studies
Yes
Yes
Yes
N/A: Not Applicable; DM: Dental Model; H: Homogeneous; EDP: Enamel, Dentine, and Pulp; CC: Cortical and Cancellous bone; O: Orthotropic bone, further details can be found in [10];
No
No
No
No
No
No
No
No
No
No
No
[10,17] with a limited number of studied subjects.
Surface Meshes
In other words, in almost all of the studies focusing on full den-

tition or single/multiple tooth analysis [2,2,6–9,11,14–16], the uti-
Availability
lized geometries, volumetric meshes, and FE models have not been

made publicly available, which makes it difficult to reproduce and
Yes
Yes
Yes
Yes
No
No
No
No
No
No
No
No
No
No
compare the results. Table 1 presents an overview of the related

studies in the literature by providing details on the studied cohort,
Mesh Conformity
discretization type, and availability of the models.

As one of the few studies with public data, the OpenMandible
Non-CM
Non-CM
Non-CM
Non-CM
Non-CM
Non-CM
Non-CM
Non-CM
[10] provides detailed geometries of one mandible structure and

all teeth obtained from a dried male skull. The study scans the
CM
CM
CM
CM
mandible in two different steps to provide detailed teeth struc-

Interface Congruency
tures, i.e., pulp and enamel. First, the mandible was scanned
using a CBCT scanner with a voxel size of 0.133 mm. Second,
Interface info
all mandibular teeth were removed from the bone before being
scanned by the micro-computed tomography (micro-CT or μCT).
†
Non-CG
Non-CG
Non-CG
Micro-CTs are high-resolution CT scans that are normally acquired

CG †
CG †
CG
CG
CG
CG
CG
CG
CG
from dead specimens or cadavers due to high x-ray exposures.

However, the detailed geometries obtained from a single mandible
CC, Semi-Adaptive
cannot cover geometrical variations across different patients in the

population.
CC, Adaptive
CC, Adaptive
CC, Adaptive
CC, Uniform
CC, Uniform
CC, Uniform
H, Adaptive
H, Adaptive
O, Adaptive
H, Uniform
Bone Type
The recently introduced OpenJaw Dataset [17] provides open-

access reconstructed geometries of three patients’ mandibles ac-
N/A
N/A
Geometry and Discretization Info
CC
quired from CBCT scans. Each patients data includes surface

meshes of the reconstructed mandible, teeth, and PDL geometries.
Uniform
Uniform
Uniform
Uniform
Uniform
Uniform
Uniform
Uniform
Uniform
Uniform
Uniform
Uniform
Uniform
PDL Type
In the geometry reconstruction step, the utilized scans in the study

(one with 0.15 mm and two with 0.3 mm voxel sizes) were upsam-
H,
H,
H,
H,
H,
H,
H,
H,
H,
H,
H,
H,
H,
pled to the same resolution of 0.15 mm [5]. Although this study

H
provides more public samples, generalizability to the population

EDP, Adaptive
EDP, Adaptive
EDP, Uniform
EDP, Uniform
EDP, Uniform
H, Adaptive
H, Adaptive
for further details., † Non-smooth top surface of PDL, see [12,13,15] for illustrations.
remains a problem. Moreover, both of the abovementioned studies

H, Uniform
H, Uniform
H, Uniform
H, Uniform
H, Uniform
Tooth Type
include manual tasks in different geometry-processing or meshing

tools making it difficult for other researchers to reproduce their
EDP
EDP
meshes using the unprocessed reconstructed meshes.

Scan Modality
Our main contributions in the Open-Full-Jaw study can be sum-

CBCT + MRI
CBCT + μCT
CBCT + μCT
marized as,
Scan Info
CBCT
CBCT
CBCT
CBCT
CBCT
CBCT
CBCT
μCT
N/A
N/A
1. We provide an open-access dataset of different patient-specific

CT
models of the human jaw, including the maxilla, mandible, full

Single Tooth
Single Tooth
3 Mandibles
1 Mandible
1 Mandible
1 Mandible
1 Mandible
dentition, and the PDL geometries obtained from CBCT scans

1 Maxilla
29 Jaws ∗
Half-arch
Half-arch
of 17 patients. It is the largest publicly available dataset with

4 Jaws
2 Jaws
2 jaws
#Jaws
validated segmented geometries and quality volumetric meshes

Cohort Info
that can directly be used in FEM studies. Furthermore, to the

it was difficult to evaluate the meshes.
#Patient
best of our knowledge, it is the first repository containing the

maxillary jaws of different patients.
17 ∗
DM
DM
1
1
1
1
1
2
1
1
2. We introduce a unique repository (https://github.com/diku-dk/

Kawamura et al. (2019) [12]
Kawamura et al. (2022) [13]
Sarrafpour et al. (2013) [15]
Open- Full- Jaw) containing (1) clinically validated segmented

Savignano et al. (2020) [6]
OpenMandible (2021) [10]

Oenning et al. (2018) [14]
Open-Full-Jaw (our study)

Benaissa et al. (2020) [2]
geometries and the resulting dense irregular surface meshes;

Boryor et al. (2009) [11]
OpenJaw (2021) [5,17]

Ortun et al. (2020) [9]
Ding et al. (2014) [8]
Lee et al. (2018) [16]
(2) the quality and adaptive volumetric/surface meshes to be

Seo et al. (2021) [7]
used in FE simulations; (3) the automatically generated FEM

files for tipping and biting scenarios used for the FE analysis
Study/Dataset
of this work; (4) the principal axes of every patient’s tooth pro-
Table 1
viding great information for the users to automatically set up

different loading conditions.
3. For reproducibility, we share our pipeline developed based
on open-source meshing tools [18,19] to generate the mod-
2
els of this study. This python-based library automates the FE contacting interface between the two adjacent/contacting domains.
model generation process, including geometry processing and The agreement of the two contacting geometries on the contact-
re-meshing tasks with minimal human intervention, by setting ing boundary or shared interface can be analyzed in two different
a few required parameters. levels; first, on the geometry level and based on the curvature of
4. This pipeline allows other researchers in the field to generate the contacting surfaces; second, on the mesh-based level focusing
quality volumetric meshes and FE models directly from dense on the agreement on the identical discretization of the contacting
and uncleaned meshes with minimal human intervention. This interface in terms of the position of vertices, edges, and faces. The
will help other researchers to easily extend their datasets with- former is called here the “interface congruency”. In the biomechan-
out spending much time and effort on manually cleaning up ical field of study, different theoretical and computational studies
the meshes and non-trivially producing conformal meshes. analyze the effect of the “ball-and-socket” joint congruencies such
5. Our pipeline ensures conformal meshes in the contacting inter- as in shoulder, hip, and temporomandibular joints to analyze their
faces without any undesired gaps or penetrations and provides instabilities and dislocations under different circumstances [20,21].
adaptive meshes that are vital for reducing the total number It should be noted that the current study only focuses on devel-
of elements while using finer meshes in specific regions, e.g., oping computational models of human mandible and maxilla for
teeth sockets, alveolar crest, and alveolar process. tooth movements and the mesh congruency of the contacting sur-
faces instead of congruency of the “ball-and-socket” joints, and the
All in all, we believe the Open-Full-Jaw dataset can be used
utilized “congruency” term needs to be distinguished from those in
for various intra- and inter-patient analyses [4,6] such as intact
biomechanical studies for analyzing the joint congruencies [20,21].
tooth movement modeling, bite force estimation, restorative pro-
As mentioned before, by using the conventional one-by-one do-
cedures modeling including cavity fillings and dental implants, just
main discretization, mesh reduction, or quality meshing processes,
to name a few. Besides, it can greatly impact the reproducibility of
it is challenging to achieve congruent surfaces due to the gener-
future studies. For the reproducibility of this study, we use open-
ation of gaps/penetrations between two contacting surfaces (see
source meshing tools. Still, the reconstructed geometries provided
Figure 1). We propose using signed distance fields in the contact-
in our dataset can be imported into any desired open-source or
ing surfaces to quantitatively evaluate the error between two con-
commercial meshing tools or FE frameworks to re-mesh and gen-
tacting regions. The signed distance function, φ (x ), for a domain
erate computational models.
⊂ R3 and an arbitrary point x ∈ R3 is defined as

2. Important traits required for a successful FEM −d (x, ∂ ), if x ∈ ,
φ ( x ) = 0 , if x ∈ ∂ , (1)
Developing a patient-specific FE model begins with segment- d (x, ∂ ), if x ∈ c ,
ing/annotating the desired regions of the medical scan obtained
from the patients. Next, the segmented regions are reconstructed where ∂ denotes the boundary of the domain , c represents
as surface meshes generally composed of irregular dense meshes the complement of , and d (x, ∂ ) indicates the Euclidean dis-
with no guarantees of manifoldness, watertightness, or absence of tance between the point x and the boundary of the domain ∂ .
self-intersection, which are crucial for developing stable and accu- We measure the congruency of the surfaces by calculating the
rate computational models. Hence, one needs to generate quality signed distance field of each contacting surface with respect to the
meshes from the exported dense and irregular meshes that are not other domain as
necessarily guaranteed to have the mentioned criteria. ε = max {|φ2 (x )|, |φ1 (y )|},
Moreover, different preprocessing steps such as geometry pro- ∀ x ∈ ∂ 1 c , y ∈ ∂ 2 c , (2)
cessing, mesh reduction, and re-meshing are essential for develop- ∂ 1 c ⊆ ∂ 1 , ∂ 2 c ⊆ ∂ 2 ,
ing FE models from image-based reconstructed geometries. When
modeling geometries with shared contacting interfaces, each of the where ∂ 1 and ∂ 2 denote the boundaries of the domains 1 and
mentioned processes can produce errors on the contacting sur- 2 , respectively, and ∂ 1c and ∂ 2c refer to the contacting sur-
faces and result in undesired gaps/penetrations between them. In faces of the domains. The computed ε values for two contacting
the cases where two adjacent segments are watertight, it is still domains are then called epsilon-congruent surfaces. The epsilon
challenging to discretize the segmented domains such that they value indicates the measurement error, where the values close to
agree on the same discretization on the shared contacting inter- zero refer to completely congruent surfaces.
faces. Therefore, the focus of this section is on essential aspects
needed to be considered for the discretization of the computa- 2.2. Conformal mesh interfaces
tional domains with shared contacting interfaces; we also discuss
potential options for developing a proper FE model of the human A special case for the congruent surfaces is “conformal meshes”
jaw. or “mesh conformity” and can be described as the identical dis-
cretization of the contacting interfaces. That is to say, for two con-
2.1. Congruent contacting interfaces tacting domains/geometries like 1 and 2 , the contacting sur-
faces ∂ 1c and ∂ 2c are assumed identical. More specifically, they
The discretization of the computational domain is an essential share identical vertices, edges, and elements on the contacting
step for developing computational models considering the biome- interfaces. In general, applying the conformal mesh criterion on
chanical behavior of tissues having sliding or bilateral contact with multi-domains with contacting surfaces is a challenging process,
other domains/tissues. It mainly affects the numerical stability of as most of the geometry and mesh processing algorithms used in
the computational models when soft structures, e.g., PDL, contact different free software produce the quality meshes per domain in-
hard tissues like bone or teeth. Besides, a coarse or poor-quality dependently. Therefore, the contacting surfaces should be identi-
discretization can itself cause a locking effect and influence the ac- fied and combined as a different step to create identical contacting
curacy of the stress/strain concentrations compared to that of the interfaces on the mesh of each domain.
values measured in vivo. Finally, the volumetric meshes can be generated by enforcing
When two contacting domains are reconstructed, discretized, the meshing algorithm to fully preserve the input surface mesh.
or re-meshed independently, undesired gaps/penetrations are in- This raises questions about the final mesh quality, as the additional
evitable, causing two different boundary definitions for an identical surface preservation constraint affects the mesh generation and the
3
Fig. 1. A visualized comparison between the ε -congruent and conformal meshes. A: the ε -congruent contacting surfaces with gaps/penetrations in the contacting interfaces.
B: the same structure with the multi-domain discretization resulting in conformal meshes. C: a close-up cross-section view of the gaps/penetrations at the root level. D: a
schematic illustration of the ε -congruent and conformal meshes.
optimization algorithms. Hence, the best approach for generating a robust filtering procedure based on either flood fill or the gener-
conformal meshes is to consider all contacting domains simultane- alized winding number [24], in both cases assuming that the user
ously while generating the volumetric meshes. is interested in a single material mesh. We propose a novel filter-
Generating multi-domain volumetric meshes would generate ing procedure for the filtering of multi-material tetrahedral meshes
zero-congruent and conformal meshes that omit any numerical er- common in medical imaging in Section 3.5.5, and we show that
rors caused by non-congruent contacting meshes. Therefore, it is our extension of TetWild is ideal for constructing our dataset, as
essential to utilize a meshing algorithm that generates volumet- it removes the expensive and tedious manual cleanup of the input
ric meshes considering multi-domain boundaries and performing surface meshes.
boolean operations on the input surface meshes by using an im-
plicit representation of the input meshes. In this study, we use 3. Controlling shared interfaces using volume mesh generation
fTetWild [18] that supports all aspects mentioned above for multi-
domain volumetric meshes. This provides a volumetric mesh of This section reviews the entire process for developing the FE
teeth-PDL-bone geometries in which there are shared points in the models of 17 jaws presented in Figure 2. To be more specific,
teeth-PDL and PDL-bone contacting surfaces. Section 3.1 describes the utilized criteria for the cohort selection
process; Sections 3.2 and 3.3 provide a detailed description for the
2.3. Proper modeling of the PDL layer CBCT segmentation and clinical validation; the geometrical varia-
tions of the reconstructed mandibles and maxillae are investigated
PDL is a thin structure that connects alveolar bone to the teeth in Section 3.4 based on widely used clinical landmarks; and finally,
cementum and acts as a shock absorber in the chewing or biting the different steps of the used pipeline for developing high-quality
process. It also plays an important role in transferring load from volumetric meshes of human jaws are studied in Section 3.5.
teeth to the bone in orthodontic treatments; when triggered with
enough orthodontic forces, it results in a bone remodeling process.
As this thin and soft structure shares interfaces with two hard tis- 3.1. Cohort selection
sues of the teeth and bone geometries, the computational model
may become unstable if the PDL geometry does not have congru- We use CBCT scans of 17 different patients with different voxel
ent surfaces with teeth and bone in the contacting interfaces. resolutions from 3Shape A/S in-house dataset. Various criteria in-
cluding the original voxel size of the scan, minimal metal filling
2.4. Quality mesh generation artifacts, and the absence of implants and severe periodontal dis-
eases are considered in selecting the mentioned cohort. Also, the
Tetrahedralization algorithms usually require a set of closed, selected cohort has no evidence of maxillofacial surgery or skeletal
self-intersection free surface meshes as input, which is a require- diseases. Sensitive information of the patients such as name, age,
ment that is challenging to guarantee in our setting without man- and gender are stripped due to the General Data Protection Regu-
ual interaction. In fact, traditional medical imaging pipeline require lation (GDPR) rules. The utilized scans were acquired from the pa-
the large majority of the manual processing in this cleanup phase, tients by their associated doctors/orthodontists as a part of treat-
as self-intersection, holes, or other imperfections normally appear ment plans, and the authors had no roles in the acquisition pro-
when processing real-world data, and are especially common (and cess.
time consuming) when dealing with thin layers, like the PDL layer.
We propose a very different approach, where we accept that 3.2. Data specifications and geometry reconstruction
these imperfection exist, and design a meshing pipeline that toler-
ates and automatically heals them. We base our approach on the To reconstruct the patient-specific geometries, first, the scans
method recently introduced in [22,23]: instead of meshing one do- are imported in 3DSlicer [30] in the standard Digital Imaging and
main at a time, TetWild meshes the entire volume of the bound- Communications in Medicine (DICOM) format. Table 2 provides de-
ing box containing the soup of triangles of all surfaces of interest. tails of the scans utilized in this study. Next, according to the pre-
The triangles are approximated with faces of the tetrahedral mesh. evaluation criteria (metal fillings or implant artifacts), we decide
This procedure does not require clean input geometry, it can toler- on which jaws are suitable to be segmented from the scan.
ate degenerate input triangles, self-intersections, and holes [22]. In The resolution of the selected CBCT scans is at most 0.3 mm
the original TetWild algorithms, the final mesh is extracted using since based on our experience, accurate tooth-bone segmentation
4
Fig. 2. Reconstructed 3D geometries from 17 CBCT scans. Note that the maxillary jaw is removed in some scans due to having several metal artifacts and/or missing teeth.
The anatomical variations in different jaws and teeth of our dataset indicate the necessity for introducing such a model repository for assessing the generalizability of the
FEM results in the related clinical/population studies.
Fig. 3. The utilized morphological landmarks for analyzing and quantification of the geometrical variation of the reconstructed mandibles and maxillae in different models.
Table 2
Specifications of the utilized CBCT scans. Different scans with various filed-of-views and slice thicknesses are used in this study. All scans are converted to an
identical slice thickness of 0.15 mm to avoid biases in the geometry reconstruction process, especially when applying smoothing filters to eliminate noises on the
segmented teeth and bone geometries. Note that the maxillary jaws of the patients with many metal artifacts and/or missing teeth are removed.
Patient’s ID Original scan Final resampled ROI Jaws Periodontal disease
Dimension Slice thickness (mm) Dimension Slice thickness (mm) Mandible Maxilla Mandible Maxilla
Patient 1 400 × 400 × 280 0.3 676 × 530 × 280 0.15 Yes - No Missing teeth
Patient 2 400 × 400 × 280 0.3 670 × 440 × 344 0.15 Yes - No Missing teeth
Patient 3 532 × 532 × 540 0.15 534 × 435 × 338 0.15 Yes Yes Moderate Moderate
Patient 4 532 × 532 × 540 0.15 534 × 435 × 338 0.15 Yes Yes Mild Mild bone loss
Patient 5 532 × 532 × 540 0.15 525 × 425 × 290 0.15 Yes Yes No No
Patient 7 400 × 400 × 280 0.3 614 × 470 × 270 0.15 Yes - No -
Patient 9 400 × 400 × 280 0.3 525 × 425 × 290 0.15 Yes - No -
Patient 10 400 × 400 × 280 0.3 525 × 425 × 290 0.15 Yes - No -
Patient 13 400 × 400 × 280 0.3 525 × 425 × 290 0.15 Yes Yes No† No∗
Patient 15 750 × 750 × 400 0.2 525 × 425 × 290 0.15 Yes Yes No† No∗

Patient 17 500 × 500 × 500 0.2 668 × 530 × 364 0.15 Yes Yes No† No
∗
The scan obtained from 3DSlicer’s “Sample Data” module, titled “CBCT-MRI Head”. † The mandible includes partially erupted wisdom tooth/teeth. The maxilla
includes an impacted wisdom tooth.
5
of the scans with slice thicknesses above 0.3 mm is very challeng-

ing. Besides, applying smoothing functions with identical kernel
sizes results in smoother segments for coarse voxel sizes, making
it difficult to remove the segmentation noise (e.g., rugged surfaces)
while preserving the fine details in desired regions (e.g., the alve-
olar crests of teeth sockets and root apexes). As a result, to avoid
biases in the geometry reconstruction process, we convert all scans
to an identical resolution of 0.15 mm, as summarized in Table 2,
Mandibular body angle: the angle between the lines connecting the left mental foramen (Ml ) to pogonion (Pg) and left gonion, inspired from [27]
Anterior nasal spine angle: the angle between the lines connecting the anterior nasal spine (S) to the left and right maxillary tuberosity (Tl , Tr )
by oversampling a cropped region of interest (ROI) containing the
jaws using linear interpolation.
The angle between the lines connecting the left canine eminence (CEl ) to the anterior nasal spine (S) and the left maxillary tuberosity (Tl )
The tooth-bone segmentation is performed based on the semi-
automatic watershed algorithm [31] provided in 3DSlicer’s Seg-
mentEditorExtraEffects extension. The result of the watershed al-
Pogonion angle: the angle between the lines connecting pogonion (Pg) to the left and right gonia (Gl , Gr ), inspired from [25,29]
gorithm is then refined to correct the misclassified teeth and bone.
Later, the geometries are smoothened using the 3DSlicer’s standard
median and joint smoothing modules [32,33]. The joint smooth-
ing method [33] smooths the adjacent segments simultaneously
Anterior maxillary basal width: the length of the line connecting the left and right canine eminences (C El , C Er ) [27]
Posterior maxillary width; the distance between the left and right maxillary tuberosity (Tl , Tr ), inspired from [27]
and enforces watertight interfaces between them. The tooth-bone
segmentation procedure proceeds until the segmentation accuracy
meets our clinical validation criteria (see Section 3.3). Finally, the
segmented regions are exported as surface meshes in the Object
file format (OBJ).
3.3. Clinical validation of the teeth-bone segmentation
The perpendicular distance between the anterior interdental (Int) and pogonion (Pg) [28]
Description of the utilized morphological landmarks obtained based on the literature [25–29], as illustrated in Figure 3.
The segmented teeth-bone geometries of all patients are val-
The perpendicular distance from the left canine eminence (CEl ) to the alveolar crest
idated by clinical experts. This is done by identifying the exis-
tence of any periodontal diseases and categorizing each patient
into one of the no, mild, or moderate periodontal disease cate-
Intergonial distance: direct distance between the left and right gonia [28]
Body height: the height of the mandible at the mental foramen (Ml ) [26]
gories. The patients with severe periodontal diseases are excluded
from the analyses based on the cohort selection criteria mentioned
Maximum breadth of mandibular body at the mental foramen [25]

in Section 3.1. Moreover, we assess the accuracy of segmentation in
areas close to teeth/bone borders, cervical regions of bone around
teeth sockets, tooth-bone interfaces, and roots.
3.4. Geometrical variations of the dataset
We evaluate the geometrical variation of the reconstructed ge-

ometries using 3D morphological and cephalometric landmarks
and measurements adapted or inspired from the literature [25–29],
as illustrated in Figure 3 and described in Table 3. Accordingly, the
measured values per patient and the mean and standard deviation
of each morphological measurement across all patients are listed
in Table 4. Note that even though the standard deviation of each
Landmark’s description
measurement may seem rather small, such small changes can lead
to significant variations in the overall shape of the jaw, indicating
the availability of high geometrical variations among different jaws
under the study.
3.5. The pipeline for generating FE model of a human jaw
The flexibility of the FE method allows it to use a wide range of

PMG Angle
MM Width
GPG Angle
MF Height
SCT Angle
GG Width
TST Angle
CE Height
Landmark
CC Width
TT Width
PI Height
spatial discretizations [34]. We opt for an unstructured tetrahedral

mesh as can be robustly generated using automatic meshing tools
[18] and can lead to similar accuracy and running time when using
high order elements as structured meshes [35].
Mandible
To eliminate the manual geometry processing, quality mesh-

Jaw type
Maxilla
Table 3
ing, or mesh decimation steps, we propose to directly use the ex-

ported surface meshes as an input to our pipeline. Figure 4 shows
a comparison between the conventional labor-intensive approach
and our method for developing FE models of a human jaw. For the
reproducibility purpose, the pipeline is implemented based on the
free open-access geometry processing library libigl [19] and the
meshing algorithm fTetWild [18]. The pipeline generates confor-
mal volumetric meshes using imperfect meshes exported from the
segmentation software with minimal human intervention. Figure 5
6
Fig. 4. Finite element models created based on irregular dense meshes (left subfigure, A) exported from the segmentation step. Right subfigure: A comparison between
the results of the conventional meshing approach (top-row) and utilized pipeline (bottom-row). The conventional approach involves time-consuming and labor-intensive
geometry and mesh processing tasks (B); this results in non-congruent contacting interfaces (F), and non-conformal meshes (G). Note that the spotted marks in F indicate
the undesired gaps/penetrations in the contacting interfaces. In contrast, the utilized method generates multi-domain volumetric meshes directly using the input irregular
meshes from the generated PDL rims (shown in red in C) and guarantees the interface congruency as well as the mesh conformity as depicted in H.
Fig. 5. Flowchart of the utilized pipeline and the characteristic of the meshes at each step. The pipeline consists of five main consecutive steps, as further described in
Section 3.5.
shows the flowchart of the utilized pipeline and the characteristics ing steps of the used pipeline, i.e., the gap and PDL rim generation
of the meshes at different steps. steps [36]. We use fTetWild [18] as a robust meshing tool to dec-
imate and “clean up” the imperfect meshes. To be more specific,
3.5.1. Preprocessing the teeth and bone geometries are tetrahedralized separately, and
As can be seen in Figure 5, the input meshes to the pipeline the boundary faces of the resulting tetrahedral meshes are then
are dense irregular meshes, which are not necessarily guaranteed extracted as the cleaned-up reduced surface meshes to be used as
to be watertight, manifold, and self-intersection-free, referred to as the input meshes to the gap generation step. Further details on the
“triangle soup” in the computer graphics [22,23,36]. Hence, we ap- surface mesh extraction process can be found in Section 3.5.6.
ply a preprocessing step to both reduce the mesh sizes and pro-
duce meshes that are manifold, watertight, and self-intersection- 3.5.2. Gap generation for the PDL tissue
free. These mesh characteristics assure meaningful values for the PDL has an average width of 0.2 mm [37]; its width can approx-
utilized signed distance functions in the next geometry process- imately be 0.15 mm around the middle third of the root and about
7
Assessment of the geometrical variations of the mandibles and maxillae based on the computed measurements illustrated in Figure 3. Note that even small changes in the values can change the overall shape of the jaw
0.21 mm [38,39] to 0.38 mm [37] near the root apex and cervical
regions.
TT Width (mm)
Reconstructing the PDL layer using CBCT scans obtained in vivo
61.15 ± 6.85
from patients in clinics is a challenging process [40] as the com-
monly used voxel dimensions for the CBCT scans, ranging from 0.2
63.66
52.62
59.87
55.36
54.78
61.01
68.12
61.51
72.85
58.99
71.66
53.33
- mm to 0.5 mm [3], are not fine enough to capture such a thin
-
-
-
structure (roughly 0.2 mm) [41,42]. Although the geometry of the
TST Angle (deg)
PDL layer can be reconstructed by segmenting it from the micro-
69.64 ± 5.43
CTs acquired in vitro, the x-ray exposure in such scans is extremely
high and harmful for the human body, and it is usually obtained
68.40
61.10
67.80
66.70
68.30
69.90
69.90
80.30
77.30
63.70
74.60
67.70
from dead specimens. Therefore, we first conduct the teeth-bone
-
-
-
-
segmentation from the CBCT scans and then apply geometry pro-
SCT Angle (deg)
cessing techniques to create a gap between the teeth and bone ge-
140.58 ± 7.96
ometries where the PDL can reside with an average thickness of
0.2 mm.
141.60
146.90
154.50
143.50
148.60
129.20
133.20
132.40
147.60
131.70
140.20
137.50
Ideally, to generate the required gaps for PDL geometries, we
shrink the bone and teeth, each by 0.1 mm. First, we shrink the
-
-
-
-
bone geometry by 0.1 mm, by moving its mesh points in the re-
CE Height (mm)
14.46 ± 2.10
verse direction of per-vertex normal with a magnitude of 0.1 mm,
which we call it explicit shrinking approach. Before performing any
explicit shrinking process, the radius curvature is locally computed
15.62
12.24
12.88
16.09
13.04
17.08
12.49
17.43
14.70
13.30
17.00
11.70
for bone vertices to identify the sharp and thin structures. This
-
-
-
-
is done by calculating the radius r of the mean curvature h at

CC Width (mm)
point t based on r (t ) = 1/h(t ), as the reciprocal of the curvature

Maxillary jaw
32.44 ± 4.55
at that point. The radius of the curvature provides useful infor-

mation about the maximum magnitude that the surface vertices
40.84
30.21
29.92
27.68
37.13
35.83
33.46
34.65
35.94
27.61
30.09
25.90
can be moved in the opposite direction of the normals of the sur-

-
-
-
-
face before a singularity occurs. This means that the movement

PI Height (mm)
of the nodes with larger values would cause self-intersection is-

21.87 ± 3.58
sues and artifacts on the surface mesh. For this reason, evaluating
the shrinking limits prior to the bone shrinking process is impor-
21.45
26.43
19.94
16.41
19.71
21.70
21.07
21.63
26.53
23.02
27.21
26.21
17.86
17.68
18.83
26.88
19.16
tant. In the case where the shrinking limit is less than 0.1 mm, the
pipeline shrinks the bone to its maximum shrinking limit, while
GG Width (mm)
more shrinking the teeth to compensate for the total desired

87.25 ± 9.79
gap.
To shrink the teeth, we use an implicit shrinking approach based
104.70
on signed distance functions presented in Equation (1). The signed

95.32
83.54
71.61
88.88
94.29
86.19
91.21
91.76
69.80
73.02
78.20
82.10
97.08
89.11
89.42
97.10
distance field contains the boundary of the geometry, i.e., zero iso-
surface [43] and the information from different offset surfaces with
GPG Angle (deg)
positive and negative values. The positive offsets/iso-surfaces rep-

considerably, indicating that high geometrical variations are available within our models.
80.14 ± 7.00
resent a dilated version of the geometry, while the negative val-

ues represent the eroded geometry. Consequently, we use an iso-
77.50
79.00
86.60
81.20
96.50
85.80
71.90
75.30
81.70
77.60
80.50
81.30
77.90
91.90
75.80
71.70
70.20
contour of φ (x ) = −0.1 to shrink the tooth with a magnitude of

0.1 mm. Since the iso-contour is still an implicit representation of
PMG Angle (deg)
the shrunk tooth, we use a contouring method based on a march-

146.91 ± 6.24
ing cubes algorithm [44] to convert it to an explicit representa-

tion including the vertices and connectivity matrix. As a result,
147.20
142.30
131.70
141.40
158.70
151.80
140.20
146.10
142.40
150.80
151.40
154.00
147.00
145.00
147.80
151.40
148.20
wherever the radius curvature might be less than the desired off-
set value, the implicit shrinking approach can prevent any singu-
larities at the root apexes. It can therefore be used as a more ro-
MF Height (mm)
bust shrinking approach, especially for thin and sharp geometries

28.64 ± 5.03
with lower radius curvature values (here, the root apex of incisors).
Note that the applied approach to the teeth creates a slightly wider
28.94
27.31
33.31
32.73
25.57
26.72
27.39
29.48
32.33
33.14
32.36
33.96
19.24
23.02
36.57
24.76
20.01
space than the desired gap in the root apexes, which is in line with
the clinical studies [39,45].
MM Width (mm)
Mandibular jaw
3.5.3. Boundary representation of PDL

45.71 ± 3.73
Instead of an explicit representation of PDL which uses closed

surface meshes, we use a boundary representation (B-Rep) ap-
46.54
46.89
41.61
42.51
38.42
47.18
41.19
49.30
50.31
48.34
51.29
46.34
45.62
49.39
45.15
46.76
40.21
proach using triangle meshes to describe the PDL domain. This

needs to be distinguished from the basis spline (B-Spline) repre-
Mean ± STD
sentation. we use B-Rep for describing the PDL domain for three
Patient’s D
Patient 11
Patient 12
Patient 13
Patient 14
Patient 15
Patient 16
Patient 17
Patient 10
Patient 1
Patient 2
Patient 3
Patient 4
Patient 5
Patient 6
Patient 7
Patient 8
Patient 9
main reasons. First, fTetwild uses the winding number informa-

Table 4
tion together with the input surface mesh to generate a volumetric

mesh with no prior assumptions such as watertight or closed sur-
faces. Therefore, a shell surface can be used as a part of the input
8
Fig. 6. The proposed tetrahedra-filtering method and its close-up view (bottom row) for proper labeling of the PDL from the raw mesh. Thetrahedra around teeth (A) and
bone (B) are used to obtain the conformal tooth-PDL and PDL-bone interfaces (C). The obtained mesh includes jagged top surfaces which can be removed by using the
positive winding numbers with respect to the PDL rims (D), resulting in smooth and clean top surfaces for the generated PDLs (E).
mesh to represent the domain boundary. Second, since fTetwild 3.5.4. Volumetric mesh generation
uses an implicit algorithm, it enables us to use B-Rep and per- A unified volumetric mesh is generated using the surface
form boolean operations on different components to describe the meshes of the teeth, bone, and produced PDL rims, It is ob-
domain boundary [18]. This is while it is not trivial to correctly de- tained by applying a union operation on the provided input sur-
fine the B-Rep models for complex geometries using the Delaunay- face meshes using fTetWild [18] with optimal envelope ( ) and
based algorithms such as TetGen [46]. Third, B-Rep models help us ideal edge length values of 2 × 10−4 and 0.01, respectively. Apply-
to achieve zero-congruent contacting interfaces and avoid numer- ing a union operation using fTetWild generates a single volumetric
ically small values for ε produced due to machine/floating-point mesh for the bounding box surrounding all input surface meshes
precision. This, in turn, assures congruency and mesh conformity called the raw mesh [18]. Note that the default filtering method in
at the contacting interfaces for modeling the PDL layer as men- fTetWild exports a labeled version of the raw mesh that assigns
tioned in Section 2.3. the tetrahedra outside the closed input surfaces as the background
The PDL geometry can be described using three main surfaces: elements.
the top surface of PDL, tooth-PDL, and PDL-bone interfaces, in
which the two last interfaces can be replaced with the teeth and 3.5.5. Proposed tetrahedra filtering method
bone geometries. To represent the free surface of the PDLs that are The default tetrahedra filtering/labeling method in fTetWild can
not in contact with the teeth and bone geometries, we generate only be used when all input meshes are closed surface meshes,
only the top surface of each PDL, called the PDL rim in this study. which is not the case for the PDL rims. Hence, a specific filter-
B-Rep model using PDL rim and teeth-bone geometries as inputs ing approach is used for assigning each tetrahedron a label asso-
to fTetWild helps to generate multi-domain volumetric mesh using ciated with one of the input domains, i.e., the teeth, PDL, or bone.
boolean operations to guarantee geometry congruency and mesh The remained tetrahedra are labeled as background. We use a dis-
conformity in the contacting interfaces. tance field-based labeling method for teeth and bone. To do so,
We utilize a gap filling method [47] to generate the rims that the barycenter of each tetrahedron is used as an average point of
connect the bone to shrunk teeth. This method detects an ini- tetrahedron vertices to decide whether the element is positioned
tial surface on the bone within the desired distance like 0.3 mm. inside or outside a surface considering the distance field of each
Next, the boundary nodes of the detected surface are smoothened barycenter with respect to each of the teeth and bone geometries.
using a Laplacian smoothing approach to provide a base surface The negative distance values indicate the points are located inside
on the tooth socket with a smooth boundary on the bone. The the surface mesh.
smoothened boundary points are then snapped to the bone surface We propose using a signed distance-based method combined
to assure the smoothened boundary is on the bone. Afterward, the with the winding number information to obtain the PDL geome-
base function with smoothened boundary is extruded with the de- tries. The PDL is expected to reside in the gap between the teeth
sired thickness (0.2 mm), and the extruded surface is snapped to and bone; hence, as shown in Figure 6, we first select all tetrahe-
the tooth surface to assure the points are located on the tooth sur- dra with positive distance values below 0.3 mm around teeth (A)
face. Finally, the PDL rim is constructed by connecting the bound- and bone (B) based on the intersection of these two sets of tetra-
ary points of the base and extruded surfaces. We only use the hedra (C). Next, the intersection of the obtained tetrahedra (C) and
PDL rims without the base/extruded surface to avoid any redun- the tetrahedra with positive winding numbers with respect to the
dant representations of tooth-PDL and PDL-bone interfaces having PDL rims (D) is utilized to produce a smooth boundary (E) for the
potential numerical errors. These errors can cause issues in the top-free surface of each PDL. The winding numbers, as also used
volumetric mesh generation process of fTetwild, when the value in fTetWild for implicit meshing, are applied to eliminate all mesh-
of the user-defined maximum deviation parameter called enve- ing artifacts introduced by the intersection operation (C). As can be
lope is less than or equal to the produced error in the duplicated seen, applying the proposed filtering method considering the infor-
surfaces. mation of the winding numbers to the raw mesh provides smooth
9
PDL top surfaces as well as conformal meshes in the tooth-PDL and In the biting scenario, a pressure load of 200 N is applied to the
PDL-bone interfaces, which include shared nodes in the contacting bottom surface of the mandible to simulate the biting force.
interfaces and provide perfect adhesions in the contacting inter- Contact definition To have a perfect adhesion in the tooth-PDL
faces. and PDL-bone interfaces, we generate a single volumetric mesh by
combining surface meshes of the shrunk teeth, PDL rim, and bone
3.5.6. Surface mesh extraction geometries using a union operation in fTetWild [18]. The nodes
After labeling all tetrahedra in the raw mesh, we select bound- at the contacting interfaces are shared between the adjacent do-
ary faces of each domain, i.e., the shared faces between the domains, thus guaranteeing a complete adhesion as well as confor-
main and background tetrahedra, and save them as reduced adap- mal meshes at the contacting interfaces. Therefore, there will be
tive surface meshes. Note that the extracted surfaces are zero- no sliding, separation, or penetration at the tooth-PDL and PDL-
congruent and conformal meshes with no meshing artifacts. bone contacting interfaces. Further details on the volumetric mesh
generation can be found in Section 3.5.4.
3.5.7. FE simulation setup Furthermore, in the both scenarios, any potential contacts be-
The used pipeline automatically sets up the FE problems based tween different teeth are modeled using the incremental potential
on two different scenarios; (1) uncontrolled tipping and (2) biting contact formulation [54].
scenarios. It automatically generates FEM files of the uncontrolled
tipping scenario for all patients and FEM files of the biting scenario 3.5.8. Technical details
for those whose scans were acquired in a natural biting position. Implicit shrinking approach The explicit representation of ge-
Note that the utilized code for defining the boundary and load- ometries only provides boundary information of them. This is
ing conditions is scenario-specific, which requires to be adapted while the implicit representation of geometries based on the
for other scenarios. The defined FE problem for each scenario is signed distance function provides useful information about the in-
exported as an XML file including the tetrahedral meshes, types of terior and exterior of the shape as well as the outline/boundary
elements, utilized material models and properties of each domain, of geometry. Therefore, we use boundary information obtained
as well as boundary and loading conditions. from the explicit representation to define the signed distance func-
In these scenarios, we use Tet4 elements for the teeth and bone tion/filed before evaluating it on the desired points in an R3 space
with negligible deformations. Furthermore, to avoid the locking and shrinking the teeth. To do so, we sample the distance field in
artifacts [48] that are common when linear tetrahedral elements the R3 space using a regular 3D grid in the bounding box of each
are used to model large deformations, we apply a quadratic ba- tooth. A fine grid sampling size of 0.1 mm is used to be smaller
sis to the elements of the PDL layer, to consider Tet10 elements than the minimum isotropic voxel size (0.15 mm) as presented in
for the PDL domain. We prefer to use an explicit tetrahedral mesh Table 2, to avoid the aliasing effect or spiky surfaces in the sam-
to model the PDL layer [16,40,49] (instead of the cheaper option pling process and have smooth geometries in the shrunk teeth.
of using shell elements) as it can more faithfully capture shear, Next, the signed distance values are obtained for each point of
bending, and buckling of the membrane, and at the same time, the grid considering the surface mesh of the tooth. Note that sam-
it simplifies both meshing and simulation, as it does not require pling using regular grids requires excessive memory for geometries
special handling for inserting shell elements in a tetrahedral mesh with large dimensions, which is not the case here as the dimen-
and coupling between the volumetric elastic model and its interac- sion of each tooth is relatively small. In general, for meshes with
tion with the shell. Further, providing a volumetric representation fine details that cover a larger space, it is suggested to use adap-
of the PDL layer allows future studies to explore the effects of uti- tive sampling approaches such as the octree-based methods [55],
lizing visco-elastic-plastic material models for the PDL layer. to densely sample voxels in specific regions near the boundary or
Material models and parameters In this study, the bone and regions with great details.
teeth tissues deform negligibly under the applied forces. There- Alternative PDL modeling An alternative approach for gener-
fore, we assume no distinctions between different structures of the ating FE models of the tooth-supporting complex can be using
bone (cortical and trabecular) and teeth (enamel, dentin, and pulp) Delaunay-based volumetric meshing tools like TetGen. Still, one
[4,50,51]. Besides, the porous fibrous periodontal ligament tissue needs to assure mesh conformity or interface congruency on the
is assumed as a homogeneous structure [49]. We use the Neo- surface meshes, and next, to enforce the explicit volumetric mesh-
Hookean material model with Poisson’s ratios of 0.3, 0.3, and 0.45, ing algorithm to preserve the surface meshes [10]. Note that pro-
and Young’s modulus values of 20 0 0 MPa, 150 0 MPa, and 68.9 MPa viding quality surface meshes with mesh conformity criteria per
to describe the mechanical behavior of the tooth, PDL, and bone se is not trivial. Besides, even though the alternative explicit ap-
domains, respectively [2,4,52]. proach can provide quality surface meshes, it may not provide op-
Note that the used simplex material models can be replaced timal quality for the tetrahedra due to the additional constraint ap-
with any complex constitutive models to properly mimic the plied to preserve the surface mesh on the boundary of the domain.
anisotropic viscoelastic behavior of the PDL [49,51,53], or or- Advantages of implicit over explicit meshing algorithms The
thotropic characteristic of the bone. This can however increase the explicit volumetric mesh generation approaches such as the
computational costs of the simulations. Delaunay-based algorithms with incremental local mesh operations
Boundary conditions In the tipping scenario, a Dirichlet bound- generate tetrahedral meshes covering the domain interior and are
ary condition is defined on all nodes located at the bottom/top highly faithful to the input mesh. This means that the algorithms
surface of the mandible/maxilla to fix displacements of the nodes cannot coarsen the triangles on the input surface meshes. There-
in all three directions. In the biting scenario, the same Dirichlet fore, they cannot produce coarse quality tetrahedra on the bound-
boundary condition is applied on the top nodes of the maxilla to ary faces of the input mesh, indicating that the element size at
fix them in all three directions. Besides, a similar Dirichlet bound- the domain boundary depends on the dense and irregular/rugged
ary condition is used to restrict the movement of the mandible surfaces exported from the segmentation step. Consequently, one
only in the anterior-posterior and medial-lateral directions with needs to prepare an adaptive surface mesh with high-quality tri-
imposing no restrictions in the third direction. angles before applying the Delaunay-based algorithms along with
Loading conditions In the tipping scenario, we apply a perpen- a sizing field function to create adaptive volumetric mesh both in
dicular force with a magnitude of 1 N at the center of each crown the domain interior and its boundary. Furthermore, since the ex-
to mimic the uncontrolled tipping motion in the lingual direction. plicit meshing algorithms fail in generating volumetric meshes in
10
Fig. 7. Comparison of the quality histograms of four different quality measurements used for the evaluation of the quality of the tetrahedra generated by OpenMandible and
Open-Full-Jaw. Left to right: The radius-edge ratio Qrl ; The volume-edge ratio Qvl ; The radius ratio Qrr ; The angle measurement Qθ . Note that Lmax , Rin , and Rout denote the
maximum edge length, the radius of the insphere, and radius of the circumsphere of a tetrahedron, respectively, and V and Smin represent the volume and minimum face
area of the tetrahedron.
the presence of any self-intersections, one needs to clean up the generated PDLs across different patients models by computing the
meshes before applying the Delaunay-based meshing algorithms to distances between the points located on the outer and inner sur-
create quality and adaptive mesh. faces of the PDL geometries. The thickness details of the generated
In contrast, implicit volumetric meshing approaches like PDLs can be seen in Table 5. As can be noticed from the table, the
fTetWild impose no assumptions on the input mesh and can han- computed minimum, maximum, and average PDL thickness values
dle imperfect input meshes with self-intersection artifacts. They are in line with those reported in the literature, i.e., 0.15, 0.3, and
can also produce an adaptive mesh that provides coarse tetrahe- 0.2mm, respectively [37].
dra at the domain boundary, slightly deviating from the input sur- Likewise, we investigate the thicknesses for the PDLs from the
face mesh, based on a user-defined input parameter. In addition, OpenMandible dataset. As it can be seen in the last row of Table 5,
fTetWild uses winding numbers along with surface meshes to gen- the obtained thicknesses for OpenMandible are about two-three
erate tetrahedral meshes based on an implicit meshing approach. times the values reported in the literature. As the PDL thickness
This enables us to use the constructive solid geometry (CSG) model plays an essential role in exerting the applied load from teeth
for describing a domain with complex boundaries by combining surfaces to the adjacent bone and teeth movements, it can be
several simpler domains using boolean operations [56]. This is deduced that under identical scenarios and loading systems the
while the explicit meshing tools cannot support CSG models [46]. OpenMandible model cannot result in stress/strain values compa-
Volumetric meshing The fTetWild algorithm uses the user- rable with the ones obtained in this study.
defined epsilon and ideal edge length parameters to control the
output mesh’s accuracy and size. The epsilon value indicates how 4.2. Mesh properties
much fTetWild can deviate from the input surface mesh for gen-
erating a quality volumetric mesh. Hence, using a smaller ep- We evaluate the total number of elements and mesh qualities
silon value preserves the input geometry’s details at the cost of of the developed models and compare the mesh qualities of our
higher computational time. In contrast, larger values can provide model with those of the OpenMandible dataset. The maximum
less accurate meshes in less computational time. Moreover, smaller mesh deviation between the input and output meshes of ftetwild
values of ideal edge length provide denser meshes, while larger is then calculated to examine whether the meshing process alters
values can produce coarser elements. Hence, these two parame- the geometries negligibly.
ters can control the output mesh size and its geometrical accu-
racy. In this work, the utilized values for epsilon and ideal edge 4.2.1. Mesh sizes
length are obtained based on a grid-search method in the intervals We develop all computational meshes of our study on a ma-
[10−6 , 10−3 ] and [0.005, 0.02], respectively. Note that the fTetWild chine with a 2.60 GHz processor and 16 GB of RAM, which takes
multiplies these parameters by the length of the diagonal of the around 40.78 ± 12.19 minutes per jaw. The details of the origi-
input meshes’ bounding box. Hence, the parameters are sensitive nal dense meshes, the generated volumetric meshes, and their ex-
to the size of the bounding box that includes all input meshes. In tracted surface meshes are summarized in Table 6. The pipeline
other words, they need to be adjusted according to the model’s produces conformal adaptive meshes from irregular dense meshes.
size; if, for instance, the model is a cut-out model with a much The utilized adaptive meshing approach helps to reduce the total
smaller bounding box dimension than the full jaw model. number of the elements, by producing coarser elements in regions
far from the teeth while generating finer elements on root apexes,
4. Results and discussion PDL layer, alveolar crest, or regions with fine details like thin walls
of the maxilla.
In this section, the generated models are assessed from differ-
ent aspects, and the results are compared with the state-of-the-art. 4.2.2. Mesh quality analysis
In Section 4.1, we measure the thickness of the generated PDL lay- The mesh quality required for an FE analysis can vary depend-
ers. In Section 4.2, different mesh properties such as mesh sizes ing on the application and utilized numerical methods [57]. In gen-
and mesh qualities are evaluated for the developed models; the re- eral, a regular tetrahedron has the highest mesh quality in com-
sults are compared to the OpenMandible model. In Section 4.3, FE putational models, and the main guideline is to avoid using low-
simulation results of different patients are presented under identi- quality/badly-shaped tetrahedra, as they can affect the accuracy of
cal tipping and biting scenarios. the numerical methods.
The OpenMandible uses TetGen to obtain volumetric meshes
4.1. Generated PDL properties from the manually generated conformal surface meshes. It sets the
upper limit of the radius-edge ratio of to-be-generated tetrahedra
We assume an average thickness of 0.2 mm for the PDL layers to 1.5. This mesh quality constraint controls the ratio between the
according to the literature [7,9] and assess the thicknesses of the radius of the circumscribed sphere and the shortest edge of each
11
Fig. 8. The obtained displacement (A) and stress (B) fields for Patient 17 under the biting scenario. Note that due to the deformation of the PDL layer caused by the exerted
biting force from the contacting mandibular teeth, the displacement field is seen on the posterior maxillary teeth (A). Likewise, stress concentrations can be seen on the
associated PDL layers, especially at the apexes and furcation regions, indicated by arrows in B.
Fig. 9. The obtained displacement (A) and stress (B) fields for the mandibular jaw of Patient 17 under the tipping scenario. Note that stress concentration is observed on
the lingual side of anterior teeth due to the labiolingual tipping load.
Table 5
The thickness of the PDL layers generated using the utilized pipeline with the statistics in line with the literature [37]. Note that there is a significant
difference between our model results and those of the OpenMandible.
Patient’s ID Mandibular jaw Maxillary jaw
Minimum Maximum Average Minimum Maximum Average
Patient 1 0.12316 0.23552 0.19530 - - -

Patient 2 0.14227 0.25042 0.19536 - - -
Patient 3 0.13779 0.33969 0.19379 0.08258 0.38502 0.20143
Patient 4 0.13779 0.21958 0.11962 0.14331 0.23754 0.19634
Patient 5 0.14537 0.24688 0.19770 0.11837 0.23631 0.19676
Patient 6 0.16350 0.24151 0.19737 0.14928 0.25741 0.19654
Patient 7 0.15258 0.25176 0.19527 - - -
Patient 8 0.13036 0.25292 0.20056 0.15210 0.32516 0.19745
Patient 9 0.13324 0.27437 0.19629 - - -
Patient 10 0.11951 0.23903 0.20100 - - -
Patient 11 0.14271 0.33307 0.19887 0.12834 0.25599 0.19465
Patient 12 0.11694 0.34803 0.20164 0.12070 0.31058 0.19219
Patient 13 0.12288 0.28818 0.20134 0.09648 0.26672 0.18296
Patient 14 0.12056 0.29998 0.19728 0.10014 0.29037 0.19222
Patient 15 0.13099 0.35793 0.20078 0.13277 0.34097 0.19148
Patient 16 0.12508 0.25914 0.19961 0.14156 0.30228 0.19535
Patient 17 0.13754 0.26163 0.19664 0.12427 0.28902 0.19610
Mean ± STD 0.13425 ± 0.0126 0.27645 ± 0.04360 0.19344 ± 0.01918 0.12416 ± 0.02193 0.29145 ± 0.044527 0.19446 ± 0.0045
OpenMandible 0.20576 1.93415 0.60772 - - -
12
T. Gholamalizadeh, F. Moshfeghifar, Z. Ferguson et al.
13
Computer Methods and Programs in Biomedicine 224 (2022) 107009

Fig. 10. A comparison of stress distributions of all utilized jaws in the tipping scenario with identical load magnitudes. The stress fields are clipped within the range [0, 1] MPa to represent the changes better using the same
color map. As can be noticed, the stress values and concentrations considerably vary across different patients, which indicates the importance of utilizing population models for multi-patient analysis and generalizability in
computational modeling studies.
tetrahedron, to prevent the production of low-quality/badly-shaped

tetrahedra. In addition to this quality constraint, OpenMandible en-
The total number of elements (×103 ) of the surface meshes (triangles) and volumetric meshes (tetrahedra). The number of elements is significantly reduced in the adaptive output surface meshes compared to that of
Total
forces TetGen to preserve the provided input surface meshes to
853
783
886
791
799
689
776
880
607
939
786
787
have conformal volumetric meshes. Preserving the surface mesh is
-
-
-
-
the only approach to produce conformal volumetric meshes when
Maxilla
Output volumetric mesh
using explicit volumetric meshing approaches. This raises the ques-
421
496
416
399
548
394
519
467
378
440
502
469
- tion of whether TetGen can achieve the specified quality constraint
-
-
-
value while enforcing another restriction to preserve the input sur-
PDLs
112
134
134
140
101
124
159
135
110
face mesh.
87
72
72
-
-
-
-
We quantitatively assess the quality of the generated volumetric
meshes of this study and those of the OpenMandible. To be more
Teeth
293
231
249
221
156
141
261
250
260
203
209
273
specific, four different quality measurements presented in [57], i.e.,
-
-
-
-
the radius-edge ratio Qrl [58], the volume-edge ratio Qvl [59,60],
the radius ratio Qrr [61,62], and the angle measurement Qθ [61],
are used to evaluate the quality of the volumetric meshes. For a
Maxilla
fair comparison, we compare the results from one of our patients

132
153
159
173
157
135
133
191
139
175
166
126
Output surface mesh
-
-
-
-
to those of the OpenMandible study, as shown in Figure 7.

In a regular tetrahedron, the values of each utilized quality
PDLs
108
96
77
93
92
85
61
77
93
70
50
50
measurement correspond to one, indicating that an ideal high-

-
-
-
-
quality mesh is expected to have a histogram peak at one. There-

Teeth
fore, the resulting distributions (normal or skewed normal) of all

100
78
84
89
75
91
67
52
67
91
50
90
Maxillary jaw’s mesh sizes
quality histograms show that our generated volumetric mesh has

-
-
-
-
higher quality elements compared to OpenMandible. The proposed

model also has narrower distributions with small discrepancies
Maxilla
Input mesh size
around their means. This holds even in the last case (Qθ ), where
1311
1894
2039
2181
1921
1531
1016
2221
2958
636
852
577
one of the peaks of the distribution (mode) for OpenMandible is

-
-
-
-
closer to one. Moreover, the OpenMandible model sees two or

Teeth
1032
more peaks in its distributions that can be modeled by mixed nor-

948
824
851
904
772
891
232
189
851
447
190
mal distributions.
-
-
-
-
4.2.3. Mesh deviation analysis

Total
517
669
788
681
664
448
459
777
430
502
681
606
497
511
487
605
630
If one does not enforce the surface preservation criterion, the

Delaunay-based volumetric meshing algorithms like TetGen include
all points of the input surface mesh and a number of additional
Mandible
Output volumetric mesh
points. Hence, the generated volumetric meshes have boundary el-

252
224
315
338
326
312
236
274
341
301
264
261
252
238
371
311
305
ements as fine as the input surface mesh. In contrast, fTetwild

provides volumetric meshes as coarse as possible on the surface
PDLs
while preserving the input geometry. This algorithm slightly devi-

111
155
127
119
149
116
110
103
87
65
71
78
74
87
72
71
99
ates from the input surface according to the user-defined maximal

deviation (envelope) value. Therefore, to generate an accurate volu-
Teeth
178
141
243
294
227
233
142
150
230
202
163
163
258
194
160
150
210
metric mesh, we must not deviate too much from the input surface
mesh.
Accordingly, we quantitatively evaluate the deviation of final
extracted surface meshes from the input surface meshes. We use
Mandible
Hausdorff distance as an error measurement between the in-

106
108
put and extracted surface meshes from the generated volumetric

Output surface mesh
75
96
98
91
73
83
99
91
79
82
78
72
97
94
70
meshes. The obtained distance values are presented in Table 7.

As it can be seen, the maximum deviation of the meshes both in
PDLs
108
101
45
77
88
82
49
54
76
71
51
49
69
60
60
50
80
mandibular and maxillary jaws and their associated teeth is negli-

gible with respect to the dimension of the entire jaw models.
Teeth
59
82
94
77
78
51
74
69
53
57
54
85
61
69
50
50
50
Mandibular jaw’s mesh sizes
4.3. FEM verification

the irregular dense input surface meshes.
As the last part of our analyses, we assess the displacement

Input irregular mesh
Mandible
and stress fields in the tipping and biting scenarios to ensure that
the stress patterns are smooth and have no unrealistic stress con-
1574
1576
1592
1553
1406
1556
1867
1969
2958
2110
2201
1045
2000
478
629
469
486
centrations. To do so, we run the simulations using the PolyFEM

[63] FE solver. PolyFEM is an FE simulation toolkit that supports
Teeth
1016
elastodynamic deformations with linear and non-linear material

198
531
768
731
526
625
869
791
708
233
188
183
877
900
770
440
models. It provides an adaptive p-refinement that allows increasing

the order of basis functions for specific domains while utilizing lin-
Patient’s ID
11
12
13
14
15
16
17
10
1
2
3
4
5
6
7
8
9
ear basis functions for the other domains. Hence, we use Tet10 el-
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
Patient
ements for the PDL layer to increase the simulations’ accuracy and
Table 6
avoid element locking issues [64]. Besides, PolyFEM uses the incre-
mental potential contact formulation [54] for contact response and
friction, which ensures valid, penetration-free meshes during the
14
Table 7
Hausdorff distances (HD) used as an error measurement to assess the deviation of the sur-
face in the output mesh from the input mesh. The HD values indicate that the output surface
meshes generated using our pipeline are close to the input surface meshes.
Patient’s ID Mandibular Jaw HD (mm) Maxillary Jaw HD (mm)
Mandible Mandibular Teeth Maxilla Maxillary Teeth
Patient 1 0.173 0.038 - -

Patient 2 0.147 0.039 - -
Patient 3 0.125 0.044 0.143 0.019
Patient 4 0.086 0.043 0.137 0.029
Patient 5 0.139 0.025 0.156 0.018
Patient 6 0.241 0.023 0.232 0.042
Patient 7 0.154 0.042 - -
Patient 8 0.454 0.029 0.135 0.039
Patient 9 0.174 0.043 - -
Patient 10 0.453 0.021 - -
Patient 11 0.152 0.027 0.298 0.023
Patient 12 0.162 0.024 0.212 0.044
Patient 13 0.158 0.022 0.320 0.024
Patient 14 0.211 0.026 0.349 0.042
Patient 15 0.205 0.033 0.231 0.020
Patient 16 0.189 0.046 0.404 0.032
Patient 17 0.180 0.028 0.257 0.016
Mean ± STD 0.200 ± 0.102 0.032 ± 0.009 0.240 ± 0.089 0.029 ± 0.011
entire simulation. The contacts are automatically detected by prox- Obtaining accurate spatial registration and high-quality volu-
imity; hence, there is no need to specify contact surfaces, which metric meshes for the human jaw can be challenging, as there are
significantly simplifies the scene setup. large variations among geometries of different patients (Figure 2),
The simulation results are visualized in ParaView [65]. such as geometrical differences in the bones and teeth, missing
Figure 8 shows the result of the biting scenario of a selected pa- teeth, and topological changes in the number of roots, e.g., the
tient, including the resulting displacement and stress fields. Due to mandibular and maxillary molars with two to four roots. There-
the deformation of the PDL layer under the applied biting load, the fore, to include different types of variations in the data for a plau-
displacement fields can be seen both on the mandibular jaw and sible deformation, one needs to have different templates covering
the posterior maxillary teeth (A). Besides, stress concentrations can missing teeth or various numbers of roots. This in turn is a time-
be observed on the associated PDL layers, especially at root apexes consuming process and can increase the complexity of the model.
and bifurcation regions indicated by arrows in (B). On the other hand, large deformations in small volumes of teeth
Figure 9 illustrates the result of the auto-generated uncon- and roots can result in distorted elements, preventing us from gen-
trolled tipping scenario for the mandibular jaw of the same pa- erating high-quality meshes, especially in the PDL layers with thin
tient. As seen in [4], the anterior teeth have higher displacement structures that need to be modeled with fine volumetric elements.
fields than the posterior teeth under an identical load (A). Hence, The current study introduces the largest-ever dataset of patient-
higher stress concentrations can be seen on the lingual side of the specific human jaws reconstructed from CBCT scans. We believe
PDLs associated with the anterior teeth (B). Additionally, simula- this unique clinically validated dataset would pave the way for fu-
tion results of different patients in the tipping scenario are shown ture population studies in the field. More specifically, data aug-
in Figure 10. As can be noticed, the stress values considerably mentation techniques using machine learning [69,70] can be ap-
change from one patient to another, indicating the importance of plied to the Open-Full-Jaw dataset to expand its size and variabil-
utilizing population models for multi-patient analysis. ity by generating plausible synthetic data. In addition, this would
It should be noted that although we have tested our models un- enable us to use deep learning methods, which require a large
der two scenarios, the developed volumetric meshes can be used amount of data for training. Still, one needs to use a dataset with
in various scenarios and different FE frameworks. Furthermore, the enough variations for sampling and assess the generated samples’
FE models can benefit from using more complex material mod- validity.
els, boundary, and loading conditions. For example, the provided
meshes can be integrated with outputs of other studies [9,10] to 6. Conclusion
consider masticatory muscles for more realistic biting scenarios.
In this work, we presented a large open-access dataset,
5. General discussion called Open-Full-Jaw (https://github.com/diku- dk/Open- Full- Jaw),
with patient-specific models of 17 human mandibles and maxil-
The utilized and conventional meshing approaches generate the lae. The dataset contains clinically validated segmented geometries
volumetric meshes using reconstructed geometries based on accu- shared as dense surface meshes and adaptive quality volumetric
rately segmented scans. However, obtaining such an accurate seg- with conformal meshes in the contacting interfaces. It also includes
mentation is inherently time-consuming and labor-intensive and, the principal axes for each patient’s teeth and the generated FEM
in some cases, could be highly challenging due to the complexity files of the uncontrolled tipping and biting scenarios for all pa-
of the problem or lack of high-resolution scans [66]. The template- tients. Finally, we share the nearly-automated pipeline used for ge-
based deformation techniques [67,68] can be used to automatically ometry processing, re-meshing, and generating volumetric meshes.
reconstruct the 3D geometries by creating a template mesh and In addition, we evaluated the generated models and quantified
deforming it according to the new data samples. Still, deforming them in terms of the mesh quality and accuracy of the models,
a template mesh using registration approaches or deep learning and compared the results with the state-of-the-art. The obtained
methods requires accurate spatial registration and high-quality vol- results indicate that the developed computational models are pre-
umetric meshes with no distorted elements. cise, considering the low error/distance from input surface meshes.
15
Moreover, the quality of the volumetric elements evaluated based [3] P.M. Cattaneo, M.A. Cornelis, Orthodontic tooth movement studied by finite el-
on different quality measurements imply that the generated volu- ement analysis: an update. what can we learn from these simulations? Current
Osteoporosis Reports (2021) 1–7.
metric meshes consist of quality elements suitable for the FEM of [4] T. Gholamalizadeh, S. Darkner, P.M. Cattaneo, P. Søndergaard, K. Erleben,
the human jaw. Hence, we believe the Open-Full-Jaw dataset can Mandibular teeth movement variations in tipping scenario: A finite ele-
be used in various FE scenarios and a wide range of intra- and ment study on several patients, in: Computational Biomechanics for Medicine,
Springer, 2021, pp. 31–43.
inter-patient analyses. [5] T. Gholamalizadeh, S. Darkner, P.L. Søndergaard, K. Erleben, A multi-patient
The shared repository includes all detailed information for analysis of the center of rotation trajectories using finite element models of
reproducing the models of this study. In addition, the utilized the human mandible, Plos one 16 (11) (2021) e0259794.
[6] R. Savignano, R.F. Viecilli, U. Oyoyo, Three-dimensional nonlinear prediction of
pipeline allows other researchers in the field to generate qual-
tooth movement from the force system and root morphology, The Angle Or-
ity volumetric meshes and FE model files directly using dense thodontist 90 (6) (2020) 811–822.
and irregular meshes with minimal human intervention. This will [7] J.-H. Seo, E. Eghan-Acquah, M.-S. Kim, J.-H. Lee, Y.-H. Jeong, T.-G. Jung,
M. Hong, W.-H. Kim, B. Kim, S.J. Lee, Comparative analysis of stress in the pe-
help other researchers easily extend their datasets without spend-
riodontal ligament and center of rotation in the tooth after orthodontic treat-
ing much time and effort on manually cleaning up the meshes ment depending on clear aligner thickness finite element analysis study, Ma-
and non-trivially producing conformal meshes. Furthermore, sim- terials 14 (2) (2021) 324.
ilar concepts as those used in this study to generate population [8] X. Ding, S.-H. Liao, X.-H. Zhu, H.M. Wang, Influence of orthotropy on biome-
chanics of peri-implant bone in complete mandible model with full dentition,
models of the tooth-supporting complex can be adapted to other BioMed research international 2014 (2014).
areas, such as pelvic girdles and hip joints [71]. [9] J. Ortún-Terrazas, J. Cegoñino, A.P.D. Palomar, In silico study of cus-
pid’periodontal ligament damage under parafunctional and traumatic con-
ditions of whole-mouth occlusions. a patient-specific evaluation, Computer
Funding methods and programs in biomedicine 184 (2020) 105107.
[10] A.M. Vukicevic, K. Zelic, D. Milasinovic, A. Sarrami-Foroushani, G. Jovicic,
P. Milovanovic, M. Djuric, N. Filipovic, A.F. Frangi, Openmandible: An open–
This project has received funding from the European Unions
source framework for highly realistic numerical modelling of lower mandible
Horizon 2020 research and innovation program under the Marie physiology, Dental Materials 37 (4) (2021) 612–624.
Sklodowska-Curie grant agreement No. 764644. This paper only [11] A. Boryor, A. Hohmann, M. Geiger, U. Wolfram, C. Sander, F.G. Sander, A down-
loadable meshed human canine tooth model with pdl and bone for finite ele-
contains the authors views, and the Research Executive Agency and
ment simulations, Dental Materials 25 (9) (2009) e57–e62.
the Commission are not responsible for any use that may be made [12] J. Kawamura, J.H. Park, Y. Kojima, Y.-A. Kook, H.-M. Kyung, J.M. Chae, Biome-
of the information it contains. 3Shape A/S provided financial sup- chanical analysis for total mesialization of the mandibular dentition: A finite
port in the form of salaries for authors TG and PS. The funders had element study, Orthodontics & craniofacial research 22 (4) (2019) 329–336.
[13] J. Kawamura, J.H. Park, N. Tamaya, J.-H. Oh, J.M. Chae, Biomechanical analysis
no role in study design, data collection/analysis, publication deci- of the maxillary molar intrusion: A finite element study, American Journal of
sion, or manuscript preparation. Orthodontics and Dentofacial Orthopedics (2022).
This work was also partially supported by the NSF CAREER [14] A.C. Oenning, A.R. Freire, A.C. Rossi, F.B. Prado, P.H.F. Caria, L. Correr-Sobrinho,
F. Haiter-Neto, Resorptive potential of impacted mandibular third molars: 3d
award under Grant No. 1652515, the NSF grants OAC-1835712, OIA- simulation by finite element analysis, Clinical oral investigations 22 (9) (2018)
1937043, CHS-1908767, CHS-1901091, NSERC DGECR-2021-00461 3195–3203.
and RGPIN-2021-03707, a Sloan Fellowship, a gift from Adobe Re- [15] B. Sarrafpour, M. Swain, Q. Li, H. Zoellner, Tooth eruption results from bone
remodelling driven by bite forces sensed by soft tissue dental follicles: a finite
search and a gift from Advanced Micro Devices, Inc. element analysis, PLoS One 8 (3) (2013) e58803.
[16] J.-S. Lee, H.-I. Choi, H. Lee, S.-J. Ahn, G. Noh, Biomechanical effect of mandibu-
lar advancement device with different protrusion positions for treatment of
Ethical approval
obstructive sleep apnoea on tooth and facial bone: A finite element study, Jour-
nal of oral rehabilitation 45 (12) (2018) 948–958.
All data used in this paper was provided, already anonymized, [17] Openjaw dataset, 2021, https://erda.ku.dk/archives/
by 3Shape A/S. The data was originally acquired for diagnostic pur- 97cd65fe80e83356f618bb9fbc7d5980/published-archive.html.
[18] Y. Hu, T. Schneider, B. Wang, D. Zorin, D. Panozzo, Fast tetrahedral meshing in
poses unrelated to this study. No other aspect of this work trig- the wild, ACM Trans. Graph. 39 (4) (2020).
gered ethical issues. [19] A. Jacobson, D. Panozzo, et al., libigl: A simple C++ geometry processing library,
2018, https://libigl.github.io/.
[20] L. Ernstbrunner, J. Werthel, T. Hatta, A. Thoreson, H. Resch, K. An, P. Moroder,
Declaration of Competing Interest Biomechanical analysis of the effect of congruence, depth and radius on the
stability ratio of a simplistic ball-and-socket joint model, Bone & joint research
5 (10) (2016) 453–460.
Authors TG and PS are affiliated with 3Shape A/S, which pro-
[21] M. Beek, J. Koolstra, L.V. Ruijven, T.V. Eijden, Three-dimensional finite element
vided support in the form of salaries and data for this study. Their analysis of the human temporomandibular joint dis, Journal of biomechanics
commercial affiliation does not alter our adherence to sharing data 33 (3) (20 0 0) 307–316.
[22] Y. Hu, Q. Zhou, X. Gao, A. Jacobson, D. Zorin, D. Panozzo, Tetrahedral meshing
and materials. Additionally, there are no patents, products in de-
in the wild, ACM Trans. Graph. 37 (4) (2018) 60–61.
velopment, or marketed products associated with this research to [23] Y. Hu, T. Schneider, B. Wang, D. Zorin, D. Panozzo, Fast tetrahedral meshing in
declare. the wild, ACM Trans. Graph. 39 (4) (2020) 117–121.
[24] A. Jacobson, L. Kavan, O. Sorkine, Robust inside-outside segmentation using
generalized winding numbers, ACM Trans. Graph. 32 (4) (2013).
Acknowledgements [25] J.E. Buikstra, Standards for data collection from human skeletal remains,
Arkansas archaeological survey research series 44 (1994).
[26] T.S. Tunis, R. Sarig, H. Cohen, B. Medlej, N. Peled, H. May, Sex estimation using
We thank NYU IT High-Performance. We also thank 3Shape A/S computed tomography of the mandible, International journal of legal medicine
for providing this study’s CBCT scans and, especially, the Dental 131 (6) (2017) 1691–1700.
CAD AI team for their support in the CBCT segmentation and com- [27] M. Bayome, J.H. Park, Y.A. Kook, New three-dimensional cephalometric anal-
yses among adults with a skeletal class i pattern and normal occlusion, The
putation of teeth axes.
Korean Journal of Orthodontics 43 (2) (2013) 62–73.
[28] R. Vallabh, J. Zhang, J. Fernandez, G. Dimitroulis, D.C. Ackland, The morphol-
ogy of the human mandible: A computational modelling study, Biomechanics
References
& Modeling in Mechanobiology 19 (4) (2020).
[29] E.A. Villanueva, A.M. Garmendia, G. Torres, G. Sánchez-Mejorada, J.A. Gómez–
[1] P. Ausiello, A. Gloria, S. Maietta, D.C. Watts, M. Martorelli, Stress distributions Valdés, Gender assessment using the mandible in the mexican population,
for hybrid composite endodontic post designs with and without a ferrule: Fea Spanish Journal of Legal Medicine 43 (4) (2017) 146–154.
study, Polymers 12 (8) (2020) 1836. [30] A. Fedorov, R. Beichel, J. Kalpathy-Cramer, J. Finet, J.-C. Fillion-Robin, S. Pu-
[2] A. Benaissa, A. Merdji, M.Z. Bendjaballah, P. Ngan, O.M. Mukdadi, Stress influ- jol, C. Bauer, D. Jennings, F. Fennessy, M. Sonka, et al., 3d slicer as an image
ence on orthodontic system components under simulated treatment loadings, computing platform for the quantitative imaging network, Magnetic resonance
Computer Methods and Programs in Biomedicine 195 (2020) 105569. imaging 30 (9) (2012) 1323–1341.
16
[31] P. Soille, Morphological image analysis: principles and applications, Springer [53] T.S. Fill, R.W. Toogood, P.W. Major, J.P. Carey, Analytically determined mechan-
Science & Business Media, 2013. ical properties of, and models for the periodontal ligament: critical review of
[32] C. Pinter, A. Lasso, G. Fichtinger, Polymorph segmentation representation for literature, Journal of Biomechanics 45 (1) (2012) 9–16.
medical image computing, Computer methods and programs in biomedicine [54] M. Li, Z. Ferguson, T. Schneider, T. Langlois, D. Zorin, D. Panozzo, C. Jiang,
171 (2019) 19–26. D.M. Kaufman, Incremental potential contact: Intersection- and inversion-free
[33] G. Taubin, T. Zhang, G. Golub, Optimal surface smoothing as filter design, in: large deformation dynamics, ACM Trans. Graph. (SIGGRAPH) 39 (4) (2020).
European Conference on Computer Vision, Springer, 1996, pp. 283–292. [55] D. Meagher, Geometric modeling using octree encoding, Computer graphics
[34] D.N. Arnold, A. Logg, M. Schlager, H. Narayanan, FEMTable, 2022, https:// and image processing 19 (2) (1982) 129–147.
www-users.cse.umn.edu/∼arnold/femtable/. [56] J.D. Foley, F.D. Van, A.V. Dam, S.K. Feiner, J.F. Hughes, J. Hughes, Computer
[35] T. Schneider, Y. Hu, X. Gao, J. Dumas, D. Zorin, D. Panozzo, A large scale com- graphics: principles and practice, Vol. 12110, Addison-Wesley Professional,
parison of tetrahedral and hexahedral elements for finite element analysis, 1996.
2019, arXiv:1903.09332 [57] J. Shewchuk, What is a good linear finite element? interpolation, conditioning,
[36] H. Xu, J. Barbic, Signed distance fields for polygon soup meshes, in: Graphics anisotropy, and quality measures (preprint), University of California at Berkeley
Interface, 2014, pp. 35–41. 73 (2002) 137.
[37] Y. Li, L.A. Jacox, S.H. Little, C.C. Ko, Orthodontic tooth movement: The biol- [58] T.J. Baker, Element quality in tetrahedral meshes, in: 7th International Confer-
ogy and clinical implications, The Kaohsiung journal of medical sciences 34 ence on Finite Element Methods in Flow Problems, 1989, p. 1018.
(4) (2018) 207–214. [59] A. Liu, B. Joe, Relationship between tetrahedron shape measures, BIT Numerical
[38] M. Baron, M. Hudson, M. Dagenais, D. Macdonald, G. Gyger, T.E. Sayegh, J. Pope, Mathematics 34 (2) (1994) 268–287.
A. Fontaine, A. Masetto, D. Matthews, et al., Relationship between disease char- [60] M.K. Misztal, K. Erleben, A. Bargteil, J. Fursund, B.B. Christensen,
acteristics and oral radiologic findings in systemic sclerosis: results from a J.A. Bærentzen, R. Bridson, Multiphase flow of immiscible fluids on un-
canadian oral health study, Arthritis care & research 68 (5) (2016) 673–680. structured moving meshes, IEEE Transactions on Visualization and Computer
[39] S.C. White, M.J. Pharoah, Oral radiology-E-Book: Principles and interpretation, Graphics 20 (1) (2013) 4–16.
Elsevier Health Sciences, 2014. [61] L.A. Freitag, C. Ollivier-Gooch, Tetrahedral mesh improvement using swapping
[40] A. Hohmann, C. Kober, P. Young, C. Dorow, M. Geiger, A. Boryor, F.M. Sander, and smoothing, International Journal for Numerical Methods in Engineering 40
C. Sander, F.G. Sander, Influence of different modeling strategies for the pe- (21) (1997) 3979–4002.
riodontal ligament on finite element simulation results, American Journal of [62] J.C. Caendish, D.A. Field, W.H. Frey, An approach to automatic three-dimen-
Orthodontics and Dentofacial Orthopedics 139 (6) (2011) 775–783. sional finite element mesh generation, International Journal for Numerical
[41] C. Dorow, J. Schneider, F.G. Sander, Finite element simulation of in vivo tooth Methods in Engineering 21 (2) (1985) 329–347.
mobility in comparison with experimental results, Journal of Mechanics in [63] T. Schneider, J. Dumas, X. Gao, D. Zorin, D. Panozzo, Polyfem, 2019, https://
Medicine and Biology 3 (01) (2003) 79–94. polyfem.github.io/.
[42] R. Savignano, R. Valentino, A. Razionale, A. Michelotti, S. Barone, V. D’Antò, [64] T. Schneider, Y. Hu, J. Dumas, X. Gao, D. Panozzo, D. Zorin, Decoupling simula-
Biomechanical effects of different auxiliary-aligner designs for the extrusion tion accuracy from mesh quality, ACM Transactions on Graphics 37 (6) (2018).
of an upper central incisor: a finite element analysis, Journal of Healthcare [65] J. Ahrens, B. Geveci, C. Law, Paraview: An end-user tool for large data visual-
Engineering 2019 (2019). ization, The visualization handbook 717 (8) (2005).
[43] S. Osher, J.A. Sethian, Fronts propagating with curvature-dependent speed: [66] Z. Zheng, H. Yan, F.C. Setzer, K.J. Shi, M. Mupparapu, J. Li, Anatomically con-
Algorithms based on hamilton-jacobi formulations, Journal of computational strained deep learning for automating dental cbct segmentation and lesion
physics 79 (1) (1988) 12–49. detection, IEEE Transactions on Automation Science and Engineering 18 (2)
[44] W.E. Lorensen, H.E. Cline, Marching cubes: A high resolution 3d surface con- (2020) 603–614.
struction algorithm, ACM siggraph computer graphics 21 (4) (1987) 163–169. [67] L. Liang, F. Kong, C. Martin, T. Pham, Q. Wang, J. Duncan, W. Sun, Machine
[45] H. Mortazavi, M. Baharvand, Review of common conditions associated with learning based 3D geometry reconstruction and modeling of aortic valve de-
periodontal ligament widening, Imaging science in dentistry 46 (4) (2016) formation using 3D computed tomography images, International journal for
229–237. numerical methods in biomedical engineering 33 (5) (2017) e2827.
[46] H. Si, Tetgen, a delaunay-based quality tetrahedral mesh generator, ACM Trans- [68] D.H. Pak, M. Liu, T. Kim, L. Liang, R. McKay, W. Sun, J.S. Duncan, Distortion
actions on Mathematical Software (TOMS) 41 (2) (2015) 1–36. energy for deep learning-based volumetric finite element mesh generation for
[47] F. Moshfeghifar, M.K. Nielsen, J.D. Tascón-Vidarte, S. Darkner, K. Erleben, A di- aortic valves, in: International Conference on Medical Image Computing and
rect geometry processing cartilage generation method using segmented bone Computer-Assisted Intervention, Springer, 2021, pp. 485–494.
models from datasets with poor cartilage visibility, 2022, arXiv:2203.10667 [69] L. Liang, M. Liu, C. Martin, J.A. Elefteriades, W. Sun, A machine learning ap-
[48] T. Schneider, Y. Hu, J. Dumas, X. Gao, D. Panozzo, D. Zorin, Decoupling simula- proach to investigate the relationship between shape features and numerically
tion accuracy from mesh quality, ACM Trans. Graph. 37 (6) (2018). predicted risk of ascending aortic aneurysm, Biomechanics and modeling in
[49] J. Ortún-Terrazas, J. Cegoñino, U. Santana-Penín, U. Santana-Mora, A.P.D. Palo- mechanobiology 16 (5) (2017) 1519–1533.
mar, Approach towards the porous fibrous structure of the periodontal liga- [70] G. Pascoletti, A. Aldieri, M. Terzini, P. Bhattacharya, M. Calì, E.M. Zanetti,
ment using micro-computerized tomography and finite element analysis, Jour- Stochastic PCA-based bone models from inverse transform sampling: Proof of
nal of the mechanical behavior of biomedical materials 79 (2018) 135–149. concept for mandibles and proximal femurs, Applied Sciences 11 (11) (2021)
[50] A. Ziegler, L. Keilig, A. Kawarizadeh, A. Jäger, C. Bourauel, Numerical simulation 5204.
of the biomechanical behaviour of multi-rooted teeth, The European Journal of [71] F. Moshfeghifar, T. Gholamalizadeh, Z. Ferguson, T. Schneider, M.B. Nielsen, D.
Orthodontics 27 (4) (2005) 333–339. Panozzo, S. Darkner, K. Erleben, Libhip: An open-access hip joint model repos-
[51] L. Qian, M. Todo, Y. Morita, Y. Matsushita, K. Koyano, Deformation analysis of itory suitable for finite element method simulation, submitted to computer
the periodontium considering the viscoelasticity of the periodontal ligament, methods and programs in biomedicine (2022)
Dental Materials 25 (10) (2009) 1285–1292.
[52] S. Benazzi, H.N. Nguyen, O. Kullmer, J.J. Hublin, Unravelling the functional
biomechanics of dental features and tooth wear, PLoS One 8 (7) (2013) e69990.
17

1 s2.0 S0169260722003911 Main

Uploaded by

Copyright:

Available Formats

1 s2.0 S0169260722003911 Main

Uploaded by

Document Information

Original Title

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

1 s2.0 S0169260722003911 Main

Uploaded by

Copyright:

Available Formats

Computer Methods and Programs in Biomedicine 224 (2022) 107009

Contents lists available at ScienceDirect

Computer Methods and Programs in Biomedicine

Open-Full-Jaw: An open-access dataset and pipeline for ﬁnite element

1. Introduction and cardiovascular surgeries [2]. More speciﬁcally, in the orthodon-

The main reason for using a single model in the literature is

In other words, in almost all of the studies focusing on full den-

lized geometries, volumetric meshes, and FE models have not been

compare the results. Table 1 presents an overview of the related

discretization type, and availability of the models.

[10] provides detailed geometries of one mandible structure and

mandible in two different steps to provide detailed teeth struc-

Micro-CTs are high-resolution CT scans that are normally acquired

from dead specimens or cadavers due to high x-ray exposures.

cannot cover geometrical variations across different patients in the

The recently introduced OpenJaw Dataset [17] provides open-

quired from CBCT scans. Each patients data includes surface

In the geometry reconstruction step, the utilized scans in the study

pled to the same resolution of 0.15 mm [5]. Although this study

provides more public samples, generalizability to the population

remains a problem. Moreover, both of the abovementioned studies

include manual tasks in different geometry-processing or meshing

meshes using the unprocessed reconstructed meshes.

Our main contributions in the Open-Full-Jaw study can be sum-

1. We provide an open-access dataset of different patient-speciﬁc

models of the human jaw, including the maxilla, mandible, full

dentition, and the PDL geometries obtained from CBCT scans

of 17 patients. It is the largest publicly available dataset with

validated segmented geometries and quality volumetric meshes

that can directly be used in FEM studies. Furthermore, to the

best of our knowledge, it is the ﬁrst repository containing the

2. We introduce a unique repository (https://github.com/diku-dk/

Sarrafpour et al. (2013) [15]

Open- Full- Jaw) containing (1) clinically validated segmented

OpenMandible (2021) [10]

Open-Full-Jaw (our study)

geometries and the resulting dense irregular surface meshes;

OpenJaw (2021) [5,17]

Lee et al. (2018) [16]

(2) the quality and adaptive volumetric/surface meshes to be

used in FE simulations; (3) the automatically generated FEM

viding great information for the users to automatically set up

Patient’s ID Original scan Final resampled ROI Jaws Periodontal disease

of the scans with slice thicknesses above 0.3 mm is very challeng-

3.3. Clinical validation of the teeth-bone segmentation

Maximum breadth of mandibular body at the mental foramen [25]

3.4. Geometrical variations of the dataset

We evaluate the geometrical variation of the reconstructed ge-

3.5. The pipeline for generating FE model of a human jaw

The ﬂexibility of the FE method allows it to use a wide range of

spatial discretizations [34]. We opt for an unstructured tetrahedral

To eliminate the manual geometry processing, quality mesh-

ing, or mesh decimation steps, we propose to directly use the ex-

PDL layer can be reconstructed by segmenting it from the micro-

is done by calculating the radius r of the mean curvature h at

point t based on r (t ) = 1/h(t ), as the reciprocal of the curvature

at that point. The radius of the curvature provides useful infor-

can be moved in the opposite direction of the normals of the sur-

face before a singularity occurs. This means that the movement