COMB3/qeq - bulk system wildly translating through simulation box


I’m observing a strange phenomenon occurring with COMB3 and the charge equilibration scheme. When performing dynamics in the NPT or NVT ensembles with a system composed of charged atoms, at one point, the whole system will start to wildly translate through the periodic boundaries of the simulation box. Here is an animation (the system is a bulk (Cu2O)250), it should be much more explicit this way: (might take a while to load) A couple of remarks about this: - LAMMPS version 1 Feb 2014c, also tried it with the svn build - I’m using the parametrization of COMB3 included in the LAMMPS distribution - polarization is off - I’ve observed this on a lot of different systems, using COMB3 - the problem doesn’t appear with fixed non-zero charges. So this is definitely related to qeq - it strongly affects the potential energy of the system, as you can see on the graph below - if you reduce the integration timestep, you will delay the phenomenon - freezing one atom of the system cancels the translation I’ve attached input files to reproduce the phenomenon. The test run is a 150ps-long NPT simulation at 300K, 0GPa, with qeq every 10 iterations, a timestep of 1fs, and default Tdamp/Pdamp values. The system is a bulk of (Cu2O)250. You can see using the outputed trajectory that the translation movement begins at t~78ps. structure.dat file is the Cu2O unit cell, a supercell is built at the beginning of the run. The run time is of course non-negligible, but I’m afraid it’s the simplest way to observe this strange translation. This is because: - COMB3 together with qeq is rather computationally expensive - you need quite large simulation cells to have converged values of energies and charges - the phenomenon occurs only once the system is equilibrated, ie after several tens of ps, depending on the system nature/size. Here are some graphs concerning the evolution of the total potential energy, the total kinetic energy, the mean charge on Cu atoms, and the sum of all charges of the system (which should ofc be very close to zero, ie charge neutrality of the system). The effect is rather impressive on the total potential energy (system translation at t>78ps): Thanks, Arthur

input.dat (2.52 KB)

structure.dat (572 Bytes)

I think Ray is the person to look at this (CCd).