Mechanics of Materials 104 (2017) 145–156
Contents lists available at ScienceDirect
Mechanics of Materials journal homepage: www.elsevier.com/locate/mechmat
Research paper
High strain rate behavior of MC2 single crystal under uniaxial compression load at high temperature: Experiments and modeling Louis Roux a,b, Hervé Chalons b, Hubert Maigre a, Daniel Nélias a,∗ a b
Univ Lyon, INSA-Lyon, CNRS UMR5259, LaMCoS, F-69621, France Safran Helicopter Engines (formerly Turbomeca), avenue Joseph Szydlowski, 64511 Bordes, France
a r t i c l e
i n f o
Article history: Received 8 December 2015 Revised 3 November 2016 Available online 10 November 2016 Keywords: Split Hopkinson pressure bars MC2 single crystal Cubic material Viscoplasticity Local Hill yield criterion Transient dynamic Signal analysis Digital images correlation High-speed impact Turbine blade
a b s t r a c t This paper presents experimental and numerical investigations on single crystal nickel base superalloy MC2. The purpose of this study is to efficiently predict the consequences of impacts on turbine blades within gas turbine engines. Dynamic compression tests on Split Hopkinson Pressure Bar (SHPB), performed at high temperature, have allowed the identification of anisotropic and viscoplastic material data. The effects of crystal orientation (1 0 0, 1 1 0 and 1 1 1), temperature (25, 750 and 10 0 0 °C) and strain rate (up to 610 s−1 ) on the compression yield strength and work hardening of MC2 single crystal have been examined. The LS-Dyna solver has been used to calibrate the parameters of a phenomenological behavior law based on a local Hill yield criterion. The results of simulations have been verified by comparison with a ballistic test on a monocrystalline plate, and on a turbine blade heated at 10 0 0 °C. The measured residual deflections and strain fields are consistent with those obtained by calculation.
1. Introduction Reliability and integrity of aeronautical gas turbine engines are key factors for flight safety. During the certification process, it must be demonstrated that any single turbine blade will be contained after failure and that no hazardous events will occur within the engine (European Aviation Safety Agency, 2010). The latter could be, for example, multiple blade failures leading to the perforation of casing and, potentially, to the release of high-energy debris. Among the acceptable means of compliance (European Aviation Safety Agency, 2012), analyses based on finite element (FE) models which faithfully represent the physics of the phenomenon are suitable to efficiently predict the consequences of impacts on turbine blades within turbo-engines operating at the nominal rotational speed. The Blade-Off phenomenon is a complex tri-dimensional problem which involves many non-linearities such as large strains and displacements, contact and plasticity. When occurring, several high velocity impacts – up to 300 m/s – are generated for which multiaxial loading, strain rates up to 103 s−1 and high plastic strain levels are involved. To efficiently handle these transient dynamic so-
∗
Corresponding author. E-mail address:
[email protected] (D. Nélias).
http://dx.doi.org/10.1016/j.mechmat.2016.11.002 0167-6636/© 2016 Elsevier Ltd. All rights reserved.
© 2016 Elsevier Ltd. All rights reserved.
licitations, explicit finite element codes are frequently used with a central difference time integration scheme to solve the motion equations (Lemaitre and Chaboche, 2009). The LS-Dyna solver (SMP double precision R7)(Hallquist, 2005), developed by Livermore Software Technology Corporation (LSTC), has been chosen because it offers many possibilities for transient calculations, like penalty based contact algorithms or boundary prescribed motions, and gives access to a large material library (LSTC, 2015b). As a matter of fact, mechanics of materials remains the starting point of explicit dynamic simulations. That explains why particular importance is typically given to viscoplastic behavior of turbine blade constitutive material which determines the main part of the material response under impacts. Because of their excellent high temperature mechanical strength, directed growth single crystal nickel-based superalloys are commonly used for high pressure turbines generally operating up to 10 0 0 °C. The aim of this work is to propose an efficient modeling strategy for a single crystal superalloy under high plastic strain rates at acceptable computational cost by using LS-Dyna solver, allowing thus to meet the industrial requirements. In the present study, the MC2 alloy, widely used to produce high pressure monocrystalline turbine blades at Safran Helicopter Engines, has been chosen to investigate its plastic dynamic behavior at high temperature.
146
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156 Table 1 Nominal composition of MC2 superalloy (Wt %). Ni
Cr
W
Ta
Co
Al
Mo
Ti
64.5
8
8
6
5
5
2
1,5
2. Background review High performance nickel-based superalloys have been designed for high temperature applications (Reed et al., 2009; Pollock and Tin, 2006) and widely used for aero-engines or industrial gas turbines. Within turbo-engines, under nominal operating conditions, the blade airfoils are subjected to elevated temperatures and high stress gradients leading mostly to corrosion, creep and fatigue damages. Thanks to their excellent fatigue and creep characteristics in the design of high-pressure turbine blades, the monocrystalline form – which means no grain boundary – is commonly used to improve the creep resistance, enabling thus to extend the blade service life (Caron et al., 2011; Caron and Lavigne, 2011). Nowadays, the research mainly focuses on the development of multi-scale crystal plasticity-based models which take into account the different damage mechanisms under normal operating conditions (creep, fatigue and corrosion) (Méric et al., 1991). Several studies on fatigue life and creep strength of single crystals (SX) can be found in the literature and provide data on the macroscopic behavior in terms of stress-strain curves, however, they are always limited to relatively low strain rates (i.e. quasistatic conditions, ε˙ < 1 s−1 ) (Segersall, 2013; Zhao et al., 2001) and low plasticity levels (ε p < 3%) (Wang et al., 2009). In these studies, the yield strength of single crystal superalloys remains high (i.e. between 800 and 1200 MPa) for temperatures up to 800 °C then decreases drastically (i.e. less than 400 MPa at 1100 °C) Caron et al. (2011). Furthermore, compression and tensile behaviors were compared at room temperature (R.T.) (Leidermark et al., 2009). The main differences are pointed for the 1 1 0 crystallographic orientation for which the yield strength is close to 930 MPa in tension and around 1070 MPa in compression (difference of 14%). For the other orientations tested in this work (i.e. 1 0 0 and 1 1 1), the differences are less than 1%. To the author’s knowledge no studies are currently available in the literature about high temperature impact tests on either monocrystalline or polycrystalline nickel-based superalloys to treat accidental situations like the blade-off event. According to the literature review, the mechanical behavior under extreme dynamic loads at high temperature of nickel-based superalloys, especially that of single crystals, is not well understood. 3. Presentation of MC2 superalloy The material presented here is a nickel-based single crystal superalloy, designated MC2, which was developed by the French aerospace laboratory ONERA in the late 1980s. As a first generation of single crystal (i.e. better creep strengthening than that of the polycrystalline form but no Rhenium addition inherent to the next generations), such as AM1 and AM3 superalloys, MC2 superalloy has been especially qualified for more powerful turbo-engines. It is now commonly used for high pressure turbine blades at Safran Helicopter Engines for operating temperatures ranging from 750 to 10 0 0 °C. The nominal composition of MC2 superalloy is presented in Table 1. Proportions of other minor elements such as C, S, Mn, P, N or O are controlled. After standard thermal treatments (i.e. solution and ageing heat treatments), the microstructure has a high volume fraction (around 70%) of cuboid precipitates (γ phase) distributed within a matrix of solution Ni (γ phase).
Because of its face centered cubic (FCC) crystalline structure, which presents 12 slip systems located on {111} octahedral planes, the macroscopic material behavior is orthotropic which means different yield strengths depending on the crystallographic orientation. Thus, three types of cylindrical sample corresponding to each main crystalline orientations – 1 0 0, 1 1 0 and 1 1 1 (according to Miller indices) – have been provided by Safran Helicopter Engines material department with an accuracy of ±10° (typical industrial criterion for crystal growth single crystal). These samples have purposely been aligned to the compression axis, see Fig. 1. The directions 1 0 0 and 1 1 1 have allowed to characterize the orthotropic behavior and the 1 1 0 direction has been used to assess the differences from modeling. Note that for the specimens which have been tested here, as for in-flight turbine blades, the material processing has ended with a standard aging heat treatment (1080 °C/6 h/air cooled + 870 °C/20 h/air cooled). 4. SHPB compression tests: typical use and experimental setup 4.1. Principle and relevant assumptions Based on unidirectional wave propagation, one of the most common techniques for dynamic characterization in the range from 102 s−1 to 104 s−1 is the Split Hopkinson Pressure Bar (SHPB) test also called Kolsky bar test (Hopkinson, 1914; Zhao and Gary, 1996). It has been widely used for defense, automotive and aeronautical applications for which materials, such as ductile metals (Liang and Khan, 1999; Ogawa, 1985), are subjected to impacts. In its genuine configuration, which has been used in the present study, this test consists in a unidirectional compression load on a cylindrical specimen located between two cylindrical bars, called respectively incident (or input) and transmitted (or output) bars, see Fig. 2. In the current setup, an initial elastic wave has been triggered by a striker impact on a first bar extremity, as shown in Fig. 2. Then, as the sample normal section, S, has been chosen smaller than that of the straight bars, the stress level is sufficient to plastically deform the specimen which is particularly hard in the case of MC2. In order to estimate the macroscopic behavior of a material under dynamic compression loads, relevant assumptions should be made. The continuity of loads and velocities is typically assumed at each interface. Local measurements are usually performed with strain gauges placed on each bar, far enough from the specimen, so as to record distinctive incident, reflected and transmitted wave signals, respectively corresponding to ε i (t), ε r (t) and ε t (t). These signals are then superimposed thanks to delay calculations based on wave celerity into the input and output bars, C1 and C2 respectively. In the following analyses, the compression stress in the specimen has been determined in the usual way from incident and reflected (i.e. traction) stress pulses, respectively σ i (t) and σ r (t), or from the transmitted stress pulse σ t (t), in accordance with the following equilibrium equation (Zhao and Gary, 1997):
σt (t ) = σi (t ) − σr (t )
(1)
Knowing physical and elastic properties of the bars (i.e. the Young modulus of the bars, EB , density ρ B , and wave celerity CB ) and strain gauge positions, the temporal evolution of velocities V1 (t) and V2 (t) can be estimated at each interface from the strain measurements, after wave transpositions, according to the following equations:
CB =
EB
ρB
V1 (t ) = −CB (εi (t ) − εr (t ))
(2) (3)
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
147
Fig. 1. Crystallographic orientation of samples showing a {111} octahedral slip plan position and compression along the z-axis.
Fig. 2. Standard setup for Hopkinson bar compression tests.
V2 (t ) = −CB εt (t )
(4)
As presented in Eqs. (5) and (6), both kinematic variables are then used to directly estimate the macroscopic strain of a L-length sample, ε (t), by integration of the strain rate, ε˙ (t ).
ε˙ (t ) = ε (t ) =
V1 (t ) − V2 (t ) L
t 0
ε˙ (t )dt
(5)
T(°C)
100 – 200 s−1
200 – 500 s−1
>500 s−1
1 0 0
25 750 10 0 0 750 10 0 0 750 10 0 0
5 5 5 3 3 3 3
5 3 3 3 3 3 3
0 1 2 2 2 2 2
1 1 0
(6)
F1 (t ) = SB EB (εi (t ) + εr (t ))
(7)
F2 (t ) = SB EB εt (t )
(8)
F1 (t ) + F2 (t ) 2S
Crystal orient.
1 1 1
According to Eq. (5), the smaller the sample length L is, the more the strain rate increases. However, according to a previous study (Iwamoto and Yokoyama, 2012), L must at least equal the sample diameter in order to limit transverse wave propagation effects (i.e. wave dispersion). Besides, in this configuration, the sample has quickly reached the equilibrium state and inertial effects can be neglected (i.e. uniform strain state has been assumed). In the same way as for strain calculation, macroscopic stress, σ (t), can be estimated from input and output loads, F1 (t) and F2 (t), on both sides of the sample as seen in Eqs. (7)–(9), where SB is the section of the bars and S that of the sample.
σ (t ) =
Table 2 Mean plastic strain rate ranges obtained from SHPB tests for various conditions.
(9)
4.2. Experimental setup A modified SHPB device, developed at LaMCoS laboratory, has been employed to investigate the dynamic response of MC2 under uniaxial compression load at temperatures of 750 and 10 0 0 °C. Note that, for modeling simplification reasons, the behavior has been assumed to be symmetrical, which means no difference in tension and in compression. Despite a simple principle, the experimental implementation of such dynamic tests is relatively complex because of extreme temperature conditions (Lennon and Ramesh, 1998) for which the nickel-based superalloy yield strength stays high. The objective is to reach relatively high plastic strain (i.e. at least 1%) on aged MC2 samples without damaging (i.e. plastically
deforming) the bars. Because of its high compression yield strength (around 500 MPa) and a mechanical impedance close to that of MC2, high strength steel has been chosen as the constitutive material of the bars. It has been observed that, in the present test rig configuration, wave damping in steel bars is negligible during the compression phase. Precaution has also been taken on gauge bridge attachment for strain measurements. Actually, it has been found that the gauge strength remains the limiting factor for the initial velocity of the striker and, consequently, to the reachable strain rates. Thus, the striker length, which is proportional to load duration, has been determined to provide a sufficient initial kinetic energy to reach the MC2 compression yield strength without damaging the gauges. The specimen has rapidly been heated from R.T. up to 10 0 0 °C by electromagnetic local induction, see Fig. 3. Thermal insulating panels, made of superwhool®fibers, have been placed all around the sample to increase heat built-up and protect steel bars from heat conduction. Thermocouple sensors have allowed to control not only the assumed homogeneous sample temperature but also the temperature at the steel bar ends. A regulation system has been used to stabilize the temperature at 750 or 10 0 0 °C. The insulating panels have been removed just before triggering the impact. A large number of preliminary tests have been carried out to adjust the described test conditions. As shown in Table 2, dynamic compression tests have been performed through temperatures ranging from 25 to 10 0 0 °C for the three crystallographic orientations. Plastic strain rate levels obtained vary from 100 s−1 to 650 s−1 . Fig. 4 shows some 1 0 0-oriented MC2 samples after dynamic compression tests for which slip lines are clearly visible on skew external surface.
148
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
Fig. 3. SHPB system for high temperature dynamic tests at LaMCoS. Specimen diameter and length are 12 mm and 15 mm, respectively; Those of high strength steel bars are 20 mm and 2500 mm, respectively.
Fig. 4. 1 0 0-oriented MC2 samples after dynamic compression tests at 10 0 0 °C for which slip lines are clearly visible on skew external surface.
Fig. 5. Post-processing example: incident (blue), reflected (red) and transmitted (green) raw tension signals. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)
4.3. Signal post-processing and analysis A 1.25 MHz frequency acquisition card has been used to record strain signals from the gauges. Examples of incident, reflected and transmitted raw signals, measured for a 1 0 0-oriented MC2 sample, are presented in Fig. 5. Data has been filtered by computing discrete Fourier transformation (i.e. FFT algorithm) and by implementing a cutoff frequency. In order to describe the transient effects at the sample location, the three required signals have been superimposed thanks to delay calculations, t1 and t2 , based on wave celerity into steel bars at R.T., CB , and gauge/sample distances which are L1 = 1330 mm and L2 = 600 mm.
CB ≈ 5100m/s
(10)
t1 =
2L1 ≈ 520μs CB
(11)
t2 =
L1 + L2 ≈ 380μs CB
(12)
From these signal operations, loads and velocities at each interface (i.e. on both sides of the sample called input and out-
put interfaces) have been calculated in accordance with Eqs. (3)– (8), as illustrated in Fig. 6. As the input/output forces are quasisimilar, it is assumed that the specimen quickly reaches the equilibrium. The compression load is maintained for about 300μs (i.e. consistent plateau). The initial peak of the input velocity means that elastic and transient effects are mixed at the very beginning of the loading (i.e. short period lasting less than 60 μ s). To clearly identify the elastic material response, a coupled experiment/simulation approach could be chosen (Gary, 2014). However, as the purpose of the current study is to characterize the viscoplastic behavior of MC2 (which mainly determines the response of the structure under impacts), these complex analyses are not useful to reach this goal. Therefore, they have not been performed. This explains why, in the following results, material response is presented in terms of compression stress and plasticity levels. The evolution of the macroscopic strain rate is calculated from Eqs. (5) and (6) and presented in Fig. 7. Eq. (13) is then used to approximate the macroscopic plastic strain by trapezoidal numerical integration of the strain rate during the identified equilibrium phase, starting at time t0 and on a period of time T. Because of the high compression yield strength of MC2 superalloy, a significant part of the initial energy has been consumed to reach a
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
149
Fig. 6. On the left-hand side: filtered and superimposed signals. On the right-hand side: evolution of estimated input/output forces and velocity.
Fig. 7. Evolution of strain rate during the dynamic compression test and identification of the effective plastic strain rate.
sufficient stress level to plastically deform the specimen. This explains why the useful period T seems to be relatively short. Regardless, the obtained residual deformation levels are significant enough (i.e. between 2 and 5%) to identity the work hardening rate.
ε¯˙ p =
1 T
t0 +T
t0
V1 (t ) − V2 (t ) dt L
Fig. 8. Effect of temperature on compression yield strength and hardening for 1 0 0-oriented samples (engineering data). Results are normalized by the 1 0 0oriented MC2 compression yield strength at 10 0 0 °C, σ 0 .
This step enables to improve readability and to provide an average behavior.
(13)
Finally, by calculating the compression stress with Eq. (9), the assumed homogeneous strain-stress state in the specimen has been estimated so as to model the macroscopic material behavior. 5. SHPB compression tests: dynamic behavior results Experimental engineering stress-strain curves have been analyzed in terms of compression yield strength and work hardening rate. The obtained macroscopic plastic strain ranges up to 5%. Only tests for which it exceeds 1% have been plotted. In what follows, σ 0 is the compression yield strength of 1 0 0-oriented MC2 at 10 0 0 °C, taken here as a reference. It should be noted that some data on the yield strength of nickel-based superalloys at high temperature can be found in the literature (Caron et al., 2011; Segersäll and Moverare, 2013). Because of significant result dispersion inherent to dynamic testing, linear functions have been fitted for each group of the stress-strain curves derived from SHPB experiments.
5.1. Effect of the temperature In quasi-static conditions, it is widely recognized that, for nickel-based superalloys, the yield strength increases with temperature reaching a peak stress temperature around 800 °C (Starenchenko et al., 2008; Simonetti and Caron, 1998). However, at high strain rates (i.e. ε˙ p > 100 s−1 ), Fig. 8 shows a significant decrease of the yield strength by about 250 MPa, as the temperature increases from R.T. to 750 °C. This observation may be explained by the movement of dislocations which is strain-rate dependent (Kear and Wilsdorf, 1962). Besides, an adiabatic heating phenomenon due to high plastic strain rate has been pointed out to explain the low work hardening rate at room temperature (2.2 GPa) compared with the one obtained at 750 and 10 0 0 °C (respectively 6.6 GPa and 3.6 GPa) (Molinari et al., 2002). In other words, thermal softening seems to be preponderant over work hardening at ambient temperature.
150
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
Fig. 10. MC2 microstructure view after impact at 10 0 0 °C obtained with SEM. γ precipitate size is less than 1μ m. View in the plan normal to 1 0 0 direction where two distinct slip lines are visible (precipitate cutting). Fig. 9. Effect of plastic strain rate on hardening at 10 0 0 °C (engineering data). Results are normalized by the 1 0 0-oriented MC2 compression yield strength at 10 0 0 °C, σ 0 .
5.2. Effect of the strain rate Fig. 9 shows the material response to different plastic strain rates. Again, it is assumed here that work hardening due to dislocation motion is linear. The work hardening rate rises from 2.5 GPa to 9.5 GPa when the strain rate increases from 153 s−1 to 610 s−1 , whereas the compression yield strength remains almost unchanged. This illustrates why the viscoplastic flow of MC2 should be taken into account for accurate modeling of impacts. 5.3. Effect of the crystal orientation For FCC single crystals, the Schmid factor for the 1 1 1 orientation is lower than the one for 1 0 0 and 1 1 0 orientations (Ramaglia, 2013; Schmid, 1924). It means that, since the shear stress is the one needed to move dislocations on the {111} plans, the yield strength for 1 1 1 oriented MC2 is higher. In accordance with the Schmid law predictions, micrographs (see Fig. 10) show sheared or cut precipitates located on the octahedral plans, reflecting a non-homogeneous deformation mechanism. Besides, a significant increase of the compression yield strength, in the order of 200 MPa, is observed for 1 1 1 oriented single crystal, see Fig. 11. This important variation must be considered for blade-off simulations because of the complex geometry of the blade and the multi-axial loading. 6. Constitutive modeling of anisotropic viscoplastic behavior 6.1. Material modeling strategy Single crystal anisotropy is usually explained by the Schmid factors (Sabnis and Mazière, 2012). On the one hand, multi-scale finite element analyses based on Schmid law, which describe several active glide systems, are ruled out for complex dynamic analyses such as blade-off (i.e. large deformation/displacement, contact, high plasticity levels) because of their expensive computational cost (Méric et al., 1991; Sabnis and Forest, 2013). On the other hand, the common isotropic von Mises yield criterion is not able to faithfully transpose yield strength variations with single
Fig. 11. Effect of crystal orientation on compression yield strength and hardening at 10 0 0 °C (Engineering data). Results are normalized by the 1 0 0-oriented MC2 compression yield strength at 10 0 0 °C, σ 0 .
crystal orientation (Han et al., 20 01a; 20 01b). Actually, a coupled orthotropic elastic and viscoplastic behavior is expected for MC2 modeling. Several interesting material cards available in LS-Dyna for anisotropic modeling are presented in Table 3 (LSTC, 2015b). Recently compatible with eight node solid elements (i.e. bricks) which are employed to mesh turbine blades, “anisotropic elastic plastic” material model (i.e. model 157 in LS-Dyna) has been selected. It offers not only anisotropic elasticity and strain rate sensitivity for anisotropic plasticity (feature also recently added) but also predicts realistic behavior for finite displacement/rotation and large strain. This material model, which is a combination of LSDyna material models 2 and 103, has already been successfully used to model composite material behavior for automotive applications (Hatt, 2014).
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
151
Fig. 12. Position of material coordinate systems attached to solid element in the global Cartesian coordinate system.
Table 3 Potential LS-Dyna material cards to model anisotropic behavior. Card #
Description
Restriction
2
Orthotropic Elastic Barlat Anisotropic Plasticity Anisotropic Viscoplastic Orthotropic Viscoelastic Orthotropic Elastic Plastic Anisotropic Elastic Plastic
No plasticity
33 103 86 108 157
hence,
f (σ ) =
(σ11 − σ22 )2 + (σ22 − σ33 )2 + (σ33 − σ11 )2 2 2 2 2 + 2L(σ12 + σ23 + σ31 ) − σ02
(19)
On a usual basis, the anisotropic parameters are determined by measuring three uni-axial yield strengths and three shear yield stresses. However, compression tests on 1 1 1 oriented MC2 provide identification data to calibrate the L parameter. Indeed, the transformation of stress tensor from the global coordinate system ,Y ,Z ) to a local (i.e. material) one B1 = (X1 , Y1 , Z1 ) can be B = (X written as follows:
Isotropic elasticity Isotropic plasticity Only for shell element Local Hill criterion
σB 1 = T T σB T
(20)
Uni-axial compression condition test leads to the following stress tensor:
6.2. Material behavior approximation First, this particular material model uses generalized Hooke’s law, reduced here to the orthotropic case, in order to estimate elastic strains as written below:
ε = e ij
with
Ci−1 jkl
σkl
(14)
εi j = εiej + εipj
(15)
where σ ij is the Cauchy stress tensor, ε kl is the strain tensor which can be decomposed as a reversible elastic part εiej and irreversible
εipj .
plastic part Cijkl is the stiffness tensor, defined here by three independent components: C11 , C12 and C44 . These three elastic parameters, being necessary and sufficient to describe orthotropic behavior (i.e. cubic elasticity), have been identified for a wide range of temperature in a previous study (Gaubert, 2013). Note that this quasi-static elastic model is useful in blade-off simulations as a stress initialization of the airfoils under centrifugal loading is typically performed before impact simulation. Then, a relative (i.e. local) Hill criterion, corresponding to a modified von Mises surface, has been employed to describe an anisotropic yield locus. Hill’s yield criterion is generally expressed as follows (Hill, 1948):
f (σ ) = F (σ11 − σ22 )2 + G(σ22 − σ33 )2 + H (σ33 − σ11 )2 2 2 2 + 2Lσ12 + 2Mσ23 + 2Nσ31
0 . 5
− σ0
(16)
where F to N are the anisotropic constants and directions 1, 2, 3 are aligned with the material axes X1 , Y1 , Z1 . The local material coordinate systems corresponding to the typical single crystal orientations are shown in Fig. 12. Thanks to the material cubic symmetry, the Hill’s orthotropic yield function can be simplified:
F =G=H= L=M=N
1 2
σ 2 0
σ11
= 0.5
σB =
σ11
0 0 0
0 0
0 0 0
(21)
For the 1 1 1 orientation, the 3D transformation matrix for the base change from B to B1 is the following:
√ 1 3 1 T = 3 1 hence,
σB 1 =
σ11 3
−1 1 1 1 −1 1
1 −1 1 −1 1 −1
(22)
1 −1 1
(23)
Stress tensor can be decomposed into the sum of deviatoric and hydrostatic tensors (diagonal matrix). In 1 1 1 material base, only shear stresses are considered by the orthotropic Hill yield criteria (Eq. (19)). Thus, it is possible to identify the L parameter thanks to 1 1 1 compression test results:
3 L= 2
σ0
σy 1 1 1
2
≈ 0.9
(24)
Material model 157 also enable to consider the viscoplastic behavior thanks to tabulated curves, f (ε p , ε˙p ). These curves have been defined by fitting and extrapolating the average experimental rate-dependent curves (see Fig. 9), assuming validity of such extrapolation to larger plastic strain levels and strain rates. First, work hardening is extrapolated to larger strains (i.e. up to 100%) by using the typical Ludwik’s equation for which the two parameters, K and n, are calibrated on the experimental quasi-static response.
(17)
σ = σ0 + K ε np
(18)
with
K = 500 MPa
(25)
et
n = 0.6
(26)
152
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
Fig. 14. Equivalent von Mises stress field obtained by tensile test simulations with different material coordinates corresponding to crystallographic orientations 1 0 0, 1 1 0 and 1 1 1. Fig. 13. Extrapolated viscoplastic model for MC2 impact simulations. Dotted lines are deriving from the experimental data. Table 4 Parameters of LS-Dyna material model 157. Local coord. system
Elastic behavior
Yield surface
X1 B1 = ⎝Y1 ⎠ Z1
C11 C12 C44
σ0
⎛ ⎞
F = 0.5 L ≈ 0.9
Work Hardening f (ε p , ε˙ p )
Secondly, in order to describe the work hardening increase under high strain rates, but with no effect on the initial yield strength σ 0 , the factors C and ε˙ 0 (the latter being a reference effective strain rate) employed in the Johnson-Cook rate-dependent factor g(ε˙ p )(Johnson and Cook, 1983), are calibrated from the experimental response of 1 0 0-oriented MC2 under dynamic loading. Thus, a logarithmic extrapolation is performed to reach higher strain rate magnitudes (i.e. up to 104 s−1 ).
g(ε˙ p ) = 1 + Cln
with
C = 2.2
ε˙ p ε˙ 0
(27)
et
ε˙ 0 = 153s−1
(28)
Fig. 13 presents the extrapolated viscoplastic model which has been used to describe the dynamic behavior of MC2 single crystal. The parameters of the current material model are also summarized in Table 4. From a practical point of view, the table referring to f (ε p , ε˙ p ), which consists in true stress versus effective plastic strain curves indexed by strain rate, is used in LS-Dyna. Intermediate values are found by linear interpolation between curves. 6.3. Verification A “dog-bone” geometry has been selected for performing tensile tests, numerically. Several transient analyses have been performed for each crystalline orientation and different strain rates (see Fig. 14). As shown in Fig. 15, the current modeling faithfully reproduces the effect of the strain rate on work hardening. A significant increase of the compression yield strength for the 1 1 1 orientation is also well predicted thanks to the Hill criterion calibration. However, for the 1 1 0 orientation, the simulation gives a stronger material response, over 150 MPa higher than σ 0 , instead of 50 MPa.
Fig. 15. Results of tensile test simulation with material model 157 showing the effect of strain rate and crystallographic orientation on the macroscopic response.
Actually, the yield strengths predicted with a quadratic criterion such as the Hill’s one cannot strictly fit that of Schmid. A comparative study has been carried out in order to estimate the error for two orientations, 1 0 0 and 1 1 0. As shown in Fig. 16 an error of up to 20% is observed for the most unfavorable case (i.e. tension/compression load along 1 1 0) when the yield strength is predicted by the local Hill criterion. However, the prediction error decreases for multi-axial loadings. Keeping in mind that the present objective is to perform blade-off simulation of a turbine stage within a turbine engine, with a good balance between complexity of the model (and identification of the material data) and computing cost, it has been decided to keep the Hill criterion to account for anisotropic plasticity. 7. Application to a normal impact on a plate By way of validation, the current model has been applied to a ballistic test for which a 1 0 0-oriented single crystal plate is impacted by a flat head cylindrical bullet made of 25CD4 steel. A
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
153
Fig. 16. Theoretical tension-shear elastic domain for 1 0 0 and 1 1 0 orientations. Local Hill function with L ≈ 0.9 is compared with isotropic von Mises yield surface and Schmid predictions for which only octahedral slips are considered (Ding et al., 2006).
Fig. 17. Numerical implementation of impact case on a single-crystal plate.
gas gun installation has been used to perform this dynamic test. Tests have been performed with plates heated to 10 0 0 °C point-set and an impact velocity of 100 m/s. With the aim of measuring residual strains, speckle patterns have been spray-painted on plate surfaces with a high-temperature resistant paint (commercial name is Ulfalux ®). Then, images before and after impact have been taken with two 4MPx cameras (i.e. stereo system). This experimental procedure has allowed to use digital images correlation (DIC)(Schreier et al., 2009), in a postmortem way, to estimate local residual strain on the impacted structure. On the one hand, the commercial DIC package, VIC-3D v7 developed by Correlated Solutions, has been used for the experimental part of the study. To ensure relevant correlation results, the speckel size, also called Subset size, has been set to 30 pixels and the virtual strain gauge size, also called step size, set to 5 pixels. Lagrange interpolation function has been used in the deformation analysis. Thus, the equivalent strain, in the meaning of von Mises using the principal plane strain formulation (i.e. only ε xx , ε yy , ε xy are not nill), is calculated in the global coordinate system as follows:
εI,II = εeq =
εxx + εyy
2
±
εxy 2 +
εI 2 − εI εII + εII 2
ε − ε 2 xx yy 2
(29) (30)
On the other hand, the numerical model for impact on MC2 clamped plate has been implemented with LS-Prepost (Fig. 17). The specimen is meshed with solid hexahedral elements with at least three elements in the thickness to retrieve good bending behavior.
Only one quarter of the model is represented and symmetrical displacement conditions are applied. The material model 157 is used with a 1 0 0 oriented coordinate system attached to each element. In order to reproduce the boundary conditions of the ballistic test, the preload is taken into account by a dynamic relaxation initialization (LSTC, 2015a). Besides, “automatic surface to surface” contact algorithm is employed to model the contact between the projectile and the plate. At the end of the simulation, a springback analysis is performed by mass damping so as to obtain the residual shape of the structure (Fig. 18). The LS-Dyna explicit solver is used to solve the following nonlinear equation of motion considering both contact and inertial effects:
MU¨ + Fint (U ) = Fext
(31)
for which M represents the lumped mass matrix of the system, U the nodal displacement vector, Fint the internal force vector, and Fext is the external force vector, including contact forces. The current simulation gives access to information that experiments cannot provide. Thus, in the vicinity of the impacted zone, the plastic strain rates are estimated around 20 0 0 s−1 and the triaxiality factor reaches ±0.66 (i.e. equi-biaxial traction/compression). As shown in Fig. 18, the plate macroscopic behavior is well predicted as measured and calculated profiles are quasi-similar, with an approximate deflection of 3 mm. In order to provide quantitative results, a comparison between experimental and simulated back face profiles, in the vicinity of the impact zone, is presented in Fig. 19. An additional test, named experiment n◦2, has been performed under the same experimental conditions as the first one (experiment n◦1) to provide more
154
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
Fig. 18. A comparison of post-mortem observation corresponding to experiment n◦1(left) and residual deflection obtained by FE simulation (right).
Fig. 19. Measured and simulated back face profiles after impact.
Fig. 21. A comparison of residual strain field after impact on back face (top: test, bottom: simulation).
Fig. 20. A comparison of residual strain field after impact on front face (top: test, bottom: simulation).
consistency. Experiment n◦1 shows a difference of 1.44% on the observed maximum deflection while experiment n◦2 presents a difference of 6.27%. According to those results, the current simulation of a normal plate impact is consistent with the experimental measures (difference of less than 10%). To supplement these observations, experimental/predicted residual strain fields are compared in Figs. 20 (front face) and 21 (back face). All around the impact zone, the calculated residual strains are close to the measured ones. At the back face, mostly subjected to tensile stress during the impact, the predicted strain at the center of the plate appears to be more intense. The maximum residual strain reaches 5% and 10% for the front and back faces, respectively. 8. Application to a bullet impact on a turbine blade We will now provide an example of application to show how the present data and methodology can be used to predict what
happens when a debris impacts a turbine blade. This example is representative of the results obtained for a series of impact tests performed at different temperatures and impact speeds, for several blade geometries and bullet mass. We have selected a test with a MC2 turbine blade heated at 10 0 0 °C and impacted by a cylindrical bullet at a speed in the range of a few hundreds of m/s. The impact has been recorded by a high speed camera. The deformation of the blade and the strain distribution on the blade surface after impact have been analyzed by DIC as in the previous section. Note that in this specific example the impact did not lead to fracture. An image of the blade at a specific time during the impact is presented in Fig. 22 (left) along with the plastic strain distribution derived from simulation (right). A comparison of the residual strain distribution is shown in Fig. 23. A quite reasonable agreement between experiment and simulation can be observed. Overall this example of results confirms the ability of the proposed methodology to tackle the mechanical behavior of a turbine blade under impact and, for a mono-crystal, the quality of the material data derived from the SHPB investigation. 9. Conclusion In order to provide experimental data for blade-off simulations, MC2 single crystal has been characterized under extreme conditions: high temperature (i.e. up to 10 0 0 °C), high plastic strain rates (i.e. >100 s−1 ) and significant plasticity levels (i.e. around 5%), reached despite of the high yield strength of MC2 superalloy. A typical SHPB installation has been modified to perform and analyze theses tests. The results presented in terms of stress-strain curves have highlighted the effects of temperature, material orientation and strain rate on the material response (i.e. compression yield strength and work hardening rate). An intermediate approach between an inaccurate material model based on the von Mises yield criterion (isotropic) and an
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
155
Fig. 22. Example of turbine blade during impact, photography (left) and simulation (right).
Fig. 23. Residual strain after impact, comparison between measurement by DIC (left) and simulation (right).
expensive crystalline plasticity model is proposed, providing industrial efficiency. Thus, it is now possible to performed blade-off simulations within a turbo-engine with the aims of a better understanding of such complex event and bring answers to airworthiness requirements. The main outcomes of the present analysis are the following: i) So as to faithfully predict the crystalline material behavior for impact applications, both anisotropic plasticity and viscoplasticity should be considered in the FE model. ii) The material model 157 available in LS-Dyna has permitted to reproduce the behavior of single crystal plate under impact with a fairly good agreement compared with ballistic experiments. This feature allows to consider first the elastic orthotropy inherent to FCC crystal, then the anisotropic plasticity by calibrating a local Hill criterion, and finally, the ratedependant work hardening. iii) This implemented efficient formulation can be used to predict macroscopic responses for industrial applications as complex as blade-off for which large displacements and deformations occur. That will be the purpose of further research. In order to provide experimental data for material modeling with strain rate dependency, a dynamic test campaign has been carried out on MC2 single crystal through a modified SHPB equipment designed for high-temperature tests up to 10 0 0 °C. Recent developments on LS-Dyna solver provide possibilities to implement an approximate material model for MC2 single crystal
under impact. Despite an error of 20% introduced for specific loading along the 1 1 0 direction, the predicted macroscopic behavior is acceptable and gives satisfactory results (i.e. predicted deflection and residual strains) in the case of high velocity impact on a monocrystalline MC2 plate. High-speed impact of a flat cylindrical bullet on a MC2 plate has been simulated with LS-Dyna and numerical results compared with those of ballistic tests. They showed a fairly good agreement. Acknowledgements The authors acknowledge the material department of Safran Helicopter Engines for technical support and providing material samples. The plate impact experiment presented here has been implemented in collaboration with Alain Veray, Safran Aircraft Engines, for whom the authors are grateful. The DIC measurements reported herein have been performed with the help of Philippe Chaudet, LaMCoS technical support, who is also acknowledged here. References Caron, P., Diologent, F., Drawin, S., 2011. Influence of chemistry on the tensile yield strength of nickel-based single crystal superalloys. Adv. Mater. Res. 278, 345–350. Caron, P., Lavigne, O., 2011. Recent studies at onera on superalloys for single crystal turbine blades. J. AerospaceLab 3, 1–14. Ding, Z.P., Liu, Y.L., Yin, Z.Y., Yang, Z.G., Cheng, X.M., 2006. Constitutive model for an fcc single-crystal material. Front. Mech. Eng. China 1, 40–47.
156
L. Roux et al. / Mechanics of Materials 104 (2017) 145–156
European Aviation Safety Agency, 2010. Certification Specifications for Engines, Amendment 3. European Aviation Safety Agency. European Aviation Safety Agency, 2012. General Acceptable Means of Compliance for Airworthiness of Products, Parts and Appliances. European Aviation Safety Agency. Gary, G., 2014. Testing with bars from dynamic to quasi-static. In: Constitutive Relations under Impact Loadings. Springer, pp. 1–58. Gaubert, A., 2013. Identification d’un modèle de viscoplasticité anisotrope pour le MC2 de 20° à 1250°. Technical Report. ONERA. Hallquist, J., 2005. LS-DYNA Theory Manual. Livermore software technology corporation (LSTC). Han, S., Li, S., Smith, D.J., 2001a. Comparison of phenomenological and crystallographic models for single crystal nickel base superalloys. i. analytical identification. Mech. Mater. 33, 251–266. Han, S., Li, S., Smith, D.J., 2001b. Comparison of phenomenological and crystallographic models for single crystal nickel base superalloys. ii. numerical simulations. Mech. Mater. 33, 251–266. Hatt, A., 2014. Anisotropic modeling of short fibers reinforced thermoplastics materials with ls-dyna. LS-DYNA Forum. Hill, R., 1948. A theory of the yielding and plastic flow of anisotropic metals. Proc. R. Soc. London A, 193–281. Hopkinson, B., 1914. A method of measuring the pressure in the deformation of high explosives by the impact of bullets. Philos. Trans. R. Soc. A213, 437–452. Iwamoto, T., Yokoyama, T., 2012. Effects of radial inertia and end friction in specimen geometry in split Hopkinson pressure bar tests: a computational study. Mech. Mater. 51, 97–109. Johnson, G.R., Cook, W.H., 1983. A constitutive model and data for metals subjected to large strains, high strain rates and high temperatures. In: Proceedings of the 7th International Symposium on Ballistics. Kear, B.H., Wilsdorf, G.F., 1962. Dislocation configurations in plastically deformed polycrystalline cu3au alloys. Trans. Metall. Soc.AIME 224, 382–386. Leidermark, D., Moverare, J.J., Simonsson, K., Sjostrom, S., Johansson, S., 2009. Room temperature yield behaviour of a single-crystal nickel-base superalloy with tension/compression asymmetry. Comput. Mater. Sci. 47, 366–372. Lemaitre, J., Chaboche, J.L., 2009. Mécanique des matériaux solides. Dunod. Lennon, A.M., Ramesh, K.T., 1998. A technique for measuring the dynamic behavior of materials at high temperatures. Int. J. Plast. 14, 1279–1292. Liang, R., Khan, A.S., 1999. A critical review of experimental results and constitutive models for bcc and fcc metals over a wide range of strain rates and temperatures. Int. J. Plast. 15, 963–980. LSTC, 2015a. LS-DYNA Keyword User’s Manual Volume I (Revision: 4686). Livermore software technology corporation (LSTC). LSTC, 2015b. LS-DYNA Keyword User’s Manual Volume II (Revision: 4686). Livermore software technology corporation (LSTC).
Méric, L., Poubanne, P., Cailletaud, G., 1991. Single crystal modeling for structural calculations: part 1 - model presentation. J. Eng. Mater. Technol. 113 (1), 162–170. Molinari, A., Musquar, C., Sutter, G., 2002. Adiabatic shear banding in high speed machining of ti-6al-4v: experiments and modeling. Int. J. Plast. 18, 443–459. Ogawa, K., 1985. Mechanical behavior of metals under tension compression loading at high strain rate. Int. J. Plast. I, 347–358. Pollock, T.M., Tin, S., 2006. Nickel-based superalloys for advanced turbine engines: chemistry, microstructure, and properties. J. Propul. Power 22, 2. Ramaglia, A.D., 2013. Application of a smooth approximation of the Schmid’s law to a single crystal gas turbine blade. J. Eng. Gas Turbines Power volume 135 Number 3 (0055), 9. Reed, R.C., Tao, T., Warnken, N., 2009. Alloys-by-design: application to nickel-based single crystal superalloys. Acta Mater. 57, 5898–5913. Sabnis, P.A., Forest, S., 2013. Crystal plasticity analysis of cylindrical indentation on a ni-base single crystal superalloy. Int. J. Plast. 51, 200–271. Sabnis, P.A., Mazière, M., 2012. Effect of secondary orientation on notch-tip plasticity in superalloy single crystals. Int. J. Plast. 28, 102–123. Schmid, E., 1924. Yield point of crystals, critical shear stress law. Proc. Int. Congr. Appl.Mech. 342. (Delft). Schreier, H., Orteu, J.J., Sutton, M.A., 2009. Image Correlation for Shape, Motion and Deformation Measurements. Springer. Segersall, M., 2013. Nickel-Based Single-Crystal Superalloys. Linköping Studies in Science and Technology Master’s thesis. Segersäll, M., Moverare, J.J., 2013. Crystallographic orientation influence on the serrated yielding behavior of a single-crystal superalloy. Materials 6, 437–444. Simonetti, M., Caron, P., 1998. Role and behaviour of m phase during deformation of a nickel-based single crystal superalloy. Mater. Sci. Eng. A254, 1–12. Starenchenko, V.A., Kozlov, E.V., Solov’eva, Y.V., Abzaev, Y.A., Koneva, N.A., 2008. Orientation dependence of the yield stress and work-hardening rate of ni3ge at different temperatures. Mater. Sci. Eng. 483–484, 602–606. Wang, Z.Q., Beyerlein, I.J., Lesar, R., 2009. Plastic anisotropy in fcc single crystals in high rate deformation. Int. J. Plast. 25, 26–48. Zhao, H., Gary, G., 1996. On the use of shpb technique to determine the dynamic behaviour of the materials in the range of small strains. Int. J. Solids Struct. 33, 3363–3375. Zhao, H., Gary, G., 1997. A new method for the separation of waves. application to the shpb technique for an unlimited duration of measurement. J. Mech. Phys. Solids 7, 1185–1202. Zhao, L.G., Tong, J., Vermeulen, B., Byrne, J., 2001. On the uniaxial mechanical behaviour of an advanced nickel base superalloy at high temperature. Mech. Mater. 33, 593–600.