From 28696ccba45b132aa6519210c4b280dcdb7e5b6a Mon Sep 17 00:00:00 2001 From: clpetix Date: Wed, 7 Feb 2024 11:32:17 -0600 Subject: [PATCH 01/13] Add volume scaling for Langevin ynamics. --- src/relentless/simulate/hoomd.py | 44 +++++++++++++++++++++++++++++++- src/relentless/simulate/md.py | 4 ++- 2 files changed, 46 insertions(+), 2 deletions(-) diff --git a/src/relentless/simulate/hoomd.py b/src/relentless/simulate/hoomd.py index 276e587e..97e59cca 100644 --- a/src/relentless/simulate/hoomd.py +++ b/src/relentless/simulate/hoomd.py @@ -656,9 +656,19 @@ class RunLangevinDynamics(_Integrator): """ - def __init__(self, steps, timestep, T, friction, seed, analyzers): + def __init__( + self, + steps, + timestep, + T, + barostat, + friction, + seed, + analyzers, + ): super().__init__(steps, timestep, analyzers) self.T = T + self.barostat = barostat self.friction = friction self.seed = seed @@ -709,6 +719,38 @@ def _call_v2(self, sim): for analyzer in self.analyzers: analyzer.pre_run(sim, self) + if isinstance(self.barostat, extent.Extent): + if sim.dimension == 2: + Lx, Ly, xy = self.barostat.as_array("HOOMD") + Lz, xz, yz = 0.0, 0.0, 0.0 + else: + Lx, Ly, Lz, xy, xz, yz = self.barostat.as_array("HOOMD") + elif all(isinstance(n, extent.Extent) for n in self.barostat): + if sim.dimension == 2: + Lx1, Ly1, xy1 = self.barostat[0].as_array("HOOMD") + Lz1, xz1, yz1 = 0.0, 0.0, 0.0 + + Lx2, Ly2, xy2 = self.barostat[1].as_array("HOOMD") + Lz2, xz2, yz2 = 0.0, 0.0, 0.0 + + Lz, xz, yz = 0.0, 0.0, 0.0 + Lx = hoomd.variant.linear_interp([(0, Lx1), (self.steps, Lx2)]) + Ly = hoomd.variant.linear_interp([(0, Ly1), (self.steps, Ly2)]) + xy = hoomd.variant.linear_interp([(0, xy1), (self.steps, xy2)]) + else: + Lx1, Ly1, Lz1, xy1, xz1, yz1 = self.barostat[0].as_array("HOOMD") + Lx2, Ly2, Lz2, xy2, xz2, yz2 = self.barostat[1].as_array("HOOMD") + + Lx = hoomd.variant.linear_interp([(0, Lx1), (self.steps, Lx2)]) + Ly = hoomd.variant.linear_interp([(0, Ly1), (self.steps, Ly2)]) + Lz = hoomd.variant.linear_interp([(0, Lz1), (self.steps, Lz2)]) + xy = hoomd.variant.linear_interp([(0, xy1), (self.steps, xy2)]) + xz = hoomd.variant.linear_interp([(0, xz1), (self.steps, xz2)]) + yz = hoomd.variant.linear_interp([(0, yz1), (self.steps, yz2)]) + + hoomd.update.box_resize( + Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, scale_particles=False + ) hoomd.run(self.steps) for analyzer in self.analyzers: diff --git a/src/relentless/simulate/md.py b/src/relentless/simulate/md.py index df660474..27f6c6a7 100644 --- a/src/relentless/simulate/md.py +++ b/src/relentless/simulate/md.py @@ -118,9 +118,10 @@ class RunLangevinDynamics(_Integrator): """ - def __init__(self, steps, timestep, T, friction, seed, analyzers=None): + def __init__(self, steps, timestep, T, barostat, friction, seed, analyzers=None): super().__init__(steps, timestep, analyzers) self.T = T + self.barostat = barostat self.friction = friction self.seed = seed @@ -130,6 +131,7 @@ def _make_delegate(self, sim): steps=self.steps, timestep=self.timestep, T=self.T, + barostat=self.barostat, friction=self.friction, seed=self.seed, analyzers=self.analyzers, From 1a846abcdd8f637ca2cc290d66542b3d0dace185 Mon Sep 17 00:00:00 2001 From: clpetix Date: Wed, 7 Feb 2024 13:52:50 -0600 Subject: [PATCH 02/13] Update implementation to match HOOMD required syntax. --- src/relentless/simulate/hoomd.py | 41 ++++++++++++-------------------- src/relentless/simulate/md.py | 4 +++- 2 files changed, 18 insertions(+), 27 deletions(-) diff --git a/src/relentless/simulate/hoomd.py b/src/relentless/simulate/hoomd.py index 97e59cca..010d9004 100644 --- a/src/relentless/simulate/hoomd.py +++ b/src/relentless/simulate/hoomd.py @@ -656,16 +656,7 @@ class RunLangevinDynamics(_Integrator): """ - def __init__( - self, - steps, - timestep, - T, - barostat, - friction, - seed, - analyzers, - ): + def __init__(self, steps, timestep, T, friction, seed, analyzers, barostat=None): super().__init__(steps, timestep, analyzers) self.T = T self.barostat = barostat @@ -722,35 +713,33 @@ def _call_v2(self, sim): if isinstance(self.barostat, extent.Extent): if sim.dimension == 2: Lx, Ly, xy = self.barostat.as_array("HOOMD") - Lz, xz, yz = 0.0, 0.0, 0.0 + Lz, xz, yz = None, None, None else: Lx, Ly, Lz, xy, xz, yz = self.barostat.as_array("HOOMD") elif all(isinstance(n, extent.Extent) for n in self.barostat): if sim.dimension == 2: Lx1, Ly1, xy1 = self.barostat[0].as_array("HOOMD") - Lz1, xz1, yz1 = 0.0, 0.0, 0.0 + Lz1, xz1, yz1 = None, None, None Lx2, Ly2, xy2 = self.barostat[1].as_array("HOOMD") - Lz2, xz2, yz2 = 0.0, 0.0, 0.0 + Lz2, xz2, yz2 = None, None, None - Lz, xz, yz = 0.0, 0.0, 0.0 - Lx = hoomd.variant.linear_interp([(0, Lx1), (self.steps, Lx2)]) - Ly = hoomd.variant.linear_interp([(0, Ly1), (self.steps, Ly2)]) - xy = hoomd.variant.linear_interp([(0, xy1), (self.steps, xy2)]) + Lz, xz, yz = None, None, None + Lx = hoomd.variant.linear_interp([(0, Lx1), (self.steps - 1, Lx2)]) + Ly = hoomd.variant.linear_interp([(0, Ly1), (self.steps - 1, Ly2)]) + xy = hoomd.variant.linear_interp([(0, xy1), (self.steps - 1, xy2)]) else: Lx1, Ly1, Lz1, xy1, xz1, yz1 = self.barostat[0].as_array("HOOMD") Lx2, Ly2, Lz2, xy2, xz2, yz2 = self.barostat[1].as_array("HOOMD") - Lx = hoomd.variant.linear_interp([(0, Lx1), (self.steps, Lx2)]) - Ly = hoomd.variant.linear_interp([(0, Ly1), (self.steps, Ly2)]) - Lz = hoomd.variant.linear_interp([(0, Lz1), (self.steps, Lz2)]) - xy = hoomd.variant.linear_interp([(0, xy1), (self.steps, xy2)]) - xz = hoomd.variant.linear_interp([(0, xz1), (self.steps, xz2)]) - yz = hoomd.variant.linear_interp([(0, yz1), (self.steps, yz2)]) + Lx = hoomd.variant.linear_interp([(0, Lx1), (self.steps - 1, Lx2)]) + Ly = hoomd.variant.linear_interp([(0, Ly1), (self.steps - 1, Ly2)]) + Lz = hoomd.variant.linear_interp([(0, Lz1), (self.steps - 1, Lz2)]) + xy = hoomd.variant.linear_interp([(0, xy1), (self.steps - 1, xy2)]) + xz = hoomd.variant.linear_interp([(0, xz1), (self.steps - 1, xz2)]) + yz = hoomd.variant.linear_interp([(0, yz1), (self.steps - 1, yz2)]) - hoomd.update.box_resize( - Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, scale_particles=False - ) + hoomd.update.box_resize(Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, phase=-1) hoomd.run(self.steps) for analyzer in self.analyzers: diff --git a/src/relentless/simulate/md.py b/src/relentless/simulate/md.py index 27f6c6a7..b761175e 100644 --- a/src/relentless/simulate/md.py +++ b/src/relentless/simulate/md.py @@ -118,7 +118,9 @@ class RunLangevinDynamics(_Integrator): """ - def __init__(self, steps, timestep, T, barostat, friction, seed, analyzers=None): + def __init__( + self, steps, timestep, T, friction, seed, barostat=None, analyzers=None + ): super().__init__(steps, timestep, analyzers) self.T = T self.barostat = barostat From 52a7edb64e1a1d392a4bc7e4c9a101e4e7290b3c Mon Sep 17 00:00:00 2001 From: clpetix Date: Wed, 7 Feb 2024 14:16:20 -0600 Subject: [PATCH 03/13] Update period to prevent triggering repeatedly. --- src/relentless/simulate/hoomd.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/src/relentless/simulate/hoomd.py b/src/relentless/simulate/hoomd.py index 010d9004..98f78001 100644 --- a/src/relentless/simulate/hoomd.py +++ b/src/relentless/simulate/hoomd.py @@ -711,12 +711,14 @@ def _call_v2(self, sim): analyzer.pre_run(sim, self) if isinstance(self.barostat, extent.Extent): + period = None if sim.dimension == 2: Lx, Ly, xy = self.barostat.as_array("HOOMD") Lz, xz, yz = None, None, None else: Lx, Ly, Lz, xy, xz, yz = self.barostat.as_array("HOOMD") elif all(isinstance(n, extent.Extent) for n in self.barostat): + period = 1 if sim.dimension == 2: Lx1, Ly1, xy1 = self.barostat[0].as_array("HOOMD") Lz1, xz1, yz1 = None, None, None @@ -739,7 +741,9 @@ def _call_v2(self, sim): xz = hoomd.variant.linear_interp([(0, xz1), (self.steps - 1, xz2)]) yz = hoomd.variant.linear_interp([(0, yz1), (self.steps - 1, yz2)]) - hoomd.update.box_resize(Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, phase=-1) + hoomd.update.box_resize( + Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, period=period + ) hoomd.run(self.steps) for analyzer in self.analyzers: From b2a16d992e6192c366fbd70a5a951f0f27666187 Mon Sep 17 00:00:00 2001 From: clpetix Date: Sat, 10 Feb 2024 17:57:48 -0600 Subject: [PATCH 04/13] Check if barostat is set before trying to scale box. --- src/relentless/simulate/hoomd.py | 91 ++++++++++++++++++++------------ 1 file changed, 57 insertions(+), 34 deletions(-) diff --git a/src/relentless/simulate/hoomd.py b/src/relentless/simulate/hoomd.py index 98f78001..113034fc 100644 --- a/src/relentless/simulate/hoomd.py +++ b/src/relentless/simulate/hoomd.py @@ -710,40 +710,63 @@ def _call_v2(self, sim): for analyzer in self.analyzers: analyzer.pre_run(sim, self) - if isinstance(self.barostat, extent.Extent): - period = None - if sim.dimension == 2: - Lx, Ly, xy = self.barostat.as_array("HOOMD") - Lz, xz, yz = None, None, None - else: - Lx, Ly, Lz, xy, xz, yz = self.barostat.as_array("HOOMD") - elif all(isinstance(n, extent.Extent) for n in self.barostat): - period = 1 - if sim.dimension == 2: - Lx1, Ly1, xy1 = self.barostat[0].as_array("HOOMD") - Lz1, xz1, yz1 = None, None, None - - Lx2, Ly2, xy2 = self.barostat[1].as_array("HOOMD") - Lz2, xz2, yz2 = None, None, None - - Lz, xz, yz = None, None, None - Lx = hoomd.variant.linear_interp([(0, Lx1), (self.steps - 1, Lx2)]) - Ly = hoomd.variant.linear_interp([(0, Ly1), (self.steps - 1, Ly2)]) - xy = hoomd.variant.linear_interp([(0, xy1), (self.steps - 1, xy2)]) - else: - Lx1, Ly1, Lz1, xy1, xz1, yz1 = self.barostat[0].as_array("HOOMD") - Lx2, Ly2, Lz2, xy2, xz2, yz2 = self.barostat[1].as_array("HOOMD") - - Lx = hoomd.variant.linear_interp([(0, Lx1), (self.steps - 1, Lx2)]) - Ly = hoomd.variant.linear_interp([(0, Ly1), (self.steps - 1, Ly2)]) - Lz = hoomd.variant.linear_interp([(0, Lz1), (self.steps - 1, Lz2)]) - xy = hoomd.variant.linear_interp([(0, xy1), (self.steps - 1, xy2)]) - xz = hoomd.variant.linear_interp([(0, xz1), (self.steps - 1, xz2)]) - yz = hoomd.variant.linear_interp([(0, yz1), (self.steps - 1, yz2)]) - - hoomd.update.box_resize( - Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, period=period - ) + if self.barostat is not None: + if isinstance(self.barostat, extent.Extent): + period = None + if sim.dimension == 2: + Lx, Ly, xy = self.barostat.as_array("HOOMD") + Lz, xz, yz = None, None, None + else: + Lx, Ly, Lz, xy, xz, yz = self.barostat.as_array("HOOMD") + elif all(isinstance(n, extent.Extent) for n in self.barostat): + period = 1 + if sim.dimension == 2: + Lx1, Ly1, xy1 = self.barostat[0].as_array("HOOMD") + Lz1, xz1, yz1 = None, None, None + + Lx2, Ly2, xy2 = self.barostat[1].as_array("HOOMD") + Lz2, xz2, yz2 = None, None, None + + Lz, xz, yz = None, None, None + Lx = hoomd.variant.linear_interp( + [(0, Lx1), (self.steps - 1, Lx2)] + ) + Ly = hoomd.variant.linear_interp( + [(0, Ly1), (self.steps - 1, Ly2)] + ) + xy = hoomd.variant.linear_interp( + [(0, xy1), (self.steps - 1, xy2)] + ) + else: + Lx1, Ly1, Lz1, xy1, xz1, yz1 = self.barostat[0].as_array( + "HOOMD" + ) + Lx2, Ly2, Lz2, xy2, xz2, yz2 = self.barostat[1].as_array( + "HOOMD" + ) + + Lx = hoomd.variant.linear_interp( + [(0, Lx1), (self.steps - 1, Lx2)] + ) + Ly = hoomd.variant.linear_interp( + [(0, Ly1), (self.steps - 1, Ly2)] + ) + Lz = hoomd.variant.linear_interp( + [(0, Lz1), (self.steps - 1, Lz2)] + ) + xy = hoomd.variant.linear_interp( + [(0, xy1), (self.steps - 1, xy2)] + ) + xz = hoomd.variant.linear_interp( + [(0, xz1), (self.steps - 1, xz2)] + ) + yz = hoomd.variant.linear_interp( + [(0, yz1), (self.steps - 1, yz2)] + ) + + hoomd.update.box_resize( + Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, period=period + ) hoomd.run(self.steps) for analyzer in self.analyzers: From 8576067f3618506c8f993205275c19f7c7e96ab4 Mon Sep 17 00:00:00 2001 From: clpetix Date: Mon, 12 Feb 2024 11:57:57 -0600 Subject: [PATCH 05/13] Add barostat to LAMMPS and Dilute simulation methods. --- src/relentless/simulate/dilute.py | 3 ++- src/relentless/simulate/lammps.py | 3 ++- 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/src/relentless/simulate/dilute.py b/src/relentless/simulate/dilute.py index f54ee9ee..540d0111 100644 --- a/src/relentless/simulate/dilute.py +++ b/src/relentless/simulate/dilute.py @@ -110,9 +110,10 @@ def _modify_ensemble(self, sim): class RunLangevinDynamics(_Integrator): - def __init__(self, steps, timestep, T, friction, seed, analyzers): + def __init__(self, steps, timestep, T, friction, seed, analyzers, barostat): super().__init__(steps, timestep, analyzers) self.T = T + self.barostat = barostat _modify_ensemble = RunBrownianDynamics._modify_ensemble diff --git a/src/relentless/simulate/lammps.py b/src/relentless/simulate/lammps.py index 036e9b34..ef5c922a 100644 --- a/src/relentless/simulate/lammps.py +++ b/src/relentless/simulate/lammps.py @@ -694,9 +694,10 @@ class RunLangevinDynamics(_Integrator): """ - def __init__(self, steps, timestep, T, friction, seed, analyzers): + def __init__(self, steps, timestep, T, friction, seed, analyzers, barostat): super().__init__(steps, timestep, analyzers) self.T = T + self.barostat = barostat self.friction = friction self.seed = seed From ac0d4c3f6499ff9c89a68c71f8dca228993a1e82 Mon Sep 17 00:00:00 2001 From: clpetix Date: Mon, 12 Feb 2024 12:00:16 -0600 Subject: [PATCH 06/13] Remove default argument from backend operation. Make barostat to be last argument. Check that barostat has at most two Extents. --- src/relentless/simulate/hoomd.py | 11 ++++------- src/relentless/simulate/md.py | 11 +++++++++-- 2 files changed, 13 insertions(+), 9 deletions(-) diff --git a/src/relentless/simulate/hoomd.py b/src/relentless/simulate/hoomd.py index 113034fc..5adc25e3 100644 --- a/src/relentless/simulate/hoomd.py +++ b/src/relentless/simulate/hoomd.py @@ -656,7 +656,7 @@ class RunLangevinDynamics(_Integrator): """ - def __init__(self, steps, timestep, T, friction, seed, analyzers, barostat=None): + def __init__(self, steps, timestep, T, friction, seed, analyzers, barostat): super().__init__(steps, timestep, analyzers) self.T = T self.barostat = barostat @@ -709,7 +709,6 @@ def _call_v2(self, sim): # run + analysis for analyzer in self.analyzers: analyzer.pre_run(sim, self) - if self.barostat is not None: if isinstance(self.barostat, extent.Extent): period = None @@ -718,15 +717,13 @@ def _call_v2(self, sim): Lz, xz, yz = None, None, None else: Lx, Ly, Lz, xy, xz, yz = self.barostat.as_array("HOOMD") - elif all(isinstance(n, extent.Extent) for n in self.barostat): + elif len(self.barostat) == 2 and all( + isinstance(n, extent.Extent) for n in self.barostat + ): period = 1 if sim.dimension == 2: Lx1, Ly1, xy1 = self.barostat[0].as_array("HOOMD") - Lz1, xz1, yz1 = None, None, None - Lx2, Ly2, xy2 = self.barostat[1].as_array("HOOMD") - Lz2, xz2, yz2 = None, None, None - Lz, xz, yz = None, None, None Lx = hoomd.variant.linear_interp( [(0, Lx1), (self.steps - 1, Lx2)] diff --git a/src/relentless/simulate/md.py b/src/relentless/simulate/md.py index b761175e..3ff1b2b3 100644 --- a/src/relentless/simulate/md.py +++ b/src/relentless/simulate/md.py @@ -119,7 +119,14 @@ class RunLangevinDynamics(_Integrator): """ def __init__( - self, steps, timestep, T, friction, seed, barostat=None, analyzers=None + self, + steps, + timestep, + T, + friction, + seed, + analyzers=None, + barostat=None, ): super().__init__(steps, timestep, analyzers) self.T = T @@ -133,10 +140,10 @@ def _make_delegate(self, sim): steps=self.steps, timestep=self.timestep, T=self.T, - barostat=self.barostat, friction=self.friction, seed=self.seed, analyzers=self.analyzers, + barostat=self.barostat, ) From 798586d109101b755479aada6e456ab083d87ba4 Mon Sep 17 00:00:00 2001 From: clpetix Date: Sat, 6 Apr 2024 10:37:04 -0500 Subject: [PATCH 07/13] Set freud version less than 3 --- requirements.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/requirements.txt b/requirements.txt index 2bb0b8d2..4d6bef1a 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,4 +1,4 @@ -freud-analysis>=2 +freud-analysis=2.13.2 gsd lammpsio>=0.4 networkx>=2.5 From b87c3c357c5e568fafe579c8e1964a190f45371f Mon Sep 17 00:00:00 2001 From: clpetix Date: Sat, 6 Apr 2024 10:37:50 -0500 Subject: [PATCH 08/13] Require freud version 2. --- requirements.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/requirements.txt b/requirements.txt index 4d6bef1a..20b6c11d 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,4 +1,4 @@ -freud-analysis=2.13.2 +freud-analysis==2.13.2 gsd lammpsio>=0.4 networkx>=2.5 From 28f55ce3168cb3429528869c8dedcceba8e76a1e Mon Sep 17 00:00:00 2001 From: clpetix Date: Sat, 6 Apr 2024 10:45:08 -0500 Subject: [PATCH 09/13] Revert freud change. --- requirements.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/requirements.txt b/requirements.txt index 20b6c11d..2bb0b8d2 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,4 +1,4 @@ -freud-analysis==2.13.2 +freud-analysis>=2 gsd lammpsio>=0.4 networkx>=2.5 From c5accde89dc1e7339076268de237a6cad913e0d8 Mon Sep 17 00:00:00 2001 From: clpetix Date: Fri, 15 May 2026 14:52:21 -0500 Subject: [PATCH 10/13] Add barostat to RunBrownianDynamics initializer in md.py --- src/relentless/simulate/md.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/src/relentless/simulate/md.py b/src/relentless/simulate/md.py index 3ff1b2b3..953bd57e 100644 --- a/src/relentless/simulate/md.py +++ b/src/relentless/simulate/md.py @@ -78,11 +78,14 @@ class RunBrownianDynamics(_Integrator): """ - def __init__(self, steps, timestep, T, friction, seed, analyzers=None): + def __init__( + self, steps, timestep, T, friction, seed, analyzers=None, barostat=None + ): super().__init__(steps, timestep, analyzers) self.T = T self.friction = friction self.seed = seed + self.barostat = barostat def _make_delegate(self, sim): return self._get_delegate( @@ -93,6 +96,7 @@ def _make_delegate(self, sim): friction=self.friction, seed=self.seed, analyzers=self.analyzers, + barostat=self.barostat, ) From 6476ddc773816f890c6aebd882ff6ac9a98e1176 Mon Sep 17 00:00:00 2001 From: clpetix Date: Fri, 15 May 2026 14:52:44 -0500 Subject: [PATCH 11/13] Add barostat to RunBrownianDynamics in dilute.py --- src/relentless/simulate/dilute.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/relentless/simulate/dilute.py b/src/relentless/simulate/dilute.py index 2de84d66..70f180e5 100644 --- a/src/relentless/simulate/dilute.py +++ b/src/relentless/simulate/dilute.py @@ -93,9 +93,10 @@ def _modify_ensemble(self, sim): class RunBrownianDynamics(_Integrator): - def __init__(self, steps, timestep, T, friction, seed, analyzers): + def __init__(self, steps, timestep, T, friction, seed, analyzers, barostat=None): super().__init__(steps, timestep, analyzers) self.T = T + self.barostat = barostat def _modify_ensemble(self, sim): thermostat = md.Thermostat(self.T) From afc47b392b43dba4a230ac8e41964af69b2e282c Mon Sep 17 00:00:00 2001 From: clpetix Date: Fri, 15 May 2026 15:11:47 -0500 Subject: [PATCH 12/13] Add box resizing functionality in HOOMD-blue. --- src/relentless/simulate/hoomd.py | 138 ++++++++++++++++++++++++++----- 1 file changed, 119 insertions(+), 19 deletions(-) diff --git a/src/relentless/simulate/hoomd.py b/src/relentless/simulate/hoomd.py index 37fc1578..6c26345e 100644 --- a/src/relentless/simulate/hoomd.py +++ b/src/relentless/simulate/hoomd.py @@ -472,6 +472,62 @@ def _make_kT(sim, thermostat, steps): else: return kB * thermostat.T + @staticmethod + def _make_box(sim, V): + """Cast an extent into a HOOMD box.""" + box_array = V.as_array("HOOMD") + if sim.dimension == 2: + return hoomd.Box(Lx=box_array[0], Ly=box_array[1], xy=box_array[2]) + else: + return hoomd.Box( + Lx=box_array[0], + Ly=box_array[1], + Lz=box_array[2], + xy=box_array[3], + xz=box_array[4], + yz=box_array[5], + ) + + def _make_box_resize(self, sim, barostat): + """Create a HOOMD box resize updater.""" + if barostat is None: + return None + elif isinstance(barostat, extent.Extent): + box = self._make_box(sim, barostat) + hoomd.update.BoxResize.update( + sim["engine"]["_hoomd"].state, box, hoomd.filter.All() + ) + return None + elif len(barostat) == 2 and all(isinstance(V, extent.Extent) for V in barostat): + box_1 = self._make_box(sim, barostat[0]) + box_2 = self._make_box(sim, barostat[1]) + if self.steps > 1: + variant = hoomd.variant.Ramp( + A=0, + B=1, + t_start=sim["engine"]["_hoomd"].timestep, + t_ramp=self.steps - 1, + ) + trigger = hoomd.trigger.Periodic(1) + if _hoomd_version >= packaging.version.Version("4.6.0"): + box = hoomd.variant.box.Interpolate(box_1, box_2, variant) + box_resize = hoomd.update.BoxResize(trigger=trigger, box=box) + else: + box_resize = hoomd.update.BoxResize( + trigger=trigger, + box1=box_1, + box2=box_2, + variant=variant, + ) + return box_resize + else: + hoomd.update.BoxResize.update( + sim["engine"]["_hoomd"].state, box_2, hoomd.filter.All() + ) + return None + else: + raise TypeError("barostat must be an Extent or a pair of Extent objects.") + class RunBrownianDynamics(_Integrator): """Perform a Brownian dynamics simulation. @@ -492,14 +548,18 @@ class RunBrownianDynamics(_Integrator): Seed used to randomly generate a uniform force. analyzers : :class:`~relentless.simulate.AnalysisOperation` or list Analysis operations to perform with run (defaults to ``None``). + barostat : :class:`~relentless.model.extent.Extent` or tuple + Target simulation extent, or pair of extents for linear box resizing. + None means no box resizing. """ - def __init__(self, steps, timestep, T, friction, seed, analyzers): + def __init__(self, steps, timestep, T, friction, seed, analyzers, barostat=None): super().__init__(steps, timestep, analyzers) self.T = T self.friction = friction self.seed = seed + self.barostat = barostat def __call__(self, sim): kT = self._make_kT(sim, self.T, self.steps) @@ -522,11 +582,16 @@ def __call__(self, sim): for analyzer in self.analyzers: analyzer.pre_run(sim, self) + box_resize = self._make_box_resize(sim, self.barostat) + if box_resize is not None: + sim["engine"]["_hoomd"].operations.updaters.append(box_resize) sim["engine"]["_hoomd"].run(self.steps, write_at_start=True) for analyzer in self.analyzers: analyzer.post_run(sim, self) + if box_resize is not None: + sim["engine"]["_hoomd"].operations.updaters.remove(box_resize) sim["engine"]["_hoomd"].operations.integrator = None del ig, bd @@ -550,15 +615,18 @@ class RunLangevinDynamics(_Integrator): Seed used to randomly generate a uniform force. analyzers : :class:`~relentless.simulate.AnalysisOperation` or list Analysis operations to perform with run (defaults to ``None``). + barostat : :class:`~relentless.model.extent.Extent` or tuple + Target simulation extent, or pair of extents for linear box resizing. + None means no box resizing. """ - def __init__(self, steps, timestep, T, friction, seed, analyzers, barostat): + def __init__(self, steps, timestep, T, friction, seed, analyzers, barostat=None): super().__init__(steps, timestep, analyzers) self.T = T - self.barostat = barostat self.friction = friction self.seed = seed + self.barostat = barostat def __call__(self, sim): kT = self._make_kT(sim, self.T, self.steps) @@ -581,11 +649,16 @@ def __call__(self, sim): for analyzer in self.analyzers: analyzer.pre_run(sim, self) + box_resize = self._make_box_resize(sim, self.barostat) + if box_resize is not None: + sim["engine"]["_hoomd"].operations.updaters.append(box_resize) sim["engine"]["_hoomd"].run(self.steps, write_at_start=True) for analyzer in self.analyzers: analyzer.post_run(sim, self) + if box_resize is not None: + sim["engine"]["_hoomd"].operations.updaters.remove(box_resize) sim["engine"]["_hoomd"].operations.integrator = None del ig, ld @@ -608,8 +681,10 @@ class RunMolecularDynamics(_Integrator): Simulation time step. thermostat : :class:`~relentless.simulate.Thermostat` Thermostat for temperature control. None means no thermostat. - barostat : :class:`~relentless.simulate.Barostat` - Barostat for pressure control. None means no barostat. + barostat : :class:`~relentless.simulate.Barostat`, \ + :class:`~relentless.model.extent.Extent`, or tuple + Barostat for pressure control, target simulation extent, or pair of + extents for linear box resizing. None means no barostat or resizing. analyzers : :class:`~relentless.simulate.AnalysisOperation` or list Analysis operations to perform with run (defaults to ``None``). @@ -631,7 +706,27 @@ def __call__(self, sim): else: kT = None - if self.thermostat is None and self.barostat is None: + # check if barostat is a box resize + if isinstance(self.barostat, extent.Extent): + self._is_box_resize = True + else: + try: + if len(self.barostat) == 2 and all( + isinstance(V, extent.Extent) for V in self.barostat + ): + self._is_box_resize = True + except TypeError: + self._is_box_resize = False + + # distinguish pressure vs box resize barostat + if self._is_box_resize: + pressure_barostat = None + resize_barostat = self.barostat + else: + pressure_barostat = self.barostat + resize_barostat = None + + if self.thermostat is None and pressure_barostat is None: if _hoomd_version.major >= 4: ig_method = hoomd.md.methods.ConstantVolume(filter=hoomd.filter.All()) else: @@ -641,7 +736,7 @@ def __call__(self, sim): ig_method = hoomd.md.methods.NVE(filter=hoomd.filter.All()) elif ( isinstance(self.thermostat, md.BerendsenThermostat) - and self.barostat is None + and pressure_barostat is None ): if _hoomd_version.major >= 4: ig_method = hoomd.md.methods.ConstantVolume( @@ -661,7 +756,7 @@ def __call__(self, sim): ) elif ( isinstance(self.thermostat, md.NoseHooverThermostat) - and self.barostat is None + and pressure_barostat is None ): if _hoomd_version.major >= 4: ig_method = hoomd.md.methods.ConstantVolume( @@ -679,12 +774,12 @@ def __call__(self, sim): kT=kT, tau=self.thermostat.tau, ) - elif self.thermostat is None and isinstance(self.barostat, md.MTKBarostat): + elif self.thermostat is None and isinstance(pressure_barostat, md.MTKBarostat): if _hoomd_version.major >= 4: ig_method = hoomd.md.methods.ConstantPressure( filter=hoomd.filter.All(), - S=self.barostat.P, - tauS=self.barostat.tau, + S=pressure_barostat.P, + tauS=pressure_barostat.tau, couple="xyz" if sim.dimension == 3 else "xy", ) else: @@ -693,18 +788,18 @@ def __call__(self, sim): warnings.simplefilter(action="ignore", category=FutureWarning) ig_method = hoomd.md.methods.NPH( filter=hoomd.filter.All(), - S=self.barostat.P, - tauS=self.barostat.tau, + S=pressure_barostat.P, + tauS=pressure_barostat.tau, couple="xyz" if sim.dimension == 3 else "xy", ) elif isinstance(self.thermostat, md.NoseHooverThermostat) and isinstance( - self.barostat, md.MTKBarostat + pressure_barostat, md.MTKBarostat ): if _hoomd_version.major >= 4: ig_method = hoomd.md.methods.ConstantPressure( filter=hoomd.filter.All(), - S=self.barostat.P, - tauS=self.barostat.tau, + S=pressure_barostat.P, + tauS=pressure_barostat.tau, couple="xyz" if sim.dimension == 3 else "xy", thermostat=hoomd.md.methods.thermostats.MTTK( kT=kT, tau=self.thermostat.tau @@ -718,8 +813,8 @@ def __call__(self, sim): filter=hoomd.filter.All(), kT=kT, tau=self.thermostat.tau, - S=self.barostat.P, - tauS=self.barostat.tau, + S=pressure_barostat.P, + tauS=pressure_barostat.tau, couple="xyz" if sim.dimension == 3 else "xy", ) else: @@ -738,11 +833,16 @@ def __call__(self, sim): for analyzer in self.analyzers: analyzer.pre_run(sim, self) + box_resize = self._make_box_resize(sim, resize_barostat) + if box_resize is not None: + sim["engine"]["_hoomd"].operations.updaters.append(box_resize) sim["engine"]["_hoomd"].run(self.steps, write_at_start=True) for analyzer in self.analyzers: analyzer.post_run(sim, self) + if box_resize is not None: + sim["engine"]["_hoomd"].operations.updaters.remove(box_resize) sim["engine"]["_hoomd"].operations.integrator = None del ig, ig_method @@ -840,7 +940,7 @@ def _get_constrained_quantities(self, sim, sim_op): if sim_op.thermostat is not None: constraints["T"] = sim_op.thermostat # conjugate pair: one or the other is set - if sim_op.barostat is not None: + if isinstance(sim_op.barostat, md.MTKBarostat): constraints["P"] = sim_op.barostat.P else: constraints["V"] = True From e0daf73dbe6f5eb25d1f7bd25e28edb4d751bb0b Mon Sep 17 00:00:00 2001 From: clpetix Date: Fri, 15 May 2026 15:13:00 -0500 Subject: [PATCH 13/13] Test box resizing in HOOMD. --- tests/simulate/test_hoomd.py | 57 ++++++++++++++++++++++++++++++++++++ 1 file changed, 57 insertions(+) diff --git a/tests/simulate/test_hoomd.py b/tests/simulate/test_hoomd.py index 9be6edcb..6f307be3 100644 --- a/tests/simulate/test_hoomd.py +++ b/tests/simulate/test_hoomd.py @@ -80,6 +80,24 @@ def ens_pot(self): return (ens, pots) + def resize_extents(self): + if self.dim == 3: + return (relentless.model.Cube(L=18.0), relentless.model.Cube(L=19.0)) + elif self.dim == 2: + return (relentless.model.Square(L=18.0), relentless.model.Square(L=19.0)) + else: + raise ValueError("HOOMD supports 2d and 3d simulations") + + def assert_box(self, sim, V): + box = sim["engine"]["_hoomd"].state.box + if self.dim == 3: + box_array = numpy.array([box.Lx, box.Ly, box.Lz, box.xy, box.xz, box.yz]) + elif self.dim == 2: + box_array = numpy.array([box.Lx, box.Ly, box.xy]) + else: + raise ValueError("HOOMD supports 2d and 3d simulations") + numpy.testing.assert_allclose(box_array, V.as_array("HOOMD")) + def ens_pot_bonds(self): if self.dim == 3: ens = relentless.model.Ensemble( @@ -705,6 +723,20 @@ def test_brownian_dynamics(self): brn.T = (ens.T, 1.5 * ens.T) h.run(pot, self.directory) + # box resizing + V1, V2 = self.resize_extents() + brn.barostat = V1 + h = relentless.simulate.HOOMD(init, brn) + sim = h.run(pot, self.directory) + self.assert_box(sim, V1) + + brn.steps = 2 + brn.barostat = (V1, V2) + sim = h.run(pot, self.directory) + self.assert_box(sim, V2) + brn.steps = 1 + brn.barostat = None + def test_langevin_dynamics(self): ens, pot = self.ens_pot() init = relentless.simulate.InitializeRandomly(seed=1, N=ens.N, V=ens.V, T=ens.T) @@ -721,6 +753,19 @@ def test_langevin_dynamics(self): lgv.T = (ens.T, 1.5 * ens.T) h.run(pot, self.directory) + # box resizing + V1, V2 = self.resize_extents() + lgv.barostat = V1 + sim = h.run(pot, self.directory) + self.assert_box(sim, V1) + + lgv.steps = 2 + lgv.barostat = (V1, V2) + sim = h.run(pot, self.directory) + self.assert_box(sim, V2) + lgv.steps = 1 + lgv.barostat = None + def test_molecular_dynamics(self): # VerletIntegrator - NVE ens, pot = self.ens_pot() @@ -746,6 +791,18 @@ def test_molecular_dynamics(self): vrl.thermostat = None h.run(pot, self.directory) + # VerletIntegrator - box resizing + V1, V2 = self.resize_extents() + vrl.barostat = V1 + sim = h.run(pot, self.directory) + self.assert_box(sim, V1) + + vrl.steps = 2 + vrl.barostat = (V1, V2) + sim = h.run(pot, self.directory) + self.assert_box(sim, V2) + vrl.steps = 1 + if relentless.mpi.world.size == 1: # VerletIntegrator - NVE (Berendsen) vrl.thermostat = relentless.simulate.BerendsenThermostat(T=1, tau=0.5)