######################################################################################
variable	t equal 353
variable 	p equal 275.7903
variable	pp equal 4000
variable	chemu equal -0.3329
variable 	seed equal 31798.0000 
######################################################################################
units		metal
atom_style      atomic
boundary 	f p p
atom_modify	map array

pair_style	lj/cut 14.0
pair_modify 	mix arithmetic
######################################################################################
read_data	test_data.txt
#read_restart 	ini_equi.restart 

group		wall1 type 1
group		wall2 type 2
group 		gas type 3

variable        rxlo equal xlo+0.0001
variable        rxhi equal xhi-0.0001
variable        rylo equal ylo+0.0001
variable        ryhi equal yhi-0.0001
variable        rzlo equal zlo+0.0001
variable        rzhi equal zhi-0.0001

region          bulk block ${rxlo} 49.999 ${rylo} ${ryhi} ${rzlo} ${rzhi} units box side in

group		wall union wall1 wall2

neighbor	1.5 bin 
neigh_modify	every 1 delay 0
compute  	mdtemp gas temp

compute         allperatom gas stress/atom NULL
compute         allp gas reduce sum c_allperatom[1] c_allperatom[2] c_allperatom[3]
variable        allpress equal -(c_allp[1]+c_allp[2]+c_allp[3])/(3*vol)

velocity        gas create $t ${seed} mom yes rot yes dist gaussian
reset_timestep	0

fix 		wallleft wall1  wall/reflect xhi 9.99999999748980e-05 units box

######################################################################################
fix		00 gas nvt temp $t $t 25
fix_modify	00 temp mdtemp

fix             record_1 gas ave/time 1 10 100 c_mdtemp v_allpress ave one file 1.realtime

dump            dump11vt all atom 50 vmd1.lammpstrj
dump_modify     dump11vt sort id

timestep 	0.001
run		5000

unfix		00
unfix		record_1

write_restart   ini_equi.restart
reset_timestep  0
######################################################################################
fix             1 gas nvt temp $t $t 25
fix_modify      1 temp mdtemp

variable	velpiston equal 1.0
velocity 	wall2 set v_velpiston 0.0 0.0 units box
fix             enepiston wall2 nve
fix             setforce wall2 setforce 0.0 0.0 0.0

variable	ttime equal time
variable	piston equal v_velpiston*(v_ttime-0.001*5000)+49.999
variable	bulkvol equal (v_piston-${rxlo})*(${ryhi}-${rylo})*(${rzhi}-${rzlo})

region          newbulk block ${rxlo} ${piston} ${rylo} ${ryhi} ${rzlo} ${rzhi} units box side in

variable	pistonad equal v_piston+0.0001
fix 		wallright wall2 wall/reflect xlo v_pistonad units box 

group		bulkgas dynamic gas region newbulk every 1

variable 	ngas equal count(gas)
variable 	nbulk equal count(bulkgas)

compute         allperatom2 bulkgas stress/atom NULL
compute         allp2 bulkgas reduce sum c_allperatom2[1] c_allperatom2[2] c_allperatom2[3]
variable        allpress2 equal -(c_allp2[1]+c_allp2[2]+c_allp2[3])/(3*v_bulkvol)

compute		locpis wall2 property/atom x
compute		locpis2 wall2 reduce ave c_locpis
variable        locpis3 equal c_locpis2

compute		velpis wall2 property/atom vx
compute		velpis2 wall2 reduce ave c_velpis 
variable        velpis3 equal c_velpis2

fix             record_3 gas ave/time 1 1 100 c_mdtemp v_allpress v_allpress2 v_ngas v_nbulk v_locpis3 v_piston ave one file press.realtime
#fix             record_4 gas ave/time 1 1000 1000 c_mdtemp v_allpress v_allpress2 v_ngas v_nbulk v_locpis3 ave one file press.final

timestep        0.001
run             5000