Dear Lammps users,
I recently simulate crosslinked epoxy resin and need to analyse its bulk property based on Dreiding forcefiled. However I encounter a difficulty at the first step that the density calculated is about 10% lower than the target value. I believe the relaxation process should be fine. Somebody in the previous discussion state there may be problems of setting parameters for the improper terms. Since the contribution of improper energy in my simulation is not significant compared to other energy component, and the improper energy I calculated using lammps is consistent with MS results (actually all the energy components I calculated are consistent with the MS restuls),I dont think it is the main reason of leading to low density. Another issue is about the Van Der Wals parameters. Mayo’s paper gives the van der Waals bond length ,R0, and I think it should corresponds to the equilibrium position in the plot of van der Waals energy. In Lammps, the van der Waals distance does not corresponds to the equilibrium position, but corresponds to the cross-section point (between the energy cure and x axis) where the van der Waals energy is zero, so the van der Waals length in Lammps should be shorter than RO in Mayo’s paper, which should be equal to R0/(2^(1/6)). The value of 2^(1/6) should be around 1.1225. Therefore I divide all R0 by 1.1225 and then input all the parameters in the pair_coeff command in the lammps input file.
One thing is, in another data file written by my friend, he used van der Waals length even shorter than I used. The length he use is about R0/(2^(1/3)). I believe that if using his parameters the density surely will be higher, but I can not see any reason why he choose R0/(2^(1/3)).
Below is my input file (please notice the atom types and help to check the pair_coeff), anyone could see any problem and teach me how to fix the low-density issue ?
units real
atom_style full
dimension 3
boundary p p p
Style
pair_style lj/cut/coul/long 12
pair_modify mix arithmetic
bond_style harmonic
angle_style harmonic
dihedral_style charmm
improper_style umbrella
kspace_style ewald 1.0e-3
special_bonds dreiding
Atom Definition
read_data in_cro_933p_0char.lmps
#Coeff
pair_coeff 1 1 0.0951 3.47299 # C_R
pair_coeff 2 2 0.0951 3.47299 # C_3
pair_coeff 3 3 0.0957 3.03315 # O_R
pair_coeff 4 4 0.0951 3.47299 # C_31
pair_coeff 5 5 0.0957 3.03315 # O_3
pair_coeff 6 6 0.0152 2.84642 # H_
pair_coeff 7 7 0.0001 2.84642 # H__A
pair_coeff 8 8 0.0774 3.2626 # N_R
Settings
dielectric 1.0
neighbor 2.0 bin
neigh_modify delay 0 every 1 check yes
timestep 1.0
Output
thermo_style custom step density vol temp press etotal pe ke evdwl ecoul ebond eangle edihed eimp
thermo 1000
1
fix 1 all nvt temp 600 600 100
velocity all create 600 8531416
compute 1 all msd
fix 11 all ave/time 5 20 1000 c_1[1] c_1[2] c_1[3] c_1[4] file msd_2ns.data
run 100000
unfix 1
fix 1 all npt temp 600 600 100 iso 1 1 100
velocity all scale 600
run 400000
unfix 1
write_restart relaxed_2ns_1.restart
2
fix 1 all nvt temp 600 600 100
velocity all scale 600
run 100000
unfix 1
fix 1 all npt temp 600 600 100 iso 1 1 100
velocity all scale 600
run 400000
unfix 1
write_restart relaxed_2ns_2.restart
3
fix 1 all nvt temp 600 600 100
velocity all scale 600
run 100000
unfix 1
fix 1 all npt temp 600 600 100 iso 1 1 100
velocity all scale 600
run 400000
unfix 1
write_restart relaxed_2ns_3.restart
4
fix 1 all nvt temp 600 600 100
velocity all scale 600
run 100000
unfix 1
fix 1 all npt temp 600 600 100 iso 1 1 100
velocity all scale 600
run 400000
unfix 1
write_restart relaxed_2ns_final.restart
write_data relaxed_2ns.data