1 Introduction

Human tooth is composed of two major calcified tissues,i.e.,dentin and enamel.Dentin occupies the majority of the tooth by weight and volume and serves as the elastic base of the hardest enamel tissue.It also protects the organic tissues in the pulp by absorbing and distributing the occlusal stresses.Thus,the mechanical behavior of dentin is of considerableimportance in understanding how stresses are potentially altered by dental restorative procedures,aging,or diseases.Dentin has a microstructure consisting of a collagen matrix reinforced by hydroxyapatite nano-platelets.One of the prominent features of dentin is the dentin tubules,which are cylindrically shaped and surrounded by the peritubular dentin^{[1, 2]}.These tubules radiate outward from the pulp chamber to the dentin-enamel junction (DEJ).The number of tubules at the DEJ is lower than that of the inner dentin,and the percentage of tubule area and the diameter of tubules in the occlusal dentin vary from about 22% and 2.5 μm near the pulp to 1% and 0.8 μm at the DEJ^[3].Nanoindentation has become a popular technique to assess the elastic mechanical properties of dentin^{[4, 5, 6]}.Young's modulus and the hardness of dentin range from 8.7 GPa to 31 GPa and from 0.12 GPa to 2.5 GPa,respectively^{[2, 7]}.In an evaluation of the unique spatial variations in properties of the constituents of dentin,Senawongse et al.^[7]demonstrated that the inner dentin possessed lower hardness and Young's modulus than the outer dentin.In the vicinity of the DEJ,Tesch et al.^[8]demonstrated that the hardness and Young's modulus of the dentin increased from the DEJ towards the coronal dentin in an area between 0 mm and 1.5 mm below the DEJ.Since the Berkovich indenter tip caused a complex stress distribution within the material,the resulted load-depth curve exhibited nonlinear responses which hided the substantial constitutive relationship of the material in plastic and damage stages.

Numerical analysis provided an alternative solution for characterizing the mechanical behavior of biological materials combined with the experimental results^{[9, 10, 11, 12, 13]}.The post-yielding behaviors of the biological hard material investigated using the numerical model and nanoindentation have been focused in recent years^{[14, 15, 16, 17, 18]}.For instance,using the elastic anisotropy,the plastic pressure-independent anisotropy,and the strain hardening,Fan and Rho^[14]performed a finite element model to characterize the nanoindentation on bone.Combining of the nanoindentation and elastic-plastic finite element analysis and using a Mohr-Coulomb cohesivefrictional strength criterion,Tai et al.^[15]investigated the ultrastructural origins of the strength of bone.With the use of anisotropic yield constitutive model,Carnelli et al.^[16]developed a finite element model to characterize the response of nanoindentation of cortical bone.The finite element model with the conical Drucker-Prager yield criterion and nanoindentation^[17]was used to investigate the plastic behavior of bone.The influence of the yield surface shape and the damage on the depth-dependent response of bone to nanoindentation using spherical and conical tips was also determined^[18].When the bone and dentin yield,the energy might dissipate from the inelastic deformation and damage^[19].The damage was then considered in several numerical models.Noticeably,Zhao et al.^[20]used a numerical model to describe the viscosity,the plastic deformation,and the damage behavior of bones based on the indentation test and discussed the relationship between the model parameters and the mechanical/physical properties.A full viscous-elastic-plastic indentation model was used to analyze the nanoindentation process for polymers and bones,which showed great promise in obtaining inelastic mechanical properties^{[21, 22]}.Lucchini et al.^[23]developed a numerical model considering damage to predict the indentation stiffness and the hardness of lamellar bone at 300 nm penetration depth along both axial and transverse directions and concluded that there was a decrease in the magnitude of indentation modulus with the penetration depth.By developing a plastic damage model accounting for the chain density in protein,the region-dependent damage behavior of enamel was demonstrated^[24].The studies provided an indication that the numerical modeling could be used to dentin to discover nonlinear mechanical behaviors in plastic and damage stages.

The primary objective of this study is to evaluate the variations in mechanical properties of dentin with aging along the thickness of coronal dentin using a combined experimental and numerical approach.Nanoindentation experiments are conducted in three regions of coronal dentin sections that are obtained from molars of patients in three age groups.The finite element analysis incorporating a plastic-damage model is used to characterize the inelastic and damage behaviors of dentin based on the experimental results.

2 Materials and methods 2.1 Specimen preparation

Non-carious human molars are obtained from participating dental hospitals in Shanghai,China according to the approved protocols issued by the Institutional Review Board of Shanghai University.The teeth have been stored in Hank's balance salt solution (HBSS) in 4 ℃ until being sliced^[25].At receipt,the teeth are divided into three age groups including the young (18 ≤ age ≤ 23;N=7),mid-aged (33 ≤ age ≤ 43;N=7),and aged (old) groups (58 ≤ age ≤ 65;N=7).A programmable slicer machine (EC400,Shenyang Kejing Instrument Co.,Ltd.,China) with a diamond abrasive slicing wheel is used to section the molar along the buccal-lingual direction with continuous water-based coolant.The dentin tubules are essentially parallel to the section surface.However,as the tubules are rarely directed entirely straight from the pulp to the DEJ,there are potential for small errors in the relative orientation between the section surface and tubules.Only one specimen with the thickness of 2 mm is machined from each tooth.The specimens are then successively polished from#600 to#5000 grit abrasive papers.All of these specimens are immersed in the HBSS before the indentation test.

2.2 Nanoindentation

The Hysitron Triboscope nanoindenter in-situ nanomechanical test system (Hysitron,Minneapolis,MN,U.S.A.) is used in this study with a Berkovich indenter (TI-039,Hysitron Inc.,Minneapolis,MN,U.S.A.).Indentations are made with a maximum load of 6 mN and performed at a rate of 0.3m·N/s,followed by 5 s dwell time,and then unloaded with the same rate as loading.The indentations are conducted on each sectioned surface along the direction of tubules from the pulp to the DEJ (as shown in Fig. 1).In summary,a total of 246 indentations are conducted in the young dentin,while 296 and 242 indentations are in the mid-aged and old dentin,respectively.Along each indentation line from the pulp to the DEJ,three regions are defined with equal distances.Specifically,the region nearest the pulp is defined as the inner dentin,the region nearest the DEJ is defined as the outer dentin,and the region between them is defined as the middle region.The elastic mechanical properties including both the hardness H and Young's modulus E are evaluated^{[26, 27]},and the corresponding load-depth curves are also obtained for further evaluation of the inelastic mechanical properties via numerical simulations.

Here,P_max is the maximum load,A is the projection contact area at the maximum load,and S is the measured stiffness.The indentation responses for each location are represented by the average load-depth curves.The average hardness and elastic modulus of dentin for the young,mid-aged,and old groups are determined.These results are statistically analyzed using the one-way ANOVA with a post-hoc tukey,and the statistical significance is set as p ≤ 0.05.

2.3 Numerical modeling for indentation

Since the indentation results in a complex stress state in the vicinity of indenter,the finite element analysis is used to simulate the nanoindentation experiments in order to characterize the inelastic mechanical behaviors of dentin.The Berkovich indenter is modelled as an analytical rigid body with a radius of 100 nm and an apex angle of 70.3°^[28].Frictionless contact between indenter and dentin is implemented^[29],and a load-control loading scheme is adopted during the simulation.According to the nanoindentation experiments,a compressive load with 6mN is applied.An axisymmetric two-dimensional finite element mesh is developed using the commercial software (ABAQUS 6.10-1) to represent the indentation loading process.A total of 24 764-node bilinear axisymmetric quadrilateral elements (CAX4R) are used to mesh the dentin (as shown in Fig. 2).The mesh is well tested for convergence and is accurate in the stress analysis in the area beneath the indenter.

Considering that dentin undergoes plastic deformation,which is evidenced by the irreversible deformation after unloading,and damage is demonstrated by the stiffness degradation in unloading,the plastic damage model proposed by Lubliner et al.^[30]and Lee and Fenves^[31]is introduced to describe the mechanical behavior of dentin.The stress can be expressed as

where d represents the damage variable characterizing the stiffness degradation during indentation,C_ijkl is the stiffness tensor,and ε_kl and ε^p_kl are the total strain and the plastic strain,respectively.In the present paper,the damage variable d is assumed to be isotropic with 0 ≤ d ≤ 1,where d=0 means that there is no damage within dentin,while d=1 indicates that the material is failed.The plastic strain follows the flow rule, where denotes the plastic multiplier,and is the effective stress defined as =σ_kl/(1 − d).f is the plastic potential which is described by the Drucker-Prager hyperbolic function, where ε is the eccentricity parameter,σ_t0 is the uniaxial tensile strength at failure, is the effective hydrostatic pressure stress, is the Misses effective equivalent stress defined in the effectivestress space,and ϕ is the dilation angle.For all the numerical simulations,the associated flow rule is adopted.Therefore,the yield function is identical to the plastic potential.

In the model,the Poisson ratio is fixed at 0.3,and the elastic modulus E is determined from the unloading slope of the load-depth curve and is listed in Table 1.Based on the similar microstructure for bone and dentin,the ratio of the equibiaxial and uniaxial compressive yield stresses is set to be 1.125,and the value of the dilation angle ϕ is chosen to be 15°^[32].The apparent yield strength is determined from the average load-depth curves by the best fitting method.In detail,an initial yield strength is assumed to calculate the deformation at a given compressive load,and the error between the simulated and experimental load-depth curves is evaluated to derive an updated yield strength.This process is repeated until the optimized yield strength,which reaches the minimum error between the simulated and experimental load-depth curves,is achieved.

Zhang et al.^[32]revealed that the damage behavior is attributed exclusively to compressive loading during the indentation test.When fitting the model to the experiment,a function of the equivalent plastic strains in compression κ_c for the damage variable d_c governed by two additional parameters a and b is expressed as follows:

The parameters a and b,which are determined by best fitting the experimental curves,are set to be 0.9 and 5,respectively.The damage variable can be used to characterize the indentation induced micro-damage,such as microcracks,breakage of macromolecular chains in protein,and failure of interface between mineral and protein.The damage behavior of dentin is evaluated in the three regions for each age group.

3 Results

The magnitudes of Young's modulus and hardness for inner,middle,and outer dentin of the young,mid-aged,and old groups are evaluated by nanoindentation,as listed in Table 1,and the spatial variations of the elastic mechanical properties of dentin are analyzed.The statistic analysis indicates that Young's modulus of inner dentin is significantly lower than those of middle dentin in the young group (p=0.003) and the mid-age group (p=0.030),and thoseof outer dentin in all age groups (p=0.002,0.016,and 0.003 for young,mid-aged,and old groups,respectively).The hardness of outer dentin is significantly higher than those of inner dentin (p=0.000 in all age groups) in all age groups and middle dentin in the young group (p=0.005) and the aged group (p=0.002).Aging has significant effects on the mechanical properties of dentin as well.Young's modulus of dentin in the old group is significantly greater than that in the young group (p=0.006,0.016,and 0.005 in inner,middle,and outer dentin,respectively) regardless of the location.The hardness of dentin in the old group is significantly higher than those in the young group (p=0.005,0.000,and 0.000 in the inner,middle,and outer regions,respectively) and the mid-aged group (p=0.002,0.034,and 0.002 in the inner,middle,and outer regions,respectively).

The nanoindentation load-depth curves are also varied with locations and aging.The representative load-depth curves for inner,middle,and outer regions of young dentin are shown in Fig. 3.Obviously,under the same indentation load,the inner dentin undergoes greater penetration depth than other regions.The average load-depth curves for the outer dentin in the young,mid-aged,and old groups are shown in Fig. 4.The penetration depth in the young group reaches the highest value.

The load-depth curves obtained from the experiments at each region in each age group are then used as the benchmark for numerical modeling,as shown in Figs.3 and 4.The plastic mechanical properties of dentin are able to be studied by adjusting parameters in the model until good agreement is achieved between the experimental and numerical results.The yield strengths of dentin σ_y are derived by fitting the experimental load-depth curve and are listed in Table 1.The apparent yield strength in the inner dentin is lower than those in the other regions in each age group,and the yield strength increases with aging.

The distribution of damage d_c for dentin corresponding to the maximum indentation depth is also plotted in Fig. 5.The maximum value of damage d_max in the element near the tip ofindenter can be found in Table 1.As indicated in Table 1,in the same age group,d_max is found in the inner dentin,while in the same location,d_max is found in the young group.Meanwhile,the damage zone can be defined as the color contours of the damage variable d_c in the vicinity of the indent in dentin.As shown in Fig. 5,the damage zone developed within the inner dentin is greater than those in the middle and outer dentin in the young group.One can also find that the size of damage zone is also affected by aging.The maximum damage zone and d_max are recognized in the young group.

4 Discussion

The region dependent elastic and plastic mechanical properties of dentin are explored using an experimental and numerical approach.The results from the nanoindentation experiments reveal that Young's modulus of the inner dentin is significantly lower than those at the outer dentin in all groups,and the hardness of the outer dentin is higher than those at the other regions in the young group.The experimental findings are in good agreement with former studies.Kinney et al.^[33]studied the spatial variation of mechanical properties of dentin and reported that the hardness and Young's modulus at the intertubular dentin decrease from the DEJ to the pulp.Angker et al.^[34]demonstrated that the hardness and Young's modulus of dentin at the inner dentin are significantly lower than those at other regions.However,the magnitude of Young's modulus in this study is slightly higher compared with Angker's findings.This can be affected by the trivial differences in indentation methods.In the former study,nanoindentation was conducted with the maximum force applied of 25 mN and 30 s delay at the maximum force,which caused a larger size of impression on the specimen surface compared with the current experiment.This was reasonable since Lucchini et al.^[23]investigated that the indentation modulus decreases with the penetration depth for the hard tissue.Despite the difference in magnitude,the variation trend is identical.The spatial variation in mechanical properties observed in this study can be accounted for the microstructure of the dentin.Kinney et al.^[33]demonstrated that the peritubular dentin possesses higher Young's modulus and hardness than the intertubular dentin,and Young's modulus and hardness of peritubular dentin are independent of locations,while both of them in the intertubular dentin near the DEJ are higher than those near the pulp.The microstructure of dentin shows that the density of dentin tubules increases with the distance from the DEJ.Accordingly,the volume of peritubular dentin also increases with the distance from the DEJ,while the volume of intertubular dentin decreases.The overall contribution of tubules and peritubular and intertubular dentin determines the spatial variation in the mechanical properties.In addition,all specimens involved in the experiment are obtained from the non-carious dentin,as the elastic modulus and hardness of the carious dentin are significantly lower^{[35, 36]}.

The spatial variation of the plastic and damage behavior of dentin is also estimated with the plastic-damage modeling.Compared with the middle and outer dentin,the inner dentin shows less resistance to penetration deformation,and the apparent yield strength of the inner dentin is lower than those in the middle and outer dentin,as listed in Table 1.The results also show that dentin is a functionally graded material with varied penetration resistances and damage variables through the thickness.This graded mechanical design is obviously in accord with the function of dentin.Dentin,as the elastic base of the hardest enamel tissue,protects the pulp tissue by absorbing and distributing the occlusal stresses.When teeth are subjected to the biting force,enamel resists the penetration force and transfers the load to dentin.In order to reduce the deformation of enamel,the outer dentin should be strong and stiff to carry the transferred load.

The mechanical properties of dentin are also affected by aging.Previous studies demonstrated that ultrastructure changed in the aged dentin.The dentin tubules are filled up with large mineral crystals^{[37, 38]},and the average intertubular crystal size decreases with aging^[39].It is also reported that there is an overall increase in the mineral content of dentin with aging^[40].The present study indicate that Young's modulus and hardness increase with aging.The same variation trend is found in its inelastic mechanical properties.The apparent yield strength of dentin and the resistance to penetration increase with aging at each location,which indicates that the aged dentin is stiffer than the young dentin.

The previous study^[41]has demonstrated that,in human dentin,there is microdamage in the form of microcracks associated with inelastic deformation,which is evidenced by the irreversible deformation after unloading.The damage associated with plastic deformation can be used to understand the nature of strength in dentin.One benefit of the numerical analysis is the ability to characterize the distributions of damage associated with plastic deformation.The plastic damage model is adopted to simulate the indentation at dentin at each location and age group.The obtained load-depth curves are fitted well with the experimental curves especially at the unloading stage.In the numerical model,the scalar damage parameter d_c,defined in (6),is used as an important parameter to characterize the gradual degradation of the modulus.The numerical simulation indicates that there are differences in damage distributions surrounding the impression in inner,middle,and outer dentin,as well as in young,mid-aged,and aged groups.The magnitude and volume of damage increase with the distance from the DEJ.This is in accord with the increase of the tubular density with the distance from the DEJ.The dentin possesses a network of tubules,which are void-like structural defects.The impression is a combination of the intertubulur and peritubular dentin and tubules.In the inner dentin,the chance of involving tubules in the impression increases and causes more damage in the area close to pulp.

Note that the damage variable d_c reaches the maximum at the center of impression and gradually decreases with the distance to the center,as shown in Fig. 5.In the real biting process,the rigid enamel transfers the distributed stress to dentin and substantially minimizes the effect of damage in dentin.The numerical simulations also indicate that the magnitude and volume of damage decrease with aging.This may be attributed to the gradual occlusion of the dentin tubules due to deposition of minerals within the lumen.According to (6),the old dentin,which possesses lower magnitude of damage variable,develops fewer plastic strains than the young dentin.This implies that little plastic deformation can be found in the old dentin,which indicates that the old dentin is more brittle compared with the young dentin.Indeed,the previous study also demonstrated that the old dentin is more brittle than the young dentin^[42].

Zhang et al.^[32]demonstrated that the measured final stiffness in nanoindentation is reduced approximate 90% compared with the initial stiffness at the onset of unloading when indenting the cortical bone.In this study for instance,in Fig. 6,the measured stiffness of the outer dentin in the young group at the initial unloading stage is 88.2μN/nm,while the stiffness is reduced to5.1 μN/nm at the final unloading.This agrees well with the former findings.The reduction mechanism can be further investigated through the numerical analysis.In fact,the contact area examined from the numerical simulation decreases during unloading,as indicated in Fig. 6.Note that the radius of the final contact area in Fig. 7(b) is approximate 1/3 of initial in Fig. 7(a).According to (2),the apparent Young's modulus in the final stage reduces 80%−85% compared with that in the initial stage of unloading.Since the damage zone is mainly distributed around the center of the indent,damage causes little influence with large contact area in the initial stage of unloading,while it causes appreciable influence with small contact area in the final stage of unloading.

5 Conclusions

The influence of location and aging on the elastic,plastic,and damage properties of human dentin is evaluated via a combination of nanoindentation and numerical simulations.Specially,the numerical simulation plays an important role in moving forward to understanding the plastic and damage mechanical behaviors of dentin based on the experimental data.Young's modulus,the hardness,and the yield strength in the inner dentin are lower than those of the middle and outer dentin in each age group.The inner dentin also undergoes greater plastic deformation and more damage compared with the middle and outer dentin.The mechanical properties of dentin are also varied with aging.Young's modulus,the hardness,and the yield strength of the old dentin are significantly greater than those of the young dentin.Under indentation loading,the aged dentin can only provide limited plastic deformation,which implies that the dentin becomes brittle and vulnerable to fracture with the increase of age.

Acknowledgements The technical supports from the Instrumental Analysis and Research Center at Shanghai University are acknowledged.