A Crystal Plasticity Model with Irradiation Effect for the Mechanical Behavior of FCC Metals

In this paper, a crystal plasticity model considering the irradiation effect based on the thermal activation theory is established. The evolutions of screw dislocations, edge dislocations, and stacking fault tetrahedrals (SFTs) (induced by irradiation) are included into the model. The interactions between dislocations and irradiation-induced SFTs are also considered. The constitutive model is numerically implemented on the ABAQUS platform through UMAT subroutine and applied to study the irradiation effect on the mechanical behavior of pure copper. The mechanical properties of single and polycrystalline copper are studied, and the simulation results show that the constitutive model can properly predict the mechanical behavior of irradiated pure copper. Especially for polycrystalline copper, the simulation results are in good agreement with the experimental data.


Introduction
In nuclear industry, the mechanical properties of metal materials are influenced by irradiation and temperature. After irradiation, metal materials with different crystal structures always show an irradiation effect, i.e., a significant increase in yield strength and a decrease in toughness [1,2]. The irradiation effect becomes more pronounced as the irradiation dose increases. When the irradiation dose increases to a certain value, an upper yield point will appear and then the stress decreases. Moreover, the intensity of irradiation effect is also affected by the irradiation temperature. The results of the post-irradiation tensile testing exhibit stronger irradiation hardening and show that materials irradiated at room temperature have less ductility compared with the high temperature's occasion [3,4].
From the microperspective, materials in the irradiation environment are impacted by energetic particles, which incurs cascade reaction subsequently and results in a large number of nanoscale irradiation defects, such as vacancies, interstitial atoms, voids, precipitates, stacking fault tetrahedrals (SFTs), small dislocation loops, and so on [5]. The type and density of defect are related to the crystal structure. For example, Singh et al. [6] investigated the microstructure and associated tensile properties of irradiated FCC (Cu, Pd and 304 stainless steel) and BCC (Fe, Mo) metals by experiments, and revealed the differences and similarities of metals with different crystal structures after irradiation.
Otherwise, many studies show that the plastic deformation of the crystal is attributed to the dislocations slipping along the corresponding slip planes. Crystal plasticity theory is a basic theory for describing the plastic deformation at mesoscale, which connects the microscopic slip mechanism with the macroscopic plastic deformation behavior [7][8][9][10]. The deformation mechanism of the irradiated FCC materials is greatly affected by the dislocations and irradiation defects. Cheong et al. [11] developed a dislocation-mechanics-based crystallographic theory to study the mechanical properties of thin polycrystalline Cu specimens. Though the irradiation effect has not been considered, the evolution model of screw and edge dislocations is useful for studying the mechanical behavior of FCC crystal. Arsenlis et al. [12] introduced an internal state variable model for the mechanical behavior of irradiated structural materials within a multi-scale framework. Using the model, the behavior of Cu tensile specimens with varying irradiation damage was simulated. De et al. [13] proposed a defect-and dislocation-density-based evolution model to capture the features of irradiation hardening as well as intra-granular softening. The Jacobian-free multi-scale method (JFMM) was further applied to improve the computational performance in the polycrystalline aggregate simulations by using a Newton-Krylov process. The mechanical responses of neutron-irradiated single and polycrystalline OFHC copper were studied and the model could capture experimentally observed grain-level phenomena. Chen et al. [14] investigated the mechanical behavior and crystallographic texture evolution of irradiated FCC metals based on a physical theoretical model. The study revealed the texture evolution along with different orientations before and after irradiation. It was concluded that irradiation-induced defects could affect the mechanical behavior and texture evolution of metals, both of which were closely related to irradiation hardening. Xiao et al. [15] established a tensorial crystal plasticity model with both irradiation and temperature effects. Using the model, the mechanical behavior and characteristics of irradiated polycrystals at different temperatures could be captured. Besides, Xiao et al. [16] also explored the effects of irradiation damage and crystal size on the mechanical behavior of FCC single crystals based on a unified size-dependent tensorial plasticity model. Krishna et al. [17] proposed a micromechanicsbased model for copper subjected to neutron irradiation. The defect evolution was considered in the model to reflect the prevention and annihilation effects of defects on dislocations.
In this paper, a crystal plasticity model considering irradiation effect is established by tracking the evolutions of screw dislocations, edge dislocations, and SFTs induced by the irradiation. The model can also reflect the interactions between dislocations and irradiation defects and has been implemented numerically on the ABAQUS finite element platform. Using the model, the mechanical behaviors of single and polycrystalline copper before and after irradiation are studied and predicted.

Crystal Plasticity Theory
According to the classical crystal plasticity theory [7][8][9][10], the total deformation gradient F is divided into two parts: F P , which denotes crystallographic slip, and F * , which represents stretching and rotation of the lattice.
where N is the number of all activated slip systems, I is the second-order identical tensor, γ α is the shear strain. m α and n α are the slip direction and normal-to-slip plane of the α slip system, respectively. The slip direction and normal-to-slip plane after deformation are given by The velocity gradient L can be decomposed into the lattice part L * and the plastic part L P : L can be also decomposed into the deformation rate tensor D and the spin rate tensor Ω.
whereγ α is the slipping rate of the α slip system. The crystal macroscopic deformation and microslipping rate are portrayed through the above formulas. It is assumed thatτ is the Jaumann derivative of Kirchhoff stress τ which takes the initial configuration as the reference state. The relationship between the shear stress and the shear strain can be derived byτ where C is the stiffness tensor. According to the work of Kocks [18], the slipping rate of the α slip system can be expressed as where Q 0 denotes the activation energy required to overcome the obstacles to dislocation motion, p and q are flow-rule-related parameters, G and G 0 represent the shear modulus at temperature T and 0 K, respectively. τ α and τ α 0 represent the lattice friction stress of each slip system at temperature T and 0 K, respectively. It is assumed that τ α 0 is the same in each sliding system, so τ α 0 is simplified to τ 0 . And the term g α represents the total slip resistance to dislocation motion, including the parts caused by dislocations and irradiation defects.

Evolution of Dislocations
The strain hardening behavior of the crystal is closely related to the dislocation density. The dislocation density increases gradually during the process of plastic deformation, making the interaction between dislocations more intense. So the dislocation motion becomes more difficult. The strain hardening mechanisms include: dislocation congestion, dislocation intersection, the formation of fixed dislocations, the reduction of Frank-Lead dislocation sources, and so on. In order to describe the strain hardening behavior of the crystal, Cheong [11] proposed a new hardening formula based on the dislocation density, here λ is a statistical coefficient, b α is the magnitude of the Burgers vector. ω 1 and ω 2 are the self and latent hardening moduli, respectively. δ αβ is the Kronecker delta. ρ α e and ρ α s are the densities of edge and screw dislocations, respectively.
It is assumed that dislocation multiplication is associated with the expansion of dislocation loops from Frank-Read-type sources and dislocation reduction is associated with the annihilation reaction when the anti-parallel dislocations meet during the motion process. According to the study of Cheong [11], the evolutions of edge dislocations and screw dislocations can be expressed aṡ The terms C e and C s characterize the contributions of edge and screw dislocations to the slip rate, respectively. d s and d e are the critical edge annihilation distances. ρ t is the total density of dislocations in a slip system. K e and K s are dimensionless constants, both of which control the motion of dislocations.

Evolution of Defects
Irradiation defects in the FCC copper crystal are mainly SFTs [4,5], with a low density of small dislocation loops. In this paper, the defect evolution model proposed by Krishna [17] is adopted where the SFT is regarded as the main irradiation defect. Besides, the annihilation of irradiation defects is due to the reactions between dislocations and SFTs. The full absorption of SFTs by screw dislocations and partial absorption of SFTs by edge dislocations are reported in many studies [19][20][21]. The absorption mechanism is simplified and then applied to the model. However, other mechanisms that lead to annihilation, such as localized plastic heating, or stress-assisted conversion from SFTs into Frank loops, have not been considered at present.
It is assumed that the defects within the critical distance are readily absorbed. The annihilation mechanism of a defect encountering a screw dislocation is shown in Fig. 1. The small triangle represents an SFT. As shown in Fig. 1, the dislocation with a length of L α s moves a distance of Y α s before arrested by forest dislocations in a time interval dt. So the annihilation area of the defects is here d q is the critical distance.
According to the similar expression obtained for edge dislocations, the additional term πd 2 q is deleted due to its inability to cross-slip. In this case, the annihilation area of the defects is here Y α is the mean free path of the dislocation segment, which is defined as the moving distance of the segment before arrested by forest dislocations (Y α e is for edge dislocation and Y α s is for screw dislocation).
l s is the mean obstacle spacing. According to the geometric considerations, l m can be expressed as Then the annihilation probabilities of SFTs are Therefore, the change rate of the SFT density iṡ here ρ α def is the areal density of the irradiation defects in the α slip system. d α def is the average defect size. It is supposed that only the overlapped or truncated SFTs are annihilated, which account for about 70%-80% of the total defects [22,23] and the proportion is represented by the term ϕ. The evolution of defects is associated with the shear strain rate through Eq. (32).
Similar to the forest dislocations, SFTs and other defects can impede the dislocation motion, resulting in the irradiation hardening. In order to simulate the mechanical behavior of the irradiated crystal Fig. 1. Schematic of reactions between dislocations and other defects material, an irradiation hardening term is introduced into Eq. (21) and the irradiated slip resistance is here i αβ is the coefficient representing the interactions between defects and dislocations in different slip systems. ω i1 and ω i2 are the self and latent hardening moduli, respectively. The areal initial density of SFTs can be calculated by where N α def is the volume density of SFTs.

Finite Element Method
The finite element platform ABAQUS provides an interface for users to define a constitutive model in a Fortran subroutine called UMAT [24]. In this paper, the tangent modulus method developed by Pierce [25] is used to implement the UMAT. It is defined that Δt is the time increment. The increment of shear strain in slip system α can be expressed as Then a linear interpolation within Δt is employed as: here θ is an integration parameter ranging from 0 to 1, while θ = 0 corresponds to the simple Eulerian time integration scheme. Peirce recommended a choice of θ between 0.5 and 1 [25]. The Taylor expansion of slipping rateγ α is written asγ Therefore, the incremental expression of shear strain is So the following formulas are obtained based on Eq. (19) The incremental expression of Δτ α and Δg α can be derived via Eqs. (18) and (19) as follows The shear strain increment Δγ α can be calculated by Eqs. (37)-(41). Once Δγ α is determined, the increments of other variables can be also derived.

Material Parameters
There are 12 crystallographic equivalent slip systems {111} 110 in a FCC metal crystal. The dislocation density of typical fully annealed FCC pure metal is about 2 × 10 12 m −2 [17]. It is assumed that the initial densities of edge and screw dislocations in each slip system are the same, and the value is about 80, 000 mm −2 . The τ 0 for unirradiated copper mainly depends on dislocations and is set as 20 MPa. For irradiated copper, besides the main radiation defect SFT, there are many other types of defects, such as voids, small clusters of self-interstitial atoms, Frank dislocation loops, precipitated phase, which would also impede the dislocation motion and result in irradiation hardening. The influence of SFT on slip resistance is considered by the evolution Eqs. (25)-(35). Besides, the influences of other defects on slip resistance also need to be considered. According to [11], the slip resistance caused by irradiation defects is assumed the same as the unirradiated lattice friction stress. Therefore, the τ 0 for irradiated copper includes the influences of dislocations and irradiation defects, and its value is set as the sum of the 20 MPa caused by dislocations and the 20 MPa caused by irradiation effects, that is, 40 MPa.
The annihilation critical distance in Eq. (25) is 1-3 nm [19][20][21], here d q = 2.4 nm. Kawasaki [26] found that the mean free path lengths of edge dislocations and screw dislocations are about 100 μm and 200 μm, respectively. The total mean free path in Eq. (26) is the average value of the above two types. The value of ϕ in Eq. (27) ranges from 0.5 to 1, here ϕ = 0.7 [27,28]. The mean size of SFTs is about 2.3-2.6 nm [5], so d def is chosen as 2.5 nm. According to the simulation and corresponding experimental results, ω i1 = 1.6, ω i2 = 1.1. Other material parameters can be found in [11,29], and all the parameters used are shown in Table 1.

Mechanical Behavior of Single Crystal
According to the tensile test sample of Takeuchi experiment [30], a Cu single-crystal model of 22.5 mm × 5 mm × 2 mm is established. In order to approach the real uniaxial tensile process, the symmetric boundary condition is applied at one end of the model, and a displacement-controlled load is applied at the other end with a true strain rate of 3 × 10 −3 s −1 . The testing temperature is set to be 295 K and the shear modulus G = 45, 000 MPa [29]. The element type is C3D8. To ensure the precision of the simulation results, the mesh independence is verified. The model is divided into 210 and 1520 grids separately, and the stress-strain curves are basically the same (Fig. 2). Consequently, the model is divided into 210 grids with less computation time.
Since the single crystal is anisotropic, the results in crystal orientations of

Mechanical Behavior of Polycrystal
Singh [5] conducted the irradiation experiment of OFHC copper and studied the change of the mechanical behavior. Since the polycrystal is composed of lots of small single crystal grains with different spatial orientations, the polycrystal displays macroscopic isotropy. With the aim to facilitate the analysis, the grain boundary is assumed as pure geometric plane, which can ensure the continuity of deformation and the compatibility of displacement. In the simplified grain boundary model, the dislocation motion is impeded by a large number of defects in boundaries, so the dislocation slipping mainly occurs at interior of grains [31]. Since the grains of OFHC Cu used by Singh [5] were relative large and the sample was fully annealed, the hypothesis is consistent with the actual structure. A plane strain model of polycrystal is established as shown in Fig. 6.  The polycrystalline model is generated by Voroni principle, and the map of 60 grains is generated in the 0.45 mm × 0.1 mm plane with an average 30-micron grain size. The average grain size is in accordance with the experimental sample.
The volume density and average size of SFTs in irradiated OFHC Cu are shown in Table 2. As can be seen, the volume density of SFTs increases accordingly with the increase of DPA. When DPA reaches a certain value, the volume density of SFTs is nearly saturated. Assuming that the initial density of defects in each system is the same and there are 12 slip systems in a FCC crystal, the N α def in each slip system is 2 × 10 22 m −3 when DPA is 0.01, and 3.75 × 10 22 m −3 when DPA is 0.1. The testing temperature is set to be 373 K, the shear modulus G decreases to 44,000 MPa and the temperature-dependent elastic constants C 11 = 163,800 MPa, C 12 = 120,400 MPa, C 44 = 73,500 MPa are calculated by Varshni's formula [32]. Based on the experimental data for irradiated copper [1,5], the coefficients h αβ and i α β are adjusted to reflect the influence of temperature on the evolution of dislocations and defects with ω 1 = 1 and ω i1 = 1.1. Then the uniaxial tensile mechanical behaviors of non-irradiated OFHC Cu as well as the OFHC Cu under the conditions of DPA = 0.01 and DPA = 0.1 are simulated respectively.
Since the deformation is not uniform in the polycrystalline uniaxial tension, the stress and strain are obtained by the volume averaging method, hereσ andε are the total average stress and strain, respectively. σ k and ε k are the average stress and strain at integration point k, respectively. υ k is the volume of integration point k.
The stress and strain contours calculated by the uniaxial tension of polycrystalline model are shown in Figs. 7 and 8. By comparing the grain distribution of the model, it is found that the stress and strain are larger at the boundary of grains, but smaller inside the grains, which reflects that the dislocations accumulate at the grain boundary, resulting in stress concentration and non-uniform deformation.
According to the above calculation, the average stress-strain curves of uniaxial tension of polycrystalline model at different conditions are simulated. The results of three sets of random crystal orientation are shown in Fig. 9. It can be seen that the stress-strain curves are very close and have the same trend, indicating that 60 grains are sufficient to reflect the tensile properties of the material. The simulation results of OFHC Cu, as shown in Fig. 10, are in good agreement with the experimental data given by Singh [5] as a whole. With the increase of irradiation dose, the strength of the material increases gradually, and the increments of the yield strength are more obvious.

Conclusion
In this paper, a crystal plasticity constitutive model with the irradiation effect is developed, which incorporates the microscopic interactions between dislocations and irradiation defects into the crystal plasticity theory. In the model, the effects of edge dislocations and screw dislocations on irradiation defects are considered, respectively. Then the model is numerically implemented and applied to simulate the uniaxial tensile deformation of single and polycrystalline copper subjected to neutron irradiation. The simulation results show that: (1) The constitutive model can reflect the anisotropy of single crystal and the mechanical behavior change of the irradiated material, such as the increment of yield strength.
(2) This constitutive model is able to reflect the effect of irradiation dose on the mechanical behavior of pure copper: with the increase of irradiation dose, the yield strength increases as well. (3) For polycrystalline copper material, the simulation results are in good agreement with the experimental data.