input.lammps

# This LAMMPS input script simulates surface tension measurement of water
# Written by Simon Gravelle (https://simongravelle.github.io/)
# Find more scripts here: https://github.com/simongravelle/lammps-input-files
# LAMMPS tutorials for beginners: https://lammpstutorials.github.io/

variable T equal 293.15 # temperature in Kelvin
variable A2m equal 1e-10 # Angstrom to meter
variable atm2Pa equal 101325 # Atmosphere to Pascal
variable N2mN equal 1e3 # Newton to milliNewton

boundary p p p
units real
atom_style full
bond_style harmonic
angle_style harmonic
pair_style lj/cut/tip4p/long 1 2 1 1 0.1546 12.0
kspace_style pppm/tip4p 1.0e-5
pair_modify shift yes mix arithmetic

lattice sc 3.4
region box block -3 3 -3 3 -7 7
create_box 2 box &
    bond/types 1 &
    angle/types 1 &
    extra/bond/per/atom 2 &
    extra/angle/per/atom 1 &
    extra/special/per/atom 2
region fluid block INF INF INF INF -3 3

include PARM.lammps

molecule h2omol ../LAMMPS-molecules/H2O-water/tip4p-2005.mol
create_atoms 0 region fluid mol h2omol 482793

# dynamics
velocity all create ${T} 4928459 mom yes rot yes dist gaussian 
fix mynve all nve
fix myber all temp/berendsen ${T} ${T} 100
fix myshk all shake 1.0e-5 200 0 b 1 a 1
timestep 2.0

# output and run
# first : equilibrate the system at temperature T
thermo 1000
run 10000

# outputs
variable xPress equal c_thermo_press[1]
variable yPress equal c_thermo_press[2]
variable zPress equal c_thermo_press[3]
variable lz equal lz
variable myst equal 0.5*${lz}*(v_zPress-0.5*(v_xPress+v_yPress))*${A2m}*${atm2Pa}*${N2mN}
dump mydmp all atom 5000 dump.lammpstrj
compute mytemp all temp
fix myat1 all ave/time 10 500 5000 c_mytemp file temperatureVStime.dat
fix myat2 all ave/time 10 500 5000 v_myst file surfacetensionVStime.dat

# longer run for surface tension measurement
run 100000