Dear all,

I am trying to implement to apply cyclic shear on a polycrystalline materials. I have used

fix NVT/SLLOD and fix deform command. The box deformation rate is time dependent.

The system is sheared by 1000 MD steps along +x direction and then under sheared condition another 1000 MD steps are simulated with the same NVT/SLLOD (Is it right to simulate the system under NVT/SLLOD when the box is not deforming anymore?). Then with the same rate the system comes back to the original position and a cycle ends with another 1000 MD steps under unsheared condition.

First question:

When MSD of the system is calculated it comes out that at a certain low strain rate (consequently lower maximum strain) MSD is becoming larger as compared to the higher strain rate. It sounds quiet strange to me. If you make some comment it would be good.

Second question:

In some cases, I found the translational drift of the centre of mass. To remove this problem I used fix momentum command. How much logical it is? What is the reason of COM movement?

Third question:

When I calculate the velocity of the top and the bottom layer (which is defined by region command) without using the fix momentum command it shows that the bottom layer has a small velocity along -x direction whereas one expects that its velocity should fluctuate around zero.

However, if I sum up the speed of the top and lower layer it matches with box deformation rate.

Why does the lower layer move with a certain velocity (-x direction)? How can I fix the problem?

The input file is shown below.

units lj

atom_style atomic

dimension 2

boundary p p p

read_data system.data

variable DT equal 0.005

timestep ${DT}

region upper block INF INF 104.409 INF INF INF

region lower block INF INF INF 1.0 INF INF

pair_style lj/cut 3.5 #LJ potential with cut off 3.5

pair_coeff 1 1 1.0 1.0 3.0 #First 1 1 says that the interaction between "1" type of atoms; Epsilon=1.0; Sigma=1.0; Cutoff=3.0

pair_coeff 2 2 1.0 1.8 5.4

pair_coeff 1 2 1.0 1.4 4.2 #First 1 2 says that the interaction between "1" and "2" type of atoms; Epsilon=1.0; Sigma=1.4; Cutoff=4.2;

pair_modify shift yes

neighbor 0.3 bin

neigh_modify delay 0 every 1 check yes

variable SS equal 1000

variable S1 equal 500

variable S2 equal 500

variable S3 equal 500

variable S4 equal 500

variable S equal "v_S1+v_S2+v_SS”

variable SR equal "v_S+v_S3"

variable f equal 0.1

reset_timestep 0

variable u equal 0.0

variable rate equal "v_f*(step*dt)" # deformation rate for the first 500 MD steps

variable displace equal "(v_f*((step*dt)*(step*dt)))/2" # Displacement of the box during the deformations

fix 2 all deform 1 xy variable v_displace v_rate remap v # Command to deform the box and v_rate and

fix 1 all nvt/sllod temp 0.01 0.01 0.1 tchain 3

compute myTemp all temp/deform

fix_modify 1 temp myTemp

compute MSD all msd com yes

fix ave_MSD all ave/time 100 1 100 c_MSD[1] c_MSD[2] c_MSD[3] c_MSD[4] file MSD.data

compute ave_vel all reduce ave vx vy vz

compute sum_vel_up all reduce/region upper ave vx vy vz

compute sum_vel_low all reduce/region lower ave vx vy vz

fix ave_vel1 all ave/time 1 1 1 c_sum_vel_up[1] c_sum_vel_up[2] c_sum_vel_low[1] c_sum_vel_low[2] v_rate c_ave_vel[1] file vel.dump

run ${S1} #by this steps box deformation rate reaches the maximum

variable u equal "v_f*v_S1*dt" #initial velocity when the box deformation rate decreases

variable rate equal "v_u-v_f*((step-v_S1)*dt)" # box deformation rate decreases

variable displace equal "(v_u*((step-v_S1)*dt)-(v_f*(((step-v_S1)*dt)*((step-v_S1)*dt)))/2)" # displacement

run ${S2} #by this MD steps period box deformation rate comes back to zero