# sample LAMMPS input script for viscosity of 2d LJ liquid
# Muller-Plathe via fix viscosity

# settings

variable	x equal 20
variable	y equal 20

variable	rho equal 0.6
variable        t equal 1.0
variable	rc equal 2.5

# problem setup

units		lj
dimension	2
atom_style	atomic
neigh_modify	delay 0 every 1

lattice         sq2 ${rho}
region          simbox block 0 $x 0 $y -0.1 0.1
create_box      1 simbox
create_atoms    1 box

pair_style      lj/cut ${rc}
pair_coeff      * * 1 1

mass            * 1.0
velocity        all create $t 97287

fix             1 all nve
fix	        2 all langevin $t $t 0.1 498094
fix	        3 all enforce2d

# equilibration run

thermo          1000
run	        5000

unfix		2

# turn on Muller-Plathe driving force and equilibrate some more

velocity	all scale $t

fix		4 all viscosity 100 x y 20

compute         layers all chunk/atom bin/1d y center 0.05 units reduced
fix		5 all ave/chunk 20 50 1000 layers vx file profile.mp.2d

# equilibration run

variable	dVx equal f_5[11][3]-f_5[1][3]

thermo          1000
thermo_style	custom step temp epair etotal press f_4 v_dVx
run		20000

# data gathering run
# reset fix viscosity to zero flux accumulation

unfix		4
fix		4 all viscosity 100 x y 20

variable	visc equal -(f_4/(2*(step*0.005-125)*lx+1.0e-10))/(v_dVx/(ly/2))
fix		vave all ave/time 1000 1 1000 v_visc ave running start 26000

thermo_style	custom step temp f_4 v_dVx v_visc f_vave

# only need to run for 10000 steps to make a good 100-frame movie

#dump	        1 all custom 100 dump.mp.2d id type x y z vx

#dump		2 all image 100 image.*.jpg vx type zoom 1.6 adiam 1.2
#dump_modify	2 pad 5 amap -0.7 0.7 ca 0.0 2 min blue max red

run	        50000
