-
Notifications
You must be signed in to change notification settings - Fork 396
Fix virial sign in NPT ensemble, wrap coordinates to unit cell, and support charge/spin configuration #1546
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: main
Are you sure you want to change the base?
Conversation
Hi @yu9824! Thank you for your pull request and welcome to our community. Action RequiredIn order to merge any pull request (code, docs, etc.), we require contributors to sign our Contributor License Agreement, and we don't seem to have one on file for you. ProcessIn order for us to review and merge your suggested changes, please sign at https://code.facebook.com/cla. If you are contributing on behalf of someone else (eg your employer), the individual CLA may not be sufficient and your employer may need to sign the corporate CLA. Once the CLA is signed, our tooling will perform checks and validations. Afterwards, the pull request will be tagged with If you have received this in error or have any questions, please contact us at [email protected]. Thanks! |
Thank you for signing our Contributor License Agreement. We can now accept your code for this (and any) Meta Open Source project. Thanks! |
It seems that the CI test stopped due to a temporary environment error (EADDRINUSE in PyTorch distributed initialization), rather than an issue with this PR’s code. |
Thanks for making this PR and contributing!, the CI tests all look good, we'll review it soon and looks mostly good! |
# virials need to be in this order: xx, yy, zz, xy, xz, yz. https://docs.lammps.org/Library_utility.html#_CPPv437lammps_fix_external_set_virial_globalPvPKcPd | ||
virial_arr = [v[0], v[4], v[8], v[1], v[2], v[5]] | ||
lmp.fix_external_set_virial_global(FIX_EXT_ID, virial_arr) | ||
def wrap_to_unit_cell( |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Can we use ASE's wrap positions instead since we already use this function in our code?
pos = wrap_positions(pos, cell, pbc=pbc, eps=0) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Agreed, we can use ASE wrap_positions inside this function. Relatedly, it would be useful to add some simple lammps box to ASE unit cell and back tests to ensure consistency between lammps and ASE representation.
# stress is defined as virial/volume in lammps | ||
assert "stress" in results, "stress must be in results to compute virial" | ||
volume = torch.det(cell).abs().item() | ||
v = (-results["stress"].detach().cpu() * volume)[0].tolist() |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks for catching this! The ASE v lammps test looks good. Can you update the comment as well: virial = - stress * volume? (in the lammps doc they didn't add the negative)
ase_kinetic, ase_pot = run_ase_npt() | ||
lammps_kinetic, lammps_pot = run_lammps("tests/lammps/lammps_npt.file") | ||
assert np.isclose(ase_kinetic, lammps_kinetic, atol=1.0) | ||
assert np.isclose(ase_pot, lammps_pot, atol=1.0) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
should we add a check for pressure as well since this is a NPT sim?
fix 1 all npt temp 300.0 300.0 0.1 iso 0.0 0.0 1.0 | ||
thermo_style custom step temp pe ke etotal press vol | ||
thermo 1 | ||
run 300 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
can we run this for 100 steps instead to reduce the test time?
Looks great, mostly just cosmetic comments, thanks for the adding the NPT example test case, very useful! |
Summary
This PR includes three main updates to improve stability and correctness in LAMMPS-based simulations.
Fix overestimated pressure in NPT ensemble
The virial term was incorrectly signed, leading to systematically higher pressure evaluations in NPT simulations.
This PR corrects the sign convention and adds a dedicated NPT test to verify the fix.
Wrap atomic positions to the unit cell
Although this normally has no effect on computed energies or forces, large atomic displacements occasionally caused PyTorch warnings due to numerical instability.
Wrapping positions mitigates such issues and improves numerical robustness.
Add charge and spin configuration support
Extended input options to allow specifying total charge and spin states, which is useful for systems that require explicit electronic parameters.
Implementation Details
Testing
Motivation
This patch improves both physical correctness (pressure evaluation) and numerical stability (coordinate wrapping).