Argon Surface tension

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.5
lx
variable yLo equal -0.5ly
variable yHi equal 0.5
ly
variable zLo equal -0.5lz
variable zHi equal 0.5
lz

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.5
lz
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}

Your message is difficult to read, since you copied your input without putting it into a block of triple backquotes (```). Like the following. Because of that, some of the special characters are interpreted as typesetting markup and thus are not shown:

# 3d Lennard-Jones melt

variable	x index 1
variable	y index 1
variable	z index 1

variable	xx equal 20*$x
variable	yy equal 20*$y
variable	zz equal 20*$z

units		lj
atom_style	atomic

lattice		fcc 0.8442
region		box block 0 ${xx} 0 ${yy} 0 ${zz}
create_box	1 box
create_atoms	1 box
mass		1 1.0

velocity	all create 1.44 87287 loop geom

pair_style	lj/cut 2.5
pair_coeff	1 1 1.0 1.0 2.5

neighbor	0.3 bin
neigh_modify	delay 0 every 20 check no

fix		1 all nve

run		100

in.ar_v7 (2.2 KB)

i uploaded it in replay thank you very much really thank you because i cant get any help

Looking at your input, I have lots of questions. I am interleaving my questions with the input below.

# Argon Surface tension
units real
atom_style  hybrid atomic sphere 

Why atom_style hybrid atomic sphere?
First off, all properties from atom style atomic are contained in atom style sphere, so there is no need for the hybrid, but also, the lj/cut pair style does not use the diameter. And nothing below uses it except fix nve/sphere, but that won’t make a difference since there is not going to be a torque on the atoms.

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 5 5 5

Why do a 1x1x1 box first and then replicate to 5x5x5? why not create a 5x5x5 box directly?

variable xLo equal -0.5*lx
variable xHi equal  0.5*lx
variable yLo equal -0.5*ly
variable yHi equal  0.5*ly
variable zLo equal -0.5*lz
variable zHi equal  0.5*lz

# Recenter system coords about the origin (0,0,0):
displace_atoms all move ${xLo} ${yLo} ${zLo} units box

Because you have periodic boundaries and a dense system, this command has effectively no impact outside of changing the image flags and the identity of atoms that are near the center or the edge of the box because they’ll be replaced by identical atoms at the same location.

# 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

Again, why not create the box you want right away, i.e from -2.5 to 2.5 and fill it with atoms.

# Increase z edge to yield "Vacuum|Argon|Vacuum" system:
variable zLoNew equal -2.0*lz
variable zHiNew equal  2.0*lz
change_box all z final ${zLoNew} ${zHiNew} units box

This does not add a vacuum but with stretch the entire system, including the atoms. If you create an empty group: group none empty and then use change_box none ... will keep the atoms in place.
But also this is a moot point, since you could have created a box from -5 to 5 lattice spacings right away and then define a second region with the same x, and y of the box but ranging from -2.5 to 2.5 in z-direction. No complex manipulation of the system needed and 2/3rd of your input commands gone.

variable Nblock equal 100000
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

This also has excessive use of variables that hinders readability like the rest of the input.

neighbor 2.0 bin
neigh_modify every 1 delay 0 check yes


fix 1 all nve/sphere

There is no need for using fix nve/sphere for this kind of system as already mentioned.

fix thermostat all langevin ${Text} ${Text} 100 1234 zero yes
fix removeMomentum all momentum 1 linear 1 1 1

Why do you need the latter if you already zero the net force? If you would also zero the net velocity once after reaching equilibrium, there should be no issue. Also for a liquid lamella, the easier, more gentle, and more physically correct way would be to apply fix spring to the system to remove a drift in z (only)


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 

Why two fix ave/time instances, you save some overhead by just using one for both since you use the same averaging settings?


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.5*lz*(v_zPress-0.5*(v_xPress+v_yPress))*${convFac}
fix st all ave/time 10 100 1000 v_st 

And another averaging with the exact same settings.

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

Using flush yes can have a negative performance impact and this is not a system that is supposed to be unstable and easily crashing so that you would have to worry about losing output.

thermo 1000

dump trj all custom 1000 Argon.trj id type x y z
dump_modify trj sort id pad 6

Same as above, sorting is extra overhead and doesn’t bring many advantages for the kind of data you write out.

timestep 1

run ${Nrun}

I really thank you for your time and advice i will take your notes and try to modify my script …thank you again

For the sake of others reading this thread later and my own curiosity, would you mind actually answering my questions? After I have given you some detailed notes, it would just be fair you show me the same courtesy, would it not?

Please also note that none of my suggestions will change the fluctuations you see. The topic of fluctuating temperatures and pressures have been discussed so many times on the mailing list, you should find the explanations you need in those archived discussions (which emphasizes the value of my first point).

i very interested in your questions and really show the courtesy but i really do not want to bother you alot

answering to your questions
1- atomic style was my choice firstly but i was trying to modify to get better results
2- this was a formal script i just modify on it to get my script so lots of commands have already been settled