RDF Normalization Mode Argument

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

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<float>(m_n_query_points - 1) / static_cast<float>(m_n_query_points);
     }

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<float> m_pcf; //!< The computed pair correlation function.
     util::ManagedArray<float>

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:

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')``

doc/source/index.rst

@@ -10,6 +10,7 @@ Table of Contents
    gettingstarted/introduction
    gettingstarted/installation
    gettingstarted/quickstart
+   gettingstarted/migration
    gettingstarted/tutorial
    gettingstarted/examples
 

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
 

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]*,

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.

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):