boundary f f p 
units real
atom_style full
bond_style harmonic
angle_style harmonic
pair_style lj/cut/coul/long 15
kspace_style pppm/electrode 1e-7
kspace_modify wire 3 



read_data ./min.data

variable unit_length equal 22.2404
variable unit_water_num equal 1000
variable unit_ion_num equal 10
variable unit_rep equal 1


variable elect_dist equal 16.28
variable Lx equal ${unit_rep}*${unit_length}
variable displace equal ${elect_dist}+${Lx}
variable Ltop equal ${elect_dist}+${Lx}
variable Lend equal ${Ltop}+${elect_dist}
variable Llimit equal ${Lend}+1


group electrolyte type 1 2 3 4
group water type 1 2
group ions type 3 4
group na type 3
group cl type 4
group bot type 5 
group top type 6

variable electrode_x0 equal bound(bot,xmax)
variable electrode_x1 equal bound(top,xmin)

region water_region block ${electrode_x0} ${electrode_x1} EDGE EDGE EDGE EDGE

compute ctemp electrolyte temp

compute layers electrolyte chunk/atom bin/1d x lower 0.04 bound x ${electrode_x0} ${electrode_x1}
compute water_layers water chunk/atom bin/1d x lower 0.04 bound x ${electrode_x0} ${electrode_x1} 
compute ion_layers ions chunk/atom bin/1d x lower 0.04 bound x ${electrode_x0} ${electrode_x1} 

compute all_dipoles electrolyte dipole mass
compute water_dipoles water dipole mass
compute ion_dipoles ions dipole mass
 
variable q atom q
variable qx atom q*(x-lx/2)
compute qtop top reduce sum v_q
compute qbot bot reduce sum v_q
compute qxtop top reduce sum v_qx
compute qxbot bot reduce sum v_qx



fix density_profile electrolyte ave/chunk 100 10 1000 layers density/mass density/number norm sample file nvt.density
fix water_profile water ave/chunk 100 10 1000 layers density/mass density/number norm sample file water.density
fix na_profile na ave/chunk 100 10 1000 layers density/mass density/number norm sample file na.density
fix cl_profile cl ave/chunk 100 10 1000 layers density/mass density/number norm sample file cl.density

fix mywall electrolyte wall/reflect ylo EDGE yhi EDGE
fix mynvt electrolyte nvt temp 298 298 100

variable bot_volt equal -1*4.0
variable top_volt equal 4.0

fix conp bot electrode/conp ${bot_volt} 1.805132 couple top ${top_volt} symm on etypes on

fix ad electrolyte ave/time 100 1 100 c_all_dipoles file all.dipole mode vector
fix wd water ave/time 100 1 100 c_water_dipoles file water.dipole mode vector
fix id ions ave/time 100 1 100 c_ion_dipoles file ion.dipole mode vector

thermo 1000
thermo_style custom step temp etotal epair ke pe evdwl ecoul elong c_qbot c_qtop c_qxbot c_qxtop
thermo_modify temp ctemp
timestep 0.5
variable etotal equal etotal
variable epair equal epair
variable ke equal ke
variable pe equal pe
variable evdwl equal evdwl
variable ecoul equal ecoul
variable elong equal elong

fix info_nvt all ave/time 100 1 100 c_ctemp v_etotal v_epair v_ke v_pe v_evdwl v_ecoul v_elong c_qbot c_qtop c_qxbot c_qxtop file nvt.info

dump nvtdump all custom 1000 nvt_298K.lammpstrj id mol type q x y z vx vy vz

run 2000000 

write_data nvt_298K.data