- Email: [email protected]
Viscoelastic finite element evaluation of transient and residual stresses in dental crowns: Design parametric study
journal of the mechanical behavior of biomedical materials 103 (2020) 103545
Contents lists available at ScienceDirect
Journal of the Mechanical Behavior of Biomedical Materials journal homepage: http://www.elsevier.com/locate/jmbbm
Viscoelastic finite element evaluation of transient and residual stresses in dental crowns: Design parametric study Sukirti Dhital a, Camila Rodrigues b, c, Yu Zhang b, Jeongho Kim a, * a
Department of Civil and Environmental Engineering, University of Connecticut, 261 Glenbrook Rd, U-3037, Storrs, CT, 06269, USA Department of Biomaterials and Biomimetics, New York University College of Dentistry, 433 First Avenue, New York, NY, 10010, USA c Post-Graduation Program in Dental Sciences, Federal University of Santa Maria, 1184 Floriano Peixoto St, Santa Maria, RS, 97015-372, Brazil b
A R T I C L E I N F O
A B S T R A C T
Keywords: Viscoelastic finite elements Dental ceramics Porcelain-veneered zirconia Residual stresses Transient stresses Coefficient of thermal contraction
Porcelain-veneered zirconia (PVZ) are one of the popular choice for crown restorations. Veneer layer of these dental restorations, however, is susceptible to chipping and delamination due to the development of transient and residual stresses during the cooling phase of veneer firing. The aim of this study is to elucidate the effect of material property mismatch, veneer to core thickness ratio, and cooling rate on these transient and residual stresses of PVZ restorations. Three-dimensional viscoelastic finite element modelling (VFEM) was performed. The VFEM model was developed using the UEXPAN subroutine in ABAQUS software and was validated for transient and residual stresses in a sandwich seal problem with experimental data available. A good agreement between the simulated VFEM results and experimental data was obtained. Using validated VFEM, two PVZ systems (PM9/ zirconia and ZirPress/zirconia), three veneer to core thickness ratios (2:1, 1:1 and 1:2), and two cooling rates controlled slow cooling at 1.74E-5 W/mm2� C (i.e. ~30 � C/min) and fast bench cooling at 1.74E-4 W/mm2� C (i.e. ~300 � C/min) were used. The results showed that PM9/zirconia has smaller thermal contraction mismatch, resulting in lesser residual stress (33.36 MPa) as compared to ZirPress/zirconia (37.94 MPa) for controlled cooling and 2:1 veneer to core ratio. In addition, in both systems with the decrease in veneer thickness, we observed a decrease in residual stresses developed. We also observed some effect of cooling rate on residual stresses. The controlled cooling resulted in lower residual stress (24.35 MPa) for PM9/zirconia with a 1:1 veneer to core thickness ratio as compared to bench cooling (28.04 MPa). The effect of cooling rate was more evident on transient stresses. For instance, in the PM9/zirconia with 1:1 thickness ratio model, the difference in transient stresses was 9.93 MPa between controlled and bench cooling. Therefore, properties such as elastic modulus and coefficient of thermal contraction (CTC), as well as the thickness ratio and cooling rate all play an important role in transient and residual stresses developed in the studied ceramic systems.
1. Introduction The use of porcelain-veneered Yttria Tetragonal Zirconia Polycrystal (Y-TZP) crowns and fixed dental prostheses has increased significantly in past years. In order to mimic the enamel characteristics, a porcelain over layer is usually applied to the zirconia core. Though they exhibit great aesthetics as well as chemical durability and biocompatibility (DeHoff et al., 2006), they are quite susceptible to chipping of the por celain veneer layer (Gherlone et al., 2014; Nicolaisen et al., 2016). These dental structures are fabricated by fusing two materials together at high temperatures. This manufacturing process results in the development of transient and residual stresses during the cooling phase. The stored
tensile residual stresses results in chipping or delamination of the veneer (Pang et al., 2015; Christensen, 2009; Denry and Kelly, 2008; Sailer et al., 2007; Zhang et al., 2012; Zhang et al., 2013. Once the veneered porcelain chips, the crown often needs to be repaired or replaced and this process can be cumbersome. Hence, researchers are driven to find an optimal combination for core and veneer layer materials to minimize the deleterious tensile residual stresses, thereby ensuring longevity. There are many factors that give rise to residual thermal stresses. Some of which are temperature-dependent material properties, geometry, layer thickness, thermal contraction mismatch and cooling rate. The residual stresses in porcelain-veneered zirconia systems have been widely studied using flat samples (Swain, 2009; Taskonak et al.,
* Corresponding author. E-mail address: [email protected] (J. Kim). https://doi.org/10.1016/j.jmbbm.2019.103545 Received 25 September 2019; Received in revised form 14 November 2019; Accepted 18 November 2019 Available online 19 November 2019 1751-6161/© 2019 Elsevier Ltd. All rights reserved.
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
analyzing strips or simpler models (Asaoka and Tesk, 1990; Dehoff and Anusavice, 1989a; Dehoff and Anusavice, 1989b; Dehoff and Anusavice, 1992; Asaoka et al., 1992; Nakatsuka and Anusavice, 1997). Several studies were carried out to understand this behavior of dental porcelains (DeHoff and Anusavice, 2004a; DeHoff and Anusavice, 2004b; Fragassa, 2016). Fragassa, 2016 investigated the viscoelastic behavior of ceramic materials using few of the available commercial finite element codes such as ANSYS, LS-DYNA, ABAQUS and MARC. Commercial software, i. e. ANSYS, was used to perform viscoelastic analysis (DeHoff and Anu savice, 2004a; DeHoff and Anusavice, 2004b; DeHoff et al., 2006, 2008). Unlike these, we developed user subroutines such as UEXPAN and UTRS for ABAQUS to better capture the viscoelastic behavior of porcelain throughout the entire firing cycle by defining shear relaxation depen dent coefficient of thermal contraction, which is not available in com mercial software like Abaqus. Knowledge of the residual stress states in porcelain-veneered struc tures are important. However, understanding the entire stress history is also very critical. Since once the transient stress exceeds the stress threshold of fracture, micro-cracks and flaws could form. Such defects can diminish the fracture resistance of ceramic restorations. Unfortu nately, in the body of dental literature, we could find only one study that analyzed transient stresses (Benetti et al., 2014). Benetti et al. have proposed that transient stresses results in the nucleation of cracks inside the porcelain material during the cooling phase, which might result in fracture under different loadings. However, Benetti et al. fell short of predicting the exact transient stresses because in their study the high temperature properties of the materials were not experimentally deter mined, rather extrapolated from the manufacturers’ data sheet or literature values of similar materials. . As alluded to in the preceding paragraph, the major advantage of our newly developed subroutines is its ability to capture the entire stress evolution during the firing cycle, especially in the cooling phase. This way, any significant transient stress can be captured and analyzed. Knowledge gained can be used as a guideline to refine the processing parameters and optimize the design of restorations. In this paper, we aim to further understand the cause of transient and residual stresses in PVZ systems. The two subroutines were first vali dated by a well-studied sandwich seal system. Such a sandwich seal system allows for real-time experimental transient stress analysis (Scherer and Rekhson, 1982). Then, the validated model was used to determine both transient and residual stresses in two PVZ systems, namely (PM9/zirconia and ZirPress/zirconia), in a stylized 3D crown. We also altered the veneer thickness with three veneer to core ratios: 2:1, 1:1 and 1:2 while keeping the total thickness constant. Further we subjected the model to two cooling rates that is controlled slow cooling (1.74E-5 W/mm2 � C or on a more familiar term 30 � C/min) and fast bench cooling (1.74E-4 W/mm2 � C or 300 � C/min). The effect of all these parameters on transient and residual stresses obtained on the veneer layer is studied.
Table 1 Materials combination and veneer to core ratio. Model Number
Veneer to Core Ratio
Core
Veneer
Model 1 Model 2 Model 3 Model 4 Model 5 Model 6
2:1 1:1 1:2 2:1 1:1 1:2
zirconia zirconia zirconia zirconia zirconia zirconia
PM9 PM9 PM9 ZirPress ZirPress ZirPress
Fig. 1. Finite element model of (a) full dental crown and b) 1/4th simulated model. Light blue represents the zirconia core and grey represents the porce lain veneer. Table 2 Material properties at the room temperature (25 � C). Ceramics PM9 ZirPress Zirconia
Young’s Modulus, E (MPa) 60,000 70,000 206,000
Poisson’s Ratio, ν 0.21 0.21 0.3
Density, ρ (kg/ mm3) 2.45E-06 2.54E-06 6.10E-06
Thermal Conductivity (W/mm ⁰C) 1.11E-03 1.15E-03 2.92E-03
Specific Heat, (J/ kg⁰C) 821 886 466
2005, 2008; Wendler et al., 2015; Choi et al., 2011, Tanaka et al., 2016). However, thermal stresses observed in flat samples do not represent the stresses in complex geometries, such as dental crowns. In addition, many experimental techniques have been reported for measuring these re sidual stresses, e.g. Vickers indentations (Baldassarri et al., 2012; Tanaka et al., 2016; Choi et al., 2011; Wendler et al., 2015), birefrin gence (Wendler et al., 2016; Belli et al., 2012; Tholey et al., 2011), hole-drilling (Mainjot et al., 2011), and fracture mechanics approaches (Taskonak et al., 2005, 2008). Such methods require highly polished surfaces, thin slices of the restoration with sufficient translucency, or complex samples preparation, which alter the residual stress states and do not reflect the actual stress profiles in an anatomically-correct restoration. To accurately determine residual stresses in anatomically-correct restorations, the numerical simulation technique finite element method has been adopted. For instance, Meira et al., 2013; Tanaka et al., 2016; and Zhang et al. determined the residual stresses in a PVZ crown by the linear elastic finite element method (LFEM), using the arbitrary values of liquid thermal contraction co efficients of porcelain. Oh et al., 2002 also used this analytical process to calculate the strength of FDPs, which they then compared with the experimental values. They observed that the calculated residual stresses were much more than those obtained experimentally and this might be due to the negligence of the viscoelastic behavior of porcelain. Hence, we can conclude that these models fail to accurately capture the stress profile in porcelain as the viscoelastic behavior of these materials above the softening temperatures has been neglected. The importance of the viscoelastic behavior has been elucidated by
2. Materials and methods 2.1. Methods Three-dimensional finite element analysis is carried out in ABAQUS. Two PVZ systems [PM9 (Vita Zahnfabrik, Bad S€ ackingen, Germany) or ZirPress (Ivoclar Vivadent, Schaan, Liechtenstein) porcelains and a zir conia core material (Zpex, Tosoh, Tokio, Japan)], three veneer to core thickness ratios (2:1, 1:1, or 1:2), and two cooling rates [controlled slow (1.74E-5 W/mm2 � C) or fast bench (1.74E-4 W/mm2 � C], were used (Table 1). The maximum thickness of the specimen was 2.1 mm near central fossa, and the veneer and core thicknesses depend upon the ratio opted. A three-dimensional (3-D) FEA model was analyzed. Owing to the symmetric geometry and boundary conditions, only 1/4th of the model (see Fig. 1) was analyzed. The simulation was carried out in two steps: 2
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 4. Specific Heat at different temperatures.
Fig. 2. Density at different temperatures.
Fig. 5. Young’s Modulus at different temperatures. Fig. 3. Conductivity at different temperatures.
heat transfer analysis to capture temperature history, and stress analysis with a thermal load. 2.2. Material properties Material properties, such as Young’s modulus, Poisson’s ratio, den sity, thermal conductivity and specific heat, at the room temperature are given in Table 2, whereas the corresponding properties at high tem peratures are plotted in Figs. 2–6, respectively. The coefficient of thermal contraction, and subsequently, the density were measured by dilatometer (L75 Platinum Series, Linseis USA, Princeton, NJ). The thermal conductivity and specific heat were deter mined using the laser flash technique (LFA 457 MicroFlash, Netzsch, Selb, Germany). The room temperature Young’s modulus and Poisson’s ratio were measured at the Edward Orton Jr Ceramic Foundation through a paid service. However, we had some difficulties to obtain the high temperature Young’s Modulus values of porcelains. This kind of tests require large specimens (typically 100 mm � 20 mm X 5 mm), which was not possible to produce with the materials we used in this study. Therefore, we assumed that all the porcelains follow similar
Fig. 6. Thermal contraction coefficient at different temperature. 3
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Table 3 Constants for thermal contraction of dental ceramics (αg ¼ a þ bT þ cT2, αl ¼ d þ eT). Materials
αe (mm/mm⁰C) a (1/⁰C) 10
PM9 ZirPress zirconia
6
8.8 8 10.066
αL (mm/mm⁰C)
b (1/⁰C)2 10
9
2.7186 5.6 1.2856
c (1/⁰C)3 10
R2
12
2 4.2 0.6697
d (1/⁰C) 10
0.997 0.978 1
5
8.716 7.184 NA
e (1/⁰C)2 10
7
1.528 1.52 NA
Fig. 8. Schematic of sandwich seal.
Fig. 7. Normalized shear relaxation function at various high temperatures (DeHoff et al., 2006).
trends at high temperatures and adopted the trend of VM9 (Kim et al., 2018) as the experimental data for this material was available. We also utilized the high temperature modulus data of zirconia (3Y-TZP) from the same reference (Kim et al., 2018). For coefficient of thermal contraction (CTCs), as shown in Fig. 6, the CTCs for veneers are divided into solid and liquid CTCs, αg(T) and αl(T) given by: (1)
αg ¼ aþbT þ cT2 and αl ¼ d þ eT
where T is in Celsius, and the five coefficients are given in Table 3. Then, for CTCs of the core, thermal contraction strain versus tem perature data were fit to a third-degree polynomial relation, from regression analysis, from which αg(T) was determined, using Eq. (1). For the veneer, αl(T) was determined by fitting the contraction data with a quadratic equation from the temperature at which creep first occurs on the heating curve to the glass transition temperature, Tg.
Fig. 9. Longitudinal stress (σyy) for a G-11/Alumina sandwich seal cooled at 3 � K/min ¼ 0.05 C/s from 618 to 20 � C.
Z
τðtÞ ¼
t
GR ðt
(2)
sÞ_γ ðsÞds
0
2.3. Viscoelastic finite element analysis
where GR(t) is the shear relaxation modulus, t is the present time, and s is the past time. The normalized shear relaxation modulus can be intro duced by (Kim et al., 2018):
The time dependent viscoelastic properties of dental ceramics, sub jected to cooling from sintering temperature to room temperature, can be introduced into ABAQUS using the viscoelastic option in the material properties which then can be used to determine residual stresses. The shear stress τ (t) in a viscoelastic material model, when subjected to small shear strain γ (t), is given by (Kim et al., 2018):
(3)
gR ðtÞ ¼ GR ðtÞ = G0
where G0 is the instantaneous shear modulus. The volumetric behavior can be written as:
Table 4 Viscoelastic coefficients for Prony series used for veneer. Temp (oC)
g1
g2
g3
g4
τ1
τ2
τ3
τ4
700
0.9960
0.0030
0.0006
0.0004
0.01316
0.10000
0.0050
0.0030
4
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
N X
gR ðtÞ ¼ 1
� gPi 1
G e t=τi
�
(5)
i¼1
where N is the number of terms (N ¼ 4 used in this study), gPi are co efficients, and τGi is shear relaxation time (Table 4) as given at the reference temperature of 700 � C (DeHoff et al., 2006). Since we have taken the reference temperature as the initial temperature and we are cooling the two different material combinations of PM9/zirconia and ZirPress/zirconia from 754 � C and 657 � C respectively, we assumed that these materials have same relaxation coefficients and relaxation time at the respective cooling temperature as those given by DeHoff at 700 � C. The viscoelastic materials show thermorheologically simple behavior in which the normalized shear relaxation function when plotted against log10t has same shape at each temperature and hence can be shifted along the x-axis, thereby resulting in a master curve at a reference temperature. Fig. 7 shows the thermorheological behavior of PM9 using the coefficients in Table 4. The leftmost solid curve is the master curve at 754 � C. The shift function, which is the shift along x-axis, that represents temperature dependence of shear relaxation functions can be defined using UTRS user subroutine in ABAQUS. This subroutine defines the Tool-Narayanaswamy approximation, which takes the form (Tool, 1946; Narayanaswamy, 1971). � � �� H 1 x 1 x ln A ¼ (6) R Tref TðtÞ Tf ðtÞ
Fig. 10. Stress for the G-11/Alumina sandwich seal when cooled at 3 K/min ¼ � 0.05 C/s, and held 4 h at 465 � C.
Z pðtÞ ¼
K0
t
kR ðt 0
sÞ_εvol ðsÞds
(4)
where p is the hydrostatic pressure, K0 is the instantaneous elastic bulk modulus, kR(t) is the normalized bulk relaxation modulus assumed to be 1.0 due to its minor influence in most stress states, and εvol is the volume strain. The normalized shear relaxation function can be described by the use of a four-term Prony series in ABAQUS using the following equation:
where Tref ¼ 1027.15 in Kelvin (754 in Celsius), H is the activation energy, R is the ideal gas constant (H/R ¼ 46400 is used), T(t) is the temperature at time t and Tf(t) is the fictive temperature at time t, and x is a material constant between 0 and 1 (x ¼ 0.3 used, e.g. (Fragassa,
Fig. 11. Temperature and Stress Profile for 3-D Slow Cooling in a) Model 1 and b) Model 4. 5
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
3. Results This section presents a validation example which shows validity of the present viscoelastic finite element procedure in predicting strains and stresses. The validated model is then used to analyze 3-D dental models in order to provide a map of stresses for different material combinations, cooling rates and veneer to core ratios. 3.1. Sandwich seal under cooling Viscoelastic finite element analyses was performed to simulate the contraction behavior of a G-11/Alumina sandwich seal cooled at 3 K/ min (Scherer and Rekhson, 1982). The authors extended the Nar ayanswamy’s model (Narayanaswamy, 1971) to accurately capture the viscoelastic behavior of the glass and then validated their formula with sandwich seal experiment. In addition, to demonstrate the applicability of our codes (UEXPAN and UTRS) in ABAQUS, we modelled the same G-11/Alumina sandwich seal (Fig. 8) subjected to cooling at 3 K/min. The following properties were used for Alumina: EAl ¼ 3.73 � 105 MPa, ν ¼ 0.3 and for G-11 glass: H/R ¼ 64500 (K) where, H/R ¼ ratio of activation enthalpy for relaxation process (H) and ideal gas constant (R) x ¼ 0.53 T ¼ 618� C ¼ Initial Temperature �
αg ðTÞ ¼ 64:7 þ 0:02T � 10 αl ðTÞ ¼ 3:43 � 10
5
7
with T in Celsius
with T in Celsius
EG11 ¼ 7.26 � 104 MPa, ν ¼ 0.3 The ratio of thickness of glass to thickness of seal was 13. The obtained FEA results, for the longitudinal stress (σyy) in this G11/Alumina model were compared with the experimental results from Scherer and Rekhson, 1982) as shown in Fig. 9. The results showed an excellent agreement with difference less than 0.7 MPa over the entire temperature history. As shown in Fig. 9, the stress increased to a maximum value upon initial cooling, but then decreased owing to the viscoelastic properties. This kind of highly nonlinear stress response cannot be captured by elastic analysis. Further, to show the effectiveness of our code in order to capture any arbitrary thermal history, including isothermal holds, we subjected the same G11-Alumina sandwich seal to a cooling rate of 3 K/min, held it for 4 h at 465 � C and then cooled it to the room temperature at the same rate. Fig. 10 shows the stress for this condition, compared with the experimental results present in Scherer and Rekhson, 1982). The agreement between the FEA predictions and experimental results was quite good with difference less than 0.6 MPa over the entire temperature history. During isothermal hold, one could observe that the stress kept increasing to a maximum value and then diminished sharply afterwards. This sandwich seal modelling showed the effectiveness of our codes in accurately capturing the viscoelastic behavior, and hence we used them to analyze ou0pr veneering ceramics.
Fig. 12. Model 1 (PM9/zirconia): Temperature; a) Temperature gradient in slow cooling b) Temperature history in Slow and Fast Cooling; Nodes repre sented are PT1-Porcelain top node, ZB1- zirconia top node and ZB2-zirconia bottom node.
2016)). The logarithmic function lnA (T, Tref, t) shifts the shear relaxa tion function along the horizontal axis for various temperatures. As the data for volume relaxation are not available and in case of glasses it has been found that the volume relaxation time is 1/4 to 1/20 of shear relaxation time, we assumed the volume relaxation time in the ceramics to be 1/10 of the shear relaxation time. Also studies have shown that volume relaxation time has almost no effect on residual stresses (DeHoff et al., 2006). As we know that fictive temperature needs to be considered during a viscoelastic analysis specially in the transition region, and Markovsky, 1984 proposed the algorithm for this fictive temperature Tf(t) which has been coded in UEXPAN. The thermal strain is given for the core and the veneer, respectively, by Z T Z T Z Tf � (7) εc ¼ αg ðTÞdT and εv ¼ αg ðTÞdT þ αl Tf dT To
Tf
3.2. Veneering ceramics
To
Viscoelastic finite element analyses were run to simulate the contraction behavior of a dilatometer specimen subjected to convective cooling (DeHoff et al., 2006). We analyzed two different material combinations, each of which were subjected to controlled cooling
where αg(T) and αl (T) are plotted in Fig. 5.
6
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 13. Residual stresses in slow cooling: a) Model 1: PM9/zirconia; b) Model 4: ZirPress/zirconia; c) Model 2: PM9/zirconia; d) Model 5: ZirPress/zirconia; e) Model 3:PM9/zirconia; f) Model 6: ZirPress/zirconia. Nodes represented are PT1-Porcelain top node, PZI- Porcelain zirconia Interface node and PB1- Porcelain Bottom Node.
(1.74E-5 W/mm2 oC) and bench cooling (1.74E-4 W/mm2 oC) on the crown surfaces exposed to air. The material combinations are as mentioned in Table 1. Fig. 11 shows the 3D temperature and stress profile for 2:1 veneer to core PM9/zirconia and ZirPress/zirconia models when subjected to control cooling. Since, this is an axisymmetric model, the temperature and stress contours are uniform through the face, as shown in the fig ures. Hence, for ease of study, only one face of the contour images are provided in Figs. 12, 13 and 15. As mentioned before, in the first step, 3-D transient heat analysis was carried out. This gave us the temperature profile of the crown when subjected to both controlled and bench cooling (see Fig. 12). As ex pected, bench cooling causes more rapid cooling than controlled cooling.
All the models follow similar temperature profile as shown in Fig. 12. After inputting these temperature profiles as load, stress analysis was performed. The maximum principal residual and transient (timedependent) stresses are plotted, respectively, in Figs. 13 and 14. Fig. 13 depicts the residual stress contour at the end of analysis, when subjected to slow cooling. It is clear that the maximum principal tensile stress in the porcelain layer, for all of the models, occur at the cusp interface of porcelain and zirconia layers. The location for maximum principal stress in the zirconia layer, however, varies from model to model. As zirconia has very high strength, we focus on tensile residual stresses developed in the porcelain layer. The maximum prin cipal residual stress contours shows that the residual stress developed in ZirPress layer is more as compared to its PM9 counterpart. For instance in PM9/zirconia models that is Model 1, Model 2 and Model 3, tensile 7
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 14. Transient stress in slow cooling: a) Model 1: PM9/zirconia; b) Model 4: ZirPress/zirconia; c) Model 2: PM9/zirconia; d) Model 5: ZirPress/zirconia; e) Model 3:PM9/zirconia; f) Model 6: ZirPress/zirconia.
residual stresses developed in the PM9 layer are 33.36 MPa, 24.35 MPa and 17.40 MPa respectively, whereas for ZirPress/zirconia models, that is Model 4, Model 5 and Model 6, tensile residual stresses developed are 37.94 MPa, 27.90 MPa and 20.02 MPa respectively in the ZirPress layer. From both the PM9 and ZirPress models, we can see that with the decrease in the thickness of veneer layer, the tensile residual stress developed in it decreases. For 2:1 veneer to core PM9/zirconia model, the residual stress developed in PM9 veneer is 33.36 MPa which is almost 10 MPa more than the residual stress in the same veneer material with 1:1 veneer to core ratio. Similar pattern can be observed in the
remaining PM9/zirconia model and all of the ZirPress/zirconia models. Besides the residual stress, the tensile transient stress developed during the cooling phase, is also one of our concern. In order to understand it better, the history of the maximum principal stress in the porcelain layer needs to be reviewed. As the abrupt rise and fall of stresses during the entire simulation period may have effect on the material behavior, one needs to look at the stress profile developed along the entire time period. Few of the key features are mentioned here. From the graphs, we could observe that for PM9/zirconia Model 1, sudden peak occurred at the node PT1 with 8
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 15. Residual stress in fast cooling: a) Model 1: PM9/zirconia; b) Model 4: ZirPress/zirconia; c) Model 2: PM9/zirconia; d) Model 5: ZirPress/zirconia; e) Model 3: PM9/zirconia; f) Model 6: ZirPress/zirconia. Nodes represented are PT1-Porcelain top node, PZI- Porcelain zirconia Interface node and PB1- Porcelain Bottom Node.
value 4.516 MPa at 37.5 s, which then gradually decreased to obtain a stable stress value. Similarly, for Model 2, at the same point and at 35 s, 2.128 MPa maximum principal tensile stress is obtained and for Model 3, at 33.5 s, 2.589 MPa of stress is observed. The maximum principal stress profile for ZirPress/zirconia is similar to that of PM9 models as shown in Fig. 14. For Model 4, at 45 s and at node PT1, the maximum principal stress obtained is 8.012 MPa, at 43 s for Model 5, the stress obtained is 4.481 MPa and 5.673 MPa for Model 6 at 42 s. The models were also subjected to fast cooling. The maximum principal residual and transient (time-dependent) stresses are plotted, respectively, in Figs. 15 and 16. Similar to slow cooling, Fig. 15 shows the residual stress contour at the end of analysis for fast cooling, in all the six models. Here too, in all cases, the maximum value of the principal tensile stress, of the porcelain layer, occurs at the cusp interface of the two layers. The tensile residual stress developed in ZirPress layer, is more when compared to PM9 models. The maximum value of the maximum principal residual stress developed in the PM9 layer are 36.36 MPa, 28.04 MPa and 20.66 MPa for Model 1, Model 2 and Model 3 respectively. Similarly, for ZirPress/zirconia model, the maximum
residual stress obtained are 39.12 MPa, 30.84 MPa and 22.66 MPa for Model 4, Model 5 and Model 6 respectively. Irrespective of the cooling rate, we can see that with decrease in the veneer thickness (as indicated by the model numbers) there is decrease in the tensile stress produced in the veneer layer, in both the materials. However, as one can conclude, tensile residual stresses obtained due to fast cooling are higher than that due to slow cooling in both the materials. For instance, for Model 1 and PM9 veneer layer, the tensile residual stress obtained due to slow cooling is 33.36 MPa and fast cooling is 36.36 MPa, which shows the increase in the residual stress. In this case, we also can observe the abrupt rise and fall in the transient stress profile for models under fast cooling, as shown in Fig. 16. However, the rise in this case is more lucid than that seen in Fig. 14. For PM9/zirconia model, we can observe the sudden peak, at node PT1, for Model 1 is 18.786 MPa at 4.5 s, for Model 2 at 4.5 s the value obtained is 12.057 MPa and for Model 3 at 4.5 s the transient tensile stress obtained is 17.811 MPa. Now, in the case of ZirPress/zirconia Model, for Model 4, at 5.5 s, the tensile transient stress obtained is 23.187 MPa, for Model 5 at 5 s we obtained 21.435 MPa and for Model 6, at 5 s, the transient 9
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 16. Transient stress in fast cooling: a) Model 1: PM9/zirconia; b) Model 4: ZirPress/zirconia; c) Model 2: PM9/zirconia; d) Model 5: ZirPress/zirconia; e) Model 3:PM9/zirconia; f) Model 6: ZirPress/zirconia.
stress obtained is 25.802 MPa. Analyzing all the 12 model combination, for different materials, thickness ratios and cooling rate, the maximum principal tensile residual stress is summarized in the bar graphs in Fig. 17.
able to accurately capture the viscoelastic behavior of materials, por celain in our case (Kim et al., 2018). This might be due to the fact that for states above glass transition temperature of porcelain, Tg, the Young’s modulus is dependent on both temperature as well as time, which means that even when the temperature is constant the modulus decreases exponentially with time and this phenomena cannot be addressed by LFEM. This results into larger residual stress obtained as compared to that from viscoelastic FEM, even when the coefficient of thermal contraction for porcelain in liquid state, used in VFEM, is much larger
4. Discussion Many researchers have used linear elastic FEM to determine residual stresses in different material systems. However, these LFEMs are not 10
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
for the sandwich seal system is compared against the test data which shows remarkably good match, thereby validating our model. This model is then used to estimate residual stresses and perform parametric study on various dental ceramic combinations with different thickness ratios, material properties and cooling rates. It has been shown through our analysis that the veneer thickness affects residual stresses developed in dental models. We analyzed three different veneer to core ratios, which are 2:1, 1:1 and 1:2. From the results, for both the material combinations, we could clearly observe that with the decrease in the veneer thickness, tensile residual stresses developed in porcelain reduced. For instance, in PM9/zirconia model, the residual stress developed in PM9 was about 33.36 MPa with 2:1 ratio, 24.35 MPa for 1:1 ratio and 17.40 MPa for 1:2 ratio, under controlled cooling. Similar pattern can be observed under bench cool ing. For ZirPress/zirconia, the results obtained were 37.94 MPa, 27.90 MPa and 20.02 MPa for 2:1, 1:1 and 1:2 ratios respectively, under controlled cooling. Thus, we could conclude that, with reduced porce lain layer thickness, the reduction in thermal gradient results in lower stresses. As we have analyzed two different material combinations, we could also see the effect of material properties in residual stresses developed. The maximum tensile stress of 33.36 MPa and 36.36 MPa occurs for controlled and bench cooling respectively, for PM9/zirconia combina tion, and similarly the maximum tensile stress of 37.94 MPa and 39.12 MPa occurs for ZirPress/zirconia model, for the same veneer core thickness ratio. To understand this difference in results, we studied the combined effect of young’s moduli and thermal contraction coefficients on the different material models analyzed. The CTC of PM9 is closer to zirconia, as compared to Zirpress. Also, the young’s modulus of PM9 is lower than that of ZirPress. Hence, residual stresses in ZirPress veneer layer are higher than those in PM9. The transient and residual stresses are a combined effect of time and temperature dependent properties such as modulus, CTC and viscoelasticity. The residual stress in zirconia is not of major concern here, as they have very high strength. However, for porcelain, with a smaller CTC mismatch and larger modulus mismatch, we obtain lower residual stresses, which is in agreement with results showed on previous studies (Kim et al., 2018). Besides the core and veneer thickness and material properties, the cooling rate also has an effect on the amount of residual and transient stresses developed. Controlled cooling produced less tensile stress as compared to bench cooling, in all combinations of materials and thick ness ratios. The fast change in temperature, in bench cooling, brings about larger temperature difference, which in turn results in higher transient and residual stresses, as verified from the models analyzed. The difference between tensile residual stresses obtained for controlled and bench cooling for Model 1, Model 2, Model 3, Model 4, Model 5 and Model 6 is 3.10 MPa, 3.69 MPa, 3.26 MPa, 1.18 MPa, 2.94 MPa and 2.64 MPa respectively. The difference in residual stresses for controlled and bench cooling might be less, which is similar to the findings of Tholey et al. (2011). However, this gap is amplified when one studies transient stresses. For the same order of six models, the difference in tensile transient stress obtained are 14.27 MPa, 9.93 MPa, 15.22 MPa, 15.18 MPa, 16.96 MPa and 20.13 MPa respectively, thereby completely elucidating the effect of cooling rate on the amount of stress produced. However, these high transient stresses occur only for a short period, but they result into high strain rates in our viscoelastic analysis, thereby increasing the susceptibility of the porcelain layer to fractures, cracking and chipping. The models analyzed tried to capture the actual conditions, however, our study had few limitations. As it is difficult to measure the Young’s Moduli of different porcelain at high temperatures, required in our study, we assumed that both PM9 and ZirPress follow exactly the same trend as that followed by VM9. However, the room temperature modulus were obtained from laboratory measurements. Besides this we also have assumed the shear relaxation coefficients of these porcelains at the required temperatures to be same as the one in DeHoff et al. (2006)
Fig. 17. Residual Stress Summary in (a) porcelain; (b) zirconia.
than that in a solid state owing to the large thermal strain. However, the thermal contraction coefficients of liquid porcelain used in LFEM are quite small, normally 3–4 times smaller than actual. VFEM in the other hand uses the temperature and time-dependent phenomenon of stress relaxation to get small stress value as compared to LFEM. The depen dence of shear relaxation modulus on the time-dependent temperature in a specimen is described by a simple shift function (Eq. (5)), supported by experimentally determined creep data (DeHoff and Anusavice, 2004a). The shear relaxation modulus data presented in Fig. 6 is uni versal for porcelains at any given temperature under cooling rates examined here. The viscoelastic finite element model that has been developed, considering the above mentioned points, is used to capture the stress profiles in two core-veneer systems, that is, PM9/zirconia and ZirPress/zirconia. To ensure that the developed VFEM models accu rately captures the obtained residual stresses, they have been validated in a sandwich seal system. The obtained stress results from VFEM model 11
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
as the data for this was also not available. Addition to shear relaxation coefficients, the volume relaxation coefficients too were assumed to be 1/10 of shear relaxation coefficients (DeHoff et al., 2006), but their influence on residual stresses is quite less. Finally, damage or crack formation is not considered during this viscoelastic analysis. Further study on the effect of crack initiation and growth on deformations and residual stresses will require damage analysis in conjunction with the present analysis.
for both the materials, residual and transient stresses developed in it also decreases. � The maximum principal stress developed in the porcelain veneer and zirconia core is less for controlled slow cooling as compared to rapid bench cooling. However, the cooling effect is more evident in the transient than residual stresses and these transient stresses contribute in producing micro-cracks, which ultimately fractures. Declaration of competing interest
5. Conclusions
All authors declare no conflict of interest.
From the results obtained, we could draw the following conclusions: � The material properties such as Young’s Modulus and Coefficient of Thermal Contraction (CTC) play a role in the amount of stress developed. If the CTC mismatch between the core and veneer is less and the porcelain modulus is lower, then the residual stress devel oped in both the layers is reduced, irrespective of the cooling rate. � The thickness of core and veneer layer play a role in the amount of stress developed. With the decrease in the thickness of veneer layer,
Acknowledgements This work was sponsored by funding from the United States National Institute of Dental & Craniofacial Research, National Institutes of Health (Grants Nos. 1R01 DE026279 and R01DE026772). The authors are also thankful to the Brazilian Federal Agency for Support and Evaluation of Graduate Education (CAPES) (finance code 001).
Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi.org/10.1016/j.jmbbm.2019.103545. Appendix B 1. Subroutine UEXPAN is called a. Density and temperature dependent elastic properties are read from the.inp file, for porcelain and coping, which are defined by the user. Specific heat for zirconia are also user defined. b. If the material is “PM9 or ZirPress” then: i. Using Narayanswami form, the ratio between relaxation time and reference (volume) relaxation time is determined as: � �� � τi ΔH 1 x 1 x ¼ exp R Tref T Tf τis where, τi ¼ structural relaxation time
τis ¼ volume relaxation time, which is taken as 1/10th of shear relaxation time H ¼ activation energy x ¼ experimentally defined constant⊲ Tref ¼ reference temperature T ¼ current temperature at end of increment Tf ¼ fictive temperature
ii. The partial fictive temperatures Tfi ðtÞ, the fictive temperature Tf ðt) and the difference between fictive temperatures (dTf) of two sequential steps are stored as state variables and calculated as: � � Tfi ðtÞ þ T � Δt τi � � Tfi ðtÞ ¼ where; Δt ¼ time increment 1 þ Δt τi �� X Tf t ¼ ci Tfi ðtÞ where; ci are the constants defined at reference temperature dTf ¼ Tf
Tf ðprevious incrementÞ
iii. The thermal contraction coefficients for glass and liquid states respectively, are obtained as: ag ¼ ag1 þ ag2 �T þ ag3 �T2 agf ¼ ag1 þ ag2 �Tf þ ag3 �T2f 12
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
alf ¼ al1 þ al2 �Tf þ al3 �T2f where, ag1 ; ag2 ; ag3 ; al1 ; al2 are the constants for thermal contraction coefficients in which subscript g is for glass and l for liquid state. iv. Finally coefficient of expansion for porcelain is calculated as: � expansion ¼ ag �temp: increment þ alf agf �fictive temp increment
2. After executing UEXPAN subroutine, UTRS is called to account for the temperature-time shift for a time domain viscoelastic analysis a. As porcelain is the only viscoelastic material, Prony series is defined by the user in the.inp file and read when required. b. The experimentally defined constant, x, is taken and activation energy Δ H, reference temperature Tref along with ideal gas constant R are defined c. Fictive temperature (Tf ) and the difference between current and previous increment fictive temperatures (dTf Þ are called from UEXPAN sub routine as Statev. d. The shift along the horizontal axis, shift factor, are determined using Tool-Narayanswamy form. The shift factor is determined at the beginning and end of the increment using: � �� � ΔH 1 x 1 x shift ¼ exp 0 0 R Tref Tðt Þ Tf ðt Þ where; Tðt Þ is the temperature at any time t and Tf ðt Þis the fictive temperature at t 0
0
0
References
0
Behav. Biomed. Mater. 82, 202–209. https://doi.org/10.1016/j. jmbbm.2018.03.020. Mainjot, A.K., Schajer, G.S., Vanheusden, A.J., Sadoun, M.J., 2011. Residual stress measurement in veneering ceramic by hole-drilling. Dent. Mater. 27, 439–444. https://doi.org/10.1016/j.dental.2010.12.002. Markovsky, A., 1984. An efficient and stable Algorithm for calculating fictive temperatures effect of variation in flame velocity profile on particle deposition rate. Mater. Sci. 56–57. Meira, J.B.C., Reis, B.R., Tanaka, C.B., Ballester, R.Y., Cesar, P.F., Versluis, A., Swain, M. V., 2013. Residual stresses in Y-TZP crowns due to changes in the thermal contraction coefficient of veneers. Dent. Mater. 29, 594–601. https://doi.org/ 10.1016/j.dental.2013.03.012. Nakatsuka, A., Anusavice, K.J., 1997. Finite element analysis of stress distribution in porcelain discs. J. Mater. Sci. 32, 3621–3627. https://doi.org/10.1023/A: 1018678313122. Narayanaswamy, O.S., 1971. A model of structural relaxation in glass. J. Am. Ceram. Soc. 54, 491–498. https://doi.org/10.1111/j.1151-2916.1971.tb12186.x. Nicolaisen, M., Bahrami, G., Schropp, L., Isidor, F., 2016. Comparison of metal-ceramic and all-ceramic three-unit posterior fixed dental prostheses: a 3-year randomized clinical trial. Int. J. Prosthodont. 29, 259–264. https://doi.org/10.11607/ijp.4504. Oh, W., G€ otzen, N., Anusavice, K.J., 2002. Influence of connector design on fracture probability of ceramic fixed-partial dentures. J. Dent. Res. 81, 623–627. https://doi. org/10.1177/154405910208100909. Pang, Z., Chughtai, A., Sailer, I., Zhang, Y., 2015. A fractographic study of clinically retrieved zirconia-ceramic and metal-ceramic fixed dental prostheses. Dent. Mater. 31, 1198–1206. https://doi.org/10.1016/j.dental.2015.07.003. Sailer, I., Feher, A., Filser, F.T., Gauckler, L.J., Lüthy, H., H€ ammerle, C.H.F., 2007. Fiveyear clinical results of zirconia frameworks for posterior fixed partial dentures. Int. J. Prosthodont. 20, 383–388. Scherer, G.W., Rekhson, S.M., 1982. Viscoelastic-elastic composites: II, sandwich seal. J. Am. Ceram. Soc. 65, 399–406. https://doi.org/10.1111/j.1151-2916.1982. tb10492.x. Swain, M.V., 2009. Unstable cracking (chipping) of veneering porcelain on all-ceramic dental crowns and fixed partial dentures. Acta Biomater. 5, 1668–1677. https://doi. org/10.1016/j.actbio.2008.12.016. Tanaka, C.B., Harisha, H., Baldassarri, M., Wolff, M.S., Tong, H., Meira, J.B.C., Zhang, Y., 2016. Experimental and finite element study of residual thermal stresses in veneered Y-TZP structures. Ceram. Int. 42, 9214–9221. https://doi.org/10.1016/j. ceramint.2016.03.018. Taskonak, B., Borges, G.A., Mecholsky, J.J., Anusavice, K.J., Moore, B.K., Yan, J., 2008. The effects of viscoelastic parameters on residual stress development in a zirconia/ glass bilayer dental ceramic. Dent. Mater. 24, 1149–1155. https://doi.org/10.1016/ j.dental.2008.01.004. Taskonak, B., Mecholsky, J.J., Anusavice, K.J., 2005. Residual stresses in bilayer dental ceramics. Biomaterials 26, 3235–3241. https://doi.org/10.1016/j. biomaterials.2004.08.025. Tholey, M.J., Swain, M.V., Thiel, N., 2011. Thermal gradients and residual stresses in veneered Y-TZP frameworks. Dent. Mater. 27, 1102–1110. https://doi.org/10.1016/ j.dental.2011.08.001. Tool, A.Q., 1946. Relation between inelastic deformability and thermal expansion of glass in its annealing range. J. Am. Ceram. Soc. 29, 240–253. https://doi.org/ 10.1111/j.1151-2916.1946.tb11592.x.
Asaoka, K., Kuwayama, N., Tesk, J.A., 1992. Influence of tempering method on residual stress in dental porcelain. J. Dent. Res. 71, 1623–1627. https://doi.org/10.1177/ 00220345920710091501. Asaoka, K., Tesk, J.A., 1990. Transient and residual stress in a porcelain-metal strip. J. Dent. Res. 69, 463–469. https://doi.org/10.1177/00220345900690020901. Baldassarri, M., Stappert, C.F.J., Wolff, M.S., Thompson, V.P., Zhang, Y., 2012. Residual stresses in porcelain-veneered zirconia prostheses. Dent. Mater. 28 (8), 873–879. Belli, R., Monteiro, S., Baratieri, L.N., Katte, H., Petschelt, A., Lohbauer, U., 2012. A photoelastic assessment of residual stresses in zirconia-veneer crowns. J. Dent. Res. 91, 316–320. https://doi.org/10.1177/0022034511435100. Benetti, P., Kelly, J.R., Sanchez, M., Della Bona, A., 2014. Influence of thermal gradients on stress state of veneered restorations. Dent. Mater. 30, 554–563. https://doi.org/ 10.1016/j.dental.2014.02.020. Choi, J.E., Waddell, J.N., Swain, M.V., 2011. Pressed ceramics onto zirconia. Part 2: indentation fracture and influence of cooling rate on residual stresses. Dent. Mater. 27, 1111–1118. https://doi.org/10.1016/j.dental.2011.08.003. Christensen, G.J., 2009. Porcelain-fused-to-metal versus zirconia-based ceramic restorations, 2009. J. Am. Dent. Assoc. 140, 1036–1039. https://doi.org/10.14219/ jada.archive.2009.0316. Dehoff, P.H., Anusavice, K.J., 1992. Analysis of tempering stresses in bilayered porcelain discs. J. Dent. Res. 71, 1139–1144. https://doi.org/10.1177/ 00220345920710050201. Dehoff, P.H., Anusavice, K.J., 1989. Effect of visco-elastic behavior on stress development in a metal-ceramic system. J. Dent. Res. 68, 1223–1230. https://doi. org/10.1177/00220345890680080201. Dehoff, P.H., Anusavice, K.J., 1989. Tempering stresses in feldspathic porcelain. J. Dent. Res. 68, 134–138. https://doi.org/10.1177/00220345890680020701. DeHoff, P.H., Anusavice, K.J., 2004. Creep functions of dental ceramics measured in a beam-bending viscometer. Dent. Mater. 20, 297–304. https://doi.org/10.1016/ S0109-5641(03)00107-6. DeHoff, P.H., Anusavice, K.J., 2004. Shear stress relaxation of dental ceramics determined from creep behavior. Dent. Mater. 20, 717–725. https://doi.org/ 10.1016/j.dental.2003.10.005. DeHoff, P.H., Anusavice, K.J., G€ otzen, N., 2006. Viscoelastic finite element analysis of an all-ceramic fixed partial denture. J. Biomech. 39, 40–48. https://doi.org/10.1016/j. jbiomech.2004.11.007. DeHoff, P.H., Barrett, A.A., Lee, R.B., Anusavice, K.J., 2008. Thermal compatibility of dental ceramic systems using cylindrical and spherical geometries. Dent. Mater. 24, 744–752. https://doi.org/10.1016/j.dental.2007.08.008. Denry, I., Kelly, J.R., 2008. State of the art of zirconia for dental applications. Dent. Mater. 24, 299–307. https://doi.org/10.1016/j.dental.2007.05.007. Fragassa, C., 2016. Modelling the viscoelastic response of ceramic materials by commercial finite elements codes. FME Trans. 44, 58–64. https://doi.org/10.5937/ fmet1601058F. Gherlone, E., Mandelli, F., Cappar� e, P., Pantaleo, G., Traini, T., Ferrini, F., 2014. A 3 years retrospective study of survival for zirconia-based single crowns fabricated from intraoral digital impressions. J. Dent. 42, 1151–1155. https://doi.org/10.1016/j. jdent.2014.06.002. Kim, J., Dhital, S., Zhivago, P., Kaizer, M.R., Zhang, Y., 2018. Viscoelastic finite element analysis of residual stresses in porcelain-veneered zirconia dental crowns. J. Mech.
13
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Wendler, M., Belli, R., Petschelt, A., Lohbauer, U., 2016. Spatial distribution of residual stresses in glass-ZrO2 sphero-cylindrical bilayers. J. Mech. Behav. Biomed. Mater. 60, 535–546. https://doi.org/10.1016/j.jmbbm.2016.03.013. Wendler, M., Belli, R., Petschelt, A., Lohbauer, U., 2015. Characterization of residual stresses in zirconia veneered bilayers assessed via sharp and blunt indentation. Dent. Mater. 31, 948–957.
Zhang, Y., Chai, H., Lee, J.J.W., Lawn, B.R., 2012. Chipping resistance of graded zirconia ceramics for dental crowns. J. Dent. Res. 91, 311–315. https://doi.org/10.1177/ 0022034511434356. Zhang, Y., Sailer, I., Lawn, B.R., 2013. Fatigue of dental ceramics. J. Dent. 41, 1135–1147. https://doi.org/10.1016/j.jdent.2013.10.007.
14
Contents lists available at ScienceDirect
Journal of the Mechanical Behavior of Biomedical Materials journal homepage: http://www.elsevier.com/locate/jmbbm
Viscoelastic finite element evaluation of transient and residual stresses in dental crowns: Design parametric study Sukirti Dhital a, Camila Rodrigues b, c, Yu Zhang b, Jeongho Kim a, * a
Department of Civil and Environmental Engineering, University of Connecticut, 261 Glenbrook Rd, U-3037, Storrs, CT, 06269, USA Department of Biomaterials and Biomimetics, New York University College of Dentistry, 433 First Avenue, New York, NY, 10010, USA c Post-Graduation Program in Dental Sciences, Federal University of Santa Maria, 1184 Floriano Peixoto St, Santa Maria, RS, 97015-372, Brazil b
A R T I C L E I N F O
A B S T R A C T
Keywords: Viscoelastic finite elements Dental ceramics Porcelain-veneered zirconia Residual stresses Transient stresses Coefficient of thermal contraction
Porcelain-veneered zirconia (PVZ) are one of the popular choice for crown restorations. Veneer layer of these dental restorations, however, is susceptible to chipping and delamination due to the development of transient and residual stresses during the cooling phase of veneer firing. The aim of this study is to elucidate the effect of material property mismatch, veneer to core thickness ratio, and cooling rate on these transient and residual stresses of PVZ restorations. Three-dimensional viscoelastic finite element modelling (VFEM) was performed. The VFEM model was developed using the UEXPAN subroutine in ABAQUS software and was validated for transient and residual stresses in a sandwich seal problem with experimental data available. A good agreement between the simulated VFEM results and experimental data was obtained. Using validated VFEM, two PVZ systems (PM9/ zirconia and ZirPress/zirconia), three veneer to core thickness ratios (2:1, 1:1 and 1:2), and two cooling rates controlled slow cooling at 1.74E-5 W/mm2� C (i.e. ~30 � C/min) and fast bench cooling at 1.74E-4 W/mm2� C (i.e. ~300 � C/min) were used. The results showed that PM9/zirconia has smaller thermal contraction mismatch, resulting in lesser residual stress (33.36 MPa) as compared to ZirPress/zirconia (37.94 MPa) for controlled cooling and 2:1 veneer to core ratio. In addition, in both systems with the decrease in veneer thickness, we observed a decrease in residual stresses developed. We also observed some effect of cooling rate on residual stresses. The controlled cooling resulted in lower residual stress (24.35 MPa) for PM9/zirconia with a 1:1 veneer to core thickness ratio as compared to bench cooling (28.04 MPa). The effect of cooling rate was more evident on transient stresses. For instance, in the PM9/zirconia with 1:1 thickness ratio model, the difference in transient stresses was 9.93 MPa between controlled and bench cooling. Therefore, properties such as elastic modulus and coefficient of thermal contraction (CTC), as well as the thickness ratio and cooling rate all play an important role in transient and residual stresses developed in the studied ceramic systems.
1. Introduction The use of porcelain-veneered Yttria Tetragonal Zirconia Polycrystal (Y-TZP) crowns and fixed dental prostheses has increased significantly in past years. In order to mimic the enamel characteristics, a porcelain over layer is usually applied to the zirconia core. Though they exhibit great aesthetics as well as chemical durability and biocompatibility (DeHoff et al., 2006), they are quite susceptible to chipping of the por celain veneer layer (Gherlone et al., 2014; Nicolaisen et al., 2016). These dental structures are fabricated by fusing two materials together at high temperatures. This manufacturing process results in the development of transient and residual stresses during the cooling phase. The stored
tensile residual stresses results in chipping or delamination of the veneer (Pang et al., 2015; Christensen, 2009; Denry and Kelly, 2008; Sailer et al., 2007; Zhang et al., 2012; Zhang et al., 2013. Once the veneered porcelain chips, the crown often needs to be repaired or replaced and this process can be cumbersome. Hence, researchers are driven to find an optimal combination for core and veneer layer materials to minimize the deleterious tensile residual stresses, thereby ensuring longevity. There are many factors that give rise to residual thermal stresses. Some of which are temperature-dependent material properties, geometry, layer thickness, thermal contraction mismatch and cooling rate. The residual stresses in porcelain-veneered zirconia systems have been widely studied using flat samples (Swain, 2009; Taskonak et al.,
* Corresponding author. E-mail address: [email protected] (J. Kim). https://doi.org/10.1016/j.jmbbm.2019.103545 Received 25 September 2019; Received in revised form 14 November 2019; Accepted 18 November 2019 Available online 19 November 2019 1751-6161/© 2019 Elsevier Ltd. All rights reserved.
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
analyzing strips or simpler models (Asaoka and Tesk, 1990; Dehoff and Anusavice, 1989a; Dehoff and Anusavice, 1989b; Dehoff and Anusavice, 1992; Asaoka et al., 1992; Nakatsuka and Anusavice, 1997). Several studies were carried out to understand this behavior of dental porcelains (DeHoff and Anusavice, 2004a; DeHoff and Anusavice, 2004b; Fragassa, 2016). Fragassa, 2016 investigated the viscoelastic behavior of ceramic materials using few of the available commercial finite element codes such as ANSYS, LS-DYNA, ABAQUS and MARC. Commercial software, i. e. ANSYS, was used to perform viscoelastic analysis (DeHoff and Anu savice, 2004a; DeHoff and Anusavice, 2004b; DeHoff et al., 2006, 2008). Unlike these, we developed user subroutines such as UEXPAN and UTRS for ABAQUS to better capture the viscoelastic behavior of porcelain throughout the entire firing cycle by defining shear relaxation depen dent coefficient of thermal contraction, which is not available in com mercial software like Abaqus. Knowledge of the residual stress states in porcelain-veneered struc tures are important. However, understanding the entire stress history is also very critical. Since once the transient stress exceeds the stress threshold of fracture, micro-cracks and flaws could form. Such defects can diminish the fracture resistance of ceramic restorations. Unfortu nately, in the body of dental literature, we could find only one study that analyzed transient stresses (Benetti et al., 2014). Benetti et al. have proposed that transient stresses results in the nucleation of cracks inside the porcelain material during the cooling phase, which might result in fracture under different loadings. However, Benetti et al. fell short of predicting the exact transient stresses because in their study the high temperature properties of the materials were not experimentally deter mined, rather extrapolated from the manufacturers’ data sheet or literature values of similar materials. . As alluded to in the preceding paragraph, the major advantage of our newly developed subroutines is its ability to capture the entire stress evolution during the firing cycle, especially in the cooling phase. This way, any significant transient stress can be captured and analyzed. Knowledge gained can be used as a guideline to refine the processing parameters and optimize the design of restorations. In this paper, we aim to further understand the cause of transient and residual stresses in PVZ systems. The two subroutines were first vali dated by a well-studied sandwich seal system. Such a sandwich seal system allows for real-time experimental transient stress analysis (Scherer and Rekhson, 1982). Then, the validated model was used to determine both transient and residual stresses in two PVZ systems, namely (PM9/zirconia and ZirPress/zirconia), in a stylized 3D crown. We also altered the veneer thickness with three veneer to core ratios: 2:1, 1:1 and 1:2 while keeping the total thickness constant. Further we subjected the model to two cooling rates that is controlled slow cooling (1.74E-5 W/mm2 � C or on a more familiar term 30 � C/min) and fast bench cooling (1.74E-4 W/mm2 � C or 300 � C/min). The effect of all these parameters on transient and residual stresses obtained on the veneer layer is studied.
Table 1 Materials combination and veneer to core ratio. Model Number
Veneer to Core Ratio
Core
Veneer
Model 1 Model 2 Model 3 Model 4 Model 5 Model 6
2:1 1:1 1:2 2:1 1:1 1:2
zirconia zirconia zirconia zirconia zirconia zirconia
PM9 PM9 PM9 ZirPress ZirPress ZirPress
Fig. 1. Finite element model of (a) full dental crown and b) 1/4th simulated model. Light blue represents the zirconia core and grey represents the porce lain veneer. Table 2 Material properties at the room temperature (25 � C). Ceramics PM9 ZirPress Zirconia
Young’s Modulus, E (MPa) 60,000 70,000 206,000
Poisson’s Ratio, ν 0.21 0.21 0.3
Density, ρ (kg/ mm3) 2.45E-06 2.54E-06 6.10E-06
Thermal Conductivity (W/mm ⁰C) 1.11E-03 1.15E-03 2.92E-03
Specific Heat, (J/ kg⁰C) 821 886 466
2005, 2008; Wendler et al., 2015; Choi et al., 2011, Tanaka et al., 2016). However, thermal stresses observed in flat samples do not represent the stresses in complex geometries, such as dental crowns. In addition, many experimental techniques have been reported for measuring these re sidual stresses, e.g. Vickers indentations (Baldassarri et al., 2012; Tanaka et al., 2016; Choi et al., 2011; Wendler et al., 2015), birefrin gence (Wendler et al., 2016; Belli et al., 2012; Tholey et al., 2011), hole-drilling (Mainjot et al., 2011), and fracture mechanics approaches (Taskonak et al., 2005, 2008). Such methods require highly polished surfaces, thin slices of the restoration with sufficient translucency, or complex samples preparation, which alter the residual stress states and do not reflect the actual stress profiles in an anatomically-correct restoration. To accurately determine residual stresses in anatomically-correct restorations, the numerical simulation technique finite element method has been adopted. For instance, Meira et al., 2013; Tanaka et al., 2016; and Zhang et al. determined the residual stresses in a PVZ crown by the linear elastic finite element method (LFEM), using the arbitrary values of liquid thermal contraction co efficients of porcelain. Oh et al., 2002 also used this analytical process to calculate the strength of FDPs, which they then compared with the experimental values. They observed that the calculated residual stresses were much more than those obtained experimentally and this might be due to the negligence of the viscoelastic behavior of porcelain. Hence, we can conclude that these models fail to accurately capture the stress profile in porcelain as the viscoelastic behavior of these materials above the softening temperatures has been neglected. The importance of the viscoelastic behavior has been elucidated by
2. Materials and methods 2.1. Methods Three-dimensional finite element analysis is carried out in ABAQUS. Two PVZ systems [PM9 (Vita Zahnfabrik, Bad S€ ackingen, Germany) or ZirPress (Ivoclar Vivadent, Schaan, Liechtenstein) porcelains and a zir conia core material (Zpex, Tosoh, Tokio, Japan)], three veneer to core thickness ratios (2:1, 1:1, or 1:2), and two cooling rates [controlled slow (1.74E-5 W/mm2 � C) or fast bench (1.74E-4 W/mm2 � C], were used (Table 1). The maximum thickness of the specimen was 2.1 mm near central fossa, and the veneer and core thicknesses depend upon the ratio opted. A three-dimensional (3-D) FEA model was analyzed. Owing to the symmetric geometry and boundary conditions, only 1/4th of the model (see Fig. 1) was analyzed. The simulation was carried out in two steps: 2
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 4. Specific Heat at different temperatures.
Fig. 2. Density at different temperatures.
Fig. 5. Young’s Modulus at different temperatures. Fig. 3. Conductivity at different temperatures.
heat transfer analysis to capture temperature history, and stress analysis with a thermal load. 2.2. Material properties Material properties, such as Young’s modulus, Poisson’s ratio, den sity, thermal conductivity and specific heat, at the room temperature are given in Table 2, whereas the corresponding properties at high tem peratures are plotted in Figs. 2–6, respectively. The coefficient of thermal contraction, and subsequently, the density were measured by dilatometer (L75 Platinum Series, Linseis USA, Princeton, NJ). The thermal conductivity and specific heat were deter mined using the laser flash technique (LFA 457 MicroFlash, Netzsch, Selb, Germany). The room temperature Young’s modulus and Poisson’s ratio were measured at the Edward Orton Jr Ceramic Foundation through a paid service. However, we had some difficulties to obtain the high temperature Young’s Modulus values of porcelains. This kind of tests require large specimens (typically 100 mm � 20 mm X 5 mm), which was not possible to produce with the materials we used in this study. Therefore, we assumed that all the porcelains follow similar
Fig. 6. Thermal contraction coefficient at different temperature. 3
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Table 3 Constants for thermal contraction of dental ceramics (αg ¼ a þ bT þ cT2, αl ¼ d þ eT). Materials
αe (mm/mm⁰C) a (1/⁰C) 10
PM9 ZirPress zirconia
6
8.8 8 10.066
αL (mm/mm⁰C)
b (1/⁰C)2 10
9
2.7186 5.6 1.2856
c (1/⁰C)3 10
R2
12
2 4.2 0.6697
d (1/⁰C) 10
0.997 0.978 1
5
8.716 7.184 NA
e (1/⁰C)2 10
7
1.528 1.52 NA
Fig. 8. Schematic of sandwich seal.
Fig. 7. Normalized shear relaxation function at various high temperatures (DeHoff et al., 2006).
trends at high temperatures and adopted the trend of VM9 (Kim et al., 2018) as the experimental data for this material was available. We also utilized the high temperature modulus data of zirconia (3Y-TZP) from the same reference (Kim et al., 2018). For coefficient of thermal contraction (CTCs), as shown in Fig. 6, the CTCs for veneers are divided into solid and liquid CTCs, αg(T) and αl(T) given by: (1)
αg ¼ aþbT þ cT2 and αl ¼ d þ eT
where T is in Celsius, and the five coefficients are given in Table 3. Then, for CTCs of the core, thermal contraction strain versus tem perature data were fit to a third-degree polynomial relation, from regression analysis, from which αg(T) was determined, using Eq. (1). For the veneer, αl(T) was determined by fitting the contraction data with a quadratic equation from the temperature at which creep first occurs on the heating curve to the glass transition temperature, Tg.
Fig. 9. Longitudinal stress (σyy) for a G-11/Alumina sandwich seal cooled at 3 � K/min ¼ 0.05 C/s from 618 to 20 � C.
Z
τðtÞ ¼
t
GR ðt
(2)
sÞ_γ ðsÞds
0
2.3. Viscoelastic finite element analysis
where GR(t) is the shear relaxation modulus, t is the present time, and s is the past time. The normalized shear relaxation modulus can be intro duced by (Kim et al., 2018):
The time dependent viscoelastic properties of dental ceramics, sub jected to cooling from sintering temperature to room temperature, can be introduced into ABAQUS using the viscoelastic option in the material properties which then can be used to determine residual stresses. The shear stress τ (t) in a viscoelastic material model, when subjected to small shear strain γ (t), is given by (Kim et al., 2018):
(3)
gR ðtÞ ¼ GR ðtÞ = G0
where G0 is the instantaneous shear modulus. The volumetric behavior can be written as:
Table 4 Viscoelastic coefficients for Prony series used for veneer. Temp (oC)
g1
g2
g3
g4
τ1
τ2
τ3
τ4
700
0.9960
0.0030
0.0006
0.0004
0.01316
0.10000
0.0050
0.0030
4
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
N X
gR ðtÞ ¼ 1
� gPi 1
G e t=τi
�
(5)
i¼1
where N is the number of terms (N ¼ 4 used in this study), gPi are co efficients, and τGi is shear relaxation time (Table 4) as given at the reference temperature of 700 � C (DeHoff et al., 2006). Since we have taken the reference temperature as the initial temperature and we are cooling the two different material combinations of PM9/zirconia and ZirPress/zirconia from 754 � C and 657 � C respectively, we assumed that these materials have same relaxation coefficients and relaxation time at the respective cooling temperature as those given by DeHoff at 700 � C. The viscoelastic materials show thermorheologically simple behavior in which the normalized shear relaxation function when plotted against log10t has same shape at each temperature and hence can be shifted along the x-axis, thereby resulting in a master curve at a reference temperature. Fig. 7 shows the thermorheological behavior of PM9 using the coefficients in Table 4. The leftmost solid curve is the master curve at 754 � C. The shift function, which is the shift along x-axis, that represents temperature dependence of shear relaxation functions can be defined using UTRS user subroutine in ABAQUS. This subroutine defines the Tool-Narayanaswamy approximation, which takes the form (Tool, 1946; Narayanaswamy, 1971). � � �� H 1 x 1 x ln A ¼ (6) R Tref TðtÞ Tf ðtÞ
Fig. 10. Stress for the G-11/Alumina sandwich seal when cooled at 3 K/min ¼ � 0.05 C/s, and held 4 h at 465 � C.
Z pðtÞ ¼
K0
t
kR ðt 0
sÞ_εvol ðsÞds
(4)
where p is the hydrostatic pressure, K0 is the instantaneous elastic bulk modulus, kR(t) is the normalized bulk relaxation modulus assumed to be 1.0 due to its minor influence in most stress states, and εvol is the volume strain. The normalized shear relaxation function can be described by the use of a four-term Prony series in ABAQUS using the following equation:
where Tref ¼ 1027.15 in Kelvin (754 in Celsius), H is the activation energy, R is the ideal gas constant (H/R ¼ 46400 is used), T(t) is the temperature at time t and Tf(t) is the fictive temperature at time t, and x is a material constant between 0 and 1 (x ¼ 0.3 used, e.g. (Fragassa,
Fig. 11. Temperature and Stress Profile for 3-D Slow Cooling in a) Model 1 and b) Model 4. 5
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
3. Results This section presents a validation example which shows validity of the present viscoelastic finite element procedure in predicting strains and stresses. The validated model is then used to analyze 3-D dental models in order to provide a map of stresses for different material combinations, cooling rates and veneer to core ratios. 3.1. Sandwich seal under cooling Viscoelastic finite element analyses was performed to simulate the contraction behavior of a G-11/Alumina sandwich seal cooled at 3 K/ min (Scherer and Rekhson, 1982). The authors extended the Nar ayanswamy’s model (Narayanaswamy, 1971) to accurately capture the viscoelastic behavior of the glass and then validated their formula with sandwich seal experiment. In addition, to demonstrate the applicability of our codes (UEXPAN and UTRS) in ABAQUS, we modelled the same G-11/Alumina sandwich seal (Fig. 8) subjected to cooling at 3 K/min. The following properties were used for Alumina: EAl ¼ 3.73 � 105 MPa, ν ¼ 0.3 and for G-11 glass: H/R ¼ 64500 (K) where, H/R ¼ ratio of activation enthalpy for relaxation process (H) and ideal gas constant (R) x ¼ 0.53 T ¼ 618� C ¼ Initial Temperature �
αg ðTÞ ¼ 64:7 þ 0:02T � 10 αl ðTÞ ¼ 3:43 � 10
5
7
with T in Celsius
with T in Celsius
EG11 ¼ 7.26 � 104 MPa, ν ¼ 0.3 The ratio of thickness of glass to thickness of seal was 13. The obtained FEA results, for the longitudinal stress (σyy) in this G11/Alumina model were compared with the experimental results from Scherer and Rekhson, 1982) as shown in Fig. 9. The results showed an excellent agreement with difference less than 0.7 MPa over the entire temperature history. As shown in Fig. 9, the stress increased to a maximum value upon initial cooling, but then decreased owing to the viscoelastic properties. This kind of highly nonlinear stress response cannot be captured by elastic analysis. Further, to show the effectiveness of our code in order to capture any arbitrary thermal history, including isothermal holds, we subjected the same G11-Alumina sandwich seal to a cooling rate of 3 K/min, held it for 4 h at 465 � C and then cooled it to the room temperature at the same rate. Fig. 10 shows the stress for this condition, compared with the experimental results present in Scherer and Rekhson, 1982). The agreement between the FEA predictions and experimental results was quite good with difference less than 0.6 MPa over the entire temperature history. During isothermal hold, one could observe that the stress kept increasing to a maximum value and then diminished sharply afterwards. This sandwich seal modelling showed the effectiveness of our codes in accurately capturing the viscoelastic behavior, and hence we used them to analyze ou0pr veneering ceramics.
Fig. 12. Model 1 (PM9/zirconia): Temperature; a) Temperature gradient in slow cooling b) Temperature history in Slow and Fast Cooling; Nodes repre sented are PT1-Porcelain top node, ZB1- zirconia top node and ZB2-zirconia bottom node.
2016)). The logarithmic function lnA (T, Tref, t) shifts the shear relaxa tion function along the horizontal axis for various temperatures. As the data for volume relaxation are not available and in case of glasses it has been found that the volume relaxation time is 1/4 to 1/20 of shear relaxation time, we assumed the volume relaxation time in the ceramics to be 1/10 of the shear relaxation time. Also studies have shown that volume relaxation time has almost no effect on residual stresses (DeHoff et al., 2006). As we know that fictive temperature needs to be considered during a viscoelastic analysis specially in the transition region, and Markovsky, 1984 proposed the algorithm for this fictive temperature Tf(t) which has been coded in UEXPAN. The thermal strain is given for the core and the veneer, respectively, by Z T Z T Z Tf � (7) εc ¼ αg ðTÞdT and εv ¼ αg ðTÞdT þ αl Tf dT To
Tf
3.2. Veneering ceramics
To
Viscoelastic finite element analyses were run to simulate the contraction behavior of a dilatometer specimen subjected to convective cooling (DeHoff et al., 2006). We analyzed two different material combinations, each of which were subjected to controlled cooling
where αg(T) and αl (T) are plotted in Fig. 5.
6
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 13. Residual stresses in slow cooling: a) Model 1: PM9/zirconia; b) Model 4: ZirPress/zirconia; c) Model 2: PM9/zirconia; d) Model 5: ZirPress/zirconia; e) Model 3:PM9/zirconia; f) Model 6: ZirPress/zirconia. Nodes represented are PT1-Porcelain top node, PZI- Porcelain zirconia Interface node and PB1- Porcelain Bottom Node.
(1.74E-5 W/mm2 oC) and bench cooling (1.74E-4 W/mm2 oC) on the crown surfaces exposed to air. The material combinations are as mentioned in Table 1. Fig. 11 shows the 3D temperature and stress profile for 2:1 veneer to core PM9/zirconia and ZirPress/zirconia models when subjected to control cooling. Since, this is an axisymmetric model, the temperature and stress contours are uniform through the face, as shown in the fig ures. Hence, for ease of study, only one face of the contour images are provided in Figs. 12, 13 and 15. As mentioned before, in the first step, 3-D transient heat analysis was carried out. This gave us the temperature profile of the crown when subjected to both controlled and bench cooling (see Fig. 12). As ex pected, bench cooling causes more rapid cooling than controlled cooling.
All the models follow similar temperature profile as shown in Fig. 12. After inputting these temperature profiles as load, stress analysis was performed. The maximum principal residual and transient (timedependent) stresses are plotted, respectively, in Figs. 13 and 14. Fig. 13 depicts the residual stress contour at the end of analysis, when subjected to slow cooling. It is clear that the maximum principal tensile stress in the porcelain layer, for all of the models, occur at the cusp interface of porcelain and zirconia layers. The location for maximum principal stress in the zirconia layer, however, varies from model to model. As zirconia has very high strength, we focus on tensile residual stresses developed in the porcelain layer. The maximum prin cipal residual stress contours shows that the residual stress developed in ZirPress layer is more as compared to its PM9 counterpart. For instance in PM9/zirconia models that is Model 1, Model 2 and Model 3, tensile 7
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 14. Transient stress in slow cooling: a) Model 1: PM9/zirconia; b) Model 4: ZirPress/zirconia; c) Model 2: PM9/zirconia; d) Model 5: ZirPress/zirconia; e) Model 3:PM9/zirconia; f) Model 6: ZirPress/zirconia.
residual stresses developed in the PM9 layer are 33.36 MPa, 24.35 MPa and 17.40 MPa respectively, whereas for ZirPress/zirconia models, that is Model 4, Model 5 and Model 6, tensile residual stresses developed are 37.94 MPa, 27.90 MPa and 20.02 MPa respectively in the ZirPress layer. From both the PM9 and ZirPress models, we can see that with the decrease in the thickness of veneer layer, the tensile residual stress developed in it decreases. For 2:1 veneer to core PM9/zirconia model, the residual stress developed in PM9 veneer is 33.36 MPa which is almost 10 MPa more than the residual stress in the same veneer material with 1:1 veneer to core ratio. Similar pattern can be observed in the
remaining PM9/zirconia model and all of the ZirPress/zirconia models. Besides the residual stress, the tensile transient stress developed during the cooling phase, is also one of our concern. In order to understand it better, the history of the maximum principal stress in the porcelain layer needs to be reviewed. As the abrupt rise and fall of stresses during the entire simulation period may have effect on the material behavior, one needs to look at the stress profile developed along the entire time period. Few of the key features are mentioned here. From the graphs, we could observe that for PM9/zirconia Model 1, sudden peak occurred at the node PT1 with 8
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 15. Residual stress in fast cooling: a) Model 1: PM9/zirconia; b) Model 4: ZirPress/zirconia; c) Model 2: PM9/zirconia; d) Model 5: ZirPress/zirconia; e) Model 3: PM9/zirconia; f) Model 6: ZirPress/zirconia. Nodes represented are PT1-Porcelain top node, PZI- Porcelain zirconia Interface node and PB1- Porcelain Bottom Node.
value 4.516 MPa at 37.5 s, which then gradually decreased to obtain a stable stress value. Similarly, for Model 2, at the same point and at 35 s, 2.128 MPa maximum principal tensile stress is obtained and for Model 3, at 33.5 s, 2.589 MPa of stress is observed. The maximum principal stress profile for ZirPress/zirconia is similar to that of PM9 models as shown in Fig. 14. For Model 4, at 45 s and at node PT1, the maximum principal stress obtained is 8.012 MPa, at 43 s for Model 5, the stress obtained is 4.481 MPa and 5.673 MPa for Model 6 at 42 s. The models were also subjected to fast cooling. The maximum principal residual and transient (time-dependent) stresses are plotted, respectively, in Figs. 15 and 16. Similar to slow cooling, Fig. 15 shows the residual stress contour at the end of analysis for fast cooling, in all the six models. Here too, in all cases, the maximum value of the principal tensile stress, of the porcelain layer, occurs at the cusp interface of the two layers. The tensile residual stress developed in ZirPress layer, is more when compared to PM9 models. The maximum value of the maximum principal residual stress developed in the PM9 layer are 36.36 MPa, 28.04 MPa and 20.66 MPa for Model 1, Model 2 and Model 3 respectively. Similarly, for ZirPress/zirconia model, the maximum
residual stress obtained are 39.12 MPa, 30.84 MPa and 22.66 MPa for Model 4, Model 5 and Model 6 respectively. Irrespective of the cooling rate, we can see that with decrease in the veneer thickness (as indicated by the model numbers) there is decrease in the tensile stress produced in the veneer layer, in both the materials. However, as one can conclude, tensile residual stresses obtained due to fast cooling are higher than that due to slow cooling in both the materials. For instance, for Model 1 and PM9 veneer layer, the tensile residual stress obtained due to slow cooling is 33.36 MPa and fast cooling is 36.36 MPa, which shows the increase in the residual stress. In this case, we also can observe the abrupt rise and fall in the transient stress profile for models under fast cooling, as shown in Fig. 16. However, the rise in this case is more lucid than that seen in Fig. 14. For PM9/zirconia model, we can observe the sudden peak, at node PT1, for Model 1 is 18.786 MPa at 4.5 s, for Model 2 at 4.5 s the value obtained is 12.057 MPa and for Model 3 at 4.5 s the transient tensile stress obtained is 17.811 MPa. Now, in the case of ZirPress/zirconia Model, for Model 4, at 5.5 s, the tensile transient stress obtained is 23.187 MPa, for Model 5 at 5 s we obtained 21.435 MPa and for Model 6, at 5 s, the transient 9
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Fig. 16. Transient stress in fast cooling: a) Model 1: PM9/zirconia; b) Model 4: ZirPress/zirconia; c) Model 2: PM9/zirconia; d) Model 5: ZirPress/zirconia; e) Model 3:PM9/zirconia; f) Model 6: ZirPress/zirconia.
stress obtained is 25.802 MPa. Analyzing all the 12 model combination, for different materials, thickness ratios and cooling rate, the maximum principal tensile residual stress is summarized in the bar graphs in Fig. 17.
able to accurately capture the viscoelastic behavior of materials, por celain in our case (Kim et al., 2018). This might be due to the fact that for states above glass transition temperature of porcelain, Tg, the Young’s modulus is dependent on both temperature as well as time, which means that even when the temperature is constant the modulus decreases exponentially with time and this phenomena cannot be addressed by LFEM. This results into larger residual stress obtained as compared to that from viscoelastic FEM, even when the coefficient of thermal contraction for porcelain in liquid state, used in VFEM, is much larger
4. Discussion Many researchers have used linear elastic FEM to determine residual stresses in different material systems. However, these LFEMs are not 10
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
for the sandwich seal system is compared against the test data which shows remarkably good match, thereby validating our model. This model is then used to estimate residual stresses and perform parametric study on various dental ceramic combinations with different thickness ratios, material properties and cooling rates. It has been shown through our analysis that the veneer thickness affects residual stresses developed in dental models. We analyzed three different veneer to core ratios, which are 2:1, 1:1 and 1:2. From the results, for both the material combinations, we could clearly observe that with the decrease in the veneer thickness, tensile residual stresses developed in porcelain reduced. For instance, in PM9/zirconia model, the residual stress developed in PM9 was about 33.36 MPa with 2:1 ratio, 24.35 MPa for 1:1 ratio and 17.40 MPa for 1:2 ratio, under controlled cooling. Similar pattern can be observed under bench cool ing. For ZirPress/zirconia, the results obtained were 37.94 MPa, 27.90 MPa and 20.02 MPa for 2:1, 1:1 and 1:2 ratios respectively, under controlled cooling. Thus, we could conclude that, with reduced porce lain layer thickness, the reduction in thermal gradient results in lower stresses. As we have analyzed two different material combinations, we could also see the effect of material properties in residual stresses developed. The maximum tensile stress of 33.36 MPa and 36.36 MPa occurs for controlled and bench cooling respectively, for PM9/zirconia combina tion, and similarly the maximum tensile stress of 37.94 MPa and 39.12 MPa occurs for ZirPress/zirconia model, for the same veneer core thickness ratio. To understand this difference in results, we studied the combined effect of young’s moduli and thermal contraction coefficients on the different material models analyzed. The CTC of PM9 is closer to zirconia, as compared to Zirpress. Also, the young’s modulus of PM9 is lower than that of ZirPress. Hence, residual stresses in ZirPress veneer layer are higher than those in PM9. The transient and residual stresses are a combined effect of time and temperature dependent properties such as modulus, CTC and viscoelasticity. The residual stress in zirconia is not of major concern here, as they have very high strength. However, for porcelain, with a smaller CTC mismatch and larger modulus mismatch, we obtain lower residual stresses, which is in agreement with results showed on previous studies (Kim et al., 2018). Besides the core and veneer thickness and material properties, the cooling rate also has an effect on the amount of residual and transient stresses developed. Controlled cooling produced less tensile stress as compared to bench cooling, in all combinations of materials and thick ness ratios. The fast change in temperature, in bench cooling, brings about larger temperature difference, which in turn results in higher transient and residual stresses, as verified from the models analyzed. The difference between tensile residual stresses obtained for controlled and bench cooling for Model 1, Model 2, Model 3, Model 4, Model 5 and Model 6 is 3.10 MPa, 3.69 MPa, 3.26 MPa, 1.18 MPa, 2.94 MPa and 2.64 MPa respectively. The difference in residual stresses for controlled and bench cooling might be less, which is similar to the findings of Tholey et al. (2011). However, this gap is amplified when one studies transient stresses. For the same order of six models, the difference in tensile transient stress obtained are 14.27 MPa, 9.93 MPa, 15.22 MPa, 15.18 MPa, 16.96 MPa and 20.13 MPa respectively, thereby completely elucidating the effect of cooling rate on the amount of stress produced. However, these high transient stresses occur only for a short period, but they result into high strain rates in our viscoelastic analysis, thereby increasing the susceptibility of the porcelain layer to fractures, cracking and chipping. The models analyzed tried to capture the actual conditions, however, our study had few limitations. As it is difficult to measure the Young’s Moduli of different porcelain at high temperatures, required in our study, we assumed that both PM9 and ZirPress follow exactly the same trend as that followed by VM9. However, the room temperature modulus were obtained from laboratory measurements. Besides this we also have assumed the shear relaxation coefficients of these porcelains at the required temperatures to be same as the one in DeHoff et al. (2006)
Fig. 17. Residual Stress Summary in (a) porcelain; (b) zirconia.
than that in a solid state owing to the large thermal strain. However, the thermal contraction coefficients of liquid porcelain used in LFEM are quite small, normally 3–4 times smaller than actual. VFEM in the other hand uses the temperature and time-dependent phenomenon of stress relaxation to get small stress value as compared to LFEM. The depen dence of shear relaxation modulus on the time-dependent temperature in a specimen is described by a simple shift function (Eq. (5)), supported by experimentally determined creep data (DeHoff and Anusavice, 2004a). The shear relaxation modulus data presented in Fig. 6 is uni versal for porcelains at any given temperature under cooling rates examined here. The viscoelastic finite element model that has been developed, considering the above mentioned points, is used to capture the stress profiles in two core-veneer systems, that is, PM9/zirconia and ZirPress/zirconia. To ensure that the developed VFEM models accu rately captures the obtained residual stresses, they have been validated in a sandwich seal system. The obtained stress results from VFEM model 11
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
as the data for this was also not available. Addition to shear relaxation coefficients, the volume relaxation coefficients too were assumed to be 1/10 of shear relaxation coefficients (DeHoff et al., 2006), but their influence on residual stresses is quite less. Finally, damage or crack formation is not considered during this viscoelastic analysis. Further study on the effect of crack initiation and growth on deformations and residual stresses will require damage analysis in conjunction with the present analysis.
for both the materials, residual and transient stresses developed in it also decreases. � The maximum principal stress developed in the porcelain veneer and zirconia core is less for controlled slow cooling as compared to rapid bench cooling. However, the cooling effect is more evident in the transient than residual stresses and these transient stresses contribute in producing micro-cracks, which ultimately fractures. Declaration of competing interest
5. Conclusions
All authors declare no conflict of interest.
From the results obtained, we could draw the following conclusions: � The material properties such as Young’s Modulus and Coefficient of Thermal Contraction (CTC) play a role in the amount of stress developed. If the CTC mismatch between the core and veneer is less and the porcelain modulus is lower, then the residual stress devel oped in both the layers is reduced, irrespective of the cooling rate. � The thickness of core and veneer layer play a role in the amount of stress developed. With the decrease in the thickness of veneer layer,
Acknowledgements This work was sponsored by funding from the United States National Institute of Dental & Craniofacial Research, National Institutes of Health (Grants Nos. 1R01 DE026279 and R01DE026772). The authors are also thankful to the Brazilian Federal Agency for Support and Evaluation of Graduate Education (CAPES) (finance code 001).
Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi.org/10.1016/j.jmbbm.2019.103545. Appendix B 1. Subroutine UEXPAN is called a. Density and temperature dependent elastic properties are read from the.inp file, for porcelain and coping, which are defined by the user. Specific heat for zirconia are also user defined. b. If the material is “PM9 or ZirPress” then: i. Using Narayanswami form, the ratio between relaxation time and reference (volume) relaxation time is determined as: � �� � τi ΔH 1 x 1 x ¼ exp R Tref T Tf τis where, τi ¼ structural relaxation time
τis ¼ volume relaxation time, which is taken as 1/10th of shear relaxation time H ¼ activation energy x ¼ experimentally defined constant⊲ Tref ¼ reference temperature T ¼ current temperature at end of increment Tf ¼ fictive temperature
ii. The partial fictive temperatures Tfi ðtÞ, the fictive temperature Tf ðt) and the difference between fictive temperatures (dTf) of two sequential steps are stored as state variables and calculated as: � � Tfi ðtÞ þ T � Δt τi � � Tfi ðtÞ ¼ where; Δt ¼ time increment 1 þ Δt τi �� X Tf t ¼ ci Tfi ðtÞ where; ci are the constants defined at reference temperature dTf ¼ Tf
Tf ðprevious incrementÞ
iii. The thermal contraction coefficients for glass and liquid states respectively, are obtained as: ag ¼ ag1 þ ag2 �T þ ag3 �T2 agf ¼ ag1 þ ag2 �Tf þ ag3 �T2f 12
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
alf ¼ al1 þ al2 �Tf þ al3 �T2f where, ag1 ; ag2 ; ag3 ; al1 ; al2 are the constants for thermal contraction coefficients in which subscript g is for glass and l for liquid state. iv. Finally coefficient of expansion for porcelain is calculated as: � expansion ¼ ag �temp: increment þ alf agf �fictive temp increment
2. After executing UEXPAN subroutine, UTRS is called to account for the temperature-time shift for a time domain viscoelastic analysis a. As porcelain is the only viscoelastic material, Prony series is defined by the user in the.inp file and read when required. b. The experimentally defined constant, x, is taken and activation energy Δ H, reference temperature Tref along with ideal gas constant R are defined c. Fictive temperature (Tf ) and the difference between current and previous increment fictive temperatures (dTf Þ are called from UEXPAN sub routine as Statev. d. The shift along the horizontal axis, shift factor, are determined using Tool-Narayanswamy form. The shift factor is determined at the beginning and end of the increment using: � �� � ΔH 1 x 1 x shift ¼ exp 0 0 R Tref Tðt Þ Tf ðt Þ where; Tðt Þ is the temperature at any time t and Tf ðt Þis the fictive temperature at t 0
0
0
References
0
Behav. Biomed. Mater. 82, 202–209. https://doi.org/10.1016/j. jmbbm.2018.03.020. Mainjot, A.K., Schajer, G.S., Vanheusden, A.J., Sadoun, M.J., 2011. Residual stress measurement in veneering ceramic by hole-drilling. Dent. Mater. 27, 439–444. https://doi.org/10.1016/j.dental.2010.12.002. Markovsky, A., 1984. An efficient and stable Algorithm for calculating fictive temperatures effect of variation in flame velocity profile on particle deposition rate. Mater. Sci. 56–57. Meira, J.B.C., Reis, B.R., Tanaka, C.B., Ballester, R.Y., Cesar, P.F., Versluis, A., Swain, M. V., 2013. Residual stresses in Y-TZP crowns due to changes in the thermal contraction coefficient of veneers. Dent. Mater. 29, 594–601. https://doi.org/ 10.1016/j.dental.2013.03.012. Nakatsuka, A., Anusavice, K.J., 1997. Finite element analysis of stress distribution in porcelain discs. J. Mater. Sci. 32, 3621–3627. https://doi.org/10.1023/A: 1018678313122. Narayanaswamy, O.S., 1971. A model of structural relaxation in glass. J. Am. Ceram. Soc. 54, 491–498. https://doi.org/10.1111/j.1151-2916.1971.tb12186.x. Nicolaisen, M., Bahrami, G., Schropp, L., Isidor, F., 2016. Comparison of metal-ceramic and all-ceramic three-unit posterior fixed dental prostheses: a 3-year randomized clinical trial. Int. J. Prosthodont. 29, 259–264. https://doi.org/10.11607/ijp.4504. Oh, W., G€ otzen, N., Anusavice, K.J., 2002. Influence of connector design on fracture probability of ceramic fixed-partial dentures. J. Dent. Res. 81, 623–627. https://doi. org/10.1177/154405910208100909. Pang, Z., Chughtai, A., Sailer, I., Zhang, Y., 2015. A fractographic study of clinically retrieved zirconia-ceramic and metal-ceramic fixed dental prostheses. Dent. Mater. 31, 1198–1206. https://doi.org/10.1016/j.dental.2015.07.003. Sailer, I., Feher, A., Filser, F.T., Gauckler, L.J., Lüthy, H., H€ ammerle, C.H.F., 2007. Fiveyear clinical results of zirconia frameworks for posterior fixed partial dentures. Int. J. Prosthodont. 20, 383–388. Scherer, G.W., Rekhson, S.M., 1982. Viscoelastic-elastic composites: II, sandwich seal. J. Am. Ceram. Soc. 65, 399–406. https://doi.org/10.1111/j.1151-2916.1982. tb10492.x. Swain, M.V., 2009. Unstable cracking (chipping) of veneering porcelain on all-ceramic dental crowns and fixed partial dentures. Acta Biomater. 5, 1668–1677. https://doi. org/10.1016/j.actbio.2008.12.016. Tanaka, C.B., Harisha, H., Baldassarri, M., Wolff, M.S., Tong, H., Meira, J.B.C., Zhang, Y., 2016. Experimental and finite element study of residual thermal stresses in veneered Y-TZP structures. Ceram. Int. 42, 9214–9221. https://doi.org/10.1016/j. ceramint.2016.03.018. Taskonak, B., Borges, G.A., Mecholsky, J.J., Anusavice, K.J., Moore, B.K., Yan, J., 2008. The effects of viscoelastic parameters on residual stress development in a zirconia/ glass bilayer dental ceramic. Dent. Mater. 24, 1149–1155. https://doi.org/10.1016/ j.dental.2008.01.004. Taskonak, B., Mecholsky, J.J., Anusavice, K.J., 2005. Residual stresses in bilayer dental ceramics. Biomaterials 26, 3235–3241. https://doi.org/10.1016/j. biomaterials.2004.08.025. Tholey, M.J., Swain, M.V., Thiel, N., 2011. Thermal gradients and residual stresses in veneered Y-TZP frameworks. Dent. Mater. 27, 1102–1110. https://doi.org/10.1016/ j.dental.2011.08.001. Tool, A.Q., 1946. Relation between inelastic deformability and thermal expansion of glass in its annealing range. J. Am. Ceram. Soc. 29, 240–253. https://doi.org/ 10.1111/j.1151-2916.1946.tb11592.x.
Asaoka, K., Kuwayama, N., Tesk, J.A., 1992. Influence of tempering method on residual stress in dental porcelain. J. Dent. Res. 71, 1623–1627. https://doi.org/10.1177/ 00220345920710091501. Asaoka, K., Tesk, J.A., 1990. Transient and residual stress in a porcelain-metal strip. J. Dent. Res. 69, 463–469. https://doi.org/10.1177/00220345900690020901. Baldassarri, M., Stappert, C.F.J., Wolff, M.S., Thompson, V.P., Zhang, Y., 2012. Residual stresses in porcelain-veneered zirconia prostheses. Dent. Mater. 28 (8), 873–879. Belli, R., Monteiro, S., Baratieri, L.N., Katte, H., Petschelt, A., Lohbauer, U., 2012. A photoelastic assessment of residual stresses in zirconia-veneer crowns. J. Dent. Res. 91, 316–320. https://doi.org/10.1177/0022034511435100. Benetti, P., Kelly, J.R., Sanchez, M., Della Bona, A., 2014. Influence of thermal gradients on stress state of veneered restorations. Dent. Mater. 30, 554–563. https://doi.org/ 10.1016/j.dental.2014.02.020. Choi, J.E., Waddell, J.N., Swain, M.V., 2011. Pressed ceramics onto zirconia. Part 2: indentation fracture and influence of cooling rate on residual stresses. Dent. Mater. 27, 1111–1118. https://doi.org/10.1016/j.dental.2011.08.003. Christensen, G.J., 2009. Porcelain-fused-to-metal versus zirconia-based ceramic restorations, 2009. J. Am. Dent. Assoc. 140, 1036–1039. https://doi.org/10.14219/ jada.archive.2009.0316. Dehoff, P.H., Anusavice, K.J., 1992. Analysis of tempering stresses in bilayered porcelain discs. J. Dent. Res. 71, 1139–1144. https://doi.org/10.1177/ 00220345920710050201. Dehoff, P.H., Anusavice, K.J., 1989. Effect of visco-elastic behavior on stress development in a metal-ceramic system. J. Dent. Res. 68, 1223–1230. https://doi. org/10.1177/00220345890680080201. Dehoff, P.H., Anusavice, K.J., 1989. Tempering stresses in feldspathic porcelain. J. Dent. Res. 68, 134–138. https://doi.org/10.1177/00220345890680020701. DeHoff, P.H., Anusavice, K.J., 2004. Creep functions of dental ceramics measured in a beam-bending viscometer. Dent. Mater. 20, 297–304. https://doi.org/10.1016/ S0109-5641(03)00107-6. DeHoff, P.H., Anusavice, K.J., 2004. Shear stress relaxation of dental ceramics determined from creep behavior. Dent. Mater. 20, 717–725. https://doi.org/ 10.1016/j.dental.2003.10.005. DeHoff, P.H., Anusavice, K.J., G€ otzen, N., 2006. Viscoelastic finite element analysis of an all-ceramic fixed partial denture. J. Biomech. 39, 40–48. https://doi.org/10.1016/j. jbiomech.2004.11.007. DeHoff, P.H., Barrett, A.A., Lee, R.B., Anusavice, K.J., 2008. Thermal compatibility of dental ceramic systems using cylindrical and spherical geometries. Dent. Mater. 24, 744–752. https://doi.org/10.1016/j.dental.2007.08.008. Denry, I., Kelly, J.R., 2008. State of the art of zirconia for dental applications. Dent. Mater. 24, 299–307. https://doi.org/10.1016/j.dental.2007.05.007. Fragassa, C., 2016. Modelling the viscoelastic response of ceramic materials by commercial finite elements codes. FME Trans. 44, 58–64. https://doi.org/10.5937/ fmet1601058F. Gherlone, E., Mandelli, F., Cappar� e, P., Pantaleo, G., Traini, T., Ferrini, F., 2014. A 3 years retrospective study of survival for zirconia-based single crowns fabricated from intraoral digital impressions. J. Dent. 42, 1151–1155. https://doi.org/10.1016/j. jdent.2014.06.002. Kim, J., Dhital, S., Zhivago, P., Kaizer, M.R., Zhang, Y., 2018. Viscoelastic finite element analysis of residual stresses in porcelain-veneered zirconia dental crowns. J. Mech.
13
S. Dhital et al.
Journal of the Mechanical Behavior of Biomedical Materials 103 (2020) 103545
Wendler, M., Belli, R., Petschelt, A., Lohbauer, U., 2016. Spatial distribution of residual stresses in glass-ZrO2 sphero-cylindrical bilayers. J. Mech. Behav. Biomed. Mater. 60, 535–546. https://doi.org/10.1016/j.jmbbm.2016.03.013. Wendler, M., Belli, R., Petschelt, A., Lohbauer, U., 2015. Characterization of residual stresses in zirconia veneered bilayers assessed via sharp and blunt indentation. Dent. Mater. 31, 948–957.
Zhang, Y., Chai, H., Lee, J.J.W., Lawn, B.R., 2012. Chipping resistance of graded zirconia ceramics for dental crowns. J. Dent. Res. 91, 311–315. https://doi.org/10.1177/ 0022034511434356. Zhang, Y., Sailer, I., Lawn, B.R., 2013. Fatigue of dental ceramics. J. Dent. 41, 1135–1147. https://doi.org/10.1016/j.jdent.2013.10.007.
14