- Email: [email protected]
An augmented hybrid constitutive model for simulation of unloading and cyclic loading behavior of conventional and highly crosslinked UHMWPE
ARTICLE IN PRESS
Biomaterials 25 (2004) 2171–2178
An augmented hybrid constitutive model for simulation of unloading and cyclic loading behavior of conventional and highly crosslinked UHMWPE . a,*, C.M. Rimnacb, S.M. Kurtza,c J.S. Bergstrom a Exponent, Inc., 21 Strathmore Road, Natick, MA 01760, USA Departments of Mechanical and Aerospace Engineering and Orthopaedics, Musculoskeletal Mechanics and Materials Laboratories, Case Western Reserve University, Cleveland, OH, USA c Implant Research Center, School of Biomedical Engineering, Science and Health Systems, Drexel University, 3141 Chestnut Street, Philadelphia, PA, USA b
Received 3 June 2003; accepted 15 August 2003
Abstract Ultra-high molecular weight polyethylene (UHMWPE) is extensively used in total joint replacements. Wear, fatigue, and fracture have limited the longevity of UHMWPE components. For this reason, significant effort has been directed towards understanding the failure and wear mechanisms of UHMWPE, both at a micro-scale and a macro-scale, within the context of joint replacements. We have previously developed, calibrated, and validated a constitutive model for predicting the loading response of conventional and highly crosslinked UHMWPE under multiaxial loading conditions (Biomaterials 24 (2003) 1365). However, to simulate in vivo changes to orthopedic components, accurate simulation of unloading behavior is of equal importance to the loading phase of the duty cycle. Consequently, in this study we have focused on understanding and predicting the mechanical response of UHMWPE during uniaxial unloading. Specifically, we have augmented our previously developed constitutive model to also allow for accurate predictions of the unloading behavior of conventional and highly crosslinked UHMWPE during cyclic loading. It is shown that our augmented hybrid model accurately captures the experimentally observed characteristics, including uniaxial cyclic loading, large strain tension, rate-effects, and multiaxial deformation histories. The augmented hybrid constitutive model will be used as a critical building block in future studies of fatigue, failure, and wear of UHMWPE. r 2003 Elsevier Ltd. All rights reserved. Keywords: Constitutive modeling; Ultra-high molecular weight polyethylene; UHMWPE; Hybrid model; FEM; Radiation crosslinking; Multiaxial mechanical behavior; Small punch test
1. Introduction Wear of the articulating surface of ultra-high molecular weight polyethylene (UHMWPE) implant components is an important problem that can significantly limit the life expectancy of total joint replacements. Wear of UHMWPE components is multifactorial and is influenced by the functional loading environment, joint kinematics, component geometry, and material properties. Recently, efforts to reduce UHMWPE wear have involved changes to the resin type, sterilization method, *Corresponding author. Tel.: +1-508-652-8514; fax: +1-508-6471899. . E-mail address: [email protected] (J.S. Bergstrom). 0142-9612/$ - see front matter r 2003 Elsevier Ltd. All rights reserved. doi:10.1016/j.biomaterials.2003.08.065
radiation crosslinking, and thermal treatments [1,2]. However, there is, at present, an incomplete understanding of the wear characteristics and mechanisms of damage evolution, complicating and impeding rapid progression and improvements in performance of UHMWPE joint replacement components. There are two complementary approaches for improving the general understanding of the wear behavior of UHMWPE components used in total joint replacements: macroscopic experimental testing and microstructural material characterization. Wear simulators and other mechanical testing techniques can provide information related to wear rates of different tribological systems and can rank the performance of materials subjected to different environments and
ARTICLE IN PRESS 2172
J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
thermomechanical histories [3]. Although useful, systematic empirical testing has thus far not enabled a priori predictions of the mechanisms causing the actual wear of UHMWPE. To predict the evolution in microscopic and macroscopic damage for new UHMWPE materials from fundamental polymer physics principles requires an understanding of the performance and response of the material on the microstructural level. Theoretical and experimental research supports the notion that the wear observed in vivo and in vitro is the result of localized high stresses and strains in the surface region of the UHMWPE component [4,5]. To better understand and predict these stresses, it is necessary to have a well-calibrated and accurate constitutive model of UHMWPE. In the orthopedic research community, the J2-plasticity model has been the most widely used approach for simulating the behavior of UHMWPE. It has been shown [6], however, that the J2-plasticity model is not an accurate general tool for predicting the largedeformation-to-failure behavior of UHMWPE. In addition, the J2-plasticity model does not accurately predict cyclic loading of UHMWPE. These are serious limitations since UHMWPE joint components undergo large deformations locally at the articulating surface and are also subject to cyclically applied loads. To address these limitations, a new constitutive model was recently developed for conventional and highly crosslinked UHMWPEs [6]. This new model, which is inspired by the physical micro-mechanisms governing the deformation resistance of polymeric materials, is an extension of specialized constitutive theories for glassy polymers that have been developed during the last 10 years. The new model, named the hybrid model (HM), has been shown to accurately predict the mechanical response of both conventional and highly crosslinked UHMWPE materials in uniaxial tension, compression, and multiaxial loading. However, to simulate in vivo changes to orthopedic components, accurate simulation of unloading behavior is also of importance. Consequently, the objective of this study was to better understand and predict the mechanical response of UHMWPE during uniaxial unloading. In this regard, we have developed an augmented hybrid constitutive model that is capable of accurately predicting the experimentally observed stress–strain response in cyclic loading for conventional and well as highly crosslinked UHMWPEs.
and uncrosslinked UHMWPE. The modification in the augmented HM specifically addresses the unloading behavior during cyclic loading. The kinematic framework used in the augmented HM is based on a decomposition of the applied deformation gradient into elastic and viscoplastic components: F=FeFp (Fig. 1). The spring and dashpot representation shown in Fig. 1a is a one-dimensional embodiment of the model framework used to capture the viscoplastic flow characteristics. With the exception of the top spring (E), all spring and dashpot elements are highly nonlinear (described in detail below). Fig. 1b depicts a map of the decomposition of a given material deformation state. This decomposition specifies how the three-dimensionality of the deformation gradients and stress tensors are connected and evolve during an applied deformation history. As in our previous modeling approaches [6,7], the deformation state is decomposed into elastic, backstress, and viscoplastic components. Compared with our previous models, we have now incorporated timedependent viscoplasticity to the backstress network to improve the predictive capabilities of the model with respect to unloading. The rationale for this approach is
2. Augmented hybrid constitutive model for predictions of UHMWPE The new augmented HM is a modification of our earlier constitutive models [6,7] aimed at predicting the large strain time-dependent behavior of both crosslinked
Fig. 1. (a) Rheological representation of the augmented HM and (b) deformation map showing the kinematics and stress tensors used in the augmented HM. These figures illustrate how the model represents the viscoplastic flow, and how the deformation state is generalized into three dimensions.
ARTICLE IN PRESS J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
as follows: the interaction between the amorphous and crystalline domains in UHMWPE is complicated by entanglements due to its very high molecular weight and also due to chemical crosslinks (when present). At large deformations, however, the underlying molecular deformation resistance, the ‘‘backstress’’ network of molecular chains, has the ability to undergo viscoplastic flow. This flow behavior is caused by the absence of an isotropic crosslinked micro-state in the material, which creates both regions with highly stretched molecular chains and regions that are less stretched. The flow behavior is a function of the highly deformed material state and the interaction between the amorphous and crystalline domains, and can be accurately captured using an energy activation representation. The kinematics of the viscoplastic flow of the backstress network is captured by decomposing the deformation gradient acting on part B of the backstress network (Fig. 1a) into elastic and viscoelastic components: Fp ¼ FeB FvB : The Cauchy stress in the system is given by the isotropic linear elastic relationship: T¼
1 ð2m Ee þ le tr½Ee 1Þ; Je e
ð1Þ
where me and le are Lame! ’s constants which can be obtained from the Young’s modulus and Poisson’s ratio by me ¼ Ee =ð2ð1 þ ne ÞÞ and le ¼ Ee ne =ðð1 þ ne Þð1 2ne ÞÞ; J e ¼ det½Fe is the relative volume change of the elastic deformation, Fe is the deformation gradient, Ee ¼ ln½Ve is the logarithmic true strain, and Ve is the left stretch tensor [8] which can be obtained from the polar decomposition of Fe : The stress acting on the equilibrium portion of the backstress network is given by the same expression as used in our earlier work [6]: TA ¼
1 ½T8chain ðFp ; mA ; llock A ; kA Þ 1 þ qA þ qA TI2 ðFp ; mA Þ;
2I p TI2 ¼ mA I1p Bp 2 1 ðBp Þ2 ; 3
ð2Þ ð3Þ
where TA is a tensor-valued function of the viscoplastic deformation gradient Fp and the material parameters lock fmA ; llock A ; kA ; qA g; where mA is the shear modulus, lA is the locking stretch, kA is the bulk modulus, and qA is a material parameter specifying the relative magnitudes of T8chain and TI2 ; and Bp is the left Cauchy–Green deformation tensor. This hyperelastic stress representation is based on the 8-chain model [9], and a term containing I2 -dependence of the strain energy density. The I2 -dependence is introduced by the crystalline domains and is manifested by the asymmetry in the response between tension and compression [6]. The stress driving the viscoplastic flow of the backstress network is obtained from the same hyper-
2173
elastic representation that was used to calculate the backstress, and has a similar framework as used in the . Bergstrom–Boyce representation of crosslinked polymers at high temperatures [10,11]: TB ¼ sB TA ðFeB Þ;
ð4Þ
where sB is a dimensionless material parameter specifying the relative stiffness of the backstress network. At small deformations, the stiffness of the backstress network is constant and the material response is linear elastic. At larger applied deformations, viscoplastic flow caused by molecular chain sliding is initiated. With increasing viscoplastic flow, the crystalline domains become distorted and provide additional molecular material to the backstress network. This is manifested by an initial reduction in the effective stiffness of the backstress network with imposed viscoplastic deformation and is captured in the model by allowing the parameter sB to evolve with the plastic deformation. The parameter sB evolves with imposed plastic deformation to capture the distributed yielding: s’B ¼ pB ðsB sBf Þ g’ C ;
ð5Þ
where pB is a material parameter specifying the transition rate of the distributed yielding event, sBf is the final value of sB reached at fully developed plastic flow, and g’ C is the magnitude of the viscoplastic flow rate (Eq. (9)): t B mB g’ vB ¼ g’ 0 base : ð6Þ tB The velocity gradient of the viscoelastic flow of the backstress network is given by dev½TB e LvB ¼ g’ vB Fe1 FB ; ð7Þ B tB where g’ vB is the rate of viscoplastic flow of the timedependent network B, tB ¼ jjdev½TB jjF ; tbase B and mB are material parameters, and g’ 0 is a constant coefficient with a value of 1=s: The yielding and plastic flow of the material is captured in the same way as in our earlier work [6,7]: p eT dev½TC e L ¼ g’ C R R ; ð8Þ tC p where Lp ¼ F’ Fp1 ; TC ¼ T ½Fe ðTA þ TB ÞFeT =J e is the stress acting on the relaxed configuration convected to the current configuration, tC ¼ jjdev½TC jjF is the effective shear stress (calculated using the Frobenius norm) driving the viscoplastic flow, mC g’ C ¼ g’ 0 ðtC =tbase C Þ
ð9Þ tbase C
is the magnitude of the viscoplastic flow, and mC are material parameters, and g’ 0 a constant coefficient with a value of 1=s: In total, the augmented HM contains 13 material parameters: two small strain elastic constants (Ee ; ne );
ARTICLE IN PRESS 2174
J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
four hyperelastic constants for the backstress network ðmA ; llock A ; kA ; qA Þ; five flow constants of the backstress network ðsBi ; sBf ; pB ; tbase B ; mB Þ; and two yield and viscoplastic flow parameters ðtbase C ; mC Þ: These parameters can readily be determined from a few select experiments, as will be discussed in Section 3.
3. Materials and methods Section 3.1 first describes the different types of UHMWPE that were examined in this study and the experimental techniques that were used to characterize the material behavior. The methods and procedures that were used to calibrate and validate the predictions from the HM are then described in Section 3.2. 3.1. Experimental In this work, we have focused on one radiation sterilized, and two highly crosslinked GUR 1050 materials. These materials have been characterized and tested in previous studies using uniaxial tension, uniaxial compression, uniaxial cyclic loading, and small punch testing. The details of the material preparations and the experimental data can found elsewhere [6,12]; however, a short summary is provided herein for clarity. 3.1.1. Previous materials and testing Ram-extruded GUR 1050 was used as the base material. All test samples were cut such that the loading direction coincided with the extrusion direction. Three groups of specimens were created. The first group was gamma radiation sterilized in nitrogen with a dose of 30 kGy (‘‘30 kGy g-N2’’), the second group was gamma irradiated with a dose of 100 kGy and then heat treated at 110 C for 2 h (‘‘100 kGy (110 C)’’), and the third group was gamma irradiated with a dose of 100 kGy and then heat treated at 150 C for 2 h (‘‘100 kGy (150 C)’’). After all material preparations, all specimens were stored in a 20 C freezer to minimize aging and oxidation effects. The micro-structure of the materials studied in this work has been extensively examined elsewhere [6,12]; e.g., the degree of crystallinity of the three materials has been determined to be 0.51 for the sterilized materials (30 kGy g-N2), 0.61 for the crosslinked material that was heat treated at 110 C, and 0.46 for the crosslinked material heat treated at 150 C. Data from three different types of room-temperature experiments were analyzed in this study. The first test type was uniaxial tension to failure at three different deformation rates (approximately corresponding to true strain rates of 0.007, 0.018, and 0.035/s). The second test type was cyclic uniaxial fully reversed tension–compression experiments. In these experiments, cylindrical specimens were cyclically loaded and unloaded to a
maximum true strain of 0.12, and a minimum true strain of 0.12. The first two load–unload cycles were analyzed. The last type of experimental data that was analyzed was from a multiaxial small punch test. In these multiaxial tests, miniaturized disc specimens with a diameter of 6.4 mm and a thickness of 0.5 mm were tested by indentation with a hemispherical head punch at a constant punch displacement rate of 0.5 mm/min. The experimental test setup recorded the punch force as a function of punch displacement.
3.2. Analytical The capability of the augmented HM to predict the response of UHMWPE was evaluated by comparing the model predictions with the aforementioned experimental data for the three materials. The first step in this effort was to calibrate the HM to the uniaxial tensile and cyclic experimental data, for each of the materials. For this purpose, the same procedure that was described in our previous work [6] was followed and is briefly summarized. The first step, the bootstrapping step, is to find an initial estimation of the material parameters. In this study, we used material parameters determined from our earlier work [6]. Then, a specialized computer program based on the Nelder–Mead simplex minimization algorithm was used to iteratively improve the correlation between the predicted data sets and the experimental data. The quality of a theoretical prediction, and therefore of the chosen material parameters, was evaluated by calculating the coefficient of determination (r2 ). The reported material parameters for each material are from the set having the highest r2 -value. After the optimal set of material parameters was found, the same parameters were then used to simulate the small punch test. This validation simulation was performed to check the capability of the augmented HM to predict a multiaxial deformation history. It is well known that many constitutive models can predict uniaxial deformation histories relatively well, but that it is significantly more difficult to accurately predict multiaxial deformation states. It has been shown, for example, that the J2 -plasticity model can accurately predict monotonic uniaxial tension or compression data for UHMWPE, but is very poor at predicting cyclic or multiaxial deformation states [6]. The small punch validation simulations were performed using the ABAQUS (HKS Inc., RI) finite element package. The simulations used an axisymmetric representation with 360 quadratic triangular elements (CAX8H) to represent the small punch geometry (see inset in Fig. 8). In the simulations the friction coefficient between the specimen and the punch, and between the specimen and the die was taken as 0.1 [6]. The quality of the validation simulation was evaluated by plotting the
ARTICLE IN PRESS J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
2175
predicted and experimental force–displacement data and by calculating the r2 -value of the predictions.
4. Results The material parameters for the three UHMWPE materials for the augmented HM are given in Table 1. As with our previous constitutive theory, nine of the parameters of the augmented HM were found to be the same for the conventional and the two highly crosslinked UHMWPEs; that is, only four material parameters are dependent on crosslinking density and thermal treatment. These four material parameters are elastic (Young’s) modulus (E), yield strength ðtbase B Þ; the effective stiffness after yield (mA ), and the limiting chain stretch ðllock A Þ; which controls the large strain behavior. A direct comparison between the experimental and the predicted data used for the calibration is shown in Figs. 2–7. Figs. 2 and 3 show the results for the sterilized GUR 1050 (30 kGy, g-N2) in monotonic large strain tension to failure, and cyclic loading with a strain amplitude of 0.12, respectively. Figs. 4 and 5 show the results for highly crosslinked GUR 1050 (100 kGy, 110 C), and Figs. 6 and 7 show the results for the highly crosslinked GUR 1050 (100 kGy, 150 C) that was heat treated at 150 C for 2 h. For all materials, the HM does a very good job of predicting both the large strain tensile data Table 1 Hybrid model (HM) material parameters for the three different types of GUR 1050 Material parameter
30 kGy g-N2
100 kGy g 110 C
100 kGy g 150 C
Ee (MPa) ne mA (MPa) llock A kA (MPa) qA sBi sBf pB tbase (MPa) B mB (MPa) tbase C mC
2020 0.46 8.22 4.40 2000 0.20 40.0 10.0 27.0 25.0 9.50 8.00 3.30
2009 0.46 10.15 2.80 2000 0.20 40.0 10.0 27.0 26.2 9.50 8.00 3.30
1270 0.46 8.14 2.52 2000 0.20 40.0 10.0 27.0 20.7 9.50 8.00 3.30
Ee is the Young’s modulus, ne is the Poisson’s ratio, mA is the shear modulus of network A, llock is the locking stretch of network A, kA is A the bulk modulus of network A, qA is a parameter specifying the asymmetry between tension and compression, sBi is a parameter that controls the initial flow resistance, sBf is a parameter that controls the final flow resistance, pB is a parameter that controls the distributed is a parameter that controls the yield strength of network yielding, tbase B B, mB is a parameter controlling the rate-dependence of network B, tbase is a parameter that controls the yield strength of network C, and C mC is a parameter that controls the rate-dependence of network B. Parameters that are unique for each material are written in bold text.
Fig. 2. Comparison between experimental uniaxial compression data and predictions from the HM for GUR 1050 (30 kGy, g-N2). The three data sets are for true strain rates of 0.007, 0.018, and 0.035/s.
Fig. 3. Comparison between experimental uniaxial cyclic tension and compression data and predictions from the HM for GUR 1050 (30 kGy, g-N2). The experimental data correspond to a true strain rate of 0.05/s.
Fig. 4. Comparison between experimental uniaxial compression data and predictions from the HM for GUR 1050 (100 kGy, 110 C). The three data sets are for true strain rates of 0.007, 0.018, and 0.035/s.
and the small strain cyclic data. The r2 -values were 0.98 or higher for all cases, except the prediction of the tensile behavior of the sterilized conventional material (30 kGy, g-N2). In this case, the r2 -value was slightly lower (0.973), mainly due to variability in the experimental data (the
ARTICLE IN PRESS 2176
J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178 Table 2 Summary of the performance of the HM to predict the response of GUR 1050
Fig. 5. Comparison between experimental uniaxial cyclic tension and compression data and predictions from the HM for GUR 1050 (100 kGy, 110 C). The experimental data correspond to a true strain rate of 0.05/s.
Fig. 6. Comparison between experimental uniaxial compression data and predictions from the HM for GUR 1050 (100 kGy, 150 C). The three data sets are for true strain rates of 0.007, 0.018, and 0.035/s.
Fig. 7. Comparison between experimental uniaxial cyclic tension and compression data and predictions from the HM for GUR 1050 (100 kGy, 150 C). The experimental data correspond to a true strain rate of 0.05/s.
stress–strain curves at different rates crossed each other at high strain levels). A summary of the predictive performance of the HM is given in Table 2.
GUR1050 material
Test mode
r2 -value
30 kGy g-N2
Uniaxial tension Uniaxial cyclic loading Small punch
0.978 0.984 0.937
100 kGy g 110 C
Uniaxial tension Uniaxial cyclic loading Small punch
0.987 0.988 0.960
100 kGy g 150 C
Uniaxial tension Uniaxial cyclic loading Small punch
0.980 0.990 0.948
Fig. 8. Comparison between experimental cyclic tension and compression data predictions from the original HM [6] for GUR 1050 (30 kGy, g-N2). The experimental data correspond to a true strain rate of 0.05/s.
The performance of the old HM [6] is illustrated in Fig. 8. The figure compares cyclic experimental data for GUR 1050 (30 kGy, g-N2) with predictions from the old HM. The material parameters that are used in this simulation are the same as those used in the original work [6]. The figure shows that the old model representation, which has been shown to work very well [6] for large strain tension, compression, and small punch loading, is not accurate at predicting cyclic loading. The new model specifically addresses this issue and enables accurate simulations of both monotonic and cyclic loading conditions using one set of material parameters. The calibrated material models were then used in a finite element model to predict the behavior in the small punch tests. The results from these validation simulations are summarized in Table 2 and in Figs. 9–11. The figures show that for all three materials the new HM does a good job of predicting also the multiaxial deformation in the small punch test, including the initial elastic slope, small strain yielding, large scale yielding, and strain localization during the biaxial stretching. The r2 -values for the small punch predictions are between 0.937 and 0.960 for the three materials.
ARTICLE IN PRESS J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
2177
5. Discussion
Fig. 9. Comparison between experimental small punch data and predictions from the HM for GUR 1050 (30 kGy, g-N2). The experimental data correspond to a punch rate of 0.5 mm/min. The figure also shows the FE mesh that was used in the small punch simulations.
Fig. 10. Comparison between experimental small punch data and predictions from the HM for GUR 1050 (100 kGy, 110 C). The experimental data correspond to a punch rate of 0.5 mm/min.
Fig. 11. Comparison between experimental small punch data and predictions from the HM for GUR 1050 (100 kGy, 150 C). The experimental data correspond to a punch rate of 0.5 mm/min.
The mechanical response of UHMWPE at large deformations is very complex, considering the nonlinear behavior during both loading and unloading. Initially, at small strains, the response is linear elastic. With increasing deformation, localized yielding is initiated at sites where the flow resistance is the lowest. The flow resistance then evolves and becomes more homogeneous in both the crystalline and the amorphous domains. Finally, at large deformations the imposed molecular chain stretching and alignment causes a stiffening in the response which continues to increase until final failure. To model these events is challenging, but necessary for developing a better understanding of the fatigue, fracture, and wear response. Despite the complexity inherent in the constitutive framework of our augmented HM, we found that only four independent material properties were needed to define the overall mechanical behavior of the conventional and highly crosslinked UHMWPE investigated in the present study for the loading and unloading histories that were considered. As in our previous constitutive model, the majority of the material parameters associated with the elastic, plastic, and backstress (recovery) behavior of UHMWPE appear to be unaffected by radiation crosslinking and thermal treatment. Although the constitutive equations used to describe our augmented HM have increased somewhat in complexity, as compared with our previous hybrid theory [6], the number of independent material properties necessary to characterize conventional and highly crosslinked UHMWPE has remained unchanged in both our previous and current theoretical frameworks. Consequently, the augmented HM outlined in our present study is proposed to be a unified constitutive theory for conventional and highly crosslinked UHMWPE materials, in the sense that it is consistent with our previous constitutive modeling approach, as well as in the sense that it appears equally applicable to conventional and highly crosslinked UHMWPE. The augmented HM has the same foundation as our previous modeling efforts [6]. The main difference is that the augmented model now also incorporates relative sliding (reptation) of the molecular chains of the backstress network that carries the main load at moderate to large deformations. The results from this study implicitly show that the relative sliding of the molecular chains in the backstress network is a unified feature of the UHMWPE, both uncrosslinked and crosslinked, mechanical behavior. Figs. 2–7 show that the HM accurately captures large strain tension and small strain cyclic loading of conventional and highly crosslinked UHMWPEs. These tests are straightforward to perform and sufficient for calibrating the model.
ARTICLE IN PRESS 2178
J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
In this study, we have focused on creating a mathematical representation of the deformation resistance and flow characteristics for conventional and highly crosslinked UHMWPE at the molecular level. This effort has focused on the physics of the deformation mechanisms by establishing the framework and equations necessary to model the behavior on the macro-scale. As already mentioned, to use the constitutive model for a given material requires a calibration step where material specific parameters are determined. A variety of numerical methods may be used to determine the material specific parameters for a constitutive theory. In our study, we chose to employ numerical optimization techniques to identify the material parameters for our constitutive theory, as opposed to graphical techniques or simple curve fitting. Of greater importance is how well the physics-inspired model framework represents the governing micromechanisms, and ultimately, how well the model can predict the behavior of a given material under different loading conditions than that for which the model was originally calibrated. The simulations of the small punch test performed in this study demonstrate that our modeling approach provides satisfactory and valid predictions of large-deformation multiaxial behavior of conventional and highly crosslinked UHMWPEs. Thus, our augmented HM yields similar consistent and valid results under large-deformation multiaxial behavior as were observed with our earlier constitutive theory [6]. However, we have now introduced a key new feature to our augmented constitutive theory, which was not incorporated in the previous HM; namely, the new ability to accurately capture the nonlinear unloading behavior of conventional and highly crosslinked UHMWPEs. In summary, the augmented HM is an accurate, validated and unified material model for simulating the loading as well as the unloading behavior of conventional and highly crosslinked UHMWPE used in joint replacements. In the present work, we have restricted our attention to cyclic uniaxial mechanical behavior at room temperature. Based on earlier testing [12], some adjustment of properties is expected for body temperature due to thermal softening. Consequently, research is ongoing to evaluate the performance of the augmented HM at body temperature during cyclic multiaxial loading. In addition, fatigue, fracture, and ultimately
wear are targeted to be studied using the augmented HM as an essential tool.
Acknowledgements This work was supported by NIH Grant 1 R01 AR 47192. Special thanks to M. Villarraga and L. Ciccarelli for assistance with the uniaxial and small punch testing.
References [1] Kurtz SM, Muratoglu OK, Evans M, Edidin AA. Advances in the processing, sterilization, and crosslinking of ultra-high molecular weight polyethylene for total joint arthroplasty. Biomaterials 1999;20:1659–88. [2] Muratoglu OK, Kurtz SM. Alternative bearing surfaces in hip replacement. In: Sinha R, editor. Hip replacement: current trends and controversies. New York: Marcel Dekker; 2002. [3] Wang A, Essner A, Polineni VK, Stark C, Dumbleton JH. Lubrication and wear of ultra-high molecular weight polyethylene in total joint replacements. Tribol Int 1998;31:17–33. [4] Pooley CM, Tabor D. Friction and molecular structure: the behavior of some thermoplastics. Proc R Soc Lond 1972;329: 251–74. [5] Edidin AA, Pruitt L, Jewett CW, Crane DJ, Roberts D, Kurtz SM. Plasticity-induced damage layer is a precursor to wear in radiation—cross-linked UHMWPE acetabular components for total hip replacement. Ultra-high-molecular-weight polyethylene. J Arthroplasty 1999;14:616–27. . JS, Rimnac CM, Kurtz SM. Prediction of multiaxial [6] Bergstrom behavior for conventional and highly crosslinked UHMWPE using a hybrid constitutive model. Biomaterials 2003;24:1365–80. . JS, Kurtz SM, Rimnac CM, Edidin AA. Constitutive [7] Bergstrom modeling of ultra-high molecular weight polyethylene under large-deformation and cyclic loading conditions. Biomaterials 2002;23:2329–43. [8] Gurtin ME. An introduction to continuum mechanics. New York: Academic Press; 1981. [9] Arruda EM, Boyce MC. A three-dimensional constitutive model for the large stretch behavior of rubber elastic materials. J Mech Phys Solids 1993;41:389–412. . JS, Boyce MC. Large strain time-dependent behavior [10] Bergstrom of filled elastomers. Mech Mater 2000;32:627–44. . JS, Boyce MC. Constitutive modelling of the large [11] Bergstrom strain time-dependent behavior of elastomers. J Mech Phys Solids 1998;46:931–54. . JS, Rimnac CM, [12] Kurtz SM, Villarraga ML, Herr MP, Bergstrom Edidin AA. Thermomechanical behavior of virgin and highly crosslinked ultra-high molecular weight polyethylene used in total joint replacements. Biomaterials 2002;23:3681–97.
Biomaterials 25 (2004) 2171–2178
An augmented hybrid constitutive model for simulation of unloading and cyclic loading behavior of conventional and highly crosslinked UHMWPE . a,*, C.M. Rimnacb, S.M. Kurtza,c J.S. Bergstrom a Exponent, Inc., 21 Strathmore Road, Natick, MA 01760, USA Departments of Mechanical and Aerospace Engineering and Orthopaedics, Musculoskeletal Mechanics and Materials Laboratories, Case Western Reserve University, Cleveland, OH, USA c Implant Research Center, School of Biomedical Engineering, Science and Health Systems, Drexel University, 3141 Chestnut Street, Philadelphia, PA, USA b
Received 3 June 2003; accepted 15 August 2003
Abstract Ultra-high molecular weight polyethylene (UHMWPE) is extensively used in total joint replacements. Wear, fatigue, and fracture have limited the longevity of UHMWPE components. For this reason, significant effort has been directed towards understanding the failure and wear mechanisms of UHMWPE, both at a micro-scale and a macro-scale, within the context of joint replacements. We have previously developed, calibrated, and validated a constitutive model for predicting the loading response of conventional and highly crosslinked UHMWPE under multiaxial loading conditions (Biomaterials 24 (2003) 1365). However, to simulate in vivo changes to orthopedic components, accurate simulation of unloading behavior is of equal importance to the loading phase of the duty cycle. Consequently, in this study we have focused on understanding and predicting the mechanical response of UHMWPE during uniaxial unloading. Specifically, we have augmented our previously developed constitutive model to also allow for accurate predictions of the unloading behavior of conventional and highly crosslinked UHMWPE during cyclic loading. It is shown that our augmented hybrid model accurately captures the experimentally observed characteristics, including uniaxial cyclic loading, large strain tension, rate-effects, and multiaxial deformation histories. The augmented hybrid constitutive model will be used as a critical building block in future studies of fatigue, failure, and wear of UHMWPE. r 2003 Elsevier Ltd. All rights reserved. Keywords: Constitutive modeling; Ultra-high molecular weight polyethylene; UHMWPE; Hybrid model; FEM; Radiation crosslinking; Multiaxial mechanical behavior; Small punch test
1. Introduction Wear of the articulating surface of ultra-high molecular weight polyethylene (UHMWPE) implant components is an important problem that can significantly limit the life expectancy of total joint replacements. Wear of UHMWPE components is multifactorial and is influenced by the functional loading environment, joint kinematics, component geometry, and material properties. Recently, efforts to reduce UHMWPE wear have involved changes to the resin type, sterilization method, *Corresponding author. Tel.: +1-508-652-8514; fax: +1-508-6471899. . E-mail address: [email protected] (J.S. Bergstrom). 0142-9612/$ - see front matter r 2003 Elsevier Ltd. All rights reserved. doi:10.1016/j.biomaterials.2003.08.065
radiation crosslinking, and thermal treatments [1,2]. However, there is, at present, an incomplete understanding of the wear characteristics and mechanisms of damage evolution, complicating and impeding rapid progression and improvements in performance of UHMWPE joint replacement components. There are two complementary approaches for improving the general understanding of the wear behavior of UHMWPE components used in total joint replacements: macroscopic experimental testing and microstructural material characterization. Wear simulators and other mechanical testing techniques can provide information related to wear rates of different tribological systems and can rank the performance of materials subjected to different environments and
ARTICLE IN PRESS 2172
J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
thermomechanical histories [3]. Although useful, systematic empirical testing has thus far not enabled a priori predictions of the mechanisms causing the actual wear of UHMWPE. To predict the evolution in microscopic and macroscopic damage for new UHMWPE materials from fundamental polymer physics principles requires an understanding of the performance and response of the material on the microstructural level. Theoretical and experimental research supports the notion that the wear observed in vivo and in vitro is the result of localized high stresses and strains in the surface region of the UHMWPE component [4,5]. To better understand and predict these stresses, it is necessary to have a well-calibrated and accurate constitutive model of UHMWPE. In the orthopedic research community, the J2-plasticity model has been the most widely used approach for simulating the behavior of UHMWPE. It has been shown [6], however, that the J2-plasticity model is not an accurate general tool for predicting the largedeformation-to-failure behavior of UHMWPE. In addition, the J2-plasticity model does not accurately predict cyclic loading of UHMWPE. These are serious limitations since UHMWPE joint components undergo large deformations locally at the articulating surface and are also subject to cyclically applied loads. To address these limitations, a new constitutive model was recently developed for conventional and highly crosslinked UHMWPEs [6]. This new model, which is inspired by the physical micro-mechanisms governing the deformation resistance of polymeric materials, is an extension of specialized constitutive theories for glassy polymers that have been developed during the last 10 years. The new model, named the hybrid model (HM), has been shown to accurately predict the mechanical response of both conventional and highly crosslinked UHMWPE materials in uniaxial tension, compression, and multiaxial loading. However, to simulate in vivo changes to orthopedic components, accurate simulation of unloading behavior is also of importance. Consequently, the objective of this study was to better understand and predict the mechanical response of UHMWPE during uniaxial unloading. In this regard, we have developed an augmented hybrid constitutive model that is capable of accurately predicting the experimentally observed stress–strain response in cyclic loading for conventional and well as highly crosslinked UHMWPEs.
and uncrosslinked UHMWPE. The modification in the augmented HM specifically addresses the unloading behavior during cyclic loading. The kinematic framework used in the augmented HM is based on a decomposition of the applied deformation gradient into elastic and viscoplastic components: F=FeFp (Fig. 1). The spring and dashpot representation shown in Fig. 1a is a one-dimensional embodiment of the model framework used to capture the viscoplastic flow characteristics. With the exception of the top spring (E), all spring and dashpot elements are highly nonlinear (described in detail below). Fig. 1b depicts a map of the decomposition of a given material deformation state. This decomposition specifies how the three-dimensionality of the deformation gradients and stress tensors are connected and evolve during an applied deformation history. As in our previous modeling approaches [6,7], the deformation state is decomposed into elastic, backstress, and viscoplastic components. Compared with our previous models, we have now incorporated timedependent viscoplasticity to the backstress network to improve the predictive capabilities of the model with respect to unloading. The rationale for this approach is
2. Augmented hybrid constitutive model for predictions of UHMWPE The new augmented HM is a modification of our earlier constitutive models [6,7] aimed at predicting the large strain time-dependent behavior of both crosslinked
Fig. 1. (a) Rheological representation of the augmented HM and (b) deformation map showing the kinematics and stress tensors used in the augmented HM. These figures illustrate how the model represents the viscoplastic flow, and how the deformation state is generalized into three dimensions.
ARTICLE IN PRESS J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
as follows: the interaction between the amorphous and crystalline domains in UHMWPE is complicated by entanglements due to its very high molecular weight and also due to chemical crosslinks (when present). At large deformations, however, the underlying molecular deformation resistance, the ‘‘backstress’’ network of molecular chains, has the ability to undergo viscoplastic flow. This flow behavior is caused by the absence of an isotropic crosslinked micro-state in the material, which creates both regions with highly stretched molecular chains and regions that are less stretched. The flow behavior is a function of the highly deformed material state and the interaction between the amorphous and crystalline domains, and can be accurately captured using an energy activation representation. The kinematics of the viscoplastic flow of the backstress network is captured by decomposing the deformation gradient acting on part B of the backstress network (Fig. 1a) into elastic and viscoelastic components: Fp ¼ FeB FvB : The Cauchy stress in the system is given by the isotropic linear elastic relationship: T¼
1 ð2m Ee þ le tr½Ee 1Þ; Je e
ð1Þ
where me and le are Lame! ’s constants which can be obtained from the Young’s modulus and Poisson’s ratio by me ¼ Ee =ð2ð1 þ ne ÞÞ and le ¼ Ee ne =ðð1 þ ne Þð1 2ne ÞÞ; J e ¼ det½Fe is the relative volume change of the elastic deformation, Fe is the deformation gradient, Ee ¼ ln½Ve is the logarithmic true strain, and Ve is the left stretch tensor [8] which can be obtained from the polar decomposition of Fe : The stress acting on the equilibrium portion of the backstress network is given by the same expression as used in our earlier work [6]: TA ¼
1 ½T8chain ðFp ; mA ; llock A ; kA Þ 1 þ qA þ qA TI2 ðFp ; mA Þ;
2I p TI2 ¼ mA I1p Bp 2 1 ðBp Þ2 ; 3
ð2Þ ð3Þ
where TA is a tensor-valued function of the viscoplastic deformation gradient Fp and the material parameters lock fmA ; llock A ; kA ; qA g; where mA is the shear modulus, lA is the locking stretch, kA is the bulk modulus, and qA is a material parameter specifying the relative magnitudes of T8chain and TI2 ; and Bp is the left Cauchy–Green deformation tensor. This hyperelastic stress representation is based on the 8-chain model [9], and a term containing I2 -dependence of the strain energy density. The I2 -dependence is introduced by the crystalline domains and is manifested by the asymmetry in the response between tension and compression [6]. The stress driving the viscoplastic flow of the backstress network is obtained from the same hyper-
2173
elastic representation that was used to calculate the backstress, and has a similar framework as used in the . Bergstrom–Boyce representation of crosslinked polymers at high temperatures [10,11]: TB ¼ sB TA ðFeB Þ;
ð4Þ
where sB is a dimensionless material parameter specifying the relative stiffness of the backstress network. At small deformations, the stiffness of the backstress network is constant and the material response is linear elastic. At larger applied deformations, viscoplastic flow caused by molecular chain sliding is initiated. With increasing viscoplastic flow, the crystalline domains become distorted and provide additional molecular material to the backstress network. This is manifested by an initial reduction in the effective stiffness of the backstress network with imposed viscoplastic deformation and is captured in the model by allowing the parameter sB to evolve with the plastic deformation. The parameter sB evolves with imposed plastic deformation to capture the distributed yielding: s’B ¼ pB ðsB sBf Þ g’ C ;
ð5Þ
where pB is a material parameter specifying the transition rate of the distributed yielding event, sBf is the final value of sB reached at fully developed plastic flow, and g’ C is the magnitude of the viscoplastic flow rate (Eq. (9)): t B mB g’ vB ¼ g’ 0 base : ð6Þ tB The velocity gradient of the viscoelastic flow of the backstress network is given by dev½TB e LvB ¼ g’ vB Fe1 FB ; ð7Þ B tB where g’ vB is the rate of viscoplastic flow of the timedependent network B, tB ¼ jjdev½TB jjF ; tbase B and mB are material parameters, and g’ 0 is a constant coefficient with a value of 1=s: The yielding and plastic flow of the material is captured in the same way as in our earlier work [6,7]: p eT dev½TC e L ¼ g’ C R R ; ð8Þ tC p where Lp ¼ F’ Fp1 ; TC ¼ T ½Fe ðTA þ TB ÞFeT =J e is the stress acting on the relaxed configuration convected to the current configuration, tC ¼ jjdev½TC jjF is the effective shear stress (calculated using the Frobenius norm) driving the viscoplastic flow, mC g’ C ¼ g’ 0 ðtC =tbase C Þ
ð9Þ tbase C
is the magnitude of the viscoplastic flow, and mC are material parameters, and g’ 0 a constant coefficient with a value of 1=s: In total, the augmented HM contains 13 material parameters: two small strain elastic constants (Ee ; ne );
ARTICLE IN PRESS 2174
J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
four hyperelastic constants for the backstress network ðmA ; llock A ; kA ; qA Þ; five flow constants of the backstress network ðsBi ; sBf ; pB ; tbase B ; mB Þ; and two yield and viscoplastic flow parameters ðtbase C ; mC Þ: These parameters can readily be determined from a few select experiments, as will be discussed in Section 3.
3. Materials and methods Section 3.1 first describes the different types of UHMWPE that were examined in this study and the experimental techniques that were used to characterize the material behavior. The methods and procedures that were used to calibrate and validate the predictions from the HM are then described in Section 3.2. 3.1. Experimental In this work, we have focused on one radiation sterilized, and two highly crosslinked GUR 1050 materials. These materials have been characterized and tested in previous studies using uniaxial tension, uniaxial compression, uniaxial cyclic loading, and small punch testing. The details of the material preparations and the experimental data can found elsewhere [6,12]; however, a short summary is provided herein for clarity. 3.1.1. Previous materials and testing Ram-extruded GUR 1050 was used as the base material. All test samples were cut such that the loading direction coincided with the extrusion direction. Three groups of specimens were created. The first group was gamma radiation sterilized in nitrogen with a dose of 30 kGy (‘‘30 kGy g-N2’’), the second group was gamma irradiated with a dose of 100 kGy and then heat treated at 110 C for 2 h (‘‘100 kGy (110 C)’’), and the third group was gamma irradiated with a dose of 100 kGy and then heat treated at 150 C for 2 h (‘‘100 kGy (150 C)’’). After all material preparations, all specimens were stored in a 20 C freezer to minimize aging and oxidation effects. The micro-structure of the materials studied in this work has been extensively examined elsewhere [6,12]; e.g., the degree of crystallinity of the three materials has been determined to be 0.51 for the sterilized materials (30 kGy g-N2), 0.61 for the crosslinked material that was heat treated at 110 C, and 0.46 for the crosslinked material heat treated at 150 C. Data from three different types of room-temperature experiments were analyzed in this study. The first test type was uniaxial tension to failure at three different deformation rates (approximately corresponding to true strain rates of 0.007, 0.018, and 0.035/s). The second test type was cyclic uniaxial fully reversed tension–compression experiments. In these experiments, cylindrical specimens were cyclically loaded and unloaded to a
maximum true strain of 0.12, and a minimum true strain of 0.12. The first two load–unload cycles were analyzed. The last type of experimental data that was analyzed was from a multiaxial small punch test. In these multiaxial tests, miniaturized disc specimens with a diameter of 6.4 mm and a thickness of 0.5 mm were tested by indentation with a hemispherical head punch at a constant punch displacement rate of 0.5 mm/min. The experimental test setup recorded the punch force as a function of punch displacement.
3.2. Analytical The capability of the augmented HM to predict the response of UHMWPE was evaluated by comparing the model predictions with the aforementioned experimental data for the three materials. The first step in this effort was to calibrate the HM to the uniaxial tensile and cyclic experimental data, for each of the materials. For this purpose, the same procedure that was described in our previous work [6] was followed and is briefly summarized. The first step, the bootstrapping step, is to find an initial estimation of the material parameters. In this study, we used material parameters determined from our earlier work [6]. Then, a specialized computer program based on the Nelder–Mead simplex minimization algorithm was used to iteratively improve the correlation between the predicted data sets and the experimental data. The quality of a theoretical prediction, and therefore of the chosen material parameters, was evaluated by calculating the coefficient of determination (r2 ). The reported material parameters for each material are from the set having the highest r2 -value. After the optimal set of material parameters was found, the same parameters were then used to simulate the small punch test. This validation simulation was performed to check the capability of the augmented HM to predict a multiaxial deformation history. It is well known that many constitutive models can predict uniaxial deformation histories relatively well, but that it is significantly more difficult to accurately predict multiaxial deformation states. It has been shown, for example, that the J2 -plasticity model can accurately predict monotonic uniaxial tension or compression data for UHMWPE, but is very poor at predicting cyclic or multiaxial deformation states [6]. The small punch validation simulations were performed using the ABAQUS (HKS Inc., RI) finite element package. The simulations used an axisymmetric representation with 360 quadratic triangular elements (CAX8H) to represent the small punch geometry (see inset in Fig. 8). In the simulations the friction coefficient between the specimen and the punch, and between the specimen and the die was taken as 0.1 [6]. The quality of the validation simulation was evaluated by plotting the
ARTICLE IN PRESS J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
2175
predicted and experimental force–displacement data and by calculating the r2 -value of the predictions.
4. Results The material parameters for the three UHMWPE materials for the augmented HM are given in Table 1. As with our previous constitutive theory, nine of the parameters of the augmented HM were found to be the same for the conventional and the two highly crosslinked UHMWPEs; that is, only four material parameters are dependent on crosslinking density and thermal treatment. These four material parameters are elastic (Young’s) modulus (E), yield strength ðtbase B Þ; the effective stiffness after yield (mA ), and the limiting chain stretch ðllock A Þ; which controls the large strain behavior. A direct comparison between the experimental and the predicted data used for the calibration is shown in Figs. 2–7. Figs. 2 and 3 show the results for the sterilized GUR 1050 (30 kGy, g-N2) in monotonic large strain tension to failure, and cyclic loading with a strain amplitude of 0.12, respectively. Figs. 4 and 5 show the results for highly crosslinked GUR 1050 (100 kGy, 110 C), and Figs. 6 and 7 show the results for the highly crosslinked GUR 1050 (100 kGy, 150 C) that was heat treated at 150 C for 2 h. For all materials, the HM does a very good job of predicting both the large strain tensile data Table 1 Hybrid model (HM) material parameters for the three different types of GUR 1050 Material parameter
30 kGy g-N2
100 kGy g 110 C
100 kGy g 150 C
Ee (MPa) ne mA (MPa) llock A kA (MPa) qA sBi sBf pB tbase (MPa) B mB (MPa) tbase C mC
2020 0.46 8.22 4.40 2000 0.20 40.0 10.0 27.0 25.0 9.50 8.00 3.30
2009 0.46 10.15 2.80 2000 0.20 40.0 10.0 27.0 26.2 9.50 8.00 3.30
1270 0.46 8.14 2.52 2000 0.20 40.0 10.0 27.0 20.7 9.50 8.00 3.30
Ee is the Young’s modulus, ne is the Poisson’s ratio, mA is the shear modulus of network A, llock is the locking stretch of network A, kA is A the bulk modulus of network A, qA is a parameter specifying the asymmetry between tension and compression, sBi is a parameter that controls the initial flow resistance, sBf is a parameter that controls the final flow resistance, pB is a parameter that controls the distributed is a parameter that controls the yield strength of network yielding, tbase B B, mB is a parameter controlling the rate-dependence of network B, tbase is a parameter that controls the yield strength of network C, and C mC is a parameter that controls the rate-dependence of network B. Parameters that are unique for each material are written in bold text.
Fig. 2. Comparison between experimental uniaxial compression data and predictions from the HM for GUR 1050 (30 kGy, g-N2). The three data sets are for true strain rates of 0.007, 0.018, and 0.035/s.
Fig. 3. Comparison between experimental uniaxial cyclic tension and compression data and predictions from the HM for GUR 1050 (30 kGy, g-N2). The experimental data correspond to a true strain rate of 0.05/s.
Fig. 4. Comparison between experimental uniaxial compression data and predictions from the HM for GUR 1050 (100 kGy, 110 C). The three data sets are for true strain rates of 0.007, 0.018, and 0.035/s.
and the small strain cyclic data. The r2 -values were 0.98 or higher for all cases, except the prediction of the tensile behavior of the sterilized conventional material (30 kGy, g-N2). In this case, the r2 -value was slightly lower (0.973), mainly due to variability in the experimental data (the
ARTICLE IN PRESS 2176
J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178 Table 2 Summary of the performance of the HM to predict the response of GUR 1050
Fig. 5. Comparison between experimental uniaxial cyclic tension and compression data and predictions from the HM for GUR 1050 (100 kGy, 110 C). The experimental data correspond to a true strain rate of 0.05/s.
Fig. 6. Comparison between experimental uniaxial compression data and predictions from the HM for GUR 1050 (100 kGy, 150 C). The three data sets are for true strain rates of 0.007, 0.018, and 0.035/s.
Fig. 7. Comparison between experimental uniaxial cyclic tension and compression data and predictions from the HM for GUR 1050 (100 kGy, 150 C). The experimental data correspond to a true strain rate of 0.05/s.
stress–strain curves at different rates crossed each other at high strain levels). A summary of the predictive performance of the HM is given in Table 2.
GUR1050 material
Test mode
r2 -value
30 kGy g-N2
Uniaxial tension Uniaxial cyclic loading Small punch
0.978 0.984 0.937
100 kGy g 110 C
Uniaxial tension Uniaxial cyclic loading Small punch
0.987 0.988 0.960
100 kGy g 150 C
Uniaxial tension Uniaxial cyclic loading Small punch
0.980 0.990 0.948
Fig. 8. Comparison between experimental cyclic tension and compression data predictions from the original HM [6] for GUR 1050 (30 kGy, g-N2). The experimental data correspond to a true strain rate of 0.05/s.
The performance of the old HM [6] is illustrated in Fig. 8. The figure compares cyclic experimental data for GUR 1050 (30 kGy, g-N2) with predictions from the old HM. The material parameters that are used in this simulation are the same as those used in the original work [6]. The figure shows that the old model representation, which has been shown to work very well [6] for large strain tension, compression, and small punch loading, is not accurate at predicting cyclic loading. The new model specifically addresses this issue and enables accurate simulations of both monotonic and cyclic loading conditions using one set of material parameters. The calibrated material models were then used in a finite element model to predict the behavior in the small punch tests. The results from these validation simulations are summarized in Table 2 and in Figs. 9–11. The figures show that for all three materials the new HM does a good job of predicting also the multiaxial deformation in the small punch test, including the initial elastic slope, small strain yielding, large scale yielding, and strain localization during the biaxial stretching. The r2 -values for the small punch predictions are between 0.937 and 0.960 for the three materials.
ARTICLE IN PRESS J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
2177
5. Discussion
Fig. 9. Comparison between experimental small punch data and predictions from the HM for GUR 1050 (30 kGy, g-N2). The experimental data correspond to a punch rate of 0.5 mm/min. The figure also shows the FE mesh that was used in the small punch simulations.
Fig. 10. Comparison between experimental small punch data and predictions from the HM for GUR 1050 (100 kGy, 110 C). The experimental data correspond to a punch rate of 0.5 mm/min.
Fig. 11. Comparison between experimental small punch data and predictions from the HM for GUR 1050 (100 kGy, 150 C). The experimental data correspond to a punch rate of 0.5 mm/min.
The mechanical response of UHMWPE at large deformations is very complex, considering the nonlinear behavior during both loading and unloading. Initially, at small strains, the response is linear elastic. With increasing deformation, localized yielding is initiated at sites where the flow resistance is the lowest. The flow resistance then evolves and becomes more homogeneous in both the crystalline and the amorphous domains. Finally, at large deformations the imposed molecular chain stretching and alignment causes a stiffening in the response which continues to increase until final failure. To model these events is challenging, but necessary for developing a better understanding of the fatigue, fracture, and wear response. Despite the complexity inherent in the constitutive framework of our augmented HM, we found that only four independent material properties were needed to define the overall mechanical behavior of the conventional and highly crosslinked UHMWPE investigated in the present study for the loading and unloading histories that were considered. As in our previous constitutive model, the majority of the material parameters associated with the elastic, plastic, and backstress (recovery) behavior of UHMWPE appear to be unaffected by radiation crosslinking and thermal treatment. Although the constitutive equations used to describe our augmented HM have increased somewhat in complexity, as compared with our previous hybrid theory [6], the number of independent material properties necessary to characterize conventional and highly crosslinked UHMWPE has remained unchanged in both our previous and current theoretical frameworks. Consequently, the augmented HM outlined in our present study is proposed to be a unified constitutive theory for conventional and highly crosslinked UHMWPE materials, in the sense that it is consistent with our previous constitutive modeling approach, as well as in the sense that it appears equally applicable to conventional and highly crosslinked UHMWPE. The augmented HM has the same foundation as our previous modeling efforts [6]. The main difference is that the augmented model now also incorporates relative sliding (reptation) of the molecular chains of the backstress network that carries the main load at moderate to large deformations. The results from this study implicitly show that the relative sliding of the molecular chains in the backstress network is a unified feature of the UHMWPE, both uncrosslinked and crosslinked, mechanical behavior. Figs. 2–7 show that the HM accurately captures large strain tension and small strain cyclic loading of conventional and highly crosslinked UHMWPEs. These tests are straightforward to perform and sufficient for calibrating the model.
ARTICLE IN PRESS 2178
J.S. Bergstrom . et al. / Biomaterials 25 (2004) 2171–2178
In this study, we have focused on creating a mathematical representation of the deformation resistance and flow characteristics for conventional and highly crosslinked UHMWPE at the molecular level. This effort has focused on the physics of the deformation mechanisms by establishing the framework and equations necessary to model the behavior on the macro-scale. As already mentioned, to use the constitutive model for a given material requires a calibration step where material specific parameters are determined. A variety of numerical methods may be used to determine the material specific parameters for a constitutive theory. In our study, we chose to employ numerical optimization techniques to identify the material parameters for our constitutive theory, as opposed to graphical techniques or simple curve fitting. Of greater importance is how well the physics-inspired model framework represents the governing micromechanisms, and ultimately, how well the model can predict the behavior of a given material under different loading conditions than that for which the model was originally calibrated. The simulations of the small punch test performed in this study demonstrate that our modeling approach provides satisfactory and valid predictions of large-deformation multiaxial behavior of conventional and highly crosslinked UHMWPEs. Thus, our augmented HM yields similar consistent and valid results under large-deformation multiaxial behavior as were observed with our earlier constitutive theory [6]. However, we have now introduced a key new feature to our augmented constitutive theory, which was not incorporated in the previous HM; namely, the new ability to accurately capture the nonlinear unloading behavior of conventional and highly crosslinked UHMWPEs. In summary, the augmented HM is an accurate, validated and unified material model for simulating the loading as well as the unloading behavior of conventional and highly crosslinked UHMWPE used in joint replacements. In the present work, we have restricted our attention to cyclic uniaxial mechanical behavior at room temperature. Based on earlier testing [12], some adjustment of properties is expected for body temperature due to thermal softening. Consequently, research is ongoing to evaluate the performance of the augmented HM at body temperature during cyclic multiaxial loading. In addition, fatigue, fracture, and ultimately
wear are targeted to be studied using the augmented HM as an essential tool.
Acknowledgements This work was supported by NIH Grant 1 R01 AR 47192. Special thanks to M. Villarraga and L. Ciccarelli for assistance with the uniaxial and small punch testing.
References [1] Kurtz SM, Muratoglu OK, Evans M, Edidin AA. Advances in the processing, sterilization, and crosslinking of ultra-high molecular weight polyethylene for total joint arthroplasty. Biomaterials 1999;20:1659–88. [2] Muratoglu OK, Kurtz SM. Alternative bearing surfaces in hip replacement. In: Sinha R, editor. Hip replacement: current trends and controversies. New York: Marcel Dekker; 2002. [3] Wang A, Essner A, Polineni VK, Stark C, Dumbleton JH. Lubrication and wear of ultra-high molecular weight polyethylene in total joint replacements. Tribol Int 1998;31:17–33. [4] Pooley CM, Tabor D. Friction and molecular structure: the behavior of some thermoplastics. Proc R Soc Lond 1972;329: 251–74. [5] Edidin AA, Pruitt L, Jewett CW, Crane DJ, Roberts D, Kurtz SM. Plasticity-induced damage layer is a precursor to wear in radiation—cross-linked UHMWPE acetabular components for total hip replacement. Ultra-high-molecular-weight polyethylene. J Arthroplasty 1999;14:616–27. . JS, Rimnac CM, Kurtz SM. Prediction of multiaxial [6] Bergstrom behavior for conventional and highly crosslinked UHMWPE using a hybrid constitutive model. Biomaterials 2003;24:1365–80. . JS, Kurtz SM, Rimnac CM, Edidin AA. Constitutive [7] Bergstrom modeling of ultra-high molecular weight polyethylene under large-deformation and cyclic loading conditions. Biomaterials 2002;23:2329–43. [8] Gurtin ME. An introduction to continuum mechanics. New York: Academic Press; 1981. [9] Arruda EM, Boyce MC. A three-dimensional constitutive model for the large stretch behavior of rubber elastic materials. J Mech Phys Solids 1993;41:389–412. . JS, Boyce MC. Large strain time-dependent behavior [10] Bergstrom of filled elastomers. Mech Mater 2000;32:627–44. . JS, Boyce MC. Constitutive modelling of the large [11] Bergstrom strain time-dependent behavior of elastomers. J Mech Phys Solids 1998;46:931–54. . JS, Rimnac CM, [12] Kurtz SM, Villarraga ML, Herr MP, Bergstrom Edidin AA. Thermomechanical behavior of virgin and highly crosslinked ultra-high molecular weight polyethylene used in total joint replacements. Biomaterials 2002;23:3681–97.