A different SHOCK method is applied here than the method used in this referenced paper and produced the same phenomenon in the high-entropy alloy (HEA) system. Here I have used the shock front absorbing boundary condition (SFABC) method. In the SFABC method, a piston with a velocity of Up is impacted on the HEA system and a shock wave is generated and propagated along the shocked direction. When the shock wave reaches the free surface a second piston with the same velocity Up is added so that the shocked state can be maintained and observed for a longer time.[1]

Simulation procedure: Here I have adopted a slightly different approach during equilibration. I have applied the Monte Carlo (MC) atom-swapping technique to simulate a more realistic system.  The system size was 356*89*89 Angstrom. The atomic interactions within the Fe-Ni-Cr-Co-Cu HEA system were defined by the embedded atom method (EAM) potential developed by Deluigi et al.[2] After a cg (conjugate gradient) minimization and NPT equilibration at 300k, MC atom swapping was applied via the fix mc command in LAMMPS. Then the simulation box changed to 600 angstrom in length and equilibrated again under NPT at 300k. After that, the SFABC was applied along the x-direction with a piston velocity of 0.8 kms-1  for 30 ps. 

Simulation results: 

1. Potential energy during MC atom swapping decreased significantly in Fig. 1 which indicates a more stable structure.

2.  When the shock wave propagates along the shocked direction the temperature increases significantly which may contribute to the phase transformation behavior of this HEA system.

3. There was a plastic wave and elastic wave zone where dislocation nucleation occurred in the plastic zone. This dislocation nucleation occurred due the the phase transition from fcc to bcc phase. The bcc phase provided favorable nucleation sites for stacking faults (SFs) to transform into hcp for the dislocation nucleation.

4. Shock stress and shear stress increased with the increase of piston velocity

For better visualization, dislocation analysis (DXA), adaptive common neighbor analysis (ACNA), and atomic shear strain analysis are shown below for 0.8 kms-1 velocity.

Figure 2: a) temperature, and b) stress distribution along the shocked direction at different time.

[1]. Alexey V. Bolesta, Lianqing Zheng, Donald L. Thompson, and Thomas D. Sewell Phys. Rev. B 76, 224108 – Published  14 December 2007

[2]. O.R. Deluigi, R.C. Pasianot, F.J. Valencia, A. Caro, D. Farkas, E.M. Bringa, Simulations of primary damage in a High Entropy Alloy: Probing enhanced radiation resistance, Acta Materialia, Volume 213, 2021, 116951

I have tried to reproduce some results that were obtained in this referenced paper and successfully I was able to produce almost similar results.  All the simulation procedures were adopted from this article. 

Atomic-layered transition metal dichalcogenides (TMDCs) have versatile applications in electronics, optoelectronics, biomedical imaging, etc. MoS2 is one of the promising TMD inorganic compounds that have various applications as transistors, memory devices, solar cells, electrocatalysts, lithium batteries, lubricants, and so on. However, due to the oxidation of MoS2, its properties especially mechanical properties might change. To understand the effect of oxidation of MoS2 on mechanical properties this simulation will be useful. 

Simulation results: 

1. From 1500k the oxidation tends to start and the MoS2 layer readily reacts with oxygen and at 1800k it becomes severe.

2.  In the oxidized MoS2, its fracture stress is lower than the pristine condition.

3. Both pristine and oxidized MoS2 layer faces phase transformation during the tensile test and the zigzag wavy peak in the stress-strain curve are the response due to the phase transformation.

Figure 1: a)front, b)front, c) isometric view of the MoS2 system with Oxygen. d) stress-strain diagram for pristine MoS2 layer and oxidized layer for different oxidation temperatures.

For better understanding of the fracture mechanism in graphene/WS2/graphene heterolayer please watch the following video.

Nanowire

•Radius* Length: 24*240 Angstrom

•The input data file was created with Atomsk and the Cu atom was replaced with the Ni atom randomly with a certain atomic percentage (%at).

Simulation Method:

• Periodic boundary condition along the Z axis.

• Time step size of 1 femtosecond (1 fs) with 4 OpenMP threads per MPI task.

• Conjugate gradient (CG) minimization.

• 30 picoseconds (ps) of microcanonical ensemble (NVE) while keeping the temperature at  300 K using a Langevin thermostat.

• Isothermal-isobaric ensemble (NPT) for 30 ps.

• Canonical (NVT) ensemble for 10 ps.

• An engineering strain rate of 109 s-1 was along the Z [0 0 1] direction applied with fix deform.

• Virial stress was calculated.

Figure: (b) Stress-Strain curve for uniaxial tensile loading, (c) and (d) are the variation of Ultimate stress and Elastic Modulus with the different Cu content in a Cu-Ni nanowire, (e) number of vacancy-interstitial pairs vs strain at different Cu content of the Cu-Ni nanowire. 

Key observations

• Cu and Ni atoms exhibit lattice misfits; therefore, when they mix in an alloy, there are a lot of point defects after equilibration compared to pure structures. The Wigner-Seitz defect analysis revealed that when Cu content increases, the number of point defects (vacancies and interstitials) also increases, resulting in lower alloy strength.

• Cu-Ni nanowire elastic modulus and ultimate strain decrease with the increase of Cu concentration.

• The stress-strain curve agrees well with the results of the functionally graded Cu-Ni Alloy nanowire. [1]

•For 10% Cu of Cu-Ni nanowire, it fails at 11.5% strain and then it plastically deforms.

 The animation of elastic strain analysis during tensile simulation is added below for a visualization of the failure mechanism of the Cu-Ni alloy nanowire.

[1] https://doi.org/10.48550/arXiv.2009.06873

Simulation Method:

• Periodic boundary condition along the X and Y axis

• Time step size of 1 femtosecond (1 fs) with 4 OpenMP threads per MPI task

• Conjugate gradient (CG) minimization

•Isothermal-isobaric ensemble (NPT) for 30 ps

• Canonical (NVT) ensemble for 20 ps

• Engineering strain rate of 10e9 per second

• Fix deform with NPT

• Virial stress was calculated.

•Potential used: Airebo for C atom, Tersoff for WSe2 and LJ for inter-layer interaction.

Figure: Initial coordinates of (a) ~21.07 nm x ~21.20 nm nanosheet of Graphene/WSe2 heterostructure showing both armchair and zigzag direction.

Failure Modes for layers:

•At armchair loading condition:

•crack propagates along the zigzag direction in monolayer graphene and h-WSe2  respectively, figure (d) and figure (e).

•At zigzag Loading condition:

•Crack propagates diagonally in the monolayer graphene and h-WSe2 respectively, figure (f) and figure (g).

Figure: (d), (e) and (f), (g) fracture modes for graphene mono layer and h-WSe2 layer while load in along armchair and zigzag directions respectively.

•In the Gr-WSE2 heterostructure, fracture first initiates at the WSe2 layer.

•Phase transformation of the WSe2 layer during the fracture for various strain levels without stress contour at 300K for zigzag loading.

•At 12.5% strain, some h-WSe2 trigonal lattice transforms to an octahedral distorted phase (t-WSe2).

Figure: However, atomic rearrangement of h-WSe2 occurs at 300K (h) along the armchair loading. and (i) along the zigzag loading, atomic rearrangement of the WSe2 layer starts at 12.5% strain, and this process continues along 45 degrees (diagonally) with the loading direction. During MD simulation, the first sign of a bond breaking is seen at 13% strain, after more stretching.

Key Findings:

•Molecular dynamics is used to study the mechanical characterization of Gr/WSe2 vertical heterostructures.

•The mechanical properties and fracture mechanisms of the Gr/WSe2 are strongly direction-dependent.

•Fracture initiates from the graphene and WSe2 layers for armchair and zigzag loading, respectively.

•The WSe2 sheet undergoes a phase transformation from h-WSe2 to t-WSe2 during zigzag loading.

Simulation Method and Results:

•Potential used: MEAM, Airebo, LJ

•Equilibration: CG minimization, NPT

•Deformation: Fix deform with NPT and strain rate is 10e9 per second

Figure:  (a) 70*66.8 Angstrom heterostructure of Graphene-Silicene-Graphene (Gr-Si-Gr) (b) mechanical characteristics of Gr-Si-Gr at 300K.

Figure: Stress-strain curve of graphene, h-WSe2, and graphene/WSe2 at 300K temperature under (c) armchair (d) zigzag loading.

Key observations:

• During the armchair loading, the failure of the heterostructure initiates in the graphene layer, and during the zigzag loading, the failure of the heterostructure initiates in the silicene layer.

• In armchair loading, the crack moves in a direction that is opposite to the direction of the load. In the zigzag sheet, the crack moves in a diagonal direction.

• The fracture stress and fracture strain of the heterostructure decrease with the increase in temperature, where the thermal vibration of the atoms leads to bond dissociation.

• Strain rate has a negligible effect on elastic modulus where the fracture stress and fracture strain increase with the increase of strain rate. At a slower strain rate, there is more time for thermal fluctuation, which eventually leads to bond dissociation.