I am trying to calculate surface tension for Argon liquid but I still have fluctuation in the surface tension results it is unstable even after 2 million timestep this is the script :
Argon Surface tension
units real
atom_style hybrid atomic sphere
lattice fcc 5.26
boundary p p p
region box block 0 1 0 1 0 1
create_box 1 box
create_atoms 1 box
forcefield
mass 1 39.948 # Ar
pair_style lj/cut 8.7675
pair_coeff 1 1 0.233 3.507 8.7675
replicate 10 10 10
variable xLo equal -0.5lx
variable xHi equal 0.5lx
variable yLo equal -0.5ly
variable yHi equal 0.5ly
variable zLo equal -0.5lz
variable zHi equal 0.5lz
Recenter system coords about the origin (0,0,0):
displace_atoms all move {xLo} {yLo} ${zLo} units box
Adjust box dimensions:
change_box all x final {xLo} {xHi} units box
change_box all y final {yLo} {yHi} units box
change_box all z final {zLo} {zHi} units box
Increase z edge to yield “Vacuum|Argon|Vacuum” system:
variable zLoNew equal -2.5lz
variable zHiNew equal 2.5lz
change_box all z final {zLoNew} {zHiNew} units box
variable Nblock equal 10000
variable Nrun equal 50*{Nblock}
variable Nr equal {Nblock}/100
variable A_in_m equal 1e-10 # Angstrom in meter
variable atm_in_Pa equal 101325 # note: 1 Pa = 1 N/m^2
variable N_in_mN equal 1e3 # Newton in milliNewton
variable convFac equal {A_in_m}*{atm_in_Pa}*${N_in_mN}
variable Text equal 100
group Argon type 1
velocity all create ${Text} 1234
neighbor 2.0 bin
neigh_modify every 1 delay 0 check yes
fix 1 all nve/sphere
fix thermostat all langevin {Text} {Text} 100 1234 zero yes
fix removeMomentum all momentum 1 linear 1 1 1
compute T all temp
fix TAve all ave/time 10 100 1000 c_T
variable P equal press
fix PAve all ave/time 10 100 1000 v_P
variable xPress equal c_thermo_press[1]
variable yPress equal c_thermo_press[2]
variable zPress equal c_thermo_press[3]
Evaluate and average surface tension in mN/m:
variable st equal 0.5lz(v_zPress-0.5*(v_xPress+v_yPress))*${convFac}
fix st all ave/time 10 100 1000 v_st
variable xyArea equal lx*ly
thermo_style custom step temp f_TAve press f_PAve etotal v_xyArea lz f_st
thermo_modify norm yes flush yes
thermo 1000
dump trj all custom 1000 Argon.trj id type x y z
dump_modify trj sort id pad 6
timestep 1
run ${Nrun}