# LAMMPS input script for standardized atomistic simulations
# Created by emc.pl v5.1, July 11, 2024 as part of EMC
# on Mon Aug 12 08:15:55 CDT 2024

# LAMMPS atomistic input script

echo		screen
units		real
atom_style	full

# Variable definitions

variable	project		index	"vlmp"		# project name
variable	source		index	../build	# data directory
variable	temperature	index	285		# system temperature
variable	tdamp		index	100		# temperature damping
variable	pressure	index	1		# system pressure
variable	pdamp		index	1000		# pressure damping
variable	dielectric	index	1		# medium dielectric
variable	kappa		index	4		# electrostatics kappa
variable	cutoff		index	9.5		# standard cutoff
variable	charge_cutoff	index	9.67		# charge cutoff
variable	precision	index	0.001		# kspace precision
variable	lseed		index	93420337	# langevin seed
variable	vseed		index	17494829	# velocity seed
variable	tequil		index	10000		# equilibration time
variable	dlimit		index	0.2		# nve/limit distance
variable	trun		index	10000000	# run time
variable	frestart	index	0		# 0: equil, 1: restart
variable	dtrestart	index	100000		# delta restart time
variable	dtdump		index	1000000		# delta dump time
variable	dtthermo	index	1000		# delta thermo time
variable	timestep	index	1		# integration time step
variable	tfreq		index	10		# profile sampling freq
variable	nsample		index	1000		# profile conf sampling
variable	dtime		equal	${tfreq}*${nsample} # profile dtime
variable	restart		index	${project}.restart

if "${frestart} != 0" then &
"variable	data		index	${restart}" &
else &
"variable	data		index	${source}/${project}.data"

# Interaction potential definition

pair_style	lj/cut/coul/long ${cutoff} ${charge_cutoff}
bond_style	harmonic
special_bonds	lj/coul 0 0 1

if "${frestart} != 0" then "read_restart ${data}" &
else "read_data ${data}"

include		${source}/${project}.params

# Integration conditions (check)

timestep	${timestep}
if "${flag_charged} != 0" then "kspace_style pppm/cg ${precision}"
dielectric	${dielectric}
fix		mom all momentum 100 linear 1 1 1 angular

# Equilibration

thermo		${dtthermo}
if "${frestart} != 0" then "jump SELF simulate"
velocity	all create ${temperature} ${vseed} &
		dist gaussian rot yes mom yes sum yes
fix		temp all langevin ${temperature} ${temperature} ${tdamp} &
		${lseed}
fix		int all nve/limit ${dlimit}
run		${tequil}
unfix		temp
unfix		int
write_restart	${project}.restart2

# Simulation

label		simulate

# Verbatim paragraph

fix		int all npt temp ${temperature} ${thigh} ${tdamp} &
		iso ${pressure} ${pressure} ${pdamp}
run		${trelax}
unfix		int

fix		int all npt temp ${thigh} ${temperature} ${tdamp} &
		iso ${pressure} ${pressure} ${pdamp}
run		${trelax}
unfix		int

write_restart	${project}.restart2

# Integrator
 
fix		int all npt temp ${temperature} ${temperature} ${tdamp} &
		iso ${pressure} ${pressure} ${pdamp}

# System sampling: energy

variable	pe equal pe
variable	ke equal ke
variable	etotal equal etotal
variable	enthalpy equal enthalpy
variable	evdwl equal evdwl
variable	ecoul equal ecoul
variable	epair equal epair
variable	ebond equal ebond
variable	eangle equal eangle
variable	edihed equal edihed
variable	eimp equal eimp
variable	emol equal emol
variable	elong equal elong
variable	etail equal etail

fix		ene all ave/time ${tfreq} ${nsample} ${dtime} &
		c_thermo_temp &
	       	v_pe v_ke v_etotal v_enthalpy v_evdwl &
		v_ecoul v_epair v_ebond v_eangle v_edihed v_eimp &
		v_emol v_elong v_etail &
		file ${project}.energy

# System sampling: pressure

fix		sample_press all ave/time ${tfreq} ${nsample} ${dtime} &
		c_thermo_temp &
		c_thermo_press[1] c_thermo_press[2] c_thermo_press[3] &
		c_thermo_press[4] c_thermo_press[5] c_thermo_press[6] &
		file ${project}.pressure

# System sampling: volume

variable	volume equal vol
variable	hxx equal lx
variable	hyy equal ly
variable	hzz equal lz
variable	hxy equal xy
variable	hxz equal xz
variable	hyz equal yz

fix		vol all ave/time ${tfreq} ${nsample} ${dtime} &
		v_volume v_hxx v_hyy v_hzz v_hyz v_hxz v_hxy &
		file ${project}.volume

# Run conditions

restart		${dtrestart} ${project}.restart1 ${project}.restart2
dump		1 all custom ${dtdump} ${project}.dump id type x y z
run		${trun}