# LAMMPS input script for water droplet on fully hydroxylated quartz in Carbon_dioxide environment slab simulation # Setup log log.fullyhydroxylated+water_droplet+CO2 echo both units real dimension 3 boundary p p f atom_style full # Potential, pair_style lj/cut/coul/long 17 bond_style harmonic angle_style harmonic pair_modify mix arithmetic kspace_style pppm 1.0e-4 kspace_modify slab 3.00 # read configuration read_data fullyhydroxylated+water_droplet+CO2_editted_data.data # Define elements charge set type 1 charge 2.100 # Si_quartz set type 2 charge -1.050 # O_quartz set type 3 charge -0.95 # O_Hydroxyl set type 4 charge 0.425 # H_Hydroxyl set type 5 charge -0.82 # O_water_droplet set type 6 charge 0.41 # H_water_droplet set type 7 charge 0.6512 # C in Carbon_dioxide set type 8 charge -0.3256 # O in Carbon_dioxide group lowerlayer type 1 2 group upperlayer type 3 4 group H2O type 5 6 group CO2 type 7 8 comm_modify cutoff 21.0 # Pair style for minimization pair_style lj/cut/coul/long 17.0 pair_coeff * * 0.0 1.0 # Minimize the system minimize 1.0e-4 1.0e-6 1000 10000 # Reset pair coefficients after minimization pair_coeff 1 1 0.0000018402 3.3019566252 # Si in quartz pair_coeff 2 2 0.1554164124 3.1655200879 # O bridging oxygen in quartz pair_coeff 3 3 0.1554164124 3.1655200879 # O in Hydroxyl pair_coeff 4 4 0.0000000000 0.0000000000 # H in Hydroxyl pair_coeff 5 5 0.1554164124 3.1655200879 # O in water droplet pair_coeff 6 6 0.0000000000 0.0000000000 # H in water droplet pair_coeff 7 7 0.0558977149 2.7570000000 # C in Carbon_dioxide pair_coeff 8 8 0.1599828400 3.0330000000 # O in Carbon_dioxide fix 1 H2O shake 0.0001 20 0 b 2 a 2 fix 2 H2O nvt temp 300.0 300.0 100.0 fix 3 upperlayer nvt temp 300.0 300.0 100.0 fix 4 lowerlayer spring/self 10.00 fix 5 lowerlayer nvt temp 300.0 300.0 100.0 fix 6 CO2 rigid/nvt molecule temp 300.00 300.00 100.0 # Apply a repulsive wall at 84.5 Å to prevent lost atoms #fix 7 all wall/harmonic zhi 84.5 1.0 84.2 7.0 #fix 7 all wall/reflect zhi 84.5 fix 7 all wall/lj126 zhi 84.5 1 2.3 2.5 # Apply a repulsive wall using wall/region #region repulsive_wall block 0.00 99.89 0.00 103.78 84.5 159.168 units box #fix 7 all wall/region repulsive_wall harmonic 1.0 84.5 6.9 velocity all create 300.00 4928459 rot yes dist gaussian # Simulation parameters timestep 2.00 thermo 1000 thermo_style custom step temp vol ke pe etotal press density lx ly lz pxx pyy pzz neighbor 2.0 bin neigh_modify every 1 delay 0 check yes # Equilibration Run dump positions1 all custom 10000 300K@fullyhydroxylated+water_droplet+CO2_equil.lammpstrj id type mol x y z dump_modify positions1 sort id run 300000 undump positions1 # Production Run dump positions2 all custom 10000 300K@fullyhydroxylated+water_droplet+CO2_prod1.lammpstrj id type mol x y z dump_modify positions2 sort id run 1500000 undump positions2 # Additional Production Run dump positions3 all custom 10000 300K@fullyhydroxylated+water_droplet+CO2_prod2.lammpstrj id type mol x y z dump_modify positions3 sort id dump dumpNVTVMD all atom 50000 dumpNVT.xyz dump nvt all atom 10000 dumpNVT.*.xyz run 50000 restart 50000 restartNVT.* undump positions3 undump dumpNVT undump dumpNVTVMD undump nvt # Density profile analysis compute density_x all chunk/atom bin/1d x 0.0 2.0 compute density_y all chunk/atom bin/1d y 0.0 2.0 compute density_z all chunk/atom bin/1d z 0.0 2.0 fix density_profile_x all ave/chunk 100 100 10000 density_x density/mass file density_profile_x.txt fix density_profile_y all ave/chunk 100 100 10000 density_y density/mass file density_profile_y.txt fix density_profile_z all ave/chunk 100 100 10000 density_z density/mass file density_profile_z.txt run 100000 write_data 300K@fullyhydroxylated+Water_droplet+CO2_simulation.data