Heterogeneous deformation in a commercially pure titanium sheet under dwell fatigue loading: Crystal plasticity modeling and experiment

ISIJ International, Vol. 61 (2021), No. 6 © 2021 ISIJ 1990 * Corresponding author: E-mail: [email protected] © 2021 The Iron and Steel Institute of Japan. This is an open access article under the terms of the Creative Commons Attribution-NonCommercial-NoDerivs license (https://creativecommons.org/licenses/by-nc-nd/4.0/). ISIJ International, Vol. 61 (2021), No. 6, pp. 1990–2001 https://doi.org/10.2355/isijinternational.ISIJINT-2020-702 1. Introduction Titanium and titanium alloys have been the mainstays of structural materials in aerospace and marine applications owing to their reduced density, excellent corrosion resistance and high fracture toughness. It has been reported that titanium alloy-manufactured compressor discs in high bypass turbofan engines suffer a significant life reduction (known as ‘dwell life debit’) under cyclic loading with a stress-dwell period in each cycle, compared with an equivalent low-cycle fatigue.1) Similar phenomenon has also been observed in commercially pure titanium (CP-Ti) at room temperature.2) The anisotropic elastoplastic response of hexagonal close-packed (HCP) α -titanium, which can cause strain incompatibility and trigger the generation of discontinuous basal facets within subsurface crack initiation sites, Heterogeneous Deformation in a Commercially Pure Titanium Sheet under Dwell Fatigue Loading: Crystal Plasticity Modeling and Experiment Liangwei YIN1)* and Osamu UMEZAWA2,3) 1) Graduate School of Engineering Science, Yokohama National University, 79-5 Tokiwadai, Hodogaya, Yokohama, 240-8501 Japan. 2) Faculty of Engineering, Yokohama National University, 79-5 Tokiwadai, Hodogaya, Yokohama, 240-8501 Japan. 3) Center of Advanced Innovation Technologies, Vysoká Škola Bánˇská - Technical University of Ostrava, 17. listopadu 15, 708 33 Ostrava-Poruba, Czech Republic. (Received on November 30, 2020; accepted on February 9, 2021; J-STAGE Advance published date: March 19, 2021) The heterogeneous deformation in a hot-rolled commercially pure titanium grade 1 sheet has been experimentally and numerically investigated under dwell fatigue loading in current paper. The residual strain fields within two regions of interest are probed by digital image correlation (DIC) after interrupted dwell fatigue test. These measurements essentially agree with predictions of a dislocation mechanismbased crystal plasticity model incorporating deformation twinning. The simulated results further indicate that axial strain localization at grain scale mainly derives from prismatic slip activity, followed by pyramidal <a>, basal slip and { 1122 } compression twinning activities. On this basis, a weighted averaged Schmid factor is proposed to correlate the axial strain accumulation with active plastic deformation modes in individual grains. Besides, the cyclic load shedding within a soft-hard-soft grain pair is captured by crystal plasticity modeling. The stress redistribution from soft to adjacent hard grain implies that the influence of deformation twinning on dwell facet nucleation is limited. The presented study confirms a robust crystal plasticity model and deepens the quantitative analysis of cold dwell fatigue in titanium and its alloys. KEY WORDS: dwell fatigue; titanium; crystal plasticity; digital image correlation; heterogeneous deformation. is recognized as a possible reason for the aforementioned premature fracture.3–5) Therefore, thoroughly understanding the deformation heterogeneity of CP-Ti under dwell fatigue loading is of essential importance. The heterogeneous deformation at grain level is closely associated with slip activities. For HCP-structured CP-Ti, there exists five active families of 30 slip systems in total, three ()0001 1120  basal, three {} 1010 1120  prismatic, six {} 1011 1120  pyramidal, twelve {} 1011 1123  1st pyramidal, and six {} 1122 1123  2nd pyramidal slip systems. A specific slip system can be activated when the resolved shear stress is higher than the critical resolved shear stress (CRSS). Many theoretical and empirical methods have been used to assess absolute CRSS values or CRSS ratios for CP-Ti.6–10) According to these works, there is little doubt that prismatic slip is the easiest to be activated,8,9) and the CRSS value for pyramidal < c+a> slip is approximately three times higher than prismatic slip.9,10) However, due to ISIJ International, Vol. 61 (2021), No. 6 © 2021 ISIJ1991 the diverse textures, grain sizes and chemical compositions (e.g. oxygen concentration) of tested materials, the CRSS ratios of basal to prismatic slip reported in previous studies scatter extensively from 1.1 to 20.0.6–11) Similarly, the CRSS ratios of pyramidal < a> to prismatic slip range from 1.2 to 12.0.6,7,11) Thus, the CRSSs for basal and pyramidal < a> slips still remain matters of debate. In addition, Zhang et al. found that the creep-induced deformation under cold dwell fatigue was intrinsically related to microscopic strain rate sensitivity (SRS).12) Thus, a careful identification of CRSS relationships and slip-dependent SRS is indispensable for understanding the dwell fatigue behaviors in CP-Ti. In contrast to slip-dominant plastic deformation in nearα and α / β titanium alloys, slip and twinning are competitive deformation mechanisms in CP-Ti with a low oxygen contents (mass% O < 0.11).13) Since the activation of highCRSS pyramidal < c+a> slip is difficult, {} 1012 tension twinning and {} 1122 compression twinning act as primary substitutes to achieve homogeneous deformation along the c-axis for hard grains (poorly oriented for slip).14) On the other hand, Luan et al. presumed that twinning nucleation and propagation in soft grains (well oriented for slip) could interrupt continuous slip, thereby reducing the dislocation pile-ups at grain boundaries and possibility of dwell facet generation in CP-Ti.15) Nevertheless, the influence of deformation twinning on dwell fatigue susceptibility of CP-Ti is still not properly understood. The conjunction of digital image correlation (DIC) and CPFE method is a commonly used approach to study the heterogeneous deformation and underlying plastic mechanisms at mesoscale.16–18) For instance, the in-situ DIC strain mapping of an oligocrystalline CP-Ti was compared to CPFE simulation results to discuss the influence of CRSS ratio on deformation behaviors during tension.16) Similar approach was extended to a coarse-grained orthorhombic α -uranium to determine the relationship between strain localization and slip activity under tension.17) Further, the micro slip traces that probed by high-resolution DIC (HRDIC) were also explicitly captured by the crystal plasticity modeling18) in a Ti-6Al-4V alloy under tension. However, most of the previous works focused on the response of oligocrystalline materials under monotonic loading conditions. Although DIC technique has been used for evaluating local deformation pattern and strain accumulation under various cyclic loading conditions,19) quantitative comparison between DIC measurements and CPFE model predictions has not been documented so far for polycrystalline titanium under dwell fatigue loading. The aim of the present study is to evaluate the heterogeneous deformation of a hot-rolled CP-Ti grade 1 sheet under dwell fatigue loading. Firstly, a series of mechanical tests and parallel simulations are conducted to calibrate the crystal plasticity parameters of CP-Ti. The residual strain fields within two regions of interest are measured by DIC after interrupted dwell fatigue test and compared to predictions of CPFE simulations. Then, the mechanisms driving axial strain localization are systematically investigated in terms of active plastic deformation modes and Schmid factors. Furthermore, the influence of deformation twinning on facet nucleation are determined based on the analysis of cyclic load shedding within a soft-hard-soft grain pair. 2. Experimental Procedure 2.1. Material The chemical compositions of tested JIS grade 1 CP-Ti are tabulated in Table 1. The material was hot-rolled into a 30-mm-thick sheet and annealed at 973 K for 1 h. The average grain size after heat treatment is 17 μ m in diameter. The flat dog-bone specimens with a thickness of 2 mm were machined parallel to the rolling direction (RD), transverse direction (TD) and 45° from RD to TD for tensile, stress relaxation and dwell fatigue tests. The configuration of specimens is shown in Fig. 1(a). Prior to mechanical tests, the specimens were first polished using SiC papers to #2400, electropolished at −40°C for 30 s in a solution of 440 ml methanol, 264 ml butanol, and 44 ml HClO4, and then etched by Kroll reagent for 60 s for microstructure characterization. 2.2. Mechanical Tests and Strain Measurement Mechanical tests were performed using a motor-driven tester (SHIMAZU AG-IS 20kNT). Tensile tests were conducted under three different strain rates (2.3 ×10−4 s−1, 2.3×10 −3 s−1, 2.3×10−2 s−1) to obtain stress-strain curves for calibrating rate-sensitive parameters. Since the tensile anisotropy of high-textured CP-Ti sheet was pointed out to be a reflection of CRSS differences, RD, TD, and 45° specimens were used to determine CRSS relationships.8) In addition, stress relaxation tests, which interrupted the tensile tests at 260 MPa by a 300 s stress relaxation period followed by unloading, were also carried out for three kinds of specimens to furtherly identify strain-hardening parameters. DIC was applied to examine the heterogeneous deformaTable 1. Chemical compositions of the grade 1 CP-Ti under investigation (mass%). O C N H Fe Ti 0.092 0.002 0.005 0.001 0.041 Bal. Fig. 1. (a) Geometry of specimens, (b) microstructure of the selected area on dwell fatigue test specimen, and (c) corresponding IPF map of RD. DIC measurements are performed in regions A and B. (Online version in color.) ISIJ International, Vol. 61 (2021), No. 6 © 2021 ISIJ 1992 tion in RD specimen under dwell fatigue loading. A large area was assigned at the center of gauge section in advance, as shown in Fig. 1(a). To combine the local strain field with initial microstructure, electron back-scattered diffraction (EBSD) measurement was conducted in this area using a JEOL JSM-6400F scanning electron microscope (SEM) equipped with a TSLTM EBSD detector. EBSD scan was performed at an accelerating voltage of 15 kV with a step size of 1.0 μ m. The data was analyzed by OIM analysis software and MTEX toolbox in MATLAB.20) The SEM image (resolution: 0.12 μ m/pixel) of the area and corresponding inverse pole figure (IPF) map of RD are shown in Figs. 1(b) and 1(c). Ex-situ DIC analysis were performed in two regions of interest (ROIs) surrounded by blue dashed lines (see Fig. 1(c)). The fabrication of speckle pattern for DIC strain mapping was referred to Ref. 21). First, 30 mg Mo powder (2.56 μ m in diameter) and 200 ml methanol were mixed and the mixture was then placed into an ultrasonic cleaner for 30 min to make a homogeneous suspension. Then ten thousand Pb balls with a diameter of 1 mm were added to shrink the Mo agglomerates. The suspension was rested for 24 h and applied to the surface of specimen using a dropper. Finally, the specimen was left for 6 h to achieve a rigid attachment of pattern. The patterns of two ROIs before and after dwell fatigue test were analyzed by an open-source MATLAB program Ncorr to provide residual strain fields22) with subset radius of 26 and subset spacing of 2. Dwell fatigue load was applied in a trapezoidal waveform with a loading/unloading time of 1 s. The stress ratio R was set to 0.01. In each cycle, the load was held at 293.55 MPa (0.95 σ y at strain rate 2.3×10 −4 s−1) for 5 s.4) It was reported that strain field established at the early stage of the dwell fatigue test,19) and strain accumulation was most pronounced at the beginning of the loading history.4) Thus, the dwell fatigue test was interrupted after 20 cycles. Meanwhile, the computation time of accompanying CPFE simulation could be reduced. 3. Crystal Plasticity Finite Simulation Framework 3.1. Crystal Plasticity Model The rate-dependent crystal plasticity framework used in present study was originally developed by Huang based on finite strain theory.23) The deformation gradient F is given by FFF ep  ................................ (1) where Fe refers to the lattice stretching and rigid body rotation, and Fp represents the deformation caused by plastic shear of material. The above formula leads to additive decomposition of velocity gradient L into elastic and plastic parts, LFFFFFFF F LFLF ee ep pe eepe               1 11 1 1 ... (2) Assuming that the plastic deformation results from all active slip and twinning systems, Lm s p          1 n ...................... (3) where n is the total number of slip and twinning systems,    is the shearing rate of the α th slip/twinning system, m α and s α are slip/twinning plane normal and slip/twinning direction, respectively. The velocity gradient L can also be written as LDWD DDWW W ep ep    ,, .......... (4) where the symmetric rate of stretching D and the antisymmetric spin tensor W can be decomposed into elastic (De, We) and plastic (Dp, Wp) parts, respectively. The plastic strain rate tensor Dp and plastic spin tensor Wp represent the symmetric and skew-symmetric parts of the plastic velocity gradient: DLL ppms sm p pp       Tn 2 2 1        ,/ .......................................... (5) WLL ms sm p pp       Tn 2 2 1       ,/ .......................................... (6) The relation between the symmetric rate of lattice stretching, De, and the corotational stress rate with respect to lattice rotation, ˆ σσ e, is given by ˆˆ :  ep pe WW ED   ............... (7) where E represents the elastic modulus tensor. The constitutive equation of a single crystal can therefore be deduced by combining Eqs. (4)–(7): ˆ::         ED Ep      1 n..... (8) Five active slip system families, {} 1012 tension twinning and {} 1122 compression twinning systems are taken into account in current model. The shearing rate    on a slip system α follows a thermally activated model considering the forward and backward dislocation un-trapped events:24)                      0 2 if c mD c bexp F kT sinh V kT               if   c .......................................... (9) where ρ m is the mobile dislocation density, b α is the Burgers vector, υ D is the Debye frequency, k is the Boltzmann constant, T is the temperature. ΔF and ΔV are the activation energy and activation volume, which significantly influence the SRS of slip systems. τ α is the resolved shear stress and   c is the effective CRSS. The {} 1012 tension and {} 1122 compression twinning systems are modeled as pseudo slips and allowed to be activated when resolved shear stresses are positive and negative, respectively. The shearing rate on twinning system obeys the following viscoplastic formula:              0 1 SSD r sgn.................... (10) where  γ 0 is the reference slipping rate and r is the rate sensitivity index,   SSD is the slip strength induced by statistically stored dislocations (SSDs). The Voce-type hardening model is applied to describe the evolution of   SSD:25) ISIJ International, Vol. 61 (2021), No. 6 © 2021 ISIJ1993     SSD nh    1........................ (11) where h αβ is the hardening moduli, which stands for the interactions among slip and twinning systems. The selfhardening moduli h αα and latent-hardening moduli h αβ are expressed as:25) hh hh dt hqh t                         0 00 0 1exp , ,     1 n .......... (12) where h0, τ ∞, τ 0 are initial hardening modulus, saturated CRSS and the initial CRSS. γ is the cumulative shear strain on all slip and twinning systems, q is the latent-hardening constant. Based on the strain-gradient crystal plasticity framework developed by Han et al.,26) the effective CRSS   c of a slip system α takes the contribution of SSDs and geometrically necessary dislocations (GNDs) into consideration:27)     c SSD GND  () () 22 ..................... (13) where   GND is the slip strength caused by GNDs. To this end, the intrinsic length scale of a material, l, is introduced to account for GNDs-induced hardening: l Gb g T    22 0 2............................... (14) where α T is an empirical coefficient that ranges from 0.3 to 0.5, G is the shear modulus, g0 represents the reference slip resistance, which is taken to be g0 = G/100. The slip strength due to the accumulation of GNDs is then computed based on Taylor dislocation model,28)    GNDGND TG ND gl Gb 0 ............. (15) where   GN D is the effective density of GNDs, which has the following form at small deformations:    GND n             xms ms 1... (16) where m α and s α are the slip plane normal vector and the slip direction, respectively. The calculation of effective density of GNDs requires the evaluation of the shear strain gradient at each element integration point, the computational execution is summarized in Ref. 29). Using Eqs. (13) and (15), the effective CRSS of a slip system α is rewritten as:     c SSD TG ND Gb()() 22 ............... (17) Twinning model is based on the conventional predominant twin reorientation (PTR) scheme.30) The volume fraction V of a twinning system is defined as: Vtwin ref    ................................. (18) where γ twin is the cumulative shear strain on a twinning system, γ ref is the reference twinning shear, which can be calculated by31)  refref a c c a c a c a ,{ },{} , 1012 1122 2 3 3 22 3                ... (19) For CP-Ti, c/a≈1.587. Thus, γ ref ,{ }1012 = 0.175, γ ref ,{ }1122 = 0.218. When the volume fraction V is greater than the threshold volume fraction Vth, lattice reorientation occurs. The orientations of all slip and twinning systems at each integration point are rotated from matrix to twinned region by the following rotation tensor Rtw: Rm mI tw tw tw 2 ...................... (20) where mtw is the normal vector of the twinning plane, I is the identity tensor. The threshold value of Vth for each twinning system is randomly generated from 0.1 to 1.0.8) The code of the crystal plasticity model is scripted in the UMAT and URDFIL subroutines incorporated within ABAQUS. 3.2. Polycrystalline Aggregate Modeling A cuboidal representative volume element (RVE) comprising 512 grains is applied to produce numerical stressstrain and stress-time curves, as shown in Fig. 2. The dimension of RVE is 100 × 100 × 100 μ m. Each grain is meshed with 1 C3D8 brick element. The x, y and z directions represent TD, RD and normal direction (ND), respectively. The boundary conditions are assigned as uX′=0, uY′=0, and uZ′=0. The uniform displacements are imposed along RD, 45°, and ND, respectively. The initial crystallographic orientations are extracted from EBSD measurements. Comparison between the (0001) pole figures of EBSD data and RVE model is illustrated in Fig. 3. For purpose of comparison with DIC strain measurements, two quasi-3D polycrystalline aggregate models are reconstructed within the abovementioned ROIs (see Fig. Fig. 2. RVE model containing 512 grains and boundary conditions. (Online version in color.) ISIJ International, Vol. 61 (2021), No. 6 © 2021 ISIJ 1994 1(c)) via an image-based approach, as shown in Fig. 4. The sketches of 2D microstructures are depicted in AutoCAD to obtain smooth grain boundaries (as illustrated in Fig. 4(a), including the number of particular grains for later discussion), extruded along z-axis and meshed by 8-node C3D8 brick elements (19 872 elements in region A and 16 542 elements in region B) in Hypermesh, as illustrated in Fig. 4(b). Two models have the same dimension of 80 × 80 × 8 μ m. Given the influence of subsurface microstructure on the free-surface deformation, each model contains a surface layer and a subsurface layer with the same thickness, mesh size and grain morphology (green area in Fig. 4(b)).16,18) Therefore, the influence of subsurface microstructure is focused on crystallography. The crystal orientations of grains on the surface layer are specified by corresponding Bunge Euler angles, whilst the orientation of grains on the subsurface layer are randomly selected from the EBSD data. The dwell fatigue load is applied on the right surface along RD and the waveform is identical to experiment. The freedoms of left, bottom and back surfaces are fixed to eliminate rigid body movement. 4. CPFE Constitutive Parameter Calibration To numerically reproduce the experimental results, CPFE constitutive parameters of tested CP-Ti should be identified. With respect to the data related to SRS, the mobile dislocation density ρ m, Burgers vector b α , and Debye frequency υ D in Eq. (9) are referred to Ref. 32), whereas the activation energy ΔF, activation volume ΔV and the rate-sensitive parameters of deformation twinning in Eq. (10) are determined via calibration. Recently, Xiong et al. have found that the SRS of prismatic slip is higher than that of basal slip in a grade 4 CP-Ti and the activation energy ΔF of basal slip is also greater than that of prismatic slip, similarly to Ti-6242 alloy.33,34) Although slip-dependent SRS variations are considered to have influence on creep deformation in cold dwell fatigue, the concrete mechanisms are still unknown. Besides, the diversity in chemical compositions should not be ignored. Thus, the activation energy and activation volume of each slip system is assumed to be the same in this study. Determination of the CRSS relationships is conducted simultaneously. The CRSSs of 1st and 2nd pyramidal <c+a> slips are taken to be the same, three times the CRSS of prismatic slip.10) Calibration procedure for the CRSS ratio between basal and prismatic slips is based on previous works. Since increased concentration of oxygen can improve the strength of CP-Ti for the solid-solution strengthening effect, and the strengthening effect of oxygen is reported to be more remarkable on prismatic slip than basal slip, the CRSS ratio between basal and prismatic slips is therefore decreases with increased amounts of oxygen.35,36) Moreover, the macroscopic strength of tested grade 1 CP-Ti in the present study is between the counterparts reported by Hama et al.8) and Xiong et al.33) In terms of these relationships, the initial CRSSbasal/CRSSprismatic in this study is optimized between 1.63 and 2.14, values reported by the abovemenFig. 3. (0001) pole figures in experiment and simulation. (Online version in color.) Fig. 4. Image-based modeling: (a) microstructure sketches of two ROIs with a dimension of 80 ×80 μ m, (b) quasi-3D models meshed by 8-node C3D8 elements and imposed boundary conditions. (Online version in color.) ISIJ International, Vol. 61 (2021), No. 6 © 2021 ISIJ1995 tioned two works. The CRSS CRSSprismatic {} / 1012 is set to 2.1 as Ref. 10) and CRSS values of pyramidal <a> slip, {} 1122 compression twinning are calibrated using curvefitting method due to a lack of systematic researches into these two systems. The anisotropic elastic stiffness matrix,37) rate-sensitive, strain gradient and strain-hardening parameters of CP-Ti are summarized in Tables 2–4. The activation energy ΔF for each slip system is 0.5447 eV, which is slightly smaller than the value of prismatic slip reported by Xiong and co-workers.33) The activation volume ΔV is 26b3, the same as the value calibrated by Gong et al.9) between 8b3–80b3, a general range for titanium and its alloys.38) The reference shearing rate  γ 0 and rate sensitive index r for twinning systems are taken to be  γ 0 = 0.001 and r = 0.05, respectively. The rank of CRSSs in ascending order is prismatic, pyramidal <a>, basal, {} 1012 tension twinning, {} 1122 compression twinning and pyramidal <c+a> slip, which resembles Ref. 8). Among these values, CRSSbasal/CRSSprismatic = 1.71 is a reasonable ratio in consideration of the biased strengthening effect caused by interstitial oxygen. CRSSpyramidal<a >/CRSSprismatic = 1.64, CRSS CRSSprismatic {} / 1122 = 2.37, are also acceptable with consideration of previous literature.10) Besides, the latent hardening constant in Eq. (12) is assumed to q=1 for all slip and twinning systems. It should be noted that the parameter determination in this study still lacks a mechanistic understanding, further investigation on some unsure parameters (e.g., the activation energy and volume, CRSSs of pyramidal <a> slip and twinning) should be carried out. 5. Results and Discussion 5.1. Macroscopic Response The data of monotonic tensile, stress relaxation tests and relevant simulation results are plotted in Figs. 5(a) and 5(b). The macroscopic SRS, strong tensile anisotropy and stress relaxation are faithfully reproduced by response of RVE model. Therefore, the rationality of CPFE constitutive parameters is verified. To study the relation between tensile anisotropy and deformation mechanisms, relative activities of different deformation modes at the strain rate of 2.3 × 10−4 s−1 are depicted in Fig. 6.8) The lower yield strength in TD can be attributed to high relative activity of prismatic slip at macroscopic yield since the CRSS of prismatic slip is the smallest. In contrast, the basal-dominant pattern at the beginning of Table 2. Anisotropic elastic constants of grade 1 CP-Ti. C11 (GPa) C22 (GPa) C12 (GPa) C13 (GPa) C44 (GPa) 162.4 180.7 92.0 69.0 46.7 Table 3. Rate-sensitive and strain gradient parameters for grade 1 CP-Ti. Slip Twinning Strain gradient Parameters Value Parameters Value Parameters Value ρ m ( μ m−2) 5  γ 0 (s−1)0.001 α T0.4 υ D (Hz) 1×1011 m0.05 G (MPa)44 230 b<a> ( μ m)0.32 b<c+a> ( μ m)0.51 ΔF (eV)0.5447 k (J/K)1.38 ×10 −23 ΔV26b3 Table 4. CRSS and strain-hardening parameters for slip and twinning systems, in MPa. Basal <a> Prism <a> Pyramidal <a> 1st Pyramidal <c+a> 2nd Pyramidal <c+a> {} 1012 Twinning {} 1122 Twinning h0850 1 750 750 2 050 2 050 1 050 1 050 τ ∞183 155 173 362 362 249 267 τ 0130 76 125 228 228 159 180 Fig. 5. Macroscopic responses obtained from (a) tensile, (b) stress relaxation tests. (Online version in color.) ISIJ International, Vol. 61 (2021), No. 6 © 2021 ISIJ 1996 plastic deformation leads to higher yield strength in RD. Additionally, a large initial hardening modulus (see Table 4) and continuous activity of the prismatic slip also lead to greater strain-hardening in TD. Although the activity of pyramidal <a> slip remains at a high level after the yield in RD, the strain-hardening effect is not as conspicuous as that of TD because of the low initial hardening modulus.39) Tensile anisotropy also manifests in twinning activities. The simulation results reveal that {} 1122 compression twinning and {} 1012 tension twinning are respectively dominant twinning systems in RD and TD, which keeps consistent with the experimental twin boundary maps after tension, as given in Figs. 7(a) and 7(b). Furthermore, relative activities of twinning at ε = 0.05 increase as the strain rate becomes higher in both RD and TD, as shown in Fig. 7(c). This observation captures the phenomenon that twinning is more active at high strain rates. Based on the calibrated parameters, the simulated cyclic responses of two ROIs are plotted in Fig. 8. The strain accumulates at a reducing rate per dwell cycle and macroscopic stress-strain behaviors are similar to the experimental trend in Ref. 2), as well as the simulated response in Ref. 4). 5.2. Strain Fields and Deformation Mechanisms In order to analyze the evolution of strain field under dwell fatigue loading, Schmid factor (SF) is employed to link the deformation mechanism with strain localization. The SF maps of basal, prismatic, pyramidal <a> slips, and {} 1122 compression twinning system within two ROIs are plotted in Fig. 9. Two ROIs are revealed to be soft-grainsdominant. Pyramidal <c+a> slip and {} 1012 tension twinning systems are not taken into account because the activations are difficult in RD, as plotted in Fig. 6(a). The evolutions of axial strain ε yy fields within two ROIs in the first dwell fatigue cycle are predicted by CPFE modelFig. 6. Evolutions of relative activities under 2.3×10 −4 s−1 strain rate during tension in (a) RD, (b) TD. (Online version in color.) Fig. 7. Twin boundary maps on (a) RD, (b) TD samples after tension. The tension and compression twin boundaries were respectively marked by red and blue lines. (c) relative activity of twinning at different strain rates when ε = 0.05. (Online version in color.) Fig. 8. Simulated cyclic responses of two ROIs. (Online version in color.) ISIJ International, Vol. 61 (2021), No. 6 © 2021 ISIJ1997 ing, as shown in Fig. 10. Due mainly to the elastic anisotropy of HCP lattice, strain heterogeneity at the beginning of dwell time (point D1) is observed even the maximum stress is under macroscopic yield strength. Due to cold creep during load-holding period, the strain accumulation at the end of dwell time (point D2) is conspicuous, which in turn aggravates the strain incompatibility between soft and hard grains at the unload state (point D3). An overall comparison between the experimental and numerical results with respect to residual strain distributions after 20 cycles are illustrated in Fig. 11. The numerical contours of residual axial strain ε yy and transverse strain ε xx are not in strict agreement with DIC measurements and manifest more heterogeneous deformation (i.e., strain incompatibility among grains), whereas most of the strain localizations are identified by CPFE modeling with fair accuracy. AdditionFig. 9. Schmid factor maps in (a) region A and (b) region B. (Online version in color.) Fig. 10. Evolutions of predicted axial strain fields within two ROIs in the first dwell fatigue cycle. (Online version in color.) ISIJ International, Vol. 61 (2021), No. 6 © 2021 ISIJ 1998 Fig. 11. Experimental and numerical results in ROIs after 20 cycles: (a) axial strain field, (b) transverse strain field. (Online version in color.) ally, in spite of the plastic strain accumulation, no significant discrepancy is discovered between the predicted axial strain distributions at the end of the 1st and 20th cycle, as illustrated in Figs. 10 and 11(a). Thus, it can be concluded that the local deformation pattern forms at the early stage of the dwell fatigue. This result is consistent with the DIC interpretation given by Littlewood et al.19) Regarding the residual axial strain ε yy field within two ROIs after 20 cycles, the experimental and numerical results exhibit similar strain localization paths, as highlighted in Fig. 11(a) by white dashed lines. In region A, the strain localizations are located at grain G1 (see Fig. 4(a)) and its neighborhood with high SFs for basal, prismatic and pyramidal <a> slip systems, as well as the center and bottom right corner including grains G3, G4, G5, and G7. In region B, grains G10-G14 and G17-G20 along the aforementioned paths demonstrate larger axial strain accumulation. In contrast, axial strain accumulation is less obvious in hard grains with low SFs for basal, prismatic and pyramidal <a> slip systems, such as grain G21. It should be noted that some deviations between DIC analysis and CPFE computation are observed at the center bottom of region A (grains G8 and G9) and the lower right corner of region B (grains G18, G19, and G20), especially the area containing grains G8 and G9. Since the spatial distribution of strain field on the free surface is likely to be affected by beneath layer, the local strain deviation between DIC strain mapping and CPFE simulation is possibly attributed to the mischaracterization of subsurface microstructure. The experimental and numerical mappings of residual transverse strain ε xx in two ROIs after 20 cycles achieve a better agreement than in the axial direction, as illustrated in Fig. 11(b). The localization band of compressive strain in region A is approximately along the top left corner to bottom right corner, whilst numerical prediction also misses the transverse strain concentration in grains G8 and G9 measured by DIC. As for the transverse strain field in region B, except for some small zones such as the top left corner, the overall compressive strain accumulation of region B is less pronounced than of region A. The strain localization directly results from slip and twinning activities. Although the operating slip and twinning systems are not incorporated explicitly as traces in crystal plasticity model, the activities can be captured by corresponding shear strain accumulation implicitly. To determine the contribution of individual plastic deformation modes to axial strain localization, the maximum shear strain accumulation on four active slip and twinning systems (basal, prismatic, pyramidal <a> slips, and {} 1122 compression twinning systems) are calculated by CPFE simulations, as given in Fig. 12. In region A, the axial strain accumulation within grains G1, G3, G4 and G7 is basically dominated by prismatic and pyramidal <a> slip activities because of high corresponding SFs. However, the contribution of basal slip activity is inconspicuous in region A. Concerning region B, the axial strain localizations in grains G11, G13, G15, G17, G18, and G19 are derived from the prismatic slip activity. Pyramidal <a> slip activity dominates the axial strain concentrations in grains G14 and G20, while basal slip activity leads to high axial strain accumulation in grains G10 and G12. This comparative analysis reveals that the activity of prismatic slip plays a primary role in axial strain localization, followed by pyramidal <a> and basal slips. In contrast, the activity of {} 1122 compression twinning demonstrates a low level in two ROIs considering the relatively high CRSS. As previously reported, the slip and twinning activations in CP-Ti under tension are observed to be prismatic slip-dominant, followed by pyramidal <a> slip, twinning and basal <a> slip.35,40,41) Besides this, previous studies have reported that the deformation mechanisms under different loading conditions (i.e., tension, fatigue and dwell fatigue loading) exhibit similar patterns with respect to frequency of different slip systems.18,42) The predicted slip activations, which are related to axial strain localization, closely match the quantitative results from previous slip trace analyses.35,40,41) However, the twinning activity is likely to be underestimated by the CPFE modeling considering only two twinning systems are modeled. Since the strain localization in an individual grain results from multiple active deformation modes, a weighted averaged