variable r equal 1 variable a equal 2.0 shell mkdir dir$r/dump # Folder for production dump files shell mkdir dir$r/restart # Folder for restart files from production runs shell mkdir dir$r/logfile # Folder for log files shell mkdir dir$r/output # Folder for production dump files variable conc equal 6 ################################################ units real atom_style full dimension 3 newton on processors * * 2 boundary p p f #read_data dir$r/input_data/KOHfinal.txt read_restart dir$r/input_data/x.r_6.8700000 ######################################## ######### force field ########## ######################################## pair_style lj/cut/coul/long 12 bond_style harmonic angle_style harmonic kspace_style pppm/electrode 1.0e-6 kspace_modify slab 3 fix zwalls all wall/reflect zlo EDGE zhi EDGE units box # Pair Coeffs # # 2 Ow # 1 Hw # 3 K # 4 O # 5 H # 6 C ## TMH parameters atom definition EPISILON SIGMA FOR OTHERS WE USED Lorentz-Berthelot rules for mixing for epsilon(Kcal per mol ) it is GM (geometric mean) and for sigma (Amstrong) it is AM (airthmatic mean) pair_coeff 2 2 0.102 3.188 # OT OT pair_coeff 1 1 0.0 0.0 # HT HT pair_coeff 3 3 0.119423 3.0 # K K pair_coeff 4 4 0.119423 2.3 # O O pair_coeff 5 5 0.0716538 1.0 # H H pair_coeff 6 6 0.0474 3.211 # C C pair_modify shift yes mix arithmetic bond_coeff 1 450.0 0.9572 bond_coeff 2 450.0 0.9572 angle_coeff 1 55.0 104.52 #Masses mass 2 15.9994 # o mass 1 1.008 # h mass 3 39 # K mass 4 15.9994 # O mass 5 1.008 # H mass 6 12.0107 # carbon group syst type 1 2 3 4 5 6 group ele type 6 group sol type 1 2 3 4 5 group pot type 3 group hydr type 4 5 group ohydr type 4 group hhydr type 5 group h1 type 1 group o1 type 2 group water type 1 2 group elell molecule 44444 group elelr molecule 44445 group elerl molecule 66665 group elerr molecule 66666 group elel molecule 44444 44445 group eler molecule 66666 66665 ## charges ### set group o1 charge -0.830 set group h1 charge 0.415 set group ele charge 0.0000 set group pot charge 0.6791 set group ohydr charge -1.1029 set group hhydr charge 0.4238 run_style verlet neighbor 2.0 bin neigh_modify exclude group ele ele check no neigh_modify delay 0 every 1 check yes #fix bal all balance 1000 1.05 shift z 20 1.05 ########### Equilibration ######################### log dir$r/logfile/log_${conc} #reset_timestep 0 timestep 0.00005 thermo 1000 compute t_water sol temp fix 1 ele setforce 0.0 0.0 0.0 fix 2 water shake 0.0001 20 0 b 1 2 a 1 t 1 2 4 5 minimize 1.0e-4 1.0e-6 2 10000 velocity sol create 0.1 86289 dump 1 syst custom 2000 dir$r/dump/dump-${conc}.lammpstrj id mol type x y z q dump_modify 1 sort id dump 2 syst custom 200000 dir$r/dump/dumpall-${conc}.lammpstrj id mol type x y z q dump_modify 2 sort id restart 200000 dir$r/output/x.r_${conc} fix 4 sol nvt temp 0.1 300.0 100 fix_modify 4 temp t_water run 50000 unfix 4 timestep 0.0001 fix 4 sol nvt temp 300 300.0 100 fix_modify 4 temp t_water run 50000 unfix 4 timestep 0.001 fix 4 sol nvt temp 300 300.0 100 fix_modify 4 temp t_water run 10000 unfix 4 timestep 0.01 fix 4 sol nvt temp 300 300.0 100 fix_modify 4 temp t_water run 10000 unfix 4 timestep 0.1 fix 4 sol nvt temp 300 300.0 100 fix_modify 4 temp t_water run 100000 write_restart dir$r/output/min.rest unfix 4 #reset_timestep 0 timestep 1 variable w equal ${a}/2 fix e eler electrode/conp -${w} 1.979 couple elel ${w} symm on etypes on fix 4 sol nvt temp 300.0 300.0 100 fix_modify 4 temp t_water ############################################################### ########################################################################################################################################### ########################################################################################################################################### # total charge on electrode calculation variable q atom q compute ql elel reduce sum v_q compute qr eler reduce sum v_q compute qll elell reduce sum v_q compute qlr elelr reduce sum v_q compute qrl elerl reduce sum v_q compute qrr elerr reduce sum v_q ########################################################################################################################################### ########################################################################################################################################### # variable vr internal 0.0 # variable vl internal 0.0 # variable qsb_right internal 0.0 # variable qsb_left internal 0.0 # fix_modify e set v eler vr # fix_modify e set v elel vl # fix_modify e set qsb eler qsb_right # fix_modify e set qsb elel qsb_left # variable mc00 internal 0.0 # variable mc01 internal 0.0 # variable mc10 internal 0.0 # variable mc11 internal 0.0 # fix_modify e set mc elel elel mc00 # fix_modify e set mc elel eler mc01 # fix_modify e set mc eler elel mc10 # fix_modify e set mc eler eler mc11 # # calculate charge as qsb + mc * voltage # # note: fix will wait until the run setup to look for its potential variables # variable qqleft equal v_qsb_left+v_mc00*v_vl+v_mc01*v_vr # variable qqright equal v_qsb_right+v_mc10*v_vl+v_mc11*v_vr ############################################################### thermo_style custom step c_t_water c_ql c_qr c_qll c_qlr c_qrl c_qrr pe etotal ecoul press spcpu thermo_modify flush yes norm no thermo_modify temp t_water run 15000000 write_restart dir$r/output/finalcoord.rest write_data dir$r/dump/finalcoord.dat