Evolution of Internal Crack in BCC Fe under Compressive Loading

A molecular dynamics model has been developed to investigate the evolution of the internal crack of nano scale during heating or compressive loading in BCC Fe. The initial configuration does not contain any pre-existing dislocations. In the case of heating, temperature shows a significant effect on crack evolution and the critical temperature at which the crack healing becomes possible is 673 K. In the case of compressive loading, the crack can be healed at 40 K at a loading rate 0.025 × 1018 Pa·m/s in 6 × 10 s. The diffusion of Fe atoms into the crack area results in the healing process. However, dislocations and voids appear during healing and their positions change continuously.


Introduction
Metals own the ability of self-healing.It is possible that low strength materials might be made more reliable at a low cost if the crack healing is applied at an intermediate stage of their service life.The self-healing effect for creep cavitation in heat resisting stainless steels has been developed by alloy design, so micro-cracks can be constrained from propagating by automatically forming compounds on the crack surfaces although they may not be annihilated [1,2].
A lot of research has been conducted on crack healing in inorganic and composite materials while less work has been performed on metals.The damage evolution in metals was analysed [3], and many kinds of heat treatment methods can be used to decrease the number of defects in metals after a long-term operation at high temperatures [4].By TEM observation, a crack in α-Fe was found to be healed completely when the temperature was increased to a critical value of 1073 K [5].It was suggested that the thermal or mechanical treatment provides the driving force to minimise the surface energy resulting in the crack healing [6].The healing of internal pre-cracks obtained using plate impact technology was studied in a vacuum at an elevated temperature [7,8], and the inhomogeneous microstructure in crack healing area was analysed [9].It has been shown that the healing of fine cracks (100 nm) begins immediately after heating then proceeds according to the vacancy diffusion mechanism in the surface layer about 10 μm thick [10].
The effect of hydrostatic pressure on the shrinkage of cavities in copper introduced by creep at 773 K and situated on grain boundaries was studied [11].The density of the specimen showed a rapid increase when a hydrostatic pressure was applied during annealing, while the change was little when annealed under atmospheric pressure.The similar result was also obtained in [12].One of the reasons why the ductility of experimental steels can be improved over the whole working temperature is to annihilate or heal micro-cracks that may occur during deformation [13].The internal crack can be healed completely and rapidly under hot plastic deformation [14].It is possible to intensify crack healing and/or refinement of the microstructure by adjusting the parameters of plastic deformation or isostatic pressure [15].
Molecular dynamics (MD) has been adopted to discover the essence of crack healing in metals.MD simulation was conducted to study crack healing in Al [16] and Cu crystals [17], both of which are FCC.Dislocation generation and motion occur during crack healing.The critical temperature for crack healing depends upon the orientation of the crack plane.BCC crystal is characterized by four close-packed directions, the <111> directions, and by the lack of a truly close-packed plane such as the octahedral plane of the FCC lattice.
The crack healing behaviour in BCC-Fe crystal was simulated, however, only heating has been considered [18].In this study, the molecular dynamics model has been further developed to investigate the evolution of nano scale crack in BCC Fe.Not only has the effect of temperature ranged 273 -1173 K been considered, but also the effect of compressive loading is taken into account.

Computational Method
The molecular dynamics method is built based on a set of Newtonian dynamics equations for a number of particles with the macro constraints and boundary condition.They are solved by a numerical method to achieve the classical trajectory and velocity.The macro physical properties and motion routine can be obtained based on the statisticcal average of the results calculated for a long duration of time.
Assuming the total number of atoms is N, which follows Newton's law.
where , U tot is the total potential energy of N atoms, is the atom mass.The position i r and velocity i of every particle in a system can be obtained by solving Equation (1).The position i and velocity of every particle at time depend on the position and velocity of that particle respectively at time t, and the leapfrog algorithm [19] is applied to calculate the positions and velocities of atoms.

 
The reliability of the molecular dynamics simulation depends on the description of the inter-atomic potential.The inter-atomic potential adopted here is the N-body potential proposed by Finnis and Sinclair according to the embedded atom method (EMA) [20], the total energy of an assembly of atoms can be written as where the second term on the right hand side of Equation ( 2) is the conventional central pair-potential summation, while the first term is the N-body term, in which where both V and Φ are functions of the inter-atomic distance r ij .Constructed by [21], the following functional forms have been adopted.
where the parameters a 1 , a   H x is the Heaviside unit-step function which gives the cut-off distance of each spline segment, r 1 and R 1 represent the cutoff radii of V and Φ respectively.

Model Description
The x-axis is along the extension line of the central crack, the y-axis along the normal of the crack plane, and the z-axis along the 1  The crack is expressed by removing some atoms in the centre of the cell [23] and the angle between the   1 10 slip plane and the crack plane is set to be 60˚ randomly.The crack is about 5.957 × 10 -9 m long and over 1.500 a wide.The width is larger than the cut-off radii of Φ, 1.300 a 0 , which means that there exists no interplaying force between the atoms on either side of the crack surface.The total number of atoms is 12582.
Because of the limit of computation scale, the number of atoms in the simulation is usually much lower than that in a real system.Boundary conditions have significant effect on the reliability of simulation results.In the case of heating, fixed boundary condition [24,25] is adopted along the x-and y-axis directions, which means that the locations of the boundary atoms are constant during the whole calculation process.In the case of loading, displacement boundary condition is adopted along the x-and y-axis directions, i.e. the boundary atoms have displacement based on the biaxial plane strain linear elastic displacement field in the solid containing crack [26].
where u x , u y are the displacement in x and y direction respectively, E is Young's modulus, γ is Poisson's ratio.
, in which  is the applied stress, a is the half length of the central crack and z = x + yi.In both cases, periodical boundary condition [17,27] is used in z direction.
The initial velocity is the Maxwell-Boltzmann distribution corresponding to a given temperature.The system temperature is maintained constant by scaling the atom velocities during the simulation.
The equilibrium of the atomic configuration is obtained after a long time relaxing at a low temperature near 0 K, as shown in Figure 1.No defects are set in this initial configuration.d y representing the minimum vertical distance between the atoms on top crack surface and bottom crack surface is adopted for assessing healing process [18].It is 1.526 a 0 in the initial configuration.Crack healing becomes possible when d y is less than the cut-off radii of V, 1.180 a 0 .

Heating
The simulated temperature in the model in this study has to be below 1185 K for keeping the parameters of Nbody potentials for BCC Fe valid.Simulations were carried out in the temperature range of 273 -1173 K and the increment is 50 K.The simulation spans 6.0 × 10 −11 s and the time step When the temperature is 573 K, the fluctuation may be over the cut-off radii of Φ, R 1 = 1.300 a 0 .When the temperature is 623 K, the fluctuation range is between R 1 = 1.300 a 0 and the cut-off radii of V, r 1 = 1.180 a 0 .When the temperature is 673 K, the fluctuation range is around r 1 = 1.180 a 0 .When the temperature reaches 773 K, the fluctuation ranges decrease rapidly to be around 0.500 a 0 .When the temperature reaches 923 K, the atoms on top and bottom surfaces move actively as shown in Figure 3, and d y becomes 0 when the simulation time is 3 × 10 −11 s, as shown in Figure 3(c).The reason why the fluctuation ranges can decrease dramatically when the temperature reaches about 773 K may be the atoms on top and bottom surfaces become active enough at this temperature and can be brought to within a close range abruptly.A threshold value might exist but need to be delivered in the future study.
The above results show that the critical temperature at which the crack healing becomes possible is 673 K.A transition temperature may exist in the range of 673 -773 K, above which the process of crack healing is accelerated dramatically.The temperature at which a complete healing may be achieved is 923 K. Figure 4 shows the process of crack healing at 1173 K.
Only after 1 × 10 −13 s, d y decreases to 0.989 a 0 , which is less than r 1 = 1.180 a 0 , as shown in Figure 4(a); after 1.1 × 10 −12 s, the crack evolves into two sections, as shown in Figure 4(b); after 1.8 × 10 −12 s, only two voids are left at the position of two crack tips, as shown in Figure 4(c); after 4.5 × 10 −12 s, the crack is healed completely, as shown in Figure 4(d); the state of crack healing remained stable after that.

Compressive Loading
The atomic configuration shown in Figure 1 is also adopted as the initial configuration for simulating crack evolution during compressive loading.The environment temperature is 40 K.The time step Δt s 1 × 10 −15 s and the loading rate is 0.025 × 10 18 Pa•m 1/2 /s. the happening of crack healing, as shown in Figure 5(a).After 6 × 10 −12 s, d y = 1.235 a 0 , which is less than R 1 = 1.300 a 0 , which means that the crack tends to be healed as shown in Figure 5(b).When the time is 3 × 10 −11 s, the crack is healed completely, as shown in Figure 5(c).Then it remains stable and the morphology when the time is 6 × 10 −11 s is shown in Figure 5(d).
During crack healing process, atoms deviates from their normal positions, which results in many defects such as dislocations and voids similar to that in Al [16] and Cu [17].Some defects in the simulation cell at 4.5 × 10 −12 , 6.0 × 10 −12 , 3.6 × 10 −11 and 6.0 × 10 −11 s after the   crack has been healed completely at 1173 K are marked in Figures 6(a)-(d) respectively.The dashed rectangle frames show original crack locations, and the circles show the positions of some defects.The distribution of the defects in the cell is inhomogenous.Defects can be found to aggregate in the region where is with the original crack, especially around the areas where are the crack tips.The positions of these defects change continuously, and the inhomogeneity of the distribution of the defects decreases with an increase of time.
The crack can be considered a very large defect in the simulation cell.The existence of nano scale crack surface results in higher system energy than that without inner crack, so there is a tendency to eliminate the crack to reduce the system energy.However, the crack healing process cannot happen when the temperature is low because an energy barrier exists between the state when the crack exists and the state when the crack is healed.Importing energy by temperature rise to activate atoms can help overcome the energy barrier.
In the process of crack healing, the activated atoms close to crack immigrate to the crack area and vacancies appear in the positions from where atoms emigrated, then the atoms adjacent to these vacancies need adjust their positions.Because the crack size is much larger than the size of any other defects appearing in the simulation cell, it may cost much more time for the atoms immigrating to the crack area to rearrange their positions than the atoms that only need adjust their positions to accommodate adjacent vacancies.This results in that the distribution of defects is inhomogeneous in the early stage of healing.
The healing of a crack means that a large defect evolves into a lot of small defects and all the defects may not be eliminated.However, the atoms adjust their positions continuously in order to reduce the system energy even after the crack has been healed.This continuous position adjustment results in more homogeneous distribution of defects with an increase of time, and the quantity of the defects in the crack area may be reduced.Some defects in crack healing region when the time is 3 × 10 −11 and 6 × 10 −11 s during compressive loading are marked in Figure 7.These defects also change their positions, but the change is not so dramatically.This may be because the atoms are not very active due to the low temperature.Most defects aggregate in the region where the original crack was.

Conclusion
The critical temperature at which nano scale crack healing of BCC Fe crystal without pre-existing dislocations becomes possible is 673 K. Temperature has a significant effect on the evolution of nano scale crack.Healing proceeds rapidly at 1173 K. Compressive loading can result in crack healing.Crack is healed completely at 40 K at loading rate 0.025 × 10 18 Pa•m 1/2 /s in 6 × 10 −12 s.Dislocations and voids appear during healing process, and their positions change continuously.In the case of heating, the distribution of the defects in the cell is inhomogenous in the early stage of healing.This inhomogeneity decreases and the quantity of the defects in the crack area may be reduced with an increase of time.In the case of compressive loading, most defects aggregate in the region where the original crack is.
12     direction.The lattice constant a for BCC Fe is 2.8665 × 10 −10 m.The length of the calculation cell along the x-axis direction is layers of atoms along the x-axis direction constitutes a cycle of 0 3 2 a long, every two layers of atoms along y-axis direction constitutes a cycle of 0 2a long, and every six layers of atoms along z-axis direction constitutes a cycle of 0 6a long.

is 1 ×
10 −16 s.Δt Relationships between d y and time at the temperatures of 573, 623, 673, 773, 873 and 923 K are shown in Fig- ure 2 for comparison.It can be seen that d y experiences a dramatic decrease in the first 1.0 -2.0 × 10 −12 s, then increases and tends to fluctuate in a narrow range in all cases.The dramatic decrease may result from the unstable calculation in this first 1.0 -2.0 × 10 −12 s.The fluctuation range moves down with an increase of temperature.