#porous electrodes
units  real
if "$(extract_setting(world_size) % 2) == 0" then "processors * * 2"
atom_style full
dimension   3
boundary p p f


#read_data initialsys5.config
read_restart restart.min

log logfileresult.txt
# 1  A
# 2  Cl
# 3  D
# 4  H
# 5  K
# 6  O


group electrolytes type 2 4 5 6
group electrode type 1 3 
group nion type 2
group pion type 5

group water type  4 6 
group bot type 1  
group top type 3 



bond_style harmonic
bond_coeff 1 100.0 1.0
angle_style harmonic
angle_coeff 1 10.0 109.47

kspace_style pppm/electrode 1e-7
kspace_modify   slab 3.0

pair_style 	lj/cut/coul/long 14
# 1  A
# 2  Cl
# 3  D
# 4  H
# 5  K
# 6  O
pair_modify	shift yes  
pair_coeff	1	1	0.11711	1
pair_coeff	1	2	0.1434	2.915
pair_coeff	1	3	0.11711	0.5
pair_coeff	1	4	0	0.5
pair_coeff	1	5	0.5258	1.92
pair_coeff	1	6	0.134881572	2.083
pair_coeff	2	2	0.01291	4.83
pair_coeff	2	3	0.1434	2.915
pair_coeff	2	4	0	2.415
pair_coeff	2	5	0.074483651	3.835
pair_coeff	2	6	0.044783574	3.998
pair_coeff	3	3	0.11711	0.5
pair_coeff	3	4	0	0.5
pair_coeff	3	5	0.5258	1.92
pair_coeff	3	6	0.134881572	2.083
pair_coeff	4	4	0	0
pair_coeff	4	5	0	1.42
pair_coeff	4	6	0	1.583
pair_coeff	5	5	0.042973	2.84
pair_coeff	5	6	0.25837677	3.003
pair_coeff	6	6	0.15535	3.166


#minimize 1.0e-4 1.0e-6 1000 1000

fix conp bot electrode/conp -1.0 1.979 couple top 1.0 symm on #algo mat_inv

fix 1 electrode setforce 0.0 0.0 0.0
fix 2 electrolytes nvt temp 300.0 300.0 100
fix 3 water shake 1e-4 100 0 b 1 a 1


neighbor	2.0 bin
neigh_modify delay 0 every 1 check yes page 500000 one 50000

velocity electrolytes create 300 4928459 rot yes dist gaussian
variable q atom q
compute qtop top reduce sum v_q
compute qbot bot reduce sum v_q
compute ctemp electrolytes temp


timestep	2
thermo		1000
thermo_style custom step temp pe ke press c_ctemp epair etotal c_qbot c_qtop v_q spcpu

thermo_modify   flush yes


dump 1 all custom 5000 fullunwrap.lammpstrj id mol type q xu yu zu 
dump_modify 1 first yes

run 5000000

#run 5000000
#write_restart  restart.min
write_restart  restartafter10M.min
print "Done"