diff --git a/ChangeLog.md b/ChangeLog.md index 2f041baaa..baf1e8a4f 100644 --- a/ChangeLog.md +++ b/ChangeLog.md @@ -4,6 +4,11 @@ The format is based on and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html). +## v3.0.0 -- YYYY-MM-DD + +### Changed +* The `normalize` argument to `freud.density.RDF` is now `normalization_mode`. + ## v2.11.0 -- 2022-08-09 ### Added diff --git a/cpp/density/RDF.cc b/cpp/density/RDF.cc index 310e18822..ab09fc5dd 100644 --- a/cpp/density/RDF.cc +++ b/cpp/density/RDF.cc @@ -11,8 +11,8 @@ namespace freud { namespace density { -RDF::RDF(unsigned int bins, float r_max, float r_min, bool normalize) - : BondHistogramCompute(), m_normalize(normalize) +RDF::RDF(unsigned int bins, float r_max, float r_min, NormalizationMode normalization_mode) + : BondHistogramCompute(), m_norm_mode(normalization_mode) { if (bins == 0) { @@ -59,7 +59,7 @@ void RDF::reduce() // Define prefactors with appropriate types to simplify and speed later code. float number_density = float(m_n_query_points) / m_box.getVolume(); - if (m_normalize) + if (m_norm_mode == NormalizationMode::finite_size) { number_density *= static_cast(m_n_query_points - 1) / static_cast(m_n_query_points); } diff --git a/cpp/density/RDF.h b/cpp/density/RDF.h index 50af6cbb9..fbf0b27a1 100644 --- a/cpp/density/RDF.h +++ b/cpp/density/RDF.h @@ -16,8 +16,16 @@ namespace freud { namespace density { class RDF : public locality::BondHistogramCompute { public: + //! Enum for each normalization mode + enum NormalizationMode + { + exact, + finite_size + }; + //! Constructor - RDF(unsigned int bins, float r_max, float r_min = 0, bool normalize = false); + RDF(unsigned int bins, float r_max, float r_min = 0, + NormalizationMode normalization_mode = NormalizationMode::exact); //! Destructor ~RDF() override = default; @@ -51,7 +59,7 @@ class RDF : public locality::BondHistogramCompute } private: - bool m_normalize; //!< Whether to enforce that the RDF should tend to 1 (instead of + NormalizationMode m_norm_mode; //!< Whether to enforce that the RDF should tend to 1 (instead of //!< num_query_points/num_points). util::ManagedArray m_pcf; //!< The computed pair correlation function. util::ManagedArray diff --git a/doc/source/conf.py b/doc/source/conf.py index 8ad6d30bb..9979531ff 100644 --- a/doc/source/conf.py +++ b/doc/source/conf.py @@ -78,7 +78,7 @@ # # This is also used if you do content translation via gettext catalogs. # Usually you set "language" from the command line for these cases. -language = None +language = "en" # There are two options for replacing |today|: either, you set today to some # non-false value, then it is used: diff --git a/doc/source/gettingstarted/migration.rst b/doc/source/gettingstarted/migration.rst new file mode 100644 index 000000000..2fadd59bc --- /dev/null +++ b/doc/source/gettingstarted/migration.rst @@ -0,0 +1,25 @@ +.. _migration: + +============================ +Migration to freud Version 3 +============================ + +Version 3 of the freud library introduces a few breaking changes to make the API +more intuitive and facilitate future development. This guide explains how to +alter scripts which use freud v2 APIs so they can be used with freud v3. + +Overview of API Changes +======================= + +.. list-table:: + :header-rows: 1 + + * - Goal + - v2 API + - v3 API + * - Use default RDF normalization. + - ``freud.density.RDF(..., normalize=False)`` + - ``freud.density.RDF(..., normalization_mode='exact')`` + * - Normalize small system RDFs to 1. + - ``freud.density.RDF(..., normalize=True)`` + - ``freud.density.RDF(..., normalization_mode='finite_size')`` diff --git a/doc/source/index.rst b/doc/source/index.rst index cc48d6c60..19e969623 100644 --- a/doc/source/index.rst +++ b/doc/source/index.rst @@ -10,6 +10,7 @@ Table of Contents gettingstarted/introduction gettingstarted/installation gettingstarted/quickstart + gettingstarted/migration gettingstarted/tutorial gettingstarted/examples diff --git a/doc/source/reference/credits.rst b/doc/source/reference/credits.rst index 81fb08b14..e63ad2bcc 100644 --- a/doc/source/reference/credits.rst +++ b/doc/source/reference/credits.rst @@ -332,6 +332,7 @@ Tommy Waltmann * ``DiffractionPattern`` now raises an error when used with non-cubic boxes. * Implement ``StaticStructureFactorDebye`` for 2D systems. * Add support for compilation with the C++17 standard. +* Update and test the ``normalization_mode`` argument to ``freud.density.RDF`` class. Maya Martirossyan diff --git a/freud/_density.pxd b/freud/_density.pxd index d2654f5d6..3eb0b55a1 100644 --- a/freud/_density.pxd +++ b/freud/_density.pxd @@ -49,7 +49,12 @@ cdef extern from "LocalDensity.h" namespace "freud::density": cdef extern from "RDF.h" namespace "freud::density": cdef cppclass RDF(BondHistogramCompute): - RDF(float, float, float, bool) except + + + ctypedef enum NormalizationMode "NormalizationMode": + exact "freud::density::RDF::NormalizationMode::exact" + finite_size "freud::density::RDF::NormalizationMode::finite_size" + + RDF(float, float, float, NormalizationMode) except + const freud._box.Box & getBox() const void accumulate(const freud._locality.NeighborQuery*, const vec3[float]*, diff --git a/freud/density.pyx b/freud/density.pyx index d3c58cdb5..120fbf139 100644 --- a/freud/density.pyx +++ b/freud/density.pyx @@ -576,10 +576,21 @@ cdef class RDF(_SpatialHistogram1D): behavior of :math:`\lim_{r \to \infty} g(r)=1` is recovered. However, for very small systems the long range behavior of the radial distribution will instead tend to :math:`\frac{N-1}{N}`. In small systems, where this - deviation is noticeable, the ``normalize`` flag may be used to rescale the - results and force the long range behavior to 1. Note that this option will - have little to no effect on larger systems (for example, for systems of 100 - particles the RDF will differ by 1%). + deviation is noticeable, the ``normalization_mode`` argument may be used to + rescale the results and force the long range behavior to 1. Note that this + option will have little to no effect on larger systems (for example, for + systems of 100 particles the RDF will differ by 1%). + + .. note:: + For correct normalization behavior, let the set of points be either: 1) + the same as the set of query points or 2) completely disjoint from the + set of query points (points shouldn't contain any particles in query + points). + + .. note:: + For correct normalization behavior when using + ``normalization_mode='finite_size'``, the ``points`` _must_ be the same + as the ``query_points`` and ``exclude_ii`` must be set to ``False``. .. note:: **2D:** :class:`freud.density.RDF` properly handles 2D boxes. @@ -593,27 +604,23 @@ cdef class RDF(_SpatialHistogram1D): r_min (float, optional): Minimum interparticle distance to include in the calculation (Default value = :code:`0`). - normalize (bool, optional): - Scale the RDF values by - :math:`\frac{N_{query\_points}}{N_{query\_points}-1}`. This - argument primarily exists to deal with standard RDF calculations - where no special ``query_points`` or ``neighbors`` are provided, - but where the number of ``query_points`` is small enough that the - long-ranged limit of :math:`g(r)` deviates significantly from - :math:`1`. It should not be used if :code:`query_points` is - provided as a different set of points, or if unusual query - arguments are provided to :meth:`~.compute`, specifically if - :code:`exclude_ii` is set to :code:`False`. This normalization is - not meaningful in such cases and will simply convolute the data. - + normalization_mode (str, optional): + There are two valid string inputs for this argument. The first + option, ``exact``, handles the normalization as shown mathematically + at the beginning of this class's docstring. The other option, + ``finite_size``, adds an extra rescaling factor of + :math:`\frac{N_{query\_points}}{N_{query\_ponts} - 1}` so the RDF + values will tend to 1 at large :math:`r` for small systems (Default + value = :code:`'exact'`). """ cdef freud._density.RDF * thisptr def __cinit__(self, unsigned int bins, float r_max, float r_min=0, - normalize=False): + normalization_mode='exact'): + norm_mode = self._validate_normalization_mode(normalization_mode) if type(self) == RDF: self.thisptr = self.histptr = new freud._density.RDF( - bins, r_max, r_min, normalize) + bins, r_max, r_min, norm_mode) # r_max is left as an attribute rather than a property for now # since that change needs to happen at the _SpatialHistogram level @@ -624,6 +631,15 @@ cdef class RDF(_SpatialHistogram1D): if type(self) == RDF: del self.thisptr + def _validate_normalization_mode(self, mode): + """Ensure the normalization mode is one of the approved values.""" + if mode == 'exact': + return freud._density.RDF.NormalizationMode.exact + elif mode == 'finite_size': + return freud._density.RDF.NormalizationMode.finite_size + else: + raise ValueError(f"invalid input {mode} for normalization_mode") + def compute(self, system, query_points=None, neighbors=None, reset=True): r"""Calculates the RDF and adds to the current RDF histogram. diff --git a/tests/test_density_RDF.py b/tests/test_density_RDF.py index f1f7bcfcc..b00e4a587 100644 --- a/tests/test_density_RDF.py +++ b/tests/test_density_RDF.py @@ -71,6 +71,19 @@ def test_invalid_rdf(self): freud.density.RDF(r_max=1, bins=10, r_min=2) with pytest.raises(ValueError): freud.density.RDF(r_max=1, bins=10, r_min=-1) + with pytest.raises(ValueError): + freud.density.RDF(r_max=1, bins=10, r_min=-1, normalization_mode="blah") + + @pytest.mark.parametrize("mode", ["exact", "finite_size"]) + def test_normalization_mode(self, mode): + """Make sure RDF can be computed with different normalization modes.""" + r_max = 10.0 + bins = 10 + num_points = 100 + box_size = r_max * 3.1 + sys = freud.data.make_random_system(box_size, num_points, is2D=True) + rdf = freud.density.RDF(r_max=r_max, bins=bins, normalization_mode=mode) + rdf.compute(sys) @pytest.mark.parametrize("r_min", [0, 0.1, 3.0]) def test_random_point(self, r_min):