From 062c0259a425782759b36ecf3a4bae877e395ff0 Mon Sep 17 00:00:00 2001 From: Ashley Knoerdel Date: Thu, 15 Aug 2024 09:03:51 -0500 Subject: [PATCH 1/4] Addition of AGCMS potential and documentation. --- doc/module-azplugins-pair.rst | 4 +- src/CMakeLists.txt | 1 + src/PairEvaluatorAGCMS.h | 142 ++++++++++++++++++++++++++++++++++ src/module.cc | 4 + src/pair.py | 103 ++++++++++++++++++++++++ src/pytest/test_pair.py | 58 ++++++++++++++ 6 files changed, 311 insertions(+), 1 deletion(-) create mode 100644 src/PairEvaluatorAGCMS.h diff --git a/doc/module-azplugins-pair.rst b/doc/module-azplugins-pair.rst index 56ecbf27..9880037e 100644 --- a/doc/module-azplugins-pair.rst +++ b/doc/module-azplugins-pair.rst @@ -17,6 +17,7 @@ azplugins.pair Hertz PerturbedLennardJones TwoPatchMorse + AGCMS .. rubric:: Details @@ -26,6 +27,7 @@ azplugins.pair DPDGeneralWeight, Hertz, PerturbedLennardJones, - TwoPatchMorse + TwoPatchMorse, + AGCMS :no-inherited-members: :show-inheritance: diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 030a8d17..ea25f270 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -32,6 +32,7 @@ set(_pair_evaluators Colloid Hertz PerturbedLennardJones + AGCMS ) # TODO: Add names of all dpd evaluators diff --git a/src/PairEvaluatorAGCMS.h b/src/PairEvaluatorAGCMS.h new file mode 100644 index 00000000..cbd04ed3 --- /dev/null +++ b/src/PairEvaluatorAGCMS.h @@ -0,0 +1,142 @@ +// Copyright (c) 2018-2020, Michael P. Howard +// Copyright (c) 2021-2024, Auburn University +// Part of azplugins, released under the BSD 3-Clause License. + +#ifndef AZPLUGINS_PAIR_EVALUATOR_AGCMS_H_ +#define AZPLUGINS_PAIR_EVALUATOR_AGCMS_H_ + +#include "PairEvaluator.h" +#include + +#ifdef __HIPCC__ +#define DEVICE __device__ +#else +#define DEVICE +#endif + +namespace hoomd + { +namespace azplugins + { +namespace detail + { + +struct PairParametersAGCMS : public PairParameters + { +#ifndef __HIPCC__ + PairParametersAGCMS() : w(0), sigma(0), a(0), q(0) { } + + PairParametersAGCMS(pybind11::dict v, bool managed = false) + { + w = v["w"].cast(); + sigma = v["sigma"].cast(); + a = v["a"].cast(); + q = v["q"].cast(); + } + + pybind11::dict asDict() + { + pybind11::dict v; + v["w"] = w; + v["sigma"] = sigma; + v["a"] = a; + v["q"] = q; + return v; + } +#endif // __HIPCC__ + + Scalar w; //!< well width [length] + Scalar sigma; //!< Hard core repulsion distance [length] + Scalar a; //!< Well depth [energy] + Scalar q; //!< steepness + } +#if HOOMD_LONGREAL_SIZE == 32 + __attribute__((aligned(16))); +#else + __attribute__((aligned(32))); +#endif + +//! Class for evaluating the adjusted generalized continuous multiple step pair potential +/*! + * This is the generalized continuous multiple step pair potential + * modified such that there is only one sigma and one minimum present + * with any combination of well depth and width causing the potential + * to have a value of zero at sigma. + */ +class PairEvaluatorAGCMS : public PairEvaluator + { + public: + typedef PairParametersAGCMS param_type; + + //! Constructor + /*! + * \param _rsq Squared distance between particles + * \param _rcutsq Cutoff radius squared + * \param _params Pair potential parameters, given by typedef above + * + * The functor initializes its members from \a _params. + */ + DEVICE PairEvaluatorAGCMS(Scalar _rsq, Scalar _rcutsq, const param_type& _params) + : PairEvaluator(_rsq, _rcutsq) + { + w = _params.w; + sigma = _params.sigma; + a = _params.a; + q = _params.q; + } + + //! Evaluate the force and energy + /*! + * \param force_divr Holds the computed force divided by r + * \param pair_eng Holds the computed pair energy + * \param energy_shift If true, the potential is shifted to zero at the cutoff + * + * \returns True if the energy calculation occurs + * + * The calculation does not occur if the pair distance is greater than the cutoff + * or if the potential is scaled to zero. + */ + + DEVICE bool evalForceAndEnergy(Scalar& force_divr, Scalar& pair_eng, bool energy_shift) + { + // compute the force divided by r in force_divr + if (rsq < rcutsq && a != Scalar(0)) + { + Scalar r = fast::sqrt(rsq); + Scalar rinv = 1 / r; + Scalar w_rinv_shifted = w / (r - sigma + w); + Scalar core_repuls = pow(w_rinv_shifted, q); + Scalar exponent = exp(q*(r - sigma - w)/w); + force_divr = -a*rinv*(q*core_repuls*w_rinv_shifted - q*exponent/(w*pow(1 + exponent, 2.0))); + pair_eng = a*(1/(1 + exponent) - core_repuls); + return true; + } + else + return false; + } + +#ifndef __HIPCC__ + //! Get the name of this potential + /*! \returns The potential name. Must be short and all lowercase, as this is the name energies + will be logged as via analyze.log. + */ + static std::string getName() + { + return std::string("agcms"); + } +#endif + + protected: + Scalar w; //!< well width + Scalar sigma; //!< Hard core repulsion distance + Scalar a; //!< well depth + Scalar q; //!< Steepness + }; + + } // end namespace detail + } // end namespace azplugins + } // end namespace hoomd + +#undef DEVICE + +#endif // AZPLUGINS_PAIR_EVALUATOR_AGCMS_H_ diff --git a/src/module.cc b/src/module.cc index f3e183d5..acb74e5d 100644 --- a/src/module.cc +++ b/src/module.cc @@ -61,6 +61,7 @@ void export_AnisoPotentialPairTwoPatchMorse(pybind11::module&); void export_PotentialPairColloid(pybind11::module&); void export_PotentialPairHertz(pybind11::module&); void export_PotentialPairPerturbedLennardJones(pybind11::module&); +void export_PotentialPairAGCMS(pybind11::module&); // dpd void export_PotentialPairDPDThermoGeneralWeight(pybind11::module&); @@ -74,6 +75,7 @@ void export_AnisoPotentialPairTwoPatchMorseGPU(pybind11::module&); void export_PotentialPairColloidGPU(pybind11::module&); void export_PotentialPairHertzGPU(pybind11::module&); void export_PotentialPairPerturbedLennardJonesGPU(pybind11::module&); +void export_PotentialPairAGCMSGPU(pybind11::module&); // dpd void export_PotentialPairDPDThermoGeneralWeightGPU(pybind11::module&); @@ -101,6 +103,7 @@ PYBIND11_MODULE(_azplugins, m) export_PotentialPairColloid(m); export_PotentialPairHertz(m); export_PotentialPairPerturbedLennardJones(m); + export_PotentialPairAGCMS(m); // dpd pair export_PotentialPairDPDThermoGeneralWeight(m); @@ -114,6 +117,7 @@ PYBIND11_MODULE(_azplugins, m) export_PotentialPairColloidGPU(m); export_PotentialPairHertzGPU(m); export_PotentialPairPerturbedLennardJonesGPU(m); + export_PotentialPairAGCMSGPU(m); // dpd pair export_PotentialPairDPDThermoGeneralWeightGPU(m); diff --git a/src/pair.py b/src/pair.py index a0a49539..af74b7b4 100644 --- a/src/pair.py +++ b/src/pair.py @@ -459,3 +459,106 @@ def __init__(self, nlist, default_r_cut=None, mode='none'): ), ) self._add_typeparam(params) + +class AGCMS(pair.Pair): + r"""Adjusted Generalized Continuous Multiple Step potential. + + Args: + nlist (hoomd.md.nlist.NeighborList): Neighbor list. + r_cut (float): Default cutoff radius :math:`[\mathrm{length}]`. + mode (str): Energy shifting/smoothing mode. + + `AGCMS` is a Generalized Continuous Mutltiple Step potential adjusted such + that all :math:`a` and :math:`w` values at a distance of :math:`\sigma` will + have a potential value of zero. + + + The original GCMS potential from |Munguia-Valadez et al.|_ is: + + .. math:: + + U(r) = U_{\rm{core}}(r, \sigma_0, w_0, q_0) + \sum^{n}_{k=1}U_{\rm{step}}(r, + \sigma_k, w_k, q_k) + + where :math:`U_{\rm{core}}` is + + .. math:: + + U_{\rm{core}} = \left(\frac{w_0}{r - \sigma_0 + w_0}\right)^{q_0} + + + and :math:`U_{\rm{step}}` is + + .. math:: + + U_{\rm{step}} = \left(\frac{a_k}{1 + \exp[q_k(r-\sigma_k-w_k)/w_k]}\right) + + The `AGCMS` potential is + + .. math:: + + U(r) = a\left[U_{\rm{step}}(r, \sigma, w, q) - + U_{\rm{core}}(r, \sigma, w, q)\right] + + where + + .. math:: + + U_{\rm{core}} = \left(\frac{w}{r - \sigma + w}\right)^q + + and :math:`U_{\rm{step}}` is + + .. math:: + + U_{\rm{step}} = \left(\frac{1}{1 + \exp[q(r-\sigma-w)/w]}\right) + + Here, :math:`w` is the well width, :math:`a` is the signed well depth, + :math:`\sigma` is the hard core repulsion diameter, and :math:`q` is the + steepness of the potential. + + Example:: + + nl = nlist.Cell() + agcms = pair.AGCMS(default_r_cut=4.0, nlist=nl) + # particle types A interacting + agcms.params[('A', 'A')] = dict(w=1.0, sigma=2.0, a=2.0, q=16.0) + + .. py:attribute:: params + + The `AGCMS` potential parameters. The dictionary has the following + keys: + + * ``w`` (`float`, **required**) - Well width :math:`w` + :math:`[\mathrm{length}]` + * ``sigma`` (`float`, **required**) - Hard core repulsion diameter + :math:`\sigma` :math:`[\mathrm{length}]` + * ``a`` (`float`, **required**) - Depth of well :math:`a` + :math:`[\mathrm{energy}]` + * ``q`` (`float`, **required**) - Steepness parameter for well + + Type: :class:`~hoomd.data.typeparam.TypeParameter` [`tuple` + [``particle_type``, ``particle_type``], `dict`] + + .. py:attribute:: mode + + Energy shifting/smoothing mode: ``"none"``. + + Type: `str` + + .. |acutei| unicode:: U+00ED .. acute i + .. |Munguia-Valadez et al.| replace:: Mungu\ |acutei|\ a-Valadez et al. + .. _Munguia-Valadez et al.: https://doi.org/10.1088/1361-648X/ac4fe8 + + """ + + _ext_module = _azplugins + _cpp_class_name = 'PotentialPairAGCMS' + + def __init__(self, nlist, default_r_cut=None, default_r_on=0, mode='none'): + super().__init__(nlist, default_r_cut, default_r_on, mode) + params = TypeParameter( + 'params', + 'particle_types', + TypeParameterDict(w=float, sigma=float, a=float, q=float, len_keys=2), + ) + self._add_typeparam(params) diff --git a/src/pytest/test_pair.py b/src/pytest/test_pair.py index 006d14ed..1838aae6 100644 --- a/src/pytest/test_pair.py +++ b/src/pytest/test_pair.py @@ -238,6 +238,62 @@ ), ] +# AGCMS +potential_tests += [ + # test the calculation of force and potential + # test goes to zero outside cutoff + PotentialTestCase( + hoomd.azplugins.pair.AGCMS, + {'w': 1.0, 'sigma': 1.0, 'a': -1.0, 'q': 16.0}, + 3.0, + False, + 4.0, + 0.0, + 0.0, + ), + # change w to check for changing + # well width + PotentialTestCase( + hoomd.azplugins.pair.AGCMS, + {'w': 1.25, 'sigma': 1.0, 'a': -1.0, 'q': 16.0}, + 3.0, + False, + 1.75, + -0.9977990978, + -0.015776422214703146, + ), + # change sigma so that now the potential will + # be shifted to the right + PotentialTestCase( + hoomd.azplugins.pair.AGCMS, + {'w': 1.0, 'sigma': 2.0, 'a': -1.0, 'q': 16.0}, + 3.0, + False, + 2.5, + -0.998142211, + 0.010875544898269139, + ), + # change well depth to increase the attractiveness + PotentialTestCase( + hoomd.azplugins.pair.AGCMS, + {'w': 1.0, 'sigma': 1.0, 'a': -5.0, 'q': 16.0}, + 3.0, + False, + 1.50, + -4.990711055, + 0.0543777724491345696, + ), + # change q to increase stiffness + PotentialTestCase( + hoomd.azplugins.pair.AGCMS, + {'w': 1.0, 'sigma': 1.0, 'a': -1.0, 'q': 50.0}, + 3.0, + False, + 1.50, + -0.9999999984, + 5.1583220989569564*10**(-8.0), + ), +] @pytest.mark.parametrize( 'potential_test', potential_tests, ids=lambda x: x.potential.__name__ @@ -284,12 +340,14 @@ def test_energy_and_force( # test that the energies match reference values, half goes to each particle energies = potential.energies e = potential_test.energy + print(energies, e) if sim.device.communicator.rank == 0: numpy.testing.assert_array_almost_equal(energies, [0.5 * e, 0.5 * e], decimal=4) # test that the forces match reference values, should be directed along x forces = potential.forces f = potential_test.force + print(forces, f) if sim.device.communicator.rank == 0: numpy.testing.assert_array_almost_equal( forces, [[-f, 0, 0], [f, 0, 0]], decimal=4 From 5f8a240292a65297d8b6c2ba033bb99f08c19c28 Mon Sep 17 00:00:00 2001 From: Ashley Knoerdel Date: Thu, 15 Aug 2024 09:22:21 -0500 Subject: [PATCH 2/4] pre-commit --- src/PairEvaluatorAGCMS.h | 24 +++++++++++++----------- src/pair.py | 1 + src/pytest/test_pair.py | 3 ++- 3 files changed, 16 insertions(+), 12 deletions(-) diff --git a/src/PairEvaluatorAGCMS.h b/src/PairEvaluatorAGCMS.h index cbd04ed3..bc6068bf 100644 --- a/src/PairEvaluatorAGCMS.h +++ b/src/PairEvaluatorAGCMS.h @@ -45,10 +45,10 @@ struct PairParametersAGCMS : public PairParameters } #endif // __HIPCC__ - Scalar w; //!< well width [length] - Scalar sigma; //!< Hard core repulsion distance [length] - Scalar a; //!< Well depth [energy] - Scalar q; //!< steepness + Scalar w; //!< well width [length] + Scalar sigma; //!< Hard core repulsion distance [length] + Scalar a; //!< Well depth [energy] + Scalar q; //!< steepness } #if HOOMD_LONGREAL_SIZE == 32 __attribute__((aligned(16))); @@ -106,9 +106,11 @@ class PairEvaluatorAGCMS : public PairEvaluator Scalar rinv = 1 / r; Scalar w_rinv_shifted = w / (r - sigma + w); Scalar core_repuls = pow(w_rinv_shifted, q); - Scalar exponent = exp(q*(r - sigma - w)/w); - force_divr = -a*rinv*(q*core_repuls*w_rinv_shifted - q*exponent/(w*pow(1 + exponent, 2.0))); - pair_eng = a*(1/(1 + exponent) - core_repuls); + Scalar exponent = exp(q * (r - sigma - w) / w); + force_divr = -a * rinv + * (q * core_repuls * w_rinv_shifted + - q * exponent / (w * pow(1 + exponent, 2.0))); + pair_eng = a * (1 / (1 + exponent) - core_repuls); return true; } else @@ -127,10 +129,10 @@ class PairEvaluatorAGCMS : public PairEvaluator #endif protected: - Scalar w; //!< well width - Scalar sigma; //!< Hard core repulsion distance - Scalar a; //!< well depth - Scalar q; //!< Steepness + Scalar w; //!< well width + Scalar sigma; //!< Hard core repulsion distance + Scalar a; //!< well depth + Scalar q; //!< Steepness }; } // end namespace detail diff --git a/src/pair.py b/src/pair.py index af74b7b4..d666055a 100644 --- a/src/pair.py +++ b/src/pair.py @@ -460,6 +460,7 @@ def __init__(self, nlist, default_r_cut=None, mode='none'): ) self._add_typeparam(params) + class AGCMS(pair.Pair): r"""Adjusted Generalized Continuous Multiple Step potential. diff --git a/src/pytest/test_pair.py b/src/pytest/test_pair.py index 1838aae6..3344d04f 100644 --- a/src/pytest/test_pair.py +++ b/src/pytest/test_pair.py @@ -291,10 +291,11 @@ False, 1.50, -0.9999999984, - 5.1583220989569564*10**(-8.0), + 5.1583220989569564 * 10 ** (-8.0), ), ] + @pytest.mark.parametrize( 'potential_test', potential_tests, ids=lambda x: x.potential.__name__ ) From 19ea2e35a32b05422cf5fc80c0bd7b0599789ba8 Mon Sep 17 00:00:00 2001 From: Ashley Knoerdel Date: Wed, 9 Oct 2024 12:48:33 -0500 Subject: [PATCH 3/4] Change name of potential and adjust documentation --- doc/module-azplugins-pair.rst | 8 +- src/CMakeLists.txt | 2 +- src/PairEvaluatorACSW.h | 144 ++++++++++++++++++++++++++++++++++ src/module.cc | 8 +- src/pair.py | 16 ++-- src/pytest/test_pair.py | 12 +-- 6 files changed, 167 insertions(+), 23 deletions(-) create mode 100644 src/PairEvaluatorACSW.h diff --git a/doc/module-azplugins-pair.rst b/doc/module-azplugins-pair.rst index 9880037e..3fd43a58 100644 --- a/doc/module-azplugins-pair.rst +++ b/doc/module-azplugins-pair.rst @@ -12,22 +12,22 @@ azplugins.pair .. autosummary:: :nosignatures: + ACSW Colloid DPDGeneralWeight Hertz PerturbedLennardJones TwoPatchMorse - AGCMS .. rubric:: Details .. automodule:: hoomd.azplugins.pair :synopsis: Pair potentials. - :members: Colloid, + :members: ACSW, + Colloid, DPDGeneralWeight, Hertz, PerturbedLennardJones, - TwoPatchMorse, - AGCMS + TwoPatchMorse :no-inherited-members: :show-inheritance: diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index ea25f270..68e51d2d 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -32,7 +32,7 @@ set(_pair_evaluators Colloid Hertz PerturbedLennardJones - AGCMS + ACSW ) # TODO: Add names of all dpd evaluators diff --git a/src/PairEvaluatorACSW.h b/src/PairEvaluatorACSW.h new file mode 100644 index 00000000..e6bcbda7 --- /dev/null +++ b/src/PairEvaluatorACSW.h @@ -0,0 +1,144 @@ +// Copyright (c) 2018-2020, Michael P. Howard +// Copyright (c) 2021-2024, Auburn University +// Part of azplugins, released under the BSD 3-Clause License. + +#ifndef AZPLUGINS_PAIR_EVALUATOR_ACSW_H_ +#define AZPLUGINS_PAIR_EVALUATOR_ACSW_H_ + +#include "PairEvaluator.h" +#include + +#ifdef __HIPCC__ +#define DEVICE __device__ +#else +#define DEVICE +#endif + +namespace hoomd + { +namespace azplugins + { +namespace detail + { + +struct PairParametersACSW : public PairParameters + { +#ifndef __HIPCC__ + PairParametersACSW() : w(0), sigma(0), a(0), q(0) { } + + PairParametersACSW(pybind11::dict v, bool managed = false) + { + w = v["w"].cast(); + sigma = v["sigma"].cast(); + a = v["a"].cast(); + q = v["q"].cast(); + } + + pybind11::dict asDict() + { + pybind11::dict v; + v["w"] = w; + v["sigma"] = sigma; + v["a"] = a; + v["q"] = q; + return v; + } +#endif // __HIPCC__ + + Scalar w; //!< well width [length] + Scalar sigma; //!< Hard core repulsion distance [length] + Scalar a; //!< Well depth [energy] + Scalar q; //!< steepness + } +#if HOOMD_LONGREAL_SIZE == 32 + __attribute__((aligned(16))); +#else + __attribute__((aligned(32))); +#endif + +//! Class for evaluating the adjusted generalized continuous multiple step pair potential +/*! + * This is the generalized continuous multiple step pair potential + * modified such that there is only one sigma and one minimum present + * with any combination of well depth and width causing the potential + * to have a value of zero at sigma. + */ +class PairEvaluatorACSW : public PairEvaluator + { + public: + typedef PairParametersACSW param_type; + + //! Constructor + /*! + * \param _rsq Squared distance between particles + * \param _rcutsq Cutoff radius squared + * \param _params Pair potential parameters, given by typedef above + * + * The functor initializes its members from \a _params. + */ + DEVICE PairEvaluatorACSW(Scalar _rsq, Scalar _rcutsq, const param_type& _params) + : PairEvaluator(_rsq, _rcutsq) + { + w = _params.w; + sigma = _params.sigma; + a = _params.a; + q = _params.q; + } + + //! Evaluate the force and energy + /*! + * \param force_divr Holds the computed force divided by r + * \param pair_eng Holds the computed pair energy + * \param energy_shift If true, the potential is shifted to zero at the cutoff + * + * \returns True if the energy calculation occurs + * + * The calculation does not occur if the pair distance is greater than the cutoff + * or if the potential is scaled to zero. + */ + + DEVICE bool evalForceAndEnergy(Scalar& force_divr, Scalar& pair_eng, bool energy_shift) + { + // compute the force divided by r in force_divr + if (rsq < rcutsq && a != Scalar(0)) + { + Scalar r = fast::sqrt(rsq); + Scalar rinv = 1 / r; + Scalar w_rinv_shifted = w / (r - sigma + w); + Scalar core_repuls = pow(w_rinv_shifted, q); + Scalar exponent = exp(q * (r - sigma - w) / w); + force_divr = -a * rinv + * (q * core_repuls * w_rinv_shifted + - q * exponent / (w * pow(1 + exponent, 2.0))); + pair_eng = a * (1 / (1 + exponent) - core_repuls); + return true; + } + else + return false; + } + +#ifndef __HIPCC__ + //! Get the name of this potential + /*! \returns The potential name. Must be short and all lowercase, as this is the name energies + will be logged as via analyze.log. + */ + static std::string getName() + { + return std::string("acsw"); + } +#endif + + protected: + Scalar w; //!< well width + Scalar sigma; //!< Hard core repulsion distance + Scalar a; //!< well depth + Scalar q; //!< Steepness + }; + + } // end namespace detail + } // end namespace azplugins + } // end namespace hoomd + +#undef DEVICE + +#endif // AZPLUGINS_PAIR_EVALUATOR_ACSW_H_ diff --git a/src/module.cc b/src/module.cc index acb74e5d..0e4bf5a7 100644 --- a/src/module.cc +++ b/src/module.cc @@ -61,7 +61,7 @@ void export_AnisoPotentialPairTwoPatchMorse(pybind11::module&); void export_PotentialPairColloid(pybind11::module&); void export_PotentialPairHertz(pybind11::module&); void export_PotentialPairPerturbedLennardJones(pybind11::module&); -void export_PotentialPairAGCMS(pybind11::module&); +void export_PotentialPairACSW(pybind11::module&); // dpd void export_PotentialPairDPDThermoGeneralWeight(pybind11::module&); @@ -75,7 +75,7 @@ void export_AnisoPotentialPairTwoPatchMorseGPU(pybind11::module&); void export_PotentialPairColloidGPU(pybind11::module&); void export_PotentialPairHertzGPU(pybind11::module&); void export_PotentialPairPerturbedLennardJonesGPU(pybind11::module&); -void export_PotentialPairAGCMSGPU(pybind11::module&); +void export_PotentialPairACSWGPU(pybind11::module&); // dpd void export_PotentialPairDPDThermoGeneralWeightGPU(pybind11::module&); @@ -103,7 +103,7 @@ PYBIND11_MODULE(_azplugins, m) export_PotentialPairColloid(m); export_PotentialPairHertz(m); export_PotentialPairPerturbedLennardJones(m); - export_PotentialPairAGCMS(m); + export_PotentialPairACSW(m); // dpd pair export_PotentialPairDPDThermoGeneralWeight(m); @@ -117,7 +117,7 @@ PYBIND11_MODULE(_azplugins, m) export_PotentialPairColloidGPU(m); export_PotentialPairHertzGPU(m); export_PotentialPairPerturbedLennardJonesGPU(m); - export_PotentialPairAGCMSGPU(m); + export_PotentialPairACSWGPU(m); // dpd pair export_PotentialPairDPDThermoGeneralWeightGPU(m); diff --git a/src/pair.py b/src/pair.py index d666055a..4fb1cb3a 100644 --- a/src/pair.py +++ b/src/pair.py @@ -461,17 +461,17 @@ def __init__(self, nlist, default_r_cut=None, mode='none'): self._add_typeparam(params) -class AGCMS(pair.Pair): - r"""Adjusted Generalized Continuous Multiple Step potential. +class ACSW(pair.Pair): + r"""Adjusted Continuous Square Well potential. Args: nlist (hoomd.md.nlist.NeighborList): Neighbor list. r_cut (float): Default cutoff radius :math:`[\mathrm{length}]`. mode (str): Energy shifting/smoothing mode. - `AGCMS` is a Generalized Continuous Mutltiple Step potential adjusted such - that all :math:`a` and :math:`w` values at a distance of :math:`\sigma` will - have a potential value of zero. + `ACSW` is a Generalized Continuous Multiple Step potential adjusted such + that it is a single continuous square well where all :math:`a` and :math:`w` + values at a distance of :math:`\sigma` will have a potential value of zero. The original GCMS potential from |Munguia-Valadez et al.|_ is: @@ -494,7 +494,7 @@ class AGCMS(pair.Pair): U_{\rm{step}} = \left(\frac{a_k}{1 + \exp[q_k(r-\sigma_k-w_k)/w_k]}\right) - The `AGCMS` potential is + The `ACSW` potential is .. math:: @@ -526,7 +526,7 @@ class AGCMS(pair.Pair): .. py:attribute:: params - The `AGCMS` potential parameters. The dictionary has the following + The `ACSW` potential parameters. The dictionary has the following keys: * ``w`` (`float`, **required**) - Well width :math:`w` @@ -553,7 +553,7 @@ class AGCMS(pair.Pair): """ _ext_module = _azplugins - _cpp_class_name = 'PotentialPairAGCMS' + _cpp_class_name = 'PotentialPairACSW' def __init__(self, nlist, default_r_cut=None, default_r_on=0, mode='none'): super().__init__(nlist, default_r_cut, default_r_on, mode) diff --git a/src/pytest/test_pair.py b/src/pytest/test_pair.py index 3344d04f..a172e8b0 100644 --- a/src/pytest/test_pair.py +++ b/src/pytest/test_pair.py @@ -238,12 +238,12 @@ ), ] -# AGCMS +# ACSW potential_tests += [ # test the calculation of force and potential # test goes to zero outside cutoff PotentialTestCase( - hoomd.azplugins.pair.AGCMS, + hoomd.azplugins.pair.ACSW, {'w': 1.0, 'sigma': 1.0, 'a': -1.0, 'q': 16.0}, 3.0, False, @@ -254,7 +254,7 @@ # change w to check for changing # well width PotentialTestCase( - hoomd.azplugins.pair.AGCMS, + hoomd.azplugins.pair.ACSW, {'w': 1.25, 'sigma': 1.0, 'a': -1.0, 'q': 16.0}, 3.0, False, @@ -265,7 +265,7 @@ # change sigma so that now the potential will # be shifted to the right PotentialTestCase( - hoomd.azplugins.pair.AGCMS, + hoomd.azplugins.pair.ACSW, {'w': 1.0, 'sigma': 2.0, 'a': -1.0, 'q': 16.0}, 3.0, False, @@ -275,7 +275,7 @@ ), # change well depth to increase the attractiveness PotentialTestCase( - hoomd.azplugins.pair.AGCMS, + hoomd.azplugins.pair.ACSW, {'w': 1.0, 'sigma': 1.0, 'a': -5.0, 'q': 16.0}, 3.0, False, @@ -285,7 +285,7 @@ ), # change q to increase stiffness PotentialTestCase( - hoomd.azplugins.pair.AGCMS, + hoomd.azplugins.pair.ACSW, {'w': 1.0, 'sigma': 1.0, 'a': -1.0, 'q': 50.0}, 3.0, False, From 7dbb32fdc33c67fa7c73a5941db0217f75022224 Mon Sep 17 00:00:00 2001 From: Ashley Knoerdel Date: Wed, 9 Oct 2024 12:50:57 -0500 Subject: [PATCH 4/4] Remove PairEvaluator file for old potential name --- src/PairEvaluatorAGCMS.h | 144 --------------------------------------- 1 file changed, 144 deletions(-) delete mode 100644 src/PairEvaluatorAGCMS.h diff --git a/src/PairEvaluatorAGCMS.h b/src/PairEvaluatorAGCMS.h deleted file mode 100644 index bc6068bf..00000000 --- a/src/PairEvaluatorAGCMS.h +++ /dev/null @@ -1,144 +0,0 @@ -// Copyright (c) 2018-2020, Michael P. Howard -// Copyright (c) 2021-2024, Auburn University -// Part of azplugins, released under the BSD 3-Clause License. - -#ifndef AZPLUGINS_PAIR_EVALUATOR_AGCMS_H_ -#define AZPLUGINS_PAIR_EVALUATOR_AGCMS_H_ - -#include "PairEvaluator.h" -#include - -#ifdef __HIPCC__ -#define DEVICE __device__ -#else -#define DEVICE -#endif - -namespace hoomd - { -namespace azplugins - { -namespace detail - { - -struct PairParametersAGCMS : public PairParameters - { -#ifndef __HIPCC__ - PairParametersAGCMS() : w(0), sigma(0), a(0), q(0) { } - - PairParametersAGCMS(pybind11::dict v, bool managed = false) - { - w = v["w"].cast(); - sigma = v["sigma"].cast(); - a = v["a"].cast(); - q = v["q"].cast(); - } - - pybind11::dict asDict() - { - pybind11::dict v; - v["w"] = w; - v["sigma"] = sigma; - v["a"] = a; - v["q"] = q; - return v; - } -#endif // __HIPCC__ - - Scalar w; //!< well width [length] - Scalar sigma; //!< Hard core repulsion distance [length] - Scalar a; //!< Well depth [energy] - Scalar q; //!< steepness - } -#if HOOMD_LONGREAL_SIZE == 32 - __attribute__((aligned(16))); -#else - __attribute__((aligned(32))); -#endif - -//! Class for evaluating the adjusted generalized continuous multiple step pair potential -/*! - * This is the generalized continuous multiple step pair potential - * modified such that there is only one sigma and one minimum present - * with any combination of well depth and width causing the potential - * to have a value of zero at sigma. - */ -class PairEvaluatorAGCMS : public PairEvaluator - { - public: - typedef PairParametersAGCMS param_type; - - //! Constructor - /*! - * \param _rsq Squared distance between particles - * \param _rcutsq Cutoff radius squared - * \param _params Pair potential parameters, given by typedef above - * - * The functor initializes its members from \a _params. - */ - DEVICE PairEvaluatorAGCMS(Scalar _rsq, Scalar _rcutsq, const param_type& _params) - : PairEvaluator(_rsq, _rcutsq) - { - w = _params.w; - sigma = _params.sigma; - a = _params.a; - q = _params.q; - } - - //! Evaluate the force and energy - /*! - * \param force_divr Holds the computed force divided by r - * \param pair_eng Holds the computed pair energy - * \param energy_shift If true, the potential is shifted to zero at the cutoff - * - * \returns True if the energy calculation occurs - * - * The calculation does not occur if the pair distance is greater than the cutoff - * or if the potential is scaled to zero. - */ - - DEVICE bool evalForceAndEnergy(Scalar& force_divr, Scalar& pair_eng, bool energy_shift) - { - // compute the force divided by r in force_divr - if (rsq < rcutsq && a != Scalar(0)) - { - Scalar r = fast::sqrt(rsq); - Scalar rinv = 1 / r; - Scalar w_rinv_shifted = w / (r - sigma + w); - Scalar core_repuls = pow(w_rinv_shifted, q); - Scalar exponent = exp(q * (r - sigma - w) / w); - force_divr = -a * rinv - * (q * core_repuls * w_rinv_shifted - - q * exponent / (w * pow(1 + exponent, 2.0))); - pair_eng = a * (1 / (1 + exponent) - core_repuls); - return true; - } - else - return false; - } - -#ifndef __HIPCC__ - //! Get the name of this potential - /*! \returns The potential name. Must be short and all lowercase, as this is the name energies - will be logged as via analyze.log. - */ - static std::string getName() - { - return std::string("agcms"); - } -#endif - - protected: - Scalar w; //!< well width - Scalar sigma; //!< Hard core repulsion distance - Scalar a; //!< well depth - Scalar q; //!< Steepness - }; - - } // end namespace detail - } // end namespace azplugins - } // end namespace hoomd - -#undef DEVICE - -#endif // AZPLUGINS_PAIR_EVALUATOR_AGCMS_H_