The Search for the Centre of Resistance of a Tooth Using a Finite Element Model and the Continuum Mechanics

In orthodontics, the CR is defined as the ideal point for the positioning of the bracket. This paper presents numerical analyses of a premolar tooth, which aim at searching the centre of resistance (CR) by applying several external torques. We built the geometry of the tooth by making a tomographic ren-derisation with the software MIMICS. The numerical simulation is performed in the software ANSYS WORKBENCH (ANSYS® Academic Research). The domain of the tooth is discretized using three-dimensional solid elements. The search for the position of CR is made through a heuristic analysis by applying external torques to the model and, afterwards, by seeking the intersection of the neutral planes (where there is no translation). The numerical results showed that the CR is not a point, but rather a region that strongly depends on the relative stiffness between the tooth and the Perio-dontium (PDL).


INTRODUCTION
In orthodontics, the CR is considered the most important point, related to dental root movement.Initially, this point was idealised in two dimensions (Burstone, 1962); however, the greatest difficulties arise when this concept is extrapolated to a spatial structure.Some authors indicated that this point does not exist in a three-dimensional and asymmetric structure, e.g. in a tooth partially restricted by a deformable body (PDL -periodontal ligament).Detailed descriptions of these difficulties have been presented by de Viecilli et al. (2013), Meyer et al. (2010), Cattaneo et al. (2008), among others.
The objective of this research is to study the centre of resistance (CR) of a tooth.The CR is defined as the ideal point for positioning the bracket, such that when a force is applied to the bracket, the dental root displaces without causing rotation.
This paper proposes the search for the position of CR by using the principles of solid mechanics.The results have shown that this point does exist in a rigid three-dimensional body, which is also subjected to small displacements and is bounded by a continuum medium with identical mechanical properties in all directions.
An additional numerical analysis is performed for determining whether a possible CR is highly dependent on both the relative flexibility between the tooth and the PDL and the variation of stiffness in the three reference directions (occlusal, buccal-lingual and mesial-distal).

STATE OF THE ART
In the literature, different works have been made for investigating on the existence and the position of the CR.In recent years, renderisation tools and numerical simulations (especially using finite element models) have played a fundamental role in accomplishing advancements in this area.The following paragraphs describe the most influential research studies related to the topic of this work.Christiansen and Burstone (1969) studied the localisation of CR experimentally by applying different forces on the maxillary central incisor.They confirmed that the CR is linked to the Moment-Force relationship.The authors showed that the CR is placed in a different position, as compared to the analytical position.This difference arose because the theoretical formulas (Burstone, 1962) were deducted in two dimensions.Tanne et al. (1991) investigated the behaviour of tooth displacements.The authors varied the length of the root and the height of alveolar bone of the tooth.A three-dimensional analysis was carried out using finite element (FE) model for discretizing a superior central incisor tooth.The tooth displacements were estimated for the different heights of the tooth and the PDL.The authors found that the position of CR is altered when these parameters are altered.The authors noted that the positions of CR are strongly influenced by the height of the dental root and by the loss of bone around the tooth.However, this work did not describe how the CR was obtained.Even so, the results indicated that two-dimensional projections of the intersection of neutral planes were used, i.e. planes in which stress and strain equal zero.Cattaneo et al. (2008) used two-dimensional, microcomputed tomographic images for generating a FE model of the tooth.The authors performed analyses with a number of loads, simulating the orthodontic movements and thus analysing the moment-force ratio.A range of values of the moment-force ratio was generated for obtaining different dental movements, thus concluding that the CR varies with the moment-force ratio.Meyer et al. (2010) compared the positions of the CRs for six different central incisors.The authors built a three-dimensional solid FE model using microcomputed tomographic images of teeth and their support structures.The CRs were estimated as being the projection of lines with zero strain, in two directions.However, according to Viecilli et al. (2013), the points obtained by Meyer et al. are not the CRs.Viecilli et al. (2013) implemented a FE model to analyse a maxilla, whose first molar is subject to several torques, varying their direction, to obtain the position of each rotation axis.A two-dimensional projection of axes of resistance was used for finding a point in which tooth displacement without rotating was feasible.Their conclusion claimed that the axes of rotation do not intersect in one single point.Therefore, the authors asserted that the CR could not exist.Dathe et al. (2013) proposed a simple example for demonstrating that, in general, the CR does not exist in three dimensions.However, the authors suggested a new concept, namely the centre of elasticity, which preserves the desired properties of CR.The authors adopted a linear relationship between loads and displacements.
Other works related to the same goal of this work are: Kojima et al. (2007), Lin et al. (2010), Nakajima et al. (2007), among others.

MATERIALS AND METHODS.
This work addresses the problem of finding the CR of a tooth as follows: A three-dimensional model was generated using a renderisation of tomographic images of two bottom premolar teeth.The images were renderisated using the commercial software MIMICS (Version 16.0, Materialise, Leuven, Belgium) 1 , as shown in Figure 1.The three-dimensional FE mesh was generated with the module 3-Matic.The mechanical analyses were performed in the commercial software (ANSYS® Academic Research), cf. Figure 2.  A heuristic analysis was carried out aiming at finding the position of CR, by applying a number of external torques to the model.The solution has been sought as the intersection of the planes where the strain is zero.
The above procedure was repeated several times by modifying the mechanical properties of materials, in an arbitrary and convenient fashion.The objectives of this analysis are: i) To verify the behaviour of the position of the CR with the variation of the relative stiffness between the tooth and the PDL.ii) To investigate the influence of the properties of PDL, the boundary conditions (fixation of alveolar bone) and the reduction of the PDL stiffness in the tooth behaviour.
Thus, the hypothesis of the existence of a CR for three-dimensional deformable bodies is tested next, especially when the elastic properties of the support tissue vary along with the direction.
In this work, all the torques were measured in .N mm , the lengths and coordinates in mm and the forces in N .Finally, a simple theoretical analysis with a rigid body and small displacements assumption was performed to determine whether such CR exists.In this model, the stiffness of PDL was assumed as constant in all directions.For these conditions, it is important to note that there is a single CR.

SIMPLIFIED ANALYTIC SOLUTION
4.1 The CR A body is considered under rotations around the orthogonal axis, in which the directions and the reference centre are arbitraries.
When a rotation, x  , is imposed around the X -axis, the displacements are as follows: , , 0, where u , v and w are the displacements in the X , Y and Z axis, respectively; and X , Y and Z are the coordinates related to an arbitrary Cartesian system.Similarly, the rotations y  and z  are applied for obtaining analogous results for the displacements.By applying the superposition principle, the total displacements of three rotations yield in: . , , The displacements can be written as a function of an arbitrary origin, X o , Y o and Z o , as follows: Latin American Journal of Solids and Structures, 2018, 15(5), e41 4/15 where K and q i represent the stiffness of the body and the reactive forces per unit of area2 , respectively.By performing the summation of forces in each direction, and by the integration of Eq. ( 3) in the respective areas, which are bounded by the elastic tissue, one obtains the following set of relationships after equalling to zero: where F i represents the force in the i direction.The formulas of Eq. ( 4) are worked out, which yields in: The conditions given by Eq. ( 5) are hold, regardless of i  , provided that the following expressions are met, , , In other words, whether arbitrary rotations are or not applied in a solid, which is at its turn connected to an elastic surface, the CR of the solid is found as being the centre of gravity of the contact surface.The contact surface is the result of the intersection between the tooth (solid) and the periodontal tissue or PDL (elastic medium).

Axes of Resistance
If a torque is applied onto an axis and the obtained structural response is a pure rotation without reactive forces, then it is said that this axis is defined as an "axis of resistance".In Solid Mechanics, the above definition matches exactly the concept of the principal moments of inertia.The procedure to calculate the principal axes of inertia in a simplified and analytical way is shown as follows.
In the centre of gravity of the contact surface, a value for rotation x  is imposed, such that y  and z  remain fixed.Hence, the displacements are given by Eq. ( 1) and the reactive forces are given by the following equations: Therefore, the acting moments in the constraints are determined by: Similarly, the value of rotation y  can be imposed, and all the other rotations remain fixed.Then the resulting moments in the constraints are as follows: In a similar manner, the following expressions are obtained, with a defined value for the rotation z  : 2 2 , , .
Latin American Journal of Solids and Structures, 2018, 15(5), e41 6/15 The moments and the products of inertia are defined as: Then Eq. ( 11) is to be replaced in Eq. ( 10) for obtaining the following expressions , , The principle of superposition is used in a form that the system of the Eq. ( 12) can be rewritten as follows: . .
According to Sadd (2009), "an important property of a tensor is that if we know its components in one coordinate system, we can find them in any other coordinate frame by using the appropriate transformation law.Considering the tensor transformation, it should be apparent that there might exist particular coordinate systems in which the components of a tensor take in maximum or minimum values".By analysing the expression (13) of the tooth problem, it is noted that the matrix represents the moments and products of inertia.With an appropriate rotation, this matrix is transformed into a diagonal system.Then, according to the aforementioned concept, the diagonal values are the principal moments of inertia, and their directions are the principal directions.At the same time, there is a point at which all directions intersect, for the principal directions of a tensor.Then, it is possible to conclude that the CR exists for the conditions used in this analytical solution.
The above results were obtained considering that the stiffness of PDL was constant in all directions.However, the stiffness of PDL is far from being constant, and this fact becomes apparent when considering that this tissue is neither isotropic, nor its thickness is uniform.Even so, the mechanical model used to simulate the PDL does have an influence on the mechanical behaviour of the whole tooth (Wood et al. 2011).Without this simplification, it Latin American Journal of Solids and Structures, 2018, 15(5), e41 7/15 would not have been possible to obtain Eq. ( 13), in other words, the tensor behaviour of the matrix would not be guaranteed, which represents the moments and products of inertia.

NUMERICAL RESULTS
This section shows the results of the numerical analyses of the model that were performed using the ANSYS software.The goal of the analyses is to find the axes of resistance, in terms of the intersection of the neutral planes.Viecilli et al. (2013) performed analyses using a FE model and questioned the existence of CR, because, to their understanding, it is not possible to find a single point in which the axes of resistance intersects between each other.However, the authors did not provide any explanation for this conclusion.
Section 4 showed that there is a CR for a rigid body, which is bounded by another body that exhibits constant stiffness in all directions.This Section further demonstrated that this CR is positioned in the GC of the surface created by the contact of the two bodies, i.e. the rigid body that is delimited by an elastic body.However, it can be proved that the latter statement is not applicable for deformable bodies.This problem is even more apparent when the elastic properties of the support body vary in the different directions.The GC of the contact surface is found using the ANSYS software (CR for a rigid body with simplified boundary conditions), cf. Figure 3, as well as the directional vectors of the main axes of this surface.The results are shown in Table 1.Due to the high stiffness of the tooth, it would be expected that the CR exists and that it is close to the GC, as shown in Table 1.Care was taken in performing this analysis, in view of the conclusions of Viecilli et al. (2013) about the existence of CR.In Figure 4, the auxiliary coordinate system is shown, which is used to reference the position of the GC, whereas the origin of the coordinates of this system is shown in Table 1.The model used in the structural analyses is composed of two premolar teeth, the PDL and bone, as shown in Figure 5 and Figure 6.The PDL coated the root of each tooth.The bone was partitioned into three parts to facilitate the application of the boundary conditions.The external faces of the bone end parts were fixed.This boundary condition ensures that the system presents non-homogenous elastic properties in the contact.The coordinates shown in Table 3 are related to the adopted initial system.However, an auxiliary system is also created to facilitate the visualization.This system is shown in Figure 4, which is located in the GC of the contact surface of the tooth and PDL.The GC's position is detailed in Table 1.The mechanical properties of materials are summarized in Table 2, which are based on a previous study (Viecilli et al., 2013).These properties were used in the initial mechanical analysis.After the initial analysis, further mechanical analyses were performed for investigating the existence of the CR.In these analyses, the values of mechanical properties were conveniently modified.The first premolar tooth was selected for the application of external forces, simulating the force exerted by the bracket on the tooth.The applied loads were torques around the GC.Thus, the applied torques were 2 .N mm in the direction of coordinates axes, and these loads are the results of the application of unitary forces (1N ) at the points shown in Table 3.In order to find the neutral planes (where the stresses and strains are zero), the model was subjected to three types of loads: • Mx is the torque around the X axis and is formed by the binary ( * ) Fy Z   and ( * Fy Z  ), where Fy is a unitary force (1 .N ) in Y direction and Z + and Z -are given in Table 3.
• My is the torque around the Y axis and is formed by the binary ( *Z+ Fx  ) and ( *Z-) Fx , where Fx is a unitary force (1 N ) in X direction.
• Mz is the torque around the Z axis and is formed by the binary ( * Fx Y  ) and ( * ) Fx Y   where Y + and Y -are given in Table 3.    Combinations of the initial mechanical properties, each with the loads described above, allowed for performing three simulations.The displacements of the tooth subjected to in the Y and Z directions are shown in Figure 7 and Figure 8, respectively.The displacements of the tooth subjected to My in the X and Z directions are shown in Figure 9 and Figure 10, respectively.Finally, the displacements of the tooth subjected to Mz in the X and Y directions are shown in Figure 11 and Figure 12, respectively.The neutral plane is depicted in black in each plot between Figure 7 and Figure 12.The resistance axes are determined as the intersection of the neutral planes, as shown in Figure 13.In Figure 13, it can be noted that the resistance axes do not intersect at one point, but they are rather close between each other.With this result, it is possible to conclude that the CR does not exist.
It is worth putting forward some few remarks, before the existence of CR to be demonstrated.A sphere of 0.15mm radius is now defined, such that the three axes pass through it.This sphere represents a small region where one might consider the CR location for practical purposes.This region has a relative difference of approximately 1.67% , related to the greatest dimension of the tooth   18mm .The methodology used to verify whether some point (test point) inside the sphere could be considered as CR consisted on applying a system of forces, with restriction of the rotations, and measuring the reaction moments.If the reaction in the test point is negligible, then this point could be considered as a good candidate for being a CR.
A trial and error approach facilitates the refinement of the determination of the CR, as follows: if in the test position, inside the sphere, which is shown in Figure 14, the rotations are restricted, and a unitary force is applied in any direction, then the reaction moment in the point must be lower than 0.09 .
N mm .Different test points that could be considered as optimal (negligible reaction moment) have been located using the above process, which are related to the torque direction (force and arm):  The points previously found were used to calculate the moments of inertia of Eq. ( 13).A unitary rotation around the axis is applied in each of these points under analysis, while keeping the remaining rotations constrained.Considering that the tensor must be as symmetrical as possible, the position was numerically found at the point

 
43.1, 57.75, 27.6   , which is closer to the point determined in section 4 that resulted from the analytical solution.Furthermore, the directions of resistance axes (cf.Table 4) are also closer to the analytically determined directions of the axes of resistance in Section 4.
Following the above results, it is possible to conclude that the GC of the contact surface is close to the test points, previously identified as candidates for being the CR.However, it is not possible at this point to conclude that the candidates are indeed the CR, and this is due to the infeasibility for determining a single point by using the numerical analysis.An additional numerical analysis was carried out with suitably modified mechanical properties, cf.Table 5.In this test, the bone and PDL stiffness were increased, and the tooth stiffness was reduced.With these modifications at hand, the deformation of the tooth is maximized, thus making it possible to verify whether this deformability influences the position of CR. Figure 15 shows the resistance axes obtained with the modified material properties.It is noted in this Figure that there is a sphere circumscribed of 0.45mm radius that tangencies all of three axes.The relative ratio, the radius of sphere related to the greatest dimension of the tooth, increases from 1.67% (initial analysis) to 5% (modified analysis).Furthermore, it is noted that the resistance axes are not orthogonal among them anymore.This conclusion about the orthogonality of axes implies in that the CR is not a geometrical property of a flexible tooth, as it is said in orthodontic technical literature.This fact confirms the initial hypotheses of Viecilli et al. (2013).Nevertheless, the 'rigid CR' can still be used for practical purposes when using the real elastic properties of the involved tissues,.

Analogous Problem
There is a difficulty to understand that the impossibility of finding the CR on a tooth is related to the fact that the tooth is a deformable body.Thus, it becomes necessary to show how the deformation of the body influences in the position of its CR.Imagining that the beam represents the tooth (dentin) and the elastic support represents the PDL and bone, one easily notices that the CR does not exist.Moreover, as the contact force is proportional to the vertical displacement, Figure 17 also indicates that the rotation centre of the tooth also depends on the relative stiffness of the tooth and the supporting tissue.When the supporting tissue presents high stiffness, the CR and the rotation axis can be used for practical applications in orthodontics.

CONCLUSIONS
The objective of this research was to study the CR of a tooth.The CR is defined as the ideal point for positioning the bracket, and this is the most important point related to dental root movement.
Three-dimensional analyses of premolar tooth were made using the software ANSYS WORKBENCH.The geometry of the models was built with a renderisation of tomographic images of a real tooth.The numerical analyses of this model are performed in section 5.
With the results of the three-dimensional model, it has been concluded that the tooth CR is not a point, but rather a region with very low values of the reactive moment.This finding can still be used for practical purposes, in a very similar manner to the CR concept.This region is strongly dependent on both the relative stiffness of the tooth and the supporting tissue.This conclusion suggests that the axes of resistance are not real geometrical properties of the tooth.In order to ensure that this conclusion is not an effect of numerical errors in the FE analyses, a simplified example has been proposed for confirming the statement.However, from both three-dimensional tooth analysis and simple analysis, the results indicated that when the support tissue is healthy, with elastic modulus near to the elastic modulus of dentin and a sufficient thin PDL, the geometrical centre (GC) could be used with a good approximation for orthodontic applications.Furthermore, a simple analysis has been seen to suffice for finding the GC of the contact surface of the tooth and PDL, as we describe in item 5.This paper used a numerical case in order to illustrate the problem.This numerical model is highly representative, since this was generated by using a real geometry.As the proposed methodology is extended for other cases, similar results to those presented in this paper should be achieved.Additional numerical models might be implemented in a next step of this research, in order to confirm the results and conclusions of this paper.
Further works are intended to improve bone and PDL material models including viscosity, plasticity and tissue regeneration, in continuing attempts for improving the understanding of orthodontic mechanics.

Figure 1 :
Figure 1: Pre-selected premolar teeth and bone structure used in this work for building the numerical model.

Figure 3 :
Figure 3: The contact surface of the tooth and PDL.

Figure 4 :
Figure 4: Auxiliary coordinate system positioned in the GC of the contact surface.

Figure 10 :
Figure 10: Z Displacement with My

Figure 15 :
Figure 15: The resistance axes (model with softened tooth and stiffened bone and PDL).

Table 6 :
Parameters of the beam used to represent the problem of finding the centroid of a deformable body.proposes a simple example that represents the problem in an analogous way.Thereby, it is proposed to analyze a beam(Coda and Greco (2004)), which is supported in part by an elastic base.A force F is exerted in the right end, cf. Figure 16.It is considered a geometrical non-linear analysis based on position description.Reissner kinematics is adopted.The parameter k that represents the stiffness of the elastic support takes values of 0.1 , 1.0 , 10.0 , 100.0 and 2 500.0/ N mm .The remaining parameters of the problem are shown in the Table6.The elastic support is modelled following the Winkler's method.

Figure 16 :
Figure 16: Beam used to represent the problem of to find the centroid of a deformable body.It is noted in Figure 17 that when the k -stiffness is increased the centroid of the reactive distributed force changes.This variation is evident when the stiffness of the elastic support becomes 2 500 / N mm .

Figure 17 :
Figure 17: The concentrated reactive force at each contact node for different values of k.

Table 1 :
Coordinates and directional vectors of the GC of the contact surface of the tooth and PDL.

Table 3 :
Application points of torque.
Latin American Journal of Solids and Structures, 2018, 15(5), e41 9/15 • Force in X and arm in Y , then the point is  

Table 4 :
Coordinates of points of resistance axes.