Skip to content

HOOMD/Topology Charge floating point error #943

@CalCraven

Description

@CalCraven

Bug summary

Create a hoomd gsd file for use in hoomd.

Code to reproduce the behavior

Please include a code snippet that can be used to reproduce this bug.

import mbuild as mb
import gmso
from gmso.parameterization import apply
import numpy as np
import unyt as u
import hoomd

mol = mb.load('CC(C)C(=O)OC', smiles=True)
top = mol.to_gmso()
ff = gmso.ForceField("oplsaa")

ptop = apply(top, ff)

print(sum([site.charge for site in ptop.sites]))
print(np.sum([site.charge for site in ptop.sites]))
print(np.sum([site.charge.in_units("elementary_charge") for site in ptop.sites]))
print(np.sum([site.charge.round(2) for site in ptop.sites]))

from gmso.external import to_hoomd_snapshot
base_units = {
    "mass": u.g / u.mol,
    "length": u.nm,
    "energy": u.kJ / u.mol,
}

snapshot, _ = to_hoomd_snapshot(ptop, base_units=base_units)
print(sum(snapshot.particles.charge))

device = hoomd.device.CPU()
sim = hoomd.Simulation(device=device, seed=1)
sim.create_state_from_snapshot(snapshot)
hoomd.write.GSD.write(
    state=sim.state, filename="charge.gsd", mode="wb"
)
with gsd.hoomd.open(name="charge.gsd", mode='r') as traj:
    # Read a specific frame (e.g., frame 0 for the initial state)
    # You can specify any frame index you need
    snapshot = traj[-1]
print(sum(snapshot.particles.charge))
8.326672684688674e-17 elementary_charge
8.326672684688674e-17
8.326672684688674e-17
8.326672684688674e-17
-2.220446049250313e-16
2.3841858e-07

As you can see the charge builds up from very small floating point issues, which I believe is compounded during the unit nondimensionalizing stage in hoomd.

Software versions

  • Which version of GMSO are you using? (python -c "import gmso; print(gmso.__version__)")
    0.14.0
  • Which version of Python (python --version)?
    3.12
  • Which operating system?
    MacOS ARM

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions