Skip to content

Hack in bonded interactions #236

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.

Sign up for GitHub

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

Draft
wants to merge 14 commits into
base: main
Choose a base branch
from
Draft

Hack in bonded interactions #236

Changes from 2 commits
Commits
Show all changes
14 commits
Select commit Hold shift + click to select a range
047da30
Filter neighbor list.
clpetix Feb 10, 2024
b1ba252
Make bonds a private variable of the initializer.
clpetix Feb 12, 2024
6e2694e
Change neighborlist method to calculate from a full neighborlist.
clpetix Apr 20, 2024
1e9b4c9
Simplify RDF calculations, including with bond exclusions
mphoward Apr 22, 2024
7fbecd5
Histogram and normalize RDFs manually
mphoward Apr 29, 2024
ccfff2d
Set number of rdf origins.
clpetix Apr 29, 2024
f8bd8cc
Merge branch 'main' into feature/bonds
mphoward May 3, 2024
3d326dc
Filter bond exclusions using Cantor pairing function.
clpetix Jun 6, 2024
fad8cd0
Merge branch 'main' of https://github.com/mphowardlab/relentless into…
clpetix Jul 15, 2024
3a0ca1b
Merge branch 'main' into feature/bonds
clpetix Aug 19, 2024
09b3ebd
Merge branch 'main' into feature/bonds
clpetix Aug 22, 2024
1d19b93
Add bond potentials and exclusions.
clpetix Oct 4, 2024
d05147b
Fix error in naming.
clpetix Oct 4, 2024
0a9ad40
Add angle exclusions.
clpetix Oct 9, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
34 changes: 31 additions & 3 deletions src/relentless/simulate/hoomd.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -70,8 +70,8 @@ def _call_v3(self, sim):
# parse masses by type
snap = sim["engine"]["_hoomd"].state.get_snapshot()
sim.masses = self._get_masses_from_snapshot(sim, snap)
sim[self]["_bonds"] = self._get_bonds_from_snapshot(sim, snap)
self._assert_dimension_safe(sim, snap)

# create the potentials, defer attaching until later
neighbor_list = hoomd.md.nlist.Tree(buffer=sim.potentials.pair.neighbor_buffer)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

One thing we do still want to do is make sure we can specify the neighbor list exclusions, and then only do the exclusions in the thermo if they are active. I think we talked about adding that as an argument / property of the PairPotentialTabulator, like the neighbor_buffer, and then we would set it here by setting the exclusions option.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

add in exclusions=('bond', 'angle')

pair_potential = hoomd.md.pair.Table(nlist=neighbor_list)
Expand Down Expand Up @@ -159,6 +159,9 @@ def _get_masses_from_snapshot(self, sim, snap):
masses.update(masses_)
return masses

def _get_bonds_from_snapshot(self, sim, snap):
return snap.bonds.group

def _assert_dimension_safe(self, sim, snap):
if sim.dimension == 3:
dim_safe = True
Expand Down Expand Up @@ -953,6 +956,7 @@ def _pre_run_v3(self, sim, sim_op):
system=sim["engine"]["_hoomd"].state,
rdf_params=self._get_rdf_params(sim),
constraints=self._get_constrained_quantities(sim, sim_op),
bonds=sim[sim.initializer]["_bonds"],
)
sim[self]["_hoomd_thermo_callback"] = hoomd.write.CustomWriter(
trigger=self.every,
Expand Down Expand Up @@ -1137,7 +1141,14 @@ class EnsembleAverageAction(Action):
flags = [Action.Flags.PRESSURE_TENSOR]

def __init__(
self, types, dimension, thermo, system, rdf_params=None, constraints=None
self,
types,
dimension,
thermo,
system,
rdf_params=None,
constraints=None,
bonds=None,
):
if dimension not in (2, 3):
raise ValueError("Only 2 or 3 dimensions supported")
Expand All @@ -1148,6 +1159,7 @@ def __init__(
self.system = system
self.rdf_params = rdf_params
self.constraints = constraints if constraints is not None else {}
self.bonds = bonds

# this method handles all the initialization
self.reset()
Expand Down Expand Up @@ -1228,11 +1240,27 @@ def act(self, timestep):
for i, j in self._rdf:
query_args = dict(self._rdf[i, j].default_query_args)
query_args.update(exclude_ii=(i == j))
neighbors = (
aabbs[j]
.query(
snap.particles.position[type_masks[i]], query_args
)
.toNeighborList()
)
if snap.bonds is not None:
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Here you would add a check if bonds are in the exclusion list as well. Probably you need to have that forwarded to this action though, since we would normally grab it from the sim.

bonds = numpy.vstack(
[self.bonds, numpy.flip(self.bonds, axis=1)]
)
# list intersect method from:
# https://stackoverflow.com/a/67113105
filter = ~(neighbors[:, None] == bonds).all(-1).any(1)
else:
filter = numpy.full((len(neighbors[:])), True)
# resetting when the samples are zero clears the RDF
self._rdf[i, j].compute(
aabbs[j],
snap.particles.position[type_masks[i]],
neighbors=query_args,
neighbors=neighbors.filter(filter),
reset=(self.num_samples == 0),
)

Expand Down