Dear developers and users of Lammps, my pace potential works poorly on stress prediction.
Potential Creation
I have created a database of 1000 DFT configurations for Cu to train an ACE potential. I’ve got 8 atoms per configurations so this test simulation is fast. After the fitting, I get a RMSE of 1.1 meV/at on energy and 3.6 meV/ang on forces which seems to indicate this ACE potential is able to represent the thermodynamical state pretty well.
The problem
Then, I’ve done a lammps calculation on the relaxed Cu system with 0 velocity on every atom and 0 step using this PACE potential. I’ve noticed the pressure is -80.45 GPa while the DFT pressure is 1.58 GPa. I find it strange that the pressure is that high since the highest DFT force on an atom is 1.6e-5 eV/ang and the pace potential predicts even lower forces. For the 3000 pressures, it varies between 80 and 88 GPa as if there was a constant to remove.
What I would expect
Although the effective force of an atom isn’t directly related to the stress, I would expect the stress to be closer to the DFT result as the pace potential seems to describe the interaction correctly.
Questions
- Is the pressure expected to works well using pace as there is only fitting on energy and forces ?
- Is there a way to output the different quantities used to compute the virial pressure for non-linear potential ?
- How could I get more insight on pace pressure ? Note that I could not find much information about pace pressure in the litterature.
- Is there a documentation about how Lammps compute virial pressure for non-linear potential ?
Additional Information
- As the pace potential is anharmonic, I cannot compute pair interaction to get dx, dy, dz, fx, fy, fz and manually compute the virial stress.
- I’ve substracted the energy of a free Cu atoms from the energy in my database
- The high pressure comes from the pressure pair term and the contribution seems equal from every atom.
- I’ve been careful with the units, as I’m using metal units which means pressure in bar.
- I don’t think the 8 atoms per configuration are the problem as I’m looking at a configuration
- I’m using the version of Lammps from the 2 Aug. 2023.
A bit of code
Here are the two line I use to get the virial pressure. This pressure is the same as the pressure in the log file (thermo) since the velocity are set to 0.
compute svirial all pressure NULL virial
fix svirial all ave/time 1 1 1 c_svirial file svirial.out mode vector
Thank you in advance