From 161eada917ae643dbceeea060500e0e51eb3ae96 Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Tue, 21 May 2024 14:57:26 -0400 Subject: [PATCH 01/25] write skeleton cmake code for building with pybind11 --- CMakeLists.txt | 3 ++ cpp/CMakeLists.txt | 40 ++++++++++++---- cpp/box/CMakeLists.txt | 23 ++++++++++ cpp/box/module-box.cc | 0 freud/CMakeLists.txt | 101 ++++++----------------------------------- 5 files changed, 70 insertions(+), 97 deletions(-) create mode 100644 cpp/box/CMakeLists.txt create mode 100644 cpp/box/module-box.cc diff --git a/CMakeLists.txt b/CMakeLists.txt index 06aca53ad..fdc51e44a 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -36,6 +36,9 @@ if(TBB_FOUND) "[${TBB_LIBRARY}][${TBB_INCLUDE_DIR}]") endif() +# go find pybind11 +find_package(pybind11 CONFIG REQUIRED) + # Fail fast if users have not cloned submodules. if(NOT WIN32) string(ASCII 27 Esc) diff --git a/cpp/CMakeLists.txt b/cpp/CMakeLists.txt index 0d8750ce5..4c4c2c1c3 100644 --- a/cpp/CMakeLists.txt +++ b/cpp/CMakeLists.txt @@ -6,16 +6,36 @@ if(WIN32) add_compile_options(/DNOMINMAX) endif() -add_subdirectory(cluster) -add_subdirectory(density) -add_subdirectory(diffraction) -add_subdirectory(environment) -add_subdirectory(locality) -add_subdirectory(order) -add_subdirectory(parallel) -add_subdirectory(pmft) -add_subdirectory(util) +# Detect when building against a conda environment set the _using_conda variable +# for use both in this file and in the parent +get_filename_component(_python_bin_dir ${PYTHON_EXECUTABLE} DIRECTORY) +if(EXISTS "${_python_bin_dir}/../conda-meta") + message("-- Detected conda environment, setting INSTALL_RPATH_USE_LINK_PATH") + set(_using_conda On) + set(_using_conda + On + PARENT_SCOPE) +else() + set(_using_conda Off) + set(_using_conda + Off + PARENT_SCOPE) +endif() + +add_subdirectory(box) +# commented out for now, uncomment them as the conversion to pybind11 progresses +#add_subdirectory(cluster) +#add_subdirectory(density) +#add_subdirectory(diffraction) +#add_subdirectory(environment) +#add_subdirectory(locality) +#add_subdirectory(order) +#add_subdirectory(parallel) +#add_subdirectory(pmft) +#add_subdirectory(util) + +#[[ add_library( libfreud SHARED $ @@ -53,3 +73,5 @@ if(CMAKE_EXPORT_COMPILE_COMMANDS) COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_BINARY_DIR}/compile_commands.json ${PROJECT_SOURCE_DIR}) endif() +]] + diff --git a/cpp/box/CMakeLists.txt b/cpp/box/CMakeLists.txt new file mode 100644 index 000000000..48927769b --- /dev/null +++ b/cpp/box/CMakeLists.txt @@ -0,0 +1,23 @@ +# create the target and add the needed properties +pybind11_add_module(box SHARED module-box.cc) +# this probably isn't needed for box, but may be needed by other modules +#target_compile_definitions( +# box +# # Avoid deprecation warnings for unsupported NumPy API versions. See +# # https://numpy.org/doc/1.19/reference/c-api/deprecations.html +# PRIVATE "NPY_NO_DEPRECATED_API=NPY_1_10_API_VERSION" +# # Default voro++ verbosity is high. +# PRIVATE "VOROPP_VERBOSE=1") + +# install +install(TARGETS box DESTINATION freud) +if(APPLE) + set_target_properties(box PROPERTIES INSTALL_RPATH "@loader_path") +else() + set_target_properties(box PROPERTIES INSTALL_RPATH "\$ORIGIN") +endif() + +if(_using_conda OR DEFINED ENV{CIBUILDWHEEL}) + set_target_properties(box + PROPERTIES INSTALL_RPATH_USE_LINK_PATH True) +endif() diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc new file mode 100644 index 000000000..e69de29bb diff --git a/freud/CMakeLists.txt b/freud/CMakeLists.txt index d1e800d84..7e258c3c9 100644 --- a/freud/CMakeLists.txt +++ b/freud/CMakeLists.txt @@ -1,21 +1,8 @@ # Need to figure out coverage before CYTHON_FLAGS are set. option(COVERAGE "Enable coverage" OFF) -# Cython flags must be set before we run find_package for Cython since the -# compiler command is created immediately. -set(CYTHON_FLAGS - "--directive binding=True,boundscheck=False,wraparound=False,embedsignature=True,always_allow_keywords=True" - CACHE STRING "The directives for Cython compilation.") - -if(COVERAGE) - set(CYTHON_FLAGS - "${CYTHON_FLAGS},linetrace=True" - CACHE STRING "The directives for Cython compilation." FORCE) -endif() - find_package(PythonLibs) find_package(PythonExtensions REQUIRED) -find_package(Cython REQUIRED) find_package(NumPy REQUIRED) include_directories(${NumPy_INCLUDE_DIRS}) @@ -33,84 +20,22 @@ if(${PYTHON_VERSION_MAJOR} EQUAL 3 add_compile_options("-Wno-deprecated-declarations") endif() -# Detect when building against a conda environment set the _using_conda variable -# for use both in this file and in the parent -get_filename_component(_python_bin_dir ${PYTHON_EXECUTABLE} DIRECTORY) -if(EXISTS "${_python_bin_dir}/../conda-meta") - message("-- Detected conda environment, setting INSTALL_RPATH_USE_LINK_PATH") - set(_using_conda On) - set(_using_conda - On - PARENT_SCOPE) -else() - set(_using_conda Off) - set(_using_conda - Off - PARENT_SCOPE) -endif() - -set(cython_modules_with_cpp - box - cluster - density - diffraction - environment - locality - order - parallel - pmft) - -set(cython_modules_without_cpp interface msd util) - -foreach(cython_module ${cython_modules_with_cpp} ${cython_modules_without_cpp}) - add_cython_target(${cython_module} PY3 CXX) - add_library(${cython_module} SHARED ${${cython_module}}) - make_python_extension_module(${cython_module}) - - target_compile_definitions( - ${cython_module} - # Avoid deprecation warnings for unsupported NumPy API versions. See - # https://numpy.org/doc/1.19/reference/c-api/deprecations.html - PRIVATE "NPY_NO_DEPRECATED_API=NPY_1_10_API_VERSION" - # Default voro++ verbosity is high. - PRIVATE "VOROPP_VERBOSE=1") - if(COVERAGE) - target_compile_definitions( - ${cython_module} # Enable line tracing for coverage purposes if requested. - PRIVATE "CYTHON_TRACE_NOGIL=1") - endif() - - target_link_libraries(${cython_module} libfreud) - - # Separate logic required for targets with C++ code. - if("${cython_module}" IN_LIST cython_modules_with_cpp) - - target_include_directories( - ${cython_module} PRIVATE ${PROJECT_SOURCE_DIR}/cpp/${cython_module}) - endif() - - install(TARGETS ${cython_module} DESTINATION freud) - - # Coverage requires the Cython-compiled C++ files for line coverage. - if(COVERAGE) - install(FILES ${${cython_module}} DESTINATION freud) - endif() - - if(APPLE) - set_target_properties(${cython_module} PROPERTIES INSTALL_RPATH - "@loader_path") - else() - set_target_properties(${cython_module} PROPERTIES INSTALL_RPATH "\$ORIGIN") - endif() +set(python_files + box.py) + #cluster.py + #density.py + #diffraction.py + #environment.py + #locality.py + #order.py + #parallel.py + #pmft.py) - if(_using_conda OR DEFINED ENV{CIBUILDWHEEL}) - set_target_properties(${cython_module} - PROPERTIES INSTALL_RPATH_USE_LINK_PATH True) - endif() -endforeach() +install(FILES python_files DESTINATION freud) +# THIS WILL NEED TO BE MOVED TO THE BUILD OF THE ORDER MODULE # The SolidLiquid class has an instance of cluster::Cluster as a member, so # including the header requires the Cluster.h header. Would prefer to inherit # this information from the _order library, but that's not possible since we're # linking to libfreud. -target_include_directories(order PUBLIC ${PROJECT_SOURCE_DIR}/cpp/cluster) +#target_include_directories(order PUBLIC ${PROJECT_SOURCE_DIR}/cpp/cluster) From a9c5eed16e483558e2de8a01f81b21991b2f7f91 Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Tue, 28 May 2024 12:24:06 -0400 Subject: [PATCH 02/25] keep working on skeleton code --- cpp/box/CMakeLists.txt | 1 + cpp/box/module-box.cc | 42 ++ freud/box.py | 990 +++++++++++++++++++++++++++++++++++++++++ 3 files changed, 1033 insertions(+) create mode 100644 freud/box.py diff --git a/cpp/box/CMakeLists.txt b/cpp/box/CMakeLists.txt index 48927769b..23cfd8ba3 100644 --- a/cpp/box/CMakeLists.txt +++ b/cpp/box/CMakeLists.txt @@ -1,5 +1,6 @@ # create the target and add the needed properties pybind11_add_module(box SHARED module-box.cc) +target_link_libraries(box PUBLIC TBB::tbb) # this probably isn't needed for box, but may be needed by other modules #target_compile_definitions( # box diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc index e69de29bb..05e6e5fe6 100644 --- a/cpp/box/module-box.cc +++ b/cpp/box/module-box.cc @@ -0,0 +1,42 @@ + +#include + +#include "Box.h" + +PYBIND11_MODULE(_box, m) +{ + pybind11::class_>(m, "Box") + // constructors + .def(pybind11::init()) + // getters and setters + .def("getLx", &Box::getLx) + .def("getLy", &Box::getLy) + .def("getLz", &Box::getLz) + .def("setL", &Box::setL) + .def("getTiltFactorXY", &Box::getTiltFactorXY) + .def("getTiltFactorXZ", &Box::getTiltFactorXZ) + .def("getTiltFactorYZ", &Box::getTiltFactorYZ) + .def("setTiltFactorXY", &Box::setTiltFactorXY) + .def("setTiltFactorXZ", &Box::setTiltFactorXZ) + .def("setTiltFactorYZ", &Box::setTiltFactorYZ) + .def("getPeriodicX", &Box::getPeriodicX) + .def("getPeriodicY", &Box::getPeriodicY) + .def("getPeriodicZ", &Box::getPeriodicZ) + .def("setPeriodic", &Box::setPeriodic) + .def("setPeriodicX", &Box::setPeriodicX) + .def("setPeriodicY", &Box::setPeriodicY) + .def("setPeriodicZ", &Box::setPeriodicZ) + .def("get2D", &Box::get2D) + .def("set2D", &Box::set2D) + .def("getVolume", &Box::getVolume) + .def("center", &Box::center) + .def("centerOfMass", &Box::centerOfMass) + // other stuff + .def("makeAbsolute", &Box::makeAbsolute) + .def("makeFractional", &Box::makeFractional) + .def("wrap", &Box::wrap) + .def("unwrap", &Box::unwrap) + .def("computeDistances", &Box::computeDistances) + .def("computeAllDistances", &Box::computeAllDistances) + .def("contains", &Box::contains); +} diff --git a/freud/box.py b/freud/box.py new file mode 100644 index 000000000..7fc88b0cc --- /dev/null +++ b/freud/box.py @@ -0,0 +1,990 @@ +# Copyright (c) 2010-2023 The Regents of the University of Michigan +# This file is from the freud project, released under the BSD 3-Clause License. + +r""" +The :class:`~.Box` class defines the geometry of a simulation box. The class +natively supports periodicity by providing the fundamental features for +wrapping vectors outside the box back into it. +""" + +import logging +import warnings + +import numpy as np + +import freud.util + +import freud._box + +logger = logging.getLogger(__name__) + + +class Box: + r"""The freud Box class for simulation boxes. + + This class defines an arbitrary triclinic geometry within which all points are confined. + By convention, the freud Box is centered at the origin (``[0, 0, 0]``), + with the extent in each dimension described by the half-open interval ``[-L/2, L/2)``. + For more information, see the `documentation + `__ + on boxes and periodic boundary conditions. + + Also available as ``freud.Box``. + + Args: + Lx (float, optional): + The x-dimension length. + Ly (float, optional): + The y-dimension length. + Lz (float, optional): + The z-dimension length (Default value = 0). + xy (float, optional): + The xy tilt factor (Default value = 0). + xz (float, optional): + The xz tilt factor (Default value = 0). + yz (float, optional): + The yz tilt factor (Default value = 0). + is2D (bool, optional): + Whether the box is 2-dimensional. Uses :code:`Lz == 0` + if :code:`None`. (Default value = :code:`None`) + """ # noqa: E501 + + def __cinit__(self, Lx, Ly, Lz=0, xy=0, xz=0, yz=0, is2D=None): + if is2D is None: + is2D = (Lz == 0) + if is2D: + if not (Lx and Ly): + raise ValueError("Lx and Ly must be nonzero for 2D boxes.") + elif Lz != 0 or xz != 0 or yz != 0: + warnings.warn( + "Specifying z-dimensions in a 2-dimensional box " + "has no effect!") + else: + if not (Lx and Ly and Lz): + raise ValueError( + "Lx, Ly, and Lz must be nonzero for 3D boxes.") + self._cpp_obj = freud._box.Box(Lx, Ly, Lz, xy, xz, yz, is2D) + + @property + def L(self): + r""":math:`\left(3, \right)` :class:`numpy.ndarray`: Get or set the + box lengths along x, y, and z.""" + return np.asarray([self._cpp_obj.getLx(), + self._cpp_obj.getLy(), + self._cpp_obj.getLz()]) + + @L.setter + def L(self, value): + try: + if len(value) != 3: + raise ValueError('setL must be called with a scalar or a list ' + 'of length 3.') + except TypeError: + # Will fail if object has no length + value = (value, value, value) + + if self.is2D and value[2] != 0: + warnings.warn( + "Specifying z-dimensions in a 2-dimensional box " + "has no effect!") + self._cpp_obj.setL(value[0], value[1], value[2]) + + @property + def Lx(self): + """float: Get or set the x-dimension length.""" + return self._cpp_obj.getLx() + + @Lx.setter + def Lx(self, value): + self.L = [value, self.Ly, self.Lz] + + @property + def Ly(self): + """float: Get or set the y-dimension length.""" + return self._cpp_obj.getLy() + + @Ly.setter + def Ly(self, value): + self.L = [self.Lx, value, self.Lz] + + @property + def Lz(self): + """float: Get or set the z-dimension length.""" + return self._cpp_obj.getLz() + + @Lz.setter + def Lz(self, value): + self.L = [self.Lx, self.Ly, value] + + @property + def xy(self): + """float: Get or set the xy tilt factor.""" + return self._cpp_obj.getTiltFactorXY() + + @xy.setter + def xy(self, value): + self._cpp_obj.setTiltFactorXY(value) + + @property + def xz(self): + """float: Get or set the xz tilt factor.""" + return self._cpp_obj.getTiltFactorXZ() + + @xz.setter + def xz(self, value): + self._cpp_obj.setTiltFactorXZ(value) + + @property + def yz(self): + """float: Get or set the yz tilt factor.""" + return self._cpp_obj.getTiltFactorYZ() + + @yz.setter + def yz(self, value): + self._cpp_obj.setTiltFactorYZ(value) + + @property + def dimensions(self): + """int: Get or set the number of dimensions (2 or 3).""" + return 2 if self.is2D else 3 + + @dimensions.setter + def dimensions(self, value): + assert value == 2 or value == 3 + self._cpp_obj.set2D(bool(value == 2)) + + @property + def is2D(self): + """bool: Whether the box is 2D.""" + return self._cpp_obj.get2D() + + @property + def L_inv(self): + r""":math:`\left(3, \right)` :class:`numpy.ndarray`: The inverse box + lengths.""" + cdef vec3[float] result = self._cpp_obj.getLinv() + return np.asarray([result.x, result.y, result.z]) + + @property + def volume(self): + """float: The box volume (area in 2D).""" + return self._cpp_obj.getVolume() + + def make_absolute(self, fractional_coordinates, out=None): + r"""Convert fractional coordinates into absolute coordinates. + + Args: + fractional_coordinates (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Fractional coordinate vector(s), between 0 and 1 within + parallelepipedal box. + out (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray` or :code:`None`): + The array in which to place the absolute coordinates. It must + be of dtype :attr:`numpy.float32`. If ``None``, this function + will return a newly allocated array (Default value = None). + + Returns: + :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: + Absolute coordinate vector(s). If ``out`` is provided, a + reference to it is returned. + """ # noqa: E501 + fractions = np.asarray(fractional_coordinates).copy() + flatten = fractions.ndim == 1 + fractions = np.atleast_2d(fractions) + fractions = freud.util._convert_array(fractions, shape=(None, 3)) + out = freud.util._convert_array( + out, shape=fractions.shape, allow_copy=False) + + cdef const float[:, ::1] l_points = fractions + cdef unsigned int Np = l_points.shape[0] + cdef float[:, ::1] l_out = out + + self._cpp_obj.makeAbsolute( + &l_points[0, 0], Np, + &l_out[0, 0]) + + return np.squeeze(out) if flatten else out + + def make_fractional(self, absolute_coordinates, out=None): + r"""Convert absolute coordinates into fractional coordinates. + + Args: + absolute_coordinates (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Absolute coordinate vector(s). + out (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray` or :code:`None`): + The array in which to place the fractional positions. It must + be of dtype :attr:`numpy.float32`. If ``None``, this function + will return a newly allocated array (Default value = None). + + Returns: + :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: + Fractional coordinate vector(s). If ``out`` is provided, a + reference to it is returned. + """ # noqa: E501 + vecs = np.asarray(absolute_coordinates).copy() + flatten = vecs.ndim == 1 + vecs = np.atleast_2d(vecs) + vecs = freud.util._convert_array(vecs, shape=(None, 3)) + out = freud.util._convert_array( + out, shape=vecs.shape, allow_copy=False) + + cdef const float[:, ::1] l_points = vecs + cdef unsigned int Np = l_points.shape[0] + cdef float[:, ::1] l_out = out + + self._cpp_obj.makeFractional( + &l_points[0, 0], Np, + &l_out[0, 0]) + + return np.squeeze(out) if flatten else out + + def get_images(self, vecs): + r"""Returns the images corresponding to unwrapped vectors. + + Args: + vecs (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Coordinates of unwrapped vector(s). + + Returns: + :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: + Image index vector(s). + """ # noqa: E501 + vecs = np.asarray(vecs) + flatten = vecs.ndim == 1 + vecs = np.atleast_2d(vecs) + vecs = freud.util._convert_array(vecs, shape=(None, 3)) + + images = np.zeros(vecs.shape, dtype=np.int32) + cdef const float[:, ::1] l_points = vecs + cdef const int[:, ::1] l_result = images + cdef unsigned int Np = l_points.shape[0] + self._cpp_obj.getImages( &l_points[0, 0], Np, + &l_result[0, 0]) + + return np.squeeze(images) if flatten else images + + def get_box_vector(self, i): + r"""Get the box vector with index :math:`i`. + + Args: + i (unsigned int): + Index (:math:`0 \leq i < d`) of the box vector, where + :math:`d` is the box dimension (2 or 3). + + Returns: + :math:`\left(3, \right)` :class:`numpy.ndarray`: + Box vector with index :math:`i`. + """ + return self.to_matrix()[:, i] + + @property + def v1(self): + """:math:`(3, )` :class:`np.ndarray`: The first box vector + :math:`(L_x, 0, 0)`.""" + return self.get_box_vector(0) + + @property + def v2(self): + r""":math:`(3, )` :class:`np.ndarray`: The second box vector + :math:`(xy \times L_y, L_y, 0)`.""" + return self.get_box_vector(1) + + @property + def v3(self): + r""":math:`(3, )` :class:`np.ndarray`: The third box vector + :math:`(xz \times L_z, yz \times L_z, L_z)`.""" + return self.get_box_vector(2) + + def wrap(self, vecs, out=None): + r"""Wrap an array of vectors into the box, using periodic boundaries. + + .. note:: Since the origin of the box is in the center, wrapping is + equivalent to applying the minimum image convention to the + input vectors. + + Args: + vecs (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Unwrapped vector(s). + out (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray` or :code:`None`): + The array in which to place the wrapped vectors. It must be of + dtype :attr:`numpy.float32`. If ``None``, this function will + return a newly allocated array (Default value = None). + + Returns: + :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: + Vector(s) wrapped into the box. If ``out`` is provided, a + reference to it is returned. + """ # noqa: E501 + vecs = np.asarray(vecs) + flatten = vecs.ndim == 1 + vecs = np.atleast_2d(vecs) + vecs = freud.util._convert_array(vecs, shape=(None, 3)) + out = freud.util._convert_array( + out, shape=vecs.shape, allow_copy=False) + + cdef const float[:, ::1] l_points = vecs + cdef unsigned int Np = l_points.shape[0] + cdef float[:, ::1] l_out = out + + self._cpp_obj.wrap( &l_points[0, 0], + Np, &l_out[0, 0]) + + return np.squeeze(out) if flatten else out + + def unwrap(self, vecs, imgs, out=None): + r"""Unwrap an array of vectors inside the box back into real space, + using an array of image indices that determine how many times to unwrap + in each dimension. + + Args: + vecs (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Vector(s) to be unwrapped. + imgs (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Image indices for vector(s). + out (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray` or :code:`None`): + The array in which to place the unwrapped vectors. It must be + of dtype :attr:`numpy.float32`. If ``None``, this function + will return a newly allocated array (Default value = None). + + Returns: + :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: + Unwrapped vector(s). If ``out`` is provided, a reference to it is + returned. + """ # noqa: E501 + vecs = np.asarray(vecs) + flatten = vecs.ndim == 1 + vecs = np.atleast_2d(vecs) + imgs = np.atleast_2d(imgs) + if vecs.shape[0] != imgs.shape[0]: + # Broadcasts (1, 3) to (N, 3) for both arrays + vecs, imgs = np.broadcast_arrays(vecs, imgs) + vecs = freud.util._convert_array(vecs, shape=(None, 3)).copy() + imgs = freud.util._convert_array( + imgs, shape=vecs.shape, dtype=np.int32) + out = freud.util._convert_array( + out, shape=vecs.shape, allow_copy=False) + + cdef const float[:, ::1] l_points = vecs + cdef const int[:, ::1] l_imgs = imgs + cdef unsigned int Np = l_points.shape[0] + cdef float[:, ::1] l_out = out + + self._cpp_obj.unwrap( &l_points[0, 0], + &l_imgs[0, 0], Np, + &l_out[0, 0]) + + return np.squeeze(out) if flatten else out + + def center_of_mass(self, vecs, masses=None): + r"""Compute center of mass of an array of vectors, using periodic boundaries. + + This calculation accounts for periodic images. `This Wikipedia page + `__ + describes the mathematics of this method. + + Example:: + + >>> import freud + >>> import numpy as np + >>> box = freud.Box.cube(10) + >>> points = [[-1, -1, 0], [-1, 1, 0], [2, 0, 0]] + >>> np.mean(points, axis=0) # Does not account for periodic images + array([0., 0., 0.]) + >>> box.center_of_mass(points) # Accounts for periodic images + array([-0.1845932, 0. , 0. ]) + + Args: + vecs (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Vectors used to find center of mass. + masses (:math:`\left(N,\right)` :class:`numpy.ndarray`): + Masses corresponding to each vector, defaulting to 1 if not + provided or :code:`None` (Default value = :code:`None`). + + Returns: + :math:`\left(3, \right)` :class:`numpy.ndarray`: + Center of mass. + """ # noqa: E501 + vecs = freud.util._convert_array(vecs, shape=(None, 3)) + cdef const float[:, ::1] l_points = vecs + + cdef float* l_masses_ptr = NULL + cdef float[::1] l_masses + if masses is not None: + l_masses = freud.util._convert_array(masses, shape=(len(vecs), )) + l_masses_ptr = &l_masses[0] + + cdef size_t Np = l_points.shape[0] + cdef vec3[float] result = self._cpp_obj.centerOfMass( + &l_points[0, 0], Np, l_masses_ptr) + return np.asarray([result.x, result.y, result.z]) + + def center(self, vecs, masses=None): + r"""Subtract center of mass from an array of vectors, using periodic boundaries. + + This calculation accounts for periodic images. `This Wikipedia page + `__ + describes the mathematics of this method. + + Example:: + + >>> import freud + >>> box = freud.Box.cube(10) + >>> points = [[-1, -1, 0], [-1, 1, 0], [2, 0, 0]] + >>> box.center(points) + array([[-0.8154068, -1.0000002, 0. ], + [-0.8154068, 1. , 0. ], + [ 2.1845937, 0. , 0. ]], dtype=float32) + + Args: + vecs (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Vectors to center. + masses (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Masses corresponding to each vector, defaulting to 1 if not + provided or :code:`None` (Default value = :code:`None`). + + Returns: + :math:`\left(N, 3\right)` :class:`numpy.ndarray`: + Vectors with center of mass subtracted. + """ # noqa: E501 + vecs = freud.util._convert_array(vecs, shape=(None, 3)).copy() + cdef const float[:, ::1] l_points = vecs + + cdef float* l_masses_ptr = NULL + cdef float[::1] l_masses + if masses is not None: + l_masses = freud.util._convert_array(masses, shape=(len(vecs), )) + l_masses_ptr = &l_masses[0] + + cdef size_t Np = l_points.shape[0] + self._cpp_obj.center( &l_points[0, 0], Np, l_masses_ptr) + return vecs + + def compute_distances(self, query_points, points): + r"""Calculate distances between two sets of points, using periodic boundaries. + + Distances are calculated row-wise, i.e. ``distances[i]`` is the + distance from ``query_points[i]`` to ``points[i]``. + + Args: + query_points (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Array of query points. + points (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Array of points. + + Returns: + :math:`\left(N, \right)` :class:`numpy.ndarray`: + Array of distances between query points and points. + """ # noqa: E501 + + query_points = freud.util._convert_array( + np.atleast_2d(query_points), shape=(None, 3)) + points = freud.util._convert_array( + np.atleast_2d(points), shape=(None, 3)) + + cdef: + const float[:, ::1] l_query_points = query_points + const float[:, ::1] l_points = points + size_t n_query_points = query_points.shape[0] + size_t n_points = points.shape[0] + float[::1] distances = np.empty( + n_query_points, dtype=np.float32) + + self._cpp_obj.computeDistances( + &l_query_points[0, 0], n_query_points, + &l_points[0, 0], n_points, + &distances[0]) + return np.asarray(distances) + + def compute_all_distances(self, query_points, points): + r"""Calculate distances between all pairs of query points and points, using periodic boundaries. + + Distances are calculated pairwise, i.e. ``distances[i, j]`` is the + distance from ``query_points[i]`` to ``points[j]``. + + Args: + query_points (:math:`\left(N_{query\_points}, 3 \right)` :class:`numpy.ndarray`): + Array of query points. + points (:math:`\left(N_{points}, 3 \right)` :class:`numpy.ndarray`): + Array of points with same length as ``query_points``. + + Returns: + :math:`\left(N_{query\_points}, N_{points}, \right)` :class:`numpy.ndarray`: + Array of distances between query points and points. + """ # noqa: E501 + query_points = freud.util._convert_array( + np.atleast_2d(query_points), shape=(None, 3)) + points = freud.util._convert_array( + np.atleast_2d(points), shape=(None, 3)) + + cdef: + const float[:, ::1] l_query_points = query_points + const float[:, ::1] l_points = points + size_t n_query_points = query_points.shape[0] + size_t n_points = points.shape[0] + float[:, ::1] distances = np.empty( + [n_query_points, n_points], dtype=np.float32) + + self._cpp_obj.computeAllDistances( + &l_query_points[0, 0], n_query_points, + &l_points[0, 0], n_points, + &distances[0, 0]) + + return np.asarray(distances) + + def contains(self, points): + r"""Compute a boolean array (mask) corresponding to point membership in a box. + + This calculation computes particle membership based on conventions + defined by :class:`Box`, ignoring periodicity. This means that in a + cubic (3D) box with dimensions ``L``, particles would be considered + inside the box if their coordinates are between ``[-L/2, L/2]``. + Particles laying at a coordinate such as ``[0, L, 0]`` would be + considered outside the box. More information about coordinate + conventions can be found in the documentation on `Using boxes + `__ + and `periodic boundary conditions + `__. + + Example:: + + >>> import freud + >>> box = freud.Box.cube(10) + >>> points = [[-4, 0, 0], [10, 0, 0], [0, -7, 0]] + >>> box.contains(points) + array([ True, False, False]) + + Args: + points (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): + Array of points. + + Returns: + :math:`\left(N, \right)` :class:`numpy.ndarray`: + Array of booleans, where ``True`` corresponds to points within + the box, and ``False`` corresponds to points outside the box. + """ # noqa: E501 + + points = freud.util._convert_array( + np.atleast_2d(points), shape=(None, 3)) + + cdef: + const float[:, ::1] l_points = points + size_t n_points = points.shape[0] + + contains_mask = freud.util._convert_array( + np.ones(n_points), dtype=bool) + cdef cpp_bool[::1] l_contains_mask = contains_mask + + self._cpp_obj.contains( + &l_points[0, 0], n_points, + &l_contains_mask[0]) + + return np.array(l_contains_mask).astype(bool) + + @property + def cubic(self): + """bool: Whether the box is a cube.""" + return ( + not self.is2D + and np.allclose( + [self.Lx, self.Lx, self.Ly, self.Ly, self.Lz, self.Lz], + [self.Ly, self.Lz, self.Lx, self.Lz, self.Lx, self.Ly], + rtol=1e-5, + atol=1e-5, + ) + and np.allclose(0, [self.xy, self.yz, self.xz], rtol=1e-5, atol=1e-5) + ) + + @property + def periodic(self): + r""":math:`\left(3, \right)` :class:`numpy.ndarray`: Get or set the + periodicity of the box in each dimension.""" + return np.asarray([self._cpp_obj.getPeriodicX(), + self._cpp_obj.getPeriodicY(), + self._cpp_obj.getPeriodicZ()]) + + @periodic.setter + def periodic(self, periodic): + # Allow passing a single value + try: + self._cpp_obj.setPeriodic(periodic[0], periodic[1], periodic[2]) + except TypeError: + # Allow single value to be passed for all directions + self._cpp_obj.setPeriodic(periodic, periodic, periodic) + + @property + def periodic_x(self): + """bool: Get or set the periodicity of the box in x.""" + return self._cpp_obj.getPeriodicX() + + @periodic_x.setter + def periodic_x(self, periodic): + self._cpp_obj.setPeriodicX(periodic) + + @property + def periodic_y(self): + """bool: Get or set the periodicity of the box in y.""" + return self._cpp_obj.getPeriodicY() + + @periodic_y.setter + def periodic_y(self, periodic): + self._cpp_obj.setPeriodicY(periodic) + + @property + def periodic_z(self): + """bool: Get or set the periodicity of the box in z.""" + return self._cpp_obj.getPeriodicZ() + + @periodic_z.setter + def periodic_z(self, periodic): + self._cpp_obj.setPeriodicZ(periodic) + + def to_dict(self): + r"""Return box as dictionary. + + Example:: + + >>> box = freud.box.Box.cube(L=10) + >>> box.to_dict() + {'Lx': 10.0, 'Ly': 10.0, 'Lz': 10.0, + 'xy': 0.0, 'xz': 0.0, 'yz': 0.0, 'dimensions': 3} + + Returns: + dict: Box parameters + """ + return { + 'Lx': self.Lx, + 'Ly': self.Ly, + 'Lz': self.Lz, + 'xy': self.xy, + 'xz': self.xz, + 'yz': self.yz, + 'dimensions': self.dimensions} + + def to_matrix(self): + r"""Returns the box matrix (3x3). + + Example:: + + >>> box = freud.box.Box.cube(L=10) + >>> box.to_matrix() + array([[10., 0., 0.], + [ 0., 10., 0.], + [ 0., 0., 10.]]) + + Returns: + :math:`\left(3, 3\right)` :class:`numpy.ndarray`: Box matrix + """ + return np.asarray([[self.Lx, self.xy * self.Ly, self.xz * self.Lz], + [0, self.Ly, self.yz * self.Lz], + [0, 0, self.Lz]]) + + def to_box_lengths_and_angles(self): + r"""Return the box lengths and angles. + + Returns: + tuple: + The box vector lengths and angles in radians + :math:`(L_1, L_2, L_3, \alpha, \beta, \gamma)`. + """ + alpha = np.arccos( + (self.xy * self.xz + self.yz) + / (np.sqrt(1 + self.xy**2) * np.sqrt(1 + self.xz**2 + self.yz**2)) + ) + beta = np.arccos(self.xz/np.sqrt(1+self.xz**2+self.yz**2)) + gamma = np.arccos(self.xy/np.sqrt(1+self.xy**2)) + L1 = self.Lx + a2 = [self.Ly*self.xy, self.Ly, 0] + a3 = [self.Lz*self.xz, self.Lz*self.yz, self.Lz] + L2 = np.linalg.norm(a2) + L3 = np.linalg.norm(a3) + return (L1, L2, L3, alpha, beta, gamma) + + def __repr__(self): + return ("freud.box.{cls}(Lx={Lx}, Ly={Ly}, Lz={Lz}, " + "xy={xy}, xz={xz}, yz={yz}, " + "is2D={is2D})").format(cls=type(self).__name__, + Lx=self.Lx, + Ly=self.Ly, + Lz=self.Lz, + xy=self.xy, + xz=self.xz, + yz=self.yz, + is2D=self.is2D) + + def __str__(self): + return repr(self) + + def __richcmp__(self, other, int op): + r"""Implement all comparisons for Cython extension classes""" + cdef Box c_other + try: + c_other = other + if op == Py_EQ: + return dereference(self._cpp_obj) == dereference(c_other._cpp_obj) + if op == Py_NE: + return dereference(self._cpp_obj) != dereference(c_other._cpp_obj) + except TypeError: + # Cython cast to Box failed + pass + return NotImplemented + + def __mul__(self, scale): + if scale > 0: + return self.__class__(Lx=self.Lx*scale, + Ly=self.Ly*scale, + Lz=self.Lz*scale, + xy=self.xy, xz=self.xz, yz=self.yz, + is2D=self.is2D) + else: + raise ValueError("Box can only be multiplied by positive values.") + + def __rmul__(self, scale): + return self * scale + + def plot(self, title=None, ax=None, image=[0, 0, 0], *args, **kwargs): + """Plot a :class:`~.box.Box` object. + + Args: + title (str): + Title of the graph. (Default value = :code:`None`). + ax (:class:`matplotlib.axes.Axes`): + Axes object to plot. If :code:`None`, make a new axes and + figure object. If plotting a 3D box, the axes must be 3D. + (Default value = :code:`None`). + image (list): + The periodic image location at which to draw the box (Default + value = :code:`[0, 0, 0]`). + *args: + Passed on to :meth:`mpl_toolkits.mplot3d.Axes3D.plot` or + :meth:`matplotlib.axes.Axes.plot`. + **kwargs: + Passed on to :meth:`mpl_toolkits.mplot3d.Axes3D.plot` or + :meth:`matplotlib.axes.Axes.plot`. + """ + import freud.plot + return freud.plot.box_plot(self, title=title, ax=ax, image=image, + *args, **kwargs) + + @classmethod + def from_box(cls, box, dimensions=None): + r"""Initialize a Box instance from a box-like object. + + Args: + box: + A box-like object + dimensions (int): + Dimensionality of the box (Default value = None) + + .. note:: Objects that can be converted to freud boxes include + lists like :code:`[Lx, Ly, Lz, xy, xz, yz]`, + dictionaries with keys + :code:`'Lx', 'Ly', 'Lz', 'xy', 'xz', 'yz', 'dimensions'`, + objects with attributes + :code:`Lx, Ly, Lz, xy, xz, yz, dimensions`, + 3x3 matrices (see :meth:`~.from_matrix`), + or existing :class:`freud.box.Box` objects. + + If any of :code:`Lz, xy, xz, yz` are not provided, they will + be set to 0. + + If all values are provided, a triclinic box will be + constructed. If only :code:`Lx, Ly, Lz` are provided, an + orthorhombic box will be constructed. If only :code:`Lx, Ly` + are provided, a rectangular (2D) box will be constructed. + + If the optional :code:`dimensions` argument is given, this + will be used as the box dimensionality. Otherwise, the box + dimensionality will be detected from the :code:`dimensions` + of the provided box. If no dimensions can be detected, the + box will be 2D if :code:`Lz == 0`, and 3D otherwise. + + Returns: + :class:`freud.box.Box`: The resulting box object. + """ + if np.asarray(box).shape == (3, 3): + # Handles 3x3 matrices + return cls.from_matrix(box, dimensions=dimensions) + try: + # Handles freud.box.Box and objects with attributes + Lx = box.Lx + Ly = box.Ly + Lz = getattr(box, 'Lz', 0) + xy = getattr(box, 'xy', 0) + xz = getattr(box, 'xz', 0) + yz = getattr(box, 'yz', 0) + if dimensions is None: + dimensions = getattr(box, 'dimensions', None) + elif dimensions != getattr(box, 'dimensions', dimensions): + raise ValueError( + "The provided dimensions argument conflicts with the " + "dimensions attribute of the provided box object.") + except AttributeError: + try: + # Handle dictionary-like + Lx = box['Lx'] + Ly = box['Ly'] + Lz = box.get('Lz', 0) + xy = box.get('xy', 0) + xz = box.get('xz', 0) + yz = box.get('yz', 0) + if dimensions is None: + dimensions = box.get('dimensions', None) + else: + if dimensions != box.get('dimensions', dimensions): + raise ValueError( + "The provided dimensions argument conflicts with " + "the dimensions attribute of the provided box " + "object.") + except (IndexError, KeyError, TypeError): + if not len(box) in [2, 3, 6]: + raise ValueError( + "List-like objects must have length 2, 3, or 6 to be " + "converted to freud.box.Box.") + # Handle list-like + Lx = box[0] + Ly = box[1] + Lz = box[2] if len(box) > 2 else 0 + xy, xz, yz = box[3:6] if len(box) == 6 else (0, 0, 0) + except: # noqa + logger.debug('Supplied box cannot be converted to type ' + 'freud.box.Box.') + raise + + # Infer dimensions if not provided. + if dimensions is None: + dimensions = 2 if Lz == 0 else 3 + is2D = (dimensions == 2) + return cls(Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, is2D=is2D) + + @classmethod + def from_matrix(cls, box_matrix, dimensions=None): + r"""Initialize a Box instance from a box matrix. + + For more information and the source for this code, see: + `HOOMD-blue's box documentation \ + `_. + + Args: + box_matrix (array-like): + A 3x3 matrix or list of lists + dimensions (int): + Number of dimensions (Default value = :code:`None`) + + Returns: + :class:`freud.box.Box`: The resulting box object. + """ + box_matrix = np.asarray(box_matrix, dtype=np.float32) + v0 = box_matrix[:, 0] + v1 = box_matrix[:, 1] + v2 = box_matrix[:, 2] + Lx = np.sqrt(np.dot(v0, v0)) + a2x = np.dot(v0, v1) / Lx + Ly = np.sqrt(np.dot(v1, v1) - a2x * a2x) + xy = a2x / Ly + v0xv1 = np.cross(v0, v1) + v0xv1mag = np.sqrt(np.dot(v0xv1, v0xv1)) + Lz = np.dot(v2, v0xv1) / v0xv1mag + if Lz != 0: + a3x = np.dot(v0, v2) / Lx + xz = a3x / Lz + yz = (np.dot(v1, v2) - a2x * a3x) / (Ly * Lz) + else: + xz = yz = 0 + if dimensions is None: + dimensions = 2 if Lz == 0 else 3 + is2D = (dimensions == 2) + return cls(Lx=Lx, Ly=Ly, Lz=Lz, + xy=xy, xz=xz, yz=yz, is2D=is2D) + + @classmethod + def cube(cls, L=None): + r"""Construct a cubic box with equal lengths. + + Args: + L (float): The edge length + + Returns: + :class:`freud.box.Box`: The resulting box object. + """ + # classmethods compiled with cython don't appear to support + # named access to positional arguments, so we keep this to + # recover the behavior + if L is None: + raise TypeError("cube() missing 1 required positional argument: L") + return cls(Lx=L, Ly=L, Lz=L, xy=0, xz=0, yz=0, is2D=False) + + @classmethod + def square(cls, L=None): + r"""Construct a 2-dimensional (square) box with equal lengths. + + Args: + L (float): The edge length + + Returns: + :class:`freud.box.Box`: The resulting box object. + """ + # classmethods compiled with cython don't appear to support + # named access to positional arguments, so we keep this to + # recover the behavior + if L is None: + raise TypeError("square() missing 1 required " + "positional argument: L") + return cls(Lx=L, Ly=L, Lz=0, xy=0, xz=0, yz=0, is2D=True) + + @classmethod + def from_box_lengths_and_angles( + cls, L1, L2, L3, alpha, beta, gamma, dimensions=None, + ): + r"""Construct a box from lengths and angles (in radians). + + All the angles provided must be between 0 and :math:`\pi`. + + Args: + L1 (float): + The length of the first lattice vector. + L2 (float): + The length of the second lattice vector. + L3 (float): + The length of the third lattice vector. + alpha (float): + The angle between second and third lattice vector in radians (must be + between 0 and :math:`\pi`). + beta (float): + The angle between first and third lattice vector in radians (must be + between 0 and :math:`\pi`). + gamma (float): + The angle between the first and second lattice vector in radians (must + be between 0 and :math:`\pi`). + dimensions (int): + The number of dimensions (Default value = :code:`None`). + + Returns: + :class:`freud.box.Box`: The resulting box object. + """ + if not 0 < alpha < np.pi: + raise ValueError("alpha must be between 0 and pi.") + if not 0 < beta < np.pi: + raise ValueError("beta must be between 0 and pi.") + if not 0 < gamma < np.pi: + raise ValueError("gamma must be between 0 and pi.") + a1 = np.array([L1, 0, 0]) + a2 = np.array([L2 * np.cos(gamma), L2 * np.sin(gamma), 0]) + a3x = np.cos(beta) + a3y = (np.cos(alpha) - np.cos(beta) * np.cos(gamma)) / np.sin(gamma) + under_sqrt = 1 - a3x**2 - a3y**2 + if under_sqrt < 0: + raise ValueError("The provided angles can not form a valid box.") + a3z = np.sqrt(under_sqrt) + a3 = np.array([L3 * a3x, L3 * a3y, L3 * a3z]) + if dimensions is None: + dimensions = 2 if L3 == 0 else 3 + return cls.from_matrix(np.array([a1, a2, a3]).T, dimensions=dimensions) + + +cdef BoxFromCPP(const freud._box.Box & cppbox): + b = Box(cppbox.getLx(), cppbox.getLy(), cppbox.getLz(), + cppbox.getTiltFactorXY(), cppbox.getTiltFactorXZ(), + cppbox.getTiltFactorYZ(), cppbox.is2D()) + b.periodic = [cppbox.getPeriodicX(), + cppbox.getPeriodicY(), + cppbox.getPeriodicZ()] + return b From 622742be460e13ffb640a1696776b67e5ea162ea Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Mon, 3 Jun 2024 12:49:56 -0400 Subject: [PATCH 03/25] start working on functions that pass array dat --- cpp/box/Box.h | 22 +++++++++++++++------- cpp/util/utils.h | 1 + freud/box.py | 22 +++++++--------------- 3 files changed, 23 insertions(+), 22 deletions(-) diff --git a/cpp/box/Box.h b/cpp/box/Box.h index 2cd72d9c8..a9a5aaca2 100644 --- a/cpp/box/Box.h +++ b/cpp/box/Box.h @@ -4,14 +4,16 @@ #ifndef BOX_H #define BOX_H -#include "utils.h" #include #include #include #include +#include "utils.h" #include "VectorMath.h" +#include + /*! \file Box.h \brief Represents simulation boxes and contains helpful wrapping functions. */ @@ -23,6 +25,8 @@ constexpr float TWO_PI = 2.0 * M_PI; namespace freud { namespace box { +using pybind11::array_t = pybind11_array; + //! Stores box dimensions and provides common routines for wrapping vectors back into the box /*! Box stores a standard HOOMD simulation box that goes from -L/2 to L/2 in each dimension, allowing Lx, Ly, Lz, and triclinic tilt factors xy, xz, and yz to be specified independently. @@ -156,9 +160,9 @@ class Box } //! Get current stored inverse of L - vec3 getLinv() const + std::array getLinv() const { - return m_Linv; + return { m_Linv.x, m_linv.y, m_Linv.z }; } //! Get tilt factor xy @@ -230,12 +234,14 @@ class Box * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void makeAbsolute(const vec3* vecs, unsigned int Nvecs, vec3* out) const + void makeAbsolute(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const { + vec3 *vecs_data = static_cast *>(vecs.data()); + vec3 *out_data = static_cast *>(out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { - out[i] = makeAbsolute(vecs[i]); + out_data[i] = makeAbsolute(vecs_data[i]); } }); } @@ -263,12 +269,14 @@ class Box * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void makeFractional(const vec3* vecs, unsigned int Nvecs, vec3* out) const + void makeFractional(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const { + vec3 *vecs_data = static_cast *>(vecs.data()); + vec3 *out_data = static_cast *>(out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { - out[i] = makeFractional(vecs[i]); + out_data[i] = makeFractional(vecs_data[i]); } }); } diff --git a/cpp/util/utils.h b/cpp/util/utils.h index 3bee6c86b..d52cf4b89 100644 --- a/cpp/util/utils.h +++ b/cpp/util/utils.h @@ -10,6 +10,7 @@ #include #include + namespace freud { namespace util { //! Clip v if it is outside the range [lo, hi]. diff --git a/freud/box.py b/freud/box.py index 7fc88b0cc..c9d811abf 100644 --- a/freud/box.py +++ b/freud/box.py @@ -49,7 +49,7 @@ class Box: if :code:`None`. (Default value = :code:`None`) """ # noqa: E501 - def __cinit__(self, Lx, Ly, Lz=0, xy=0, xz=0, yz=0, is2D=None): + def __init__(self, Lx, Ly, Lz=0, xy=0, xz=0, yz=0, is2D=None): if is2D is None: is2D = (Lz == 0) if is2D: @@ -162,8 +162,8 @@ def is2D(self): def L_inv(self): r""":math:`\left(3, \right)` :class:`numpy.ndarray`: The inverse box lengths.""" - cdef vec3[float] result = self._cpp_obj.getLinv() - return np.asarray([result.x, result.y, result.z]) + result = self._cpp_obj.getLinv() + return np.asarray(result) @property def volume(self): @@ -194,13 +194,9 @@ def make_absolute(self, fractional_coordinates, out=None): out = freud.util._convert_array( out, shape=fractions.shape, allow_copy=False) - cdef const float[:, ::1] l_points = fractions - cdef unsigned int Np = l_points.shape[0] - cdef float[:, ::1] l_out = out + Np = fractions.shape[0] - self._cpp_obj.makeAbsolute( - &l_points[0, 0], Np, - &l_out[0, 0]) + self._cpp_obj.makeAbsolute(fractions, Np, out) return np.squeeze(out) if flatten else out @@ -227,13 +223,9 @@ def make_fractional(self, absolute_coordinates, out=None): out = freud.util._convert_array( out, shape=vecs.shape, allow_copy=False) - cdef const float[:, ::1] l_points = vecs - cdef unsigned int Np = l_points.shape[0] - cdef float[:, ::1] l_out = out + Np = vecs.shape[0] - self._cpp_obj.makeFractional( - &l_points[0, 0], Np, - &l_out[0, 0]) + self._cpp_obj.makeFractional(vecs, Np, out) return np.squeeze(out) if flatten else out From b79e7c2961534f7aa01926f7f1e3eef3b2bf262f Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Mon, 3 Jun 2024 17:47:06 -0400 Subject: [PATCH 04/25] changes so tests pass with pybind11 --- CMakeLists.txt | 3 - cpp/box/Box.h | 113 ++++++++++++++++++----------- cpp/box/CMakeLists.txt | 15 ++-- cpp/box/module-box.cc | 15 ++-- freud/CMakeLists.txt | 5 +- freud/__init__.py | 62 ++++++++-------- freud/box.py | 123 +++++++++++--------------------- freud/util.py | 158 +++++++++++++++++++++++++++++++++++++++++ 8 files changed, 322 insertions(+), 172 deletions(-) create mode 100644 freud/util.py diff --git a/CMakeLists.txt b/CMakeLists.txt index fdc51e44a..836f8f711 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -80,6 +80,3 @@ set(ignoreMe "${SKBUILD}") add_subdirectory(cpp) add_subdirectory(freud) -if(_using_conda OR DEFINED ENV{CIBUILDWHEEL}) - set_target_properties(libfreud PROPERTIES INSTALL_RPATH_USE_LINK_PATH True) -endif() diff --git a/cpp/box/Box.h b/cpp/box/Box.h index a9a5aaca2..c87db25fe 100644 --- a/cpp/box/Box.h +++ b/cpp/box/Box.h @@ -8,6 +8,7 @@ #include #include #include +#include #include "utils.h" #include "VectorMath.h" @@ -25,7 +26,8 @@ constexpr float TWO_PI = 2.0 * M_PI; namespace freud { namespace box { -using pybind11::array_t = pybind11_array; +template +using pybind11_array = pybind11::array_t; //! Stores box dimensions and provides common routines for wrapping vectors back into the box /*! Box stores a standard HOOMD simulation box that goes from -L/2 to L/2 in each dimension, allowing Lx, Ly, @@ -97,12 +99,6 @@ class Box return !(*this == b); } - //! Set L, box lengths, inverses. Box is also centered at zero. - void setL(const vec3& L) - { - setL(L.x, L.y, L.z); - } - //! Set L, box lengths, inverses. Box is also centered at zero. void setL(const float Lx, const float Ly, const float Lz) { @@ -162,7 +158,7 @@ class Box //! Get current stored inverse of L std::array getLinv() const { - return { m_Linv.x, m_linv.y, m_Linv.z }; + return { m_Linv.x, m_Linv.y, m_Linv.z }; } //! Get tilt factor xy @@ -234,10 +230,10 @@ class Box * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void makeAbsolute(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const + void makeAbsolutePython(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const { - vec3 *vecs_data = static_cast *>(vecs.data()); - vec3 *out_data = static_cast *>(out.data()); + vec3 *vecs_data = (vec3 *)(vecs.data()); + vec3 *out_data = (vec3 *)(out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { @@ -269,10 +265,10 @@ class Box * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void makeFractional(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const + void makeFractionalPython(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const { - vec3 *vecs_data = static_cast *>(vecs.data()); - vec3 *out_data = static_cast *>(out.data()); + vec3 *vecs_data = (vec3 *)(vecs.data()); + vec3 *out_data = (vec3 *)(out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { @@ -302,12 +298,15 @@ class Box * \param Nvecs Number of vectors \param res Array to save the images */ - void getImages(vec3* vecs, unsigned int Nvecs, vec3* res) const + void getImages(pybind11_array vecs, unsigned int Nvecs, + pybind11::array_t res) const { + vec3 *vecs_data = (vec3 *)(vecs.data()); + vec3 *out_data = (vec3 *)(res.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { - getImage(vecs[i], res[i]); + getImage(vecs_data[i], out_data[i]); } }); } @@ -345,12 +344,14 @@ class Box * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void wrap(const vec3* vecs, unsigned int Nvecs, vec3* out) const + void wrapPython(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const { + vec3 *vecs_data = (vec3 *)(vecs.data()); + vec3 *out_data = (vec3 *)(out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { - out[i] = wrap(vecs[i]); + out_data[i] = wrap(vecs_data[i]); } }); } @@ -361,16 +362,20 @@ class Box \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void unwrap(const vec3* vecs, const vec3* images, unsigned int Nvecs, vec3* out) const + void unwrap(pybind11_array vecs, pybind11_array images, unsigned int Nvecs, + pybind11_array out) const { + vec3 *vecs_data = (vec3 *)(vecs.data()); + vec3 *images_data = (vec3 *)(images.data()); + vec3 *out_data = (vec3 *)(out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { - out[i] = vecs[i] + getLatticeVector(0) * float(images[i].x) - + getLatticeVector(1) * float(images[i].y); + out_data[i] = vecs_data[i] + getLatticeVector(0) * float(images_data[i].x) + + getLatticeVector(1) * float(images_data[i].y); if (!m_2d) { - out[i] += getLatticeVector(2) * float(images[i].z); + out_data[i] += getLatticeVector(2) * float(images_data[i].z); } } }); @@ -382,7 +387,7 @@ class Box * \param masses Optional array of masses, of length Nvecs * \return Center of mass as a vec3 */ - vec3 centerOfMass(vec3* vecs, size_t Nvecs, const float* masses = nullptr) const + vec3 centerOfMass(vec3 *vecs, size_t Nvecs, const float *masses) const { // This roughly follows the implementation in // https://en.wikipedia.org/wiki/Center_of_mass#Systems_with_periodic_boundary_conditions @@ -392,30 +397,44 @@ class Box for (size_t i = 0; i < Nvecs; ++i) { vec3 phase(constants::TWO_PI * makeFractional(vecs[i])); - vec3> xi(std::polar(float(1.0), phase.x), std::polar(float(1.0), phase.y), + vec3> xi(std::polar(float(1.0), phase.x), + std::polar(float(1.0), phase.y), std::polar(float(1.0), phase.z)); - float mass = (masses != nullptr) ? masses[i] : float(1.0); + float mass = masses[i]; total_mass += mass; xi_mean += std::complex(mass, 0) * xi; } xi_mean /= std::complex(total_mass, 0); - return wrap(makeAbsolute(vec3(std::arg(xi_mean.x), std::arg(xi_mean.y), std::arg(xi_mean.z)) + return wrap(makeAbsolute(vec3(std::arg(xi_mean.x), + std::arg(xi_mean.y), + std::arg(xi_mean.z)) / constants::TWO_PI)); } + std::array centerOfMassPython(pybind11_array vecs, size_t Nvecs, pybind11_array masses) const + { + vec3 *vecs_data = (vec3 *)(vecs.data()); + float *masses_data = (float *)(masses.data()); + auto com = centerOfMass(vecs_data, Nvecs, masses_data); + return { com.x, com.y, com.z }; + } + //! Subtract center of mass from vectors /*! \param vecs Vectors to center * \param Nvecs Number of vectors * \param masses Optional array of masses, of length Nvecs */ - void center(vec3* vecs, unsigned int Nvecs, const float* masses = nullptr) const + void center(pybind11_array vecs, unsigned int Nvecs, pybind11_array masses) const { - vec3 com(centerOfMass(vecs, Nvecs, masses)); + vec3 *vecs_data = (vec3 *)(vecs.data()); + float *masses_data = (float *)(masses.data()); + + vec3 com(centerOfMass(vecs_data, Nvecs, masses_data)); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { - vecs[i] = wrap(vecs[i] - com); + vecs_data[i] = wrap(vecs_data[i] - com); } }); } @@ -438,9 +457,14 @@ class Box \param distances Pointer to array of length n_query_points containing distances between each point and query_point (overwritten in place). */ - void computeDistances(const vec3* query_points, const unsigned int n_query_points, - const vec3* points, const unsigned int n_points, float* distances) const + void computeDistances(pybind11_array query_points, const unsigned int n_query_points, + pybind11_array points, const unsigned int n_points, + pybind11_array distances) const { + vec3 *query_points_data = (vec3 *)(query_points.data()); + vec3 *points_data = (vec3 *)(points.data()); + float *distances_data = (float *)(distances.data()); + if (n_query_points != n_points) { throw std::invalid_argument("The number of query points and points must match."); @@ -448,7 +472,7 @@ class Box util::forLoopWrapper(0, n_query_points, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { - distances[i] = computeDistance(query_points[i], points[i]); + distances_data[i] = computeDistance(query_points_data[i], points_data[i]); } }); } @@ -461,16 +485,22 @@ class Box \param distances Pointer to array of length n_query_points*n_points containing distances between points and query_points (overwritten in place). */ - void computeAllDistances(const vec3* query_points, const unsigned int n_query_points, - const vec3* points, const unsigned int n_points, float* distances) const + void computeAllDistances(pybind11_array query_points, const unsigned int n_query_points, + pybind11_array points, const unsigned int n_points, + pybind11_array distances) const { + vec3 *query_points_data = (vec3 *)(query_points.data()); + vec3 *points_data = (vec3 *)(points.data()); + float *distances_data = (float *)(distances.data()); + util::forLoopWrapper2D( 0, n_query_points, 0, n_points, [&](size_t begin_n, size_t end_n, size_t begin_m, size_t end_m) { for (size_t i = begin_n; i < end_n; ++i) { for (size_t j = begin_m; j < end_m; ++j) { - distances[i * n_points + j] = computeDistance(query_points[i], points[j]); + distances_data[i * n_points + j] = computeDistance( + query_points_data[i], points_data[j]); } } }); @@ -481,10 +511,14 @@ class Box \param n_points The number of points. \param contains_mask Mask of points inside the box. */ - void contains(const vec3* points, const unsigned int n_points, bool* contains_mask) const + void contains(pybind11_array points, const unsigned int n_points, + pybind11_array contains_mask) const { + vec3 *points_data = (vec3 *)(points.data()); + bool *contains_mask_data = (bool *)(contains_mask.data()); + util::forLoopWrapper(0, n_points, [&](size_t begin, size_t end) { - std::transform(&points[begin], &points[end], &contains_mask[begin], + std::transform(&points_data[begin], &points_data[end], &contains_mask_data[begin], [this](const vec3& point) -> bool { vec3 image(0, 0, 0); getImage(point, image); @@ -536,11 +570,6 @@ class Box /*! \param periodic Flags to set * \post Period flags are set to \a periodic */ - void setPeriodic(vec3 periodic) - { - m_periodic = periodic; - } - void setPeriodic(bool x, bool y, bool z) { m_periodic = vec3(x, y, z); diff --git a/cpp/box/CMakeLists.txt b/cpp/box/CMakeLists.txt index 23cfd8ba3..223e6c47e 100644 --- a/cpp/box/CMakeLists.txt +++ b/cpp/box/CMakeLists.txt @@ -1,6 +1,6 @@ # create the target and add the needed properties -pybind11_add_module(box SHARED module-box.cc) -target_link_libraries(box PUBLIC TBB::tbb) +pybind11_add_module(_box SHARED module-box.cc) +target_link_libraries(_box PUBLIC TBB::tbb) # this probably isn't needed for box, but may be needed by other modules #target_compile_definitions( # box @@ -10,15 +10,16 @@ target_link_libraries(box PUBLIC TBB::tbb) # # Default voro++ verbosity is high. # PRIVATE "VOROPP_VERBOSE=1") -# install -install(TARGETS box DESTINATION freud) if(APPLE) - set_target_properties(box PROPERTIES INSTALL_RPATH "@loader_path") + set_target_properties(_box PROPERTIES INSTALL_RPATH "@loader_path") else() - set_target_properties(box PROPERTIES INSTALL_RPATH "\$ORIGIN") + set_target_properties(_box PROPERTIES INSTALL_RPATH "\$ORIGIN") endif() if(_using_conda OR DEFINED ENV{CIBUILDWHEEL}) - set_target_properties(box + set_target_properties(_box PROPERTIES INSTALL_RPATH_USE_LINK_PATH True) endif() + +# install +install(TARGETS _box DESTINATION freud) diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc index 05e6e5fe6..37441af84 100644 --- a/cpp/box/module-box.cc +++ b/cpp/box/module-box.cc @@ -1,8 +1,11 @@ #include +#include #include "Box.h" +using namespace freud::box; + PYBIND11_MODULE(_box, m) { pybind11::class_>(m, "Box") @@ -13,6 +16,7 @@ PYBIND11_MODULE(_box, m) .def("getLy", &Box::getLy) .def("getLz", &Box::getLz) .def("setL", &Box::setL) + .def("getLinv", &Box::getLinv) .def("getTiltFactorXY", &Box::getTiltFactorXY) .def("getTiltFactorXZ", &Box::getTiltFactorXZ) .def("getTiltFactorYZ", &Box::getTiltFactorYZ) @@ -26,16 +30,17 @@ PYBIND11_MODULE(_box, m) .def("setPeriodicX", &Box::setPeriodicX) .def("setPeriodicY", &Box::setPeriodicY) .def("setPeriodicZ", &Box::setPeriodicZ) - .def("get2D", &Box::get2D) + .def("is2D", &Box::is2D) .def("set2D", &Box::set2D) .def("getVolume", &Box::getVolume) .def("center", &Box::center) - .def("centerOfMass", &Box::centerOfMass) + .def("centerOfMass", &Box::centerOfMassPython) // other stuff - .def("makeAbsolute", &Box::makeAbsolute) - .def("makeFractional", &Box::makeFractional) - .def("wrap", &Box::wrap) + .def("makeAbsolute", &Box::makeAbsolutePython) + .def("makeFractional", &Box::makeFractionalPython) + .def("wrap", &Box::wrapPython) .def("unwrap", &Box::unwrap) + .def("getImages", &Box::getImages) .def("computeDistances", &Box::computeDistances) .def("computeAllDistances", &Box::computeAllDistances) .def("contains", &Box::contains); diff --git a/freud/CMakeLists.txt b/freud/CMakeLists.txt index 7e258c3c9..f8fc35192 100644 --- a/freud/CMakeLists.txt +++ b/freud/CMakeLists.txt @@ -21,7 +21,8 @@ if(${PYTHON_VERSION_MAJOR} EQUAL 3 endif() set(python_files - box.py) + box.py + util.py) #cluster.py #density.py #diffraction.py @@ -31,7 +32,7 @@ set(python_files #parallel.py #pmft.py) -install(FILES python_files DESTINATION freud) +install(FILES ${python_files} DESTINATION freud) # THIS WILL NEED TO BE MOVED TO THE BUILD OF THE ORDER MODULE # The SolidLiquid class has an instance of cluster::Cluster as a member, so diff --git a/freud/__init__.py b/freud/__init__.py index 6cdd0ad71..aadf17ec1 100644 --- a/freud/__init__.py +++ b/freud/__init__.py @@ -3,49 +3,49 @@ from . import ( box, - cluster, - data, - density, - diffraction, - environment, - interface, - locality, - msd, - order, - parallel, - pmft, + #cluster, + #data, + #density, + #diffraction, + #environment, + #interface, + #locality, + #msd, + #order, + #parallel, + #pmft, ) from .box import Box -from .locality import AABBQuery, LinkCell, NeighborList -from .parallel import NumThreads, get_num_threads, set_num_threads +#from .locality import AABBQuery, LinkCell, NeighborList +#from .parallel import NumThreads, get_num_threads, set_num_threads # Override TBB's default autoselection. This is necessary because once the # automatic selection runs, the user cannot change it. -set_num_threads(0) +#set_num_threads(0) __version__ = "3.0.0" __all__ = [ "__version__", "box", - "cluster", - "data", - "density", - "diffraction", - "environment", - "interface", - "locality", - "msd", - "order", - "parallel", - "pmft", + #"cluster", + #"data", + #"density", + #"diffraction", + #"environment", + #"interface", + #"locality", + #"msd", + #"order", + #"parallel", + #"pmft", "Box", - "AABBQuery", - "LinkCell", - "NeighborList", - "get_num_threads", - "set_num_threads", - "NumThreads", + #"AABBQuery", + #"LinkCell", + #"NeighborList", + #"get_num_threads", + #"set_num_threads", + #"NumThreads", ] __citation__ = """@article{freud2020, diff --git a/freud/box.py b/freud/box.py index c9d811abf..acec10eac 100644 --- a/freud/box.py +++ b/freud/box.py @@ -143,6 +143,14 @@ def yz(self): def yz(self, value): self._cpp_obj.setTiltFactorYZ(value) + def __eq__(self, other): + if type(other) != freud.box.Box: + return False + return self.Lx == other.Lx and self.Ly == other.Ly and self.Lz == other.Lz \ + and self.xy == other.xy and self.xz == other.xz and self.yz == other.yz \ + and self.is2D == other.is2D and self.periodic_x == other.periodic_x and \ + self.periodic_y == other.periodic_y and self.periodic_z == other.periodic_z + @property def dimensions(self): """int: Get or set the number of dimensions (2 or 3).""" @@ -156,7 +164,7 @@ def dimensions(self, value): @property def is2D(self): """bool: Whether the box is 2D.""" - return self._cpp_obj.get2D() + return self._cpp_obj.is2D() @property def L_inv(self): @@ -246,11 +254,8 @@ def get_images(self, vecs): vecs = freud.util._convert_array(vecs, shape=(None, 3)) images = np.zeros(vecs.shape, dtype=np.int32) - cdef const float[:, ::1] l_points = vecs - cdef const int[:, ::1] l_result = images - cdef unsigned int Np = l_points.shape[0] - self._cpp_obj.getImages( &l_points[0, 0], Np, - &l_result[0, 0]) + Np = vecs.shape[0] + self._cpp_obj.getImages(vecs, Np, images) return np.squeeze(images) if flatten else images @@ -313,12 +318,8 @@ def wrap(self, vecs, out=None): out = freud.util._convert_array( out, shape=vecs.shape, allow_copy=False) - cdef const float[:, ::1] l_points = vecs - cdef unsigned int Np = l_points.shape[0] - cdef float[:, ::1] l_out = out - - self._cpp_obj.wrap( &l_points[0, 0], - Np, &l_out[0, 0]) + Np = vecs.shape[0] + self._cpp_obj.wrap(vecs, Np, out) return np.squeeze(out) if flatten else out @@ -355,14 +356,8 @@ def unwrap(self, vecs, imgs, out=None): out = freud.util._convert_array( out, shape=vecs.shape, allow_copy=False) - cdef const float[:, ::1] l_points = vecs - cdef const int[:, ::1] l_imgs = imgs - cdef unsigned int Np = l_points.shape[0] - cdef float[:, ::1] l_out = out - - self._cpp_obj.unwrap( &l_points[0, 0], - &l_imgs[0, 0], Np, - &l_out[0, 0]) + Np = vecs.shape[0] + self._cpp_obj.unwrap(vecs, imgs, Np, out) return np.squeeze(out) if flatten else out @@ -396,18 +391,15 @@ def center_of_mass(self, vecs, masses=None): Center of mass. """ # noqa: E501 vecs = freud.util._convert_array(vecs, shape=(None, 3)) - cdef const float[:, ::1] l_points = vecs + Np = vecs.shape[0] - cdef float* l_masses_ptr = NULL - cdef float[::1] l_masses if masses is not None: - l_masses = freud.util._convert_array(masses, shape=(len(vecs), )) - l_masses_ptr = &l_masses[0] + masses = freud.util._convert_array(masses, shape=(len(vecs), )) + else: + masses = np.ones(Np) - cdef size_t Np = l_points.shape[0] - cdef vec3[float] result = self._cpp_obj.centerOfMass( - &l_points[0, 0], Np, l_masses_ptr) - return np.asarray([result.x, result.y, result.z]) + result = self._cpp_obj.centerOfMass(vecs, Np, masses) + return np.asarray(result) def center(self, vecs, masses=None): r"""Subtract center of mass from an array of vectors, using periodic boundaries. @@ -438,16 +430,14 @@ def center(self, vecs, masses=None): Vectors with center of mass subtracted. """ # noqa: E501 vecs = freud.util._convert_array(vecs, shape=(None, 3)).copy() - cdef const float[:, ::1] l_points = vecs + Np = vecs.shape[0] - cdef float* l_masses_ptr = NULL - cdef float[::1] l_masses if masses is not None: - l_masses = freud.util._convert_array(masses, shape=(len(vecs), )) - l_masses_ptr = &l_masses[0] + masses = freud.util._convert_array(masses, shape=(len(vecs), )) + else: + masses = np.ones(Np) - cdef size_t Np = l_points.shape[0] - self._cpp_obj.center( &l_points[0, 0], Np, l_masses_ptr) + self._cpp_obj.center(vecs, Np, masses) return vecs def compute_distances(self, query_points, points): @@ -472,18 +462,12 @@ def compute_distances(self, query_points, points): points = freud.util._convert_array( np.atleast_2d(points), shape=(None, 3)) - cdef: - const float[:, ::1] l_query_points = query_points - const float[:, ::1] l_points = points - size_t n_query_points = query_points.shape[0] - size_t n_points = points.shape[0] - float[::1] distances = np.empty( - n_query_points, dtype=np.float32) - - self._cpp_obj.computeDistances( - &l_query_points[0, 0], n_query_points, - &l_points[0, 0], n_points, - &distances[0]) + n_query_points = query_points.shape[0] + n_points = points.shape[0] + distances = np.empty(n_query_points, dtype=np.float32) + + self._cpp_obj.computeDistances(query_points, n_query_points, points, + n_points, distances) return np.asarray(distances) def compute_all_distances(self, query_points, points): @@ -507,18 +491,12 @@ def compute_all_distances(self, query_points, points): points = freud.util._convert_array( np.atleast_2d(points), shape=(None, 3)) - cdef: - const float[:, ::1] l_query_points = query_points - const float[:, ::1] l_points = points - size_t n_query_points = query_points.shape[0] - size_t n_points = points.shape[0] - float[:, ::1] distances = np.empty( - [n_query_points, n_points], dtype=np.float32) + n_query_points = query_points.shape[0] + n_points = points.shape[0] + distances = np.empty([n_query_points, n_points], dtype=np.float32) - self._cpp_obj.computeAllDistances( - &l_query_points[0, 0], n_query_points, - &l_points[0, 0], n_points, - &distances[0, 0]) + self._cpp_obj.computeAllDistances(query_points, n_query_points, + points, n_points, distances) return np.asarray(distances) @@ -557,19 +535,14 @@ def contains(self, points): points = freud.util._convert_array( np.atleast_2d(points), shape=(None, 3)) - cdef: - const float[:, ::1] l_points = points - size_t n_points = points.shape[0] + n_points = points.shape[0] contains_mask = freud.util._convert_array( np.ones(n_points), dtype=bool) - cdef cpp_bool[::1] l_contains_mask = contains_mask - self._cpp_obj.contains( - &l_points[0, 0], n_points, - &l_contains_mask[0]) + self._cpp_obj.contains(points, n_points, contains_mask) - return np.array(l_contains_mask).astype(bool) + return np.array(contains_mask).astype(bool) @property def cubic(self): @@ -705,20 +678,6 @@ def __repr__(self): def __str__(self): return repr(self) - def __richcmp__(self, other, int op): - r"""Implement all comparisons for Cython extension classes""" - cdef Box c_other - try: - c_other = other - if op == Py_EQ: - return dereference(self._cpp_obj) == dereference(c_other._cpp_obj) - if op == Py_NE: - return dereference(self._cpp_obj) != dereference(c_other._cpp_obj) - except TypeError: - # Cython cast to Box failed - pass - return NotImplemented - def __mul__(self, scale): if scale > 0: return self.__class__(Lx=self.Lx*scale, @@ -972,7 +931,7 @@ def from_box_lengths_and_angles( return cls.from_matrix(np.array([a1, a2, a3]).T, dimensions=dimensions) -cdef BoxFromCPP(const freud._box.Box & cppbox): +def BoxFromCPP(cppbox): b = Box(cppbox.getLx(), cppbox.getLy(), cppbox.getLz(), cppbox.getTiltFactorXY(), cppbox.getTiltFactorXZ(), cppbox.getTiltFactorYZ(), cppbox.is2D()) diff --git a/freud/util.py b/freud/util.py new file mode 100644 index 000000000..11a1d2866 --- /dev/null +++ b/freud/util.py @@ -0,0 +1,158 @@ +# Copyright (c) 2010-2023 The Regents of the University of Michigan +# This file is from the freud project, released under the BSD 3-Clause License. + +from functools import wraps + +import numpy as np + +import freud.box + + +class _Compute(object): + r"""Parent class for all compute classes in freud. + + The primary purpose of this class is to prevent access of uncomputed + values. This is accomplished by maintaining a boolean flag to track whether + the compute method in a class has been called and decorating class + properties that rely on compute having been called. + + To use this class, one would write, for example, + + .. code-block:: python + class Cluster(_Compute): + + def compute(...) + ... + + @_Compute._computed_property + def cluster_idx(self): + return ... + + Attributes: + _called_compute (bool): + Flag representing whether the compute method has been called. + """ + + def __init__(self): + self._called_compute = False + + def __getattribute__(self, attr): + """Compute methods set a flag to indicate that quantities have been + computed. Compute must be called before plotting.""" + attribute = object.__getattribute__(self, attr) + if attr == 'compute': + # Set the attribute *after* computing. This enables + # self._called_compute to be used in the compute method itself. + compute = attribute + + @wraps(compute) + def compute_wrapper(*args, **kwargs): + return_value = compute(*args, **kwargs) + self._called_compute = True + return return_value + return compute_wrapper + elif attr == 'plot': + if not self._called_compute: + raise AttributeError( + "The compute method must be called before calling plot.") + return attribute + + @staticmethod + def _computed_property(prop): + r"""Decorator that makes a class method to be a property with limited access. + + Args: + prop (callable): The property function. + + Returns: + Decorator decorating appropriate property method. + """ + + @property + @wraps(prop) + def wrapper(self, *args, **kwargs): + if not self._called_compute: + raise AttributeError( + "Property not computed. Call compute first.") + return prop(self, *args, **kwargs) + return wrapper + + def __str__(self): + return repr(self) + + +def _convert_array(array, shape=None, dtype=np.float32, requirements=("C", ), + allow_copy=True): + """Function which takes a given array, checks the dimensions and shape, + and converts to a supplied dtype. + + Args: + array (:class:`numpy.ndarray` or :code:`None`): Array to check and convert. + If :code:`None`, an empty array of given shape and type will be initialized + (Default value: :code:`None`). + shape: (tuple of int and :code:`None`): Expected shape of the array. + Only the dimensions that are not :code:`None` are checked. + (Default value = :code:`None`). + dtype: :code:`dtype` to convert the array to if :code:`array.dtype` + is different. If :code:`None`, :code:`dtype` will not be changed + (Default value = :attr:`numpy.float32`). + requirements (Sequence[str]): A sequence of string flags to be passed to + :func:`numpy.require`. + allow_copy (bool): If :code:`False` and the input array does not already + conform to the required dtype and other requirements, this function + will raise an error rather than coercing the array into a copy that + does satisfy the requirements (Default value = :code:`True`). + + Returns: + :class:`numpy.ndarray`: Array. + """ + if array is None: + return np.empty(shape, dtype=dtype) + + array = np.asarray(array) + return_arr = np.require(array, dtype=dtype, requirements=requirements) + + if not allow_copy and return_arr is not array: + raise ValueError("The provided output array must have dtype " + f"{dtype}, and have the following array flags: " + f"{', '.join(requirements)}.") + + if shape is not None: + if return_arr.ndim != len(shape): + raise ValueError("array.ndim = {}; expected ndim = {}".format( + return_arr.ndim, len(shape))) + + for i, s in enumerate(shape): + if s is not None and return_arr.shape[i] != s: + shape_str = "(" + ", ".join(str(i) if i is not None + else "..." for i in shape) + ")" + raise ValueError('array.shape= {}; expected shape = {}'.format( + return_arr.shape, shape_str)) + + return return_arr + + +def _convert_box(box, dimensions=None): + """Function which takes a box-like object and attempts to convert it to + :class:`freud.box.Box`. Existing :class:`freud.box.Box` objects are + used directly. + + Args: + box (box-like object (see :meth:`freud.box.Box.from_box`)): Box to + check and convert if needed. + dimensions (int): Number of dimensions the box should be. If not None, + used to verify the box dimensions (Default value = :code:`None`). + + Returns: + :class:`freud.box.Box`: freud box. + """ + if not isinstance(box, freud.box.Box): + try: + box = freud.box.Box.from_box(box) + except ValueError: + raise + + if dimensions is not None and box.dimensions != dimensions: + raise ValueError("The box must be {}-dimensional.".format(dimensions)) + + return box From 88eed2c11181f3ed91e4d8cc0d124d4e60593855 Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Mon, 3 Jun 2024 17:48:52 -0400 Subject: [PATCH 05/25] remove files we don't need anymore --- freud/_box.pxd | 63 --- freud/box.pxd | 10 - freud/box.pyx | 1002 ------------------------------------------------ 3 files changed, 1075 deletions(-) delete mode 100644 freud/_box.pxd delete mode 100644 freud/box.pxd delete mode 100644 freud/box.pyx diff --git a/freud/_box.pxd b/freud/_box.pxd deleted file mode 100644 index 1d8843e81..000000000 --- a/freud/_box.pxd +++ /dev/null @@ -1,63 +0,0 @@ -# Copyright (c) 2010-2023 The Regents of the University of Michigan -# This file is from the freud project, released under the BSD 3-Clause License. - -from libcpp cimport bool - -from freud.util cimport vec3 - -ctypedef unsigned int uint - -cdef extern from "Box.h" namespace "freud::box": - cdef cppclass Box: - Box() - Box(float, bool) - Box(float, float, float, bool) - Box(float, float, float, float, float, float, bool) - - bool operator==(const Box &) const - bool operator!=(const Box &) const - - void setL(vec3[float]) - void setL(float, float, float) - - void set2D(bool) - bool is2D() const - - float getLx() const - float getLy() const - float getLz() const - - vec3[float] getL() const - vec3[float] getLinv() const - - float getTiltFactorXY() const - float getTiltFactorXZ() const - float getTiltFactorYZ() const - - void setTiltFactorXY(float) - void setTiltFactorXZ(float) - void setTiltFactorYZ(float) - - float getVolume() const - void makeAbsolute(const vec3[float]*, unsigned int, vec3[float]*) const - void makeFractional(const vec3[float]*, unsigned int, vec3[float]*) const - void getImages(vec3[float]*, unsigned int, vec3[int]*) const - void wrap(const vec3[float]*, unsigned int, vec3[float]*) const - void unwrap(const vec3[float]*, const vec3[int]*, - unsigned int, vec3[float]*) const - vec3[float] centerOfMass(vec3[float]*, size_t, float*) const - void center(vec3[float]*, size_t, float*) const - void computeDistances(vec3[float]*, unsigned int, - vec3[float]*, unsigned int, float* - ) except + - void computeAllDistances(vec3[float]*, unsigned int, - vec3[float]*, unsigned int, float*) - void contains(vec3[float]*, unsigned int, bool*) const - vec3[bool] getPeriodic() const - bool getPeriodicX() const - bool getPeriodicY() const - bool getPeriodicZ() const - void setPeriodic(bool, bool, bool) - void setPeriodicX(bool) - void setPeriodicY(bool) - void setPeriodicZ(bool) diff --git a/freud/box.pxd b/freud/box.pxd deleted file mode 100644 index 7a760f7d6..000000000 --- a/freud/box.pxd +++ /dev/null @@ -1,10 +0,0 @@ -# Copyright (c) 2010-2023 The Regents of the University of Michigan -# This file is from the freud project, released under the BSD 3-Clause License. - -cimport freud._box - - -cdef class Box: - cdef freud._box.Box * thisptr - -cdef BoxFromCPP(const freud._box.Box & cppbox) diff --git a/freud/box.pyx b/freud/box.pyx deleted file mode 100644 index f72151bbf..000000000 --- a/freud/box.pyx +++ /dev/null @@ -1,1002 +0,0 @@ -# Copyright (c) 2010-2023 The Regents of the University of Michigan -# This file is from the freud project, released under the BSD 3-Clause License. - -r""" -The :class:`~.Box` class defines the geometry of a simulation box. The class -natively supports periodicity by providing the fundamental features for -wrapping vectors outside the box back into it. -""" - -from cpython.object cimport Py_EQ, Py_NE -from cython.operator cimport dereference -from libcpp cimport bool as cpp_bool - -from freud.util cimport vec3 - -import logging -import warnings - -import numpy as np - -import freud.util - -cimport numpy as np - -cimport freud._box - -logger = logging.getLogger(__name__) - -# numpy must be initialized. When using numpy from C or Cython you must -# _always_ do that, or you will have segfaults -np.import_array() - -cdef class Box: - r"""The freud Box class for simulation boxes. - - This class defines an arbitrary triclinic geometry within which all points are confined. - By convention, the freud Box is centered at the origin (``[0, 0, 0]``), - with the extent in each dimension described by the half-open interval ``[-L/2, L/2)``. - For more information, see the `documentation - `__ - on boxes and periodic boundary conditions. - - Also available as ``freud.Box``. - - Args: - Lx (float, optional): - The x-dimension length. - Ly (float, optional): - The y-dimension length. - Lz (float, optional): - The z-dimension length (Default value = 0). - xy (float, optional): - The xy tilt factor (Default value = 0). - xz (float, optional): - The xz tilt factor (Default value = 0). - yz (float, optional): - The yz tilt factor (Default value = 0). - is2D (bool, optional): - Whether the box is 2-dimensional. Uses :code:`Lz == 0` - if :code:`None`. (Default value = :code:`None`) - """ # noqa: E501 - - def __cinit__(self, Lx, Ly, Lz=0, xy=0, xz=0, yz=0, is2D=None): - if is2D is None: - is2D = (Lz == 0) - if is2D: - if not (Lx and Ly): - raise ValueError("Lx and Ly must be nonzero for 2D boxes.") - elif Lz != 0 or xz != 0 or yz != 0: - warnings.warn( - "Specifying z-dimensions in a 2-dimensional box " - "has no effect!") - else: - if not (Lx and Ly and Lz): - raise ValueError( - "Lx, Ly, and Lz must be nonzero for 3D boxes.") - self.thisptr = new freud._box.Box(Lx, Ly, Lz, xy, xz, yz, is2D) - - def __dealloc__(self): - del self.thisptr - - @property - def L(self): - r""":math:`\left(3, \right)` :class:`numpy.ndarray`: Get or set the - box lengths along x, y, and z.""" - cdef vec3[float] result = self.thisptr.getL() - return np.asarray([result.x, result.y, result.z]) - - @L.setter - def L(self, value): - try: - if len(value) != 3: - raise ValueError('setL must be called with a scalar or a list ' - 'of length 3.') - except TypeError: - # Will fail if object has no length - value = (value, value, value) - - if self.is2D and value[2] != 0: - warnings.warn( - "Specifying z-dimensions in a 2-dimensional box " - "has no effect!") - self.thisptr.setL(value[0], value[1], value[2]) - - @property - def Lx(self): - """float: Get or set the x-dimension length.""" - return self.thisptr.getLx() - - @Lx.setter - def Lx(self, value): - self.L = [value, self.Ly, self.Lz] - - @property - def Ly(self): - """float: Get or set the y-dimension length.""" - return self.thisptr.getLy() - - @Ly.setter - def Ly(self, value): - self.L = [self.Lx, value, self.Lz] - - @property - def Lz(self): - """float: Get or set the z-dimension length.""" - return self.thisptr.getLz() - - @Lz.setter - def Lz(self, value): - self.L = [self.Lx, self.Ly, value] - - @property - def xy(self): - """float: Get or set the xy tilt factor.""" - return self.thisptr.getTiltFactorXY() - - @xy.setter - def xy(self, value): - self.thisptr.setTiltFactorXY(value) - - @property - def xz(self): - """float: Get or set the xz tilt factor.""" - return self.thisptr.getTiltFactorXZ() - - @xz.setter - def xz(self, value): - self.thisptr.setTiltFactorXZ(value) - - @property - def yz(self): - """float: Get or set the yz tilt factor.""" - return self.thisptr.getTiltFactorYZ() - - @yz.setter - def yz(self, value): - self.thisptr.setTiltFactorYZ(value) - - @property - def dimensions(self): - """int: Get or set the number of dimensions (2 or 3).""" - return 2 if self.is2D else 3 - - @dimensions.setter - def dimensions(self, value): - assert value == 2 or value == 3 - self.thisptr.set2D(bool(value == 2)) - - @property - def is2D(self): - """bool: Whether the box is 2D.""" - return self.thisptr.is2D() - - @property - def L_inv(self): - r""":math:`\left(3, \right)` :class:`numpy.ndarray`: The inverse box - lengths.""" - cdef vec3[float] result = self.thisptr.getLinv() - return np.asarray([result.x, result.y, result.z]) - - @property - def volume(self): - """float: The box volume (area in 2D).""" - return self.thisptr.getVolume() - - def make_absolute(self, fractional_coordinates, out=None): - r"""Convert fractional coordinates into absolute coordinates. - - Args: - fractional_coordinates (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Fractional coordinate vector(s), between 0 and 1 within - parallelepipedal box. - out (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray` or :code:`None`): - The array in which to place the absolute coordinates. It must - be of dtype :attr:`numpy.float32`. If ``None``, this function - will return a newly allocated array (Default value = None). - - Returns: - :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: - Absolute coordinate vector(s). If ``out`` is provided, a - reference to it is returned. - """ # noqa: E501 - fractions = np.asarray(fractional_coordinates).copy() - flatten = fractions.ndim == 1 - fractions = np.atleast_2d(fractions) - fractions = freud.util._convert_array(fractions, shape=(None, 3)) - out = freud.util._convert_array( - out, shape=fractions.shape, allow_copy=False) - - cdef const float[:, ::1] l_points = fractions - cdef unsigned int Np = l_points.shape[0] - cdef float[:, ::1] l_out = out - - self.thisptr.makeAbsolute( - &l_points[0, 0], Np, - &l_out[0, 0]) - - return np.squeeze(out) if flatten else out - - def make_fractional(self, absolute_coordinates, out=None): - r"""Convert absolute coordinates into fractional coordinates. - - Args: - absolute_coordinates (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Absolute coordinate vector(s). - out (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray` or :code:`None`): - The array in which to place the fractional positions. It must - be of dtype :attr:`numpy.float32`. If ``None``, this function - will return a newly allocated array (Default value = None). - - Returns: - :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: - Fractional coordinate vector(s). If ``out`` is provided, a - reference to it is returned. - """ # noqa: E501 - vecs = np.asarray(absolute_coordinates).copy() - flatten = vecs.ndim == 1 - vecs = np.atleast_2d(vecs) - vecs = freud.util._convert_array(vecs, shape=(None, 3)) - out = freud.util._convert_array( - out, shape=vecs.shape, allow_copy=False) - - cdef const float[:, ::1] l_points = vecs - cdef unsigned int Np = l_points.shape[0] - cdef float[:, ::1] l_out = out - - self.thisptr.makeFractional( - &l_points[0, 0], Np, - &l_out[0, 0]) - - return np.squeeze(out) if flatten else out - - def get_images(self, vecs): - r"""Returns the images corresponding to unwrapped vectors. - - Args: - vecs (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Coordinates of unwrapped vector(s). - - Returns: - :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: - Image index vector(s). - """ # noqa: E501 - vecs = np.asarray(vecs) - flatten = vecs.ndim == 1 - vecs = np.atleast_2d(vecs) - vecs = freud.util._convert_array(vecs, shape=(None, 3)) - - images = np.zeros(vecs.shape, dtype=np.int32) - cdef const float[:, ::1] l_points = vecs - cdef const int[:, ::1] l_result = images - cdef unsigned int Np = l_points.shape[0] - self.thisptr.getImages( &l_points[0, 0], Np, - &l_result[0, 0]) - - return np.squeeze(images) if flatten else images - - def get_box_vector(self, i): - r"""Get the box vector with index :math:`i`. - - Args: - i (unsigned int): - Index (:math:`0 \leq i < d`) of the box vector, where - :math:`d` is the box dimension (2 or 3). - - Returns: - :math:`\left(3, \right)` :class:`numpy.ndarray`: - Box vector with index :math:`i`. - """ - return self.to_matrix()[:, i] - - @property - def v1(self): - """:math:`(3, )` :class:`np.ndarray`: The first box vector - :math:`(L_x, 0, 0)`.""" - return self.get_box_vector(0) - - @property - def v2(self): - r""":math:`(3, )` :class:`np.ndarray`: The second box vector - :math:`(xy \times L_y, L_y, 0)`.""" - return self.get_box_vector(1) - - @property - def v3(self): - r""":math:`(3, )` :class:`np.ndarray`: The third box vector - :math:`(xz \times L_z, yz \times L_z, L_z)`.""" - return self.get_box_vector(2) - - def wrap(self, vecs, out=None): - r"""Wrap an array of vectors into the box, using periodic boundaries. - - .. note:: Since the origin of the box is in the center, wrapping is - equivalent to applying the minimum image convention to the - input vectors. - - Args: - vecs (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Unwrapped vector(s). - out (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray` or :code:`None`): - The array in which to place the wrapped vectors. It must be of - dtype :attr:`numpy.float32`. If ``None``, this function will - return a newly allocated array (Default value = None). - - Returns: - :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: - Vector(s) wrapped into the box. If ``out`` is provided, a - reference to it is returned. - """ # noqa: E501 - vecs = np.asarray(vecs) - flatten = vecs.ndim == 1 - vecs = np.atleast_2d(vecs) - vecs = freud.util._convert_array(vecs, shape=(None, 3)) - out = freud.util._convert_array( - out, shape=vecs.shape, allow_copy=False) - - cdef const float[:, ::1] l_points = vecs - cdef unsigned int Np = l_points.shape[0] - cdef float[:, ::1] l_out = out - - self.thisptr.wrap( &l_points[0, 0], - Np, &l_out[0, 0]) - - return np.squeeze(out) if flatten else out - - def unwrap(self, vecs, imgs, out=None): - r"""Unwrap an array of vectors inside the box back into real space, - using an array of image indices that determine how many times to unwrap - in each dimension. - - Args: - vecs (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Vector(s) to be unwrapped. - imgs (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Image indices for vector(s). - out (:math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray` or :code:`None`): - The array in which to place the unwrapped vectors. It must be - of dtype :attr:`numpy.float32`. If ``None``, this function - will return a newly allocated array (Default value = None). - - Returns: - :math:`\left(3, \right)` or :math:`\left(N, 3\right)` :class:`numpy.ndarray`: - Unwrapped vector(s). If ``out`` is provided, a reference to it is - returned. - """ # noqa: E501 - vecs = np.asarray(vecs) - flatten = vecs.ndim == 1 - vecs = np.atleast_2d(vecs) - imgs = np.atleast_2d(imgs) - if vecs.shape[0] != imgs.shape[0]: - # Broadcasts (1, 3) to (N, 3) for both arrays - vecs, imgs = np.broadcast_arrays(vecs, imgs) - vecs = freud.util._convert_array(vecs, shape=(None, 3)).copy() - imgs = freud.util._convert_array( - imgs, shape=vecs.shape, dtype=np.int32) - out = freud.util._convert_array( - out, shape=vecs.shape, allow_copy=False) - - cdef const float[:, ::1] l_points = vecs - cdef const int[:, ::1] l_imgs = imgs - cdef unsigned int Np = l_points.shape[0] - cdef float[:, ::1] l_out = out - - self.thisptr.unwrap( &l_points[0, 0], - &l_imgs[0, 0], Np, - &l_out[0, 0]) - - return np.squeeze(out) if flatten else out - - def center_of_mass(self, vecs, masses=None): - r"""Compute center of mass of an array of vectors, using periodic boundaries. - - This calculation accounts for periodic images. `This Wikipedia page - `__ - describes the mathematics of this method. - - Example:: - - >>> import freud - >>> import numpy as np - >>> box = freud.Box.cube(10) - >>> points = [[-1, -1, 0], [-1, 1, 0], [2, 0, 0]] - >>> np.mean(points, axis=0) # Does not account for periodic images - array([0., 0., 0.]) - >>> box.center_of_mass(points) # Accounts for periodic images - array([-0.1845932, 0. , 0. ]) - - Args: - vecs (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Vectors used to find center of mass. - masses (:math:`\left(N,\right)` :class:`numpy.ndarray`): - Masses corresponding to each vector, defaulting to 1 if not - provided or :code:`None` (Default value = :code:`None`). - - Returns: - :math:`\left(3, \right)` :class:`numpy.ndarray`: - Center of mass. - """ # noqa: E501 - vecs = freud.util._convert_array(vecs, shape=(None, 3)) - cdef const float[:, ::1] l_points = vecs - - cdef float* l_masses_ptr = NULL - cdef float[::1] l_masses - if masses is not None: - l_masses = freud.util._convert_array(masses, shape=(len(vecs), )) - l_masses_ptr = &l_masses[0] - - cdef size_t Np = l_points.shape[0] - cdef vec3[float] result = self.thisptr.centerOfMass( - &l_points[0, 0], Np, l_masses_ptr) - return np.asarray([result.x, result.y, result.z]) - - def center(self, vecs, masses=None): - r"""Subtract center of mass from an array of vectors, using periodic boundaries. - - This calculation accounts for periodic images. `This Wikipedia page - `__ - describes the mathematics of this method. - - Example:: - - >>> import freud - >>> box = freud.Box.cube(10) - >>> points = [[-1, -1, 0], [-1, 1, 0], [2, 0, 0]] - >>> box.center(points) - array([[-0.8154068, -1.0000002, 0. ], - [-0.8154068, 1. , 0. ], - [ 2.1845937, 0. , 0. ]], dtype=float32) - - Args: - vecs (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Vectors to center. - masses (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Masses corresponding to each vector, defaulting to 1 if not - provided or :code:`None` (Default value = :code:`None`). - - Returns: - :math:`\left(N, 3\right)` :class:`numpy.ndarray`: - Vectors with center of mass subtracted. - """ # noqa: E501 - vecs = freud.util._convert_array(vecs, shape=(None, 3)).copy() - cdef const float[:, ::1] l_points = vecs - - cdef float* l_masses_ptr = NULL - cdef float[::1] l_masses - if masses is not None: - l_masses = freud.util._convert_array(masses, shape=(len(vecs), )) - l_masses_ptr = &l_masses[0] - - cdef size_t Np = l_points.shape[0] - self.thisptr.center( &l_points[0, 0], Np, l_masses_ptr) - return vecs - - def compute_distances(self, query_points, points): - r"""Calculate distances between two sets of points, using periodic boundaries. - - Distances are calculated row-wise, i.e. ``distances[i]`` is the - distance from ``query_points[i]`` to ``points[i]``. - - Args: - query_points (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Array of query points. - points (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Array of points. - - Returns: - :math:`\left(N, \right)` :class:`numpy.ndarray`: - Array of distances between query points and points. - """ # noqa: E501 - - query_points = freud.util._convert_array( - np.atleast_2d(query_points), shape=(None, 3)) - points = freud.util._convert_array( - np.atleast_2d(points), shape=(None, 3)) - - cdef: - const float[:, ::1] l_query_points = query_points - const float[:, ::1] l_points = points - size_t n_query_points = query_points.shape[0] - size_t n_points = points.shape[0] - float[::1] distances = np.empty( - n_query_points, dtype=np.float32) - - self.thisptr.computeDistances( - &l_query_points[0, 0], n_query_points, - &l_points[0, 0], n_points, - &distances[0]) - return np.asarray(distances) - - def compute_all_distances(self, query_points, points): - r"""Calculate distances between all pairs of query points and points, using periodic boundaries. - - Distances are calculated pairwise, i.e. ``distances[i, j]`` is the - distance from ``query_points[i]`` to ``points[j]``. - - Args: - query_points (:math:`\left(N_{query\_points}, 3 \right)` :class:`numpy.ndarray`): - Array of query points. - points (:math:`\left(N_{points}, 3 \right)` :class:`numpy.ndarray`): - Array of points with same length as ``query_points``. - - Returns: - :math:`\left(N_{query\_points}, N_{points}, \right)` :class:`numpy.ndarray`: - Array of distances between query points and points. - """ # noqa: E501 - query_points = freud.util._convert_array( - np.atleast_2d(query_points), shape=(None, 3)) - points = freud.util._convert_array( - np.atleast_2d(points), shape=(None, 3)) - - cdef: - const float[:, ::1] l_query_points = query_points - const float[:, ::1] l_points = points - size_t n_query_points = query_points.shape[0] - size_t n_points = points.shape[0] - float[:, ::1] distances = np.empty( - [n_query_points, n_points], dtype=np.float32) - - self.thisptr.computeAllDistances( - &l_query_points[0, 0], n_query_points, - &l_points[0, 0], n_points, - &distances[0, 0]) - - return np.asarray(distances) - - def contains(self, points): - r"""Compute a boolean array (mask) corresponding to point membership in a box. - - This calculation computes particle membership based on conventions - defined by :class:`Box`, ignoring periodicity. This means that in a - cubic (3D) box with dimensions ``L``, particles would be considered - inside the box if their coordinates are between ``[-L/2, L/2]``. - Particles laying at a coordinate such as ``[0, L, 0]`` would be - considered outside the box. More information about coordinate - conventions can be found in the documentation on `Using boxes - `__ - and `periodic boundary conditions - `__. - - Example:: - - >>> import freud - >>> box = freud.Box.cube(10) - >>> points = [[-4, 0, 0], [10, 0, 0], [0, -7, 0]] - >>> box.contains(points) - array([ True, False, False]) - - Args: - points (:math:`\left(N, 3\right)` :class:`numpy.ndarray`): - Array of points. - - Returns: - :math:`\left(N, \right)` :class:`numpy.ndarray`: - Array of booleans, where ``True`` corresponds to points within - the box, and ``False`` corresponds to points outside the box. - """ # noqa: E501 - - points = freud.util._convert_array( - np.atleast_2d(points), shape=(None, 3)) - - cdef: - const float[:, ::1] l_points = points - size_t n_points = points.shape[0] - - contains_mask = freud.util._convert_array( - np.ones(n_points), dtype=bool) - cdef cpp_bool[::1] l_contains_mask = contains_mask - - self.thisptr.contains( - &l_points[0, 0], n_points, - &l_contains_mask[0]) - - return np.array(l_contains_mask).astype(bool) - - @property - def cubic(self): - """bool: Whether the box is a cube.""" - return ( - not self.is2D - and np.allclose( - [self.Lx, self.Lx, self.Ly, self.Ly, self.Lz, self.Lz], - [self.Ly, self.Lz, self.Lx, self.Lz, self.Lx, self.Ly], - rtol=1e-5, - atol=1e-5, - ) - and np.allclose(0, [self.xy, self.yz, self.xz], rtol=1e-5, atol=1e-5) - ) - - @property - def periodic(self): - r""":math:`\left(3, \right)` :class:`numpy.ndarray`: Get or set the - periodicity of the box in each dimension.""" - periodic = self.thisptr.getPeriodic() - return np.asarray([periodic.x, periodic.y, periodic.z]) - - @periodic.setter - def periodic(self, periodic): - # Allow passing a single value - try: - self.thisptr.setPeriodic(periodic[0], periodic[1], periodic[2]) - except TypeError: - # Allow single value to be passed for all directions - self.thisptr.setPeriodic(periodic, periodic, periodic) - - @property - def periodic_x(self): - """bool: Get or set the periodicity of the box in x.""" - return self.thisptr.getPeriodicX() - - @periodic_x.setter - def periodic_x(self, periodic): - self.thisptr.setPeriodicX(periodic) - - @property - def periodic_y(self): - """bool: Get or set the periodicity of the box in y.""" - return self.thisptr.getPeriodicY() - - @periodic_y.setter - def periodic_y(self, periodic): - self.thisptr.setPeriodicY(periodic) - - @property - def periodic_z(self): - """bool: Get or set the periodicity of the box in z.""" - return self.thisptr.getPeriodicZ() - - @periodic_z.setter - def periodic_z(self, periodic): - self.thisptr.setPeriodicZ(periodic) - - def to_dict(self): - r"""Return box as dictionary. - - Example:: - - >>> box = freud.box.Box.cube(L=10) - >>> box.to_dict() - {'Lx': 10.0, 'Ly': 10.0, 'Lz': 10.0, - 'xy': 0.0, 'xz': 0.0, 'yz': 0.0, 'dimensions': 3} - - Returns: - dict: Box parameters - """ - return { - 'Lx': self.Lx, - 'Ly': self.Ly, - 'Lz': self.Lz, - 'xy': self.xy, - 'xz': self.xz, - 'yz': self.yz, - 'dimensions': self.dimensions} - - def to_matrix(self): - r"""Returns the box matrix (3x3). - - Example:: - - >>> box = freud.box.Box.cube(L=10) - >>> box.to_matrix() - array([[10., 0., 0.], - [ 0., 10., 0.], - [ 0., 0., 10.]]) - - Returns: - :math:`\left(3, 3\right)` :class:`numpy.ndarray`: Box matrix - """ - return np.asarray([[self.Lx, self.xy * self.Ly, self.xz * self.Lz], - [0, self.Ly, self.yz * self.Lz], - [0, 0, self.Lz]]) - - def to_box_lengths_and_angles(self): - r"""Return the box lengths and angles. - - Returns: - tuple: - The box vector lengths and angles in radians - :math:`(L_1, L_2, L_3, \alpha, \beta, \gamma)`. - """ - alpha = np.arccos( - (self.xy * self.xz + self.yz) - / (np.sqrt(1 + self.xy**2) * np.sqrt(1 + self.xz**2 + self.yz**2)) - ) - beta = np.arccos(self.xz/np.sqrt(1+self.xz**2+self.yz**2)) - gamma = np.arccos(self.xy/np.sqrt(1+self.xy**2)) - L1 = self.Lx - a2 = [self.Ly*self.xy, self.Ly, 0] - a3 = [self.Lz*self.xz, self.Lz*self.yz, self.Lz] - L2 = np.linalg.norm(a2) - L3 = np.linalg.norm(a3) - return (L1, L2, L3, alpha, beta, gamma) - - def __repr__(self): - return ("freud.box.{cls}(Lx={Lx}, Ly={Ly}, Lz={Lz}, " - "xy={xy}, xz={xz}, yz={yz}, " - "is2D={is2D})").format(cls=type(self).__name__, - Lx=self.Lx, - Ly=self.Ly, - Lz=self.Lz, - xy=self.xy, - xz=self.xz, - yz=self.yz, - is2D=self.is2D) - - def __str__(self): - return repr(self) - - def __richcmp__(self, other, int op): - r"""Implement all comparisons for Cython extension classes""" - cdef Box c_other - try: - c_other = other - if op == Py_EQ: - return dereference(self.thisptr) == dereference(c_other.thisptr) - if op == Py_NE: - return dereference(self.thisptr) != dereference(c_other.thisptr) - except TypeError: - # Cython cast to Box failed - pass - return NotImplemented - - def __mul__(self, scale): - if scale > 0: - return self.__class__(Lx=self.Lx*scale, - Ly=self.Ly*scale, - Lz=self.Lz*scale, - xy=self.xy, xz=self.xz, yz=self.yz, - is2D=self.is2D) - else: - raise ValueError("Box can only be multiplied by positive values.") - - def __rmul__(self, scale): - return self * scale - - def plot(self, title=None, ax=None, image=[0, 0, 0], *args, **kwargs): - """Plot a :class:`~.box.Box` object. - - Args: - title (str): - Title of the graph. (Default value = :code:`None`). - ax (:class:`matplotlib.axes.Axes`): - Axes object to plot. If :code:`None`, make a new axes and - figure object. If plotting a 3D box, the axes must be 3D. - (Default value = :code:`None`). - image (list): - The periodic image location at which to draw the box (Default - value = :code:`[0, 0, 0]`). - *args: - Passed on to :meth:`mpl_toolkits.mplot3d.Axes3D.plot` or - :meth:`matplotlib.axes.Axes.plot`. - **kwargs: - Passed on to :meth:`mpl_toolkits.mplot3d.Axes3D.plot` or - :meth:`matplotlib.axes.Axes.plot`. - """ - import freud.plot - return freud.plot.box_plot(self, title=title, ax=ax, image=image, - *args, **kwargs) - - @classmethod - def from_box(cls, box, dimensions=None): - r"""Initialize a Box instance from a box-like object. - - Args: - box: - A box-like object - dimensions (int): - Dimensionality of the box (Default value = None) - - .. note:: Objects that can be converted to freud boxes include - lists like :code:`[Lx, Ly, Lz, xy, xz, yz]`, - dictionaries with keys - :code:`'Lx', 'Ly', 'Lz', 'xy', 'xz', 'yz', 'dimensions'`, - objects with attributes - :code:`Lx, Ly, Lz, xy, xz, yz, dimensions`, - 3x3 matrices (see :meth:`~.from_matrix`), - or existing :class:`freud.box.Box` objects. - - If any of :code:`Lz, xy, xz, yz` are not provided, they will - be set to 0. - - If all values are provided, a triclinic box will be - constructed. If only :code:`Lx, Ly, Lz` are provided, an - orthorhombic box will be constructed. If only :code:`Lx, Ly` - are provided, a rectangular (2D) box will be constructed. - - If the optional :code:`dimensions` argument is given, this - will be used as the box dimensionality. Otherwise, the box - dimensionality will be detected from the :code:`dimensions` - of the provided box. If no dimensions can be detected, the - box will be 2D if :code:`Lz == 0`, and 3D otherwise. - - Returns: - :class:`freud.box.Box`: The resulting box object. - """ - if np.asarray(box).shape == (3, 3): - # Handles 3x3 matrices - return cls.from_matrix(box, dimensions=dimensions) - try: - # Handles freud.box.Box and objects with attributes - Lx = box.Lx - Ly = box.Ly - Lz = getattr(box, 'Lz', 0) - xy = getattr(box, 'xy', 0) - xz = getattr(box, 'xz', 0) - yz = getattr(box, 'yz', 0) - if dimensions is None: - dimensions = getattr(box, 'dimensions', None) - elif dimensions != getattr(box, 'dimensions', dimensions): - raise ValueError( - "The provided dimensions argument conflicts with the " - "dimensions attribute of the provided box object.") - except AttributeError: - try: - # Handle dictionary-like - Lx = box['Lx'] - Ly = box['Ly'] - Lz = box.get('Lz', 0) - xy = box.get('xy', 0) - xz = box.get('xz', 0) - yz = box.get('yz', 0) - if dimensions is None: - dimensions = box.get('dimensions', None) - else: - if dimensions != box.get('dimensions', dimensions): - raise ValueError( - "The provided dimensions argument conflicts with " - "the dimensions attribute of the provided box " - "object.") - except (IndexError, KeyError, TypeError): - if not len(box) in [2, 3, 6]: - raise ValueError( - "List-like objects must have length 2, 3, or 6 to be " - "converted to freud.box.Box.") - # Handle list-like - Lx = box[0] - Ly = box[1] - Lz = box[2] if len(box) > 2 else 0 - xy, xz, yz = box[3:6] if len(box) == 6 else (0, 0, 0) - except: # noqa - logger.debug('Supplied box cannot be converted to type ' - 'freud.box.Box.') - raise - - # Infer dimensions if not provided. - if dimensions is None: - dimensions = 2 if Lz == 0 else 3 - is2D = (dimensions == 2) - return cls(Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, is2D=is2D) - - @classmethod - def from_matrix(cls, box_matrix, dimensions=None): - r"""Initialize a Box instance from a box matrix. - - For more information and the source for this code, see: - `HOOMD-blue's box documentation \ - `_. - - Args: - box_matrix (array-like): - A 3x3 matrix or list of lists - dimensions (int): - Number of dimensions (Default value = :code:`None`) - - Returns: - :class:`freud.box.Box`: The resulting box object. - """ - box_matrix = np.asarray(box_matrix, dtype=np.float32) - v0 = box_matrix[:, 0] - v1 = box_matrix[:, 1] - v2 = box_matrix[:, 2] - Lx = np.sqrt(np.dot(v0, v0)) - a2x = np.dot(v0, v1) / Lx - Ly = np.sqrt(np.dot(v1, v1) - a2x * a2x) - xy = a2x / Ly - v0xv1 = np.cross(v0, v1) - v0xv1mag = np.sqrt(np.dot(v0xv1, v0xv1)) - Lz = np.dot(v2, v0xv1) / v0xv1mag - if Lz != 0: - a3x = np.dot(v0, v2) / Lx - xz = a3x / Lz - yz = (np.dot(v1, v2) - a2x * a3x) / (Ly * Lz) - else: - xz = yz = 0 - if dimensions is None: - dimensions = 2 if Lz == 0 else 3 - is2D = (dimensions == 2) - return cls(Lx=Lx, Ly=Ly, Lz=Lz, - xy=xy, xz=xz, yz=yz, is2D=is2D) - - @classmethod - def cube(cls, L=None): - r"""Construct a cubic box with equal lengths. - - Args: - L (float): The edge length - - Returns: - :class:`freud.box.Box`: The resulting box object. - """ - # classmethods compiled with cython don't appear to support - # named access to positional arguments, so we keep this to - # recover the behavior - if L is None: - raise TypeError("cube() missing 1 required positional argument: L") - return cls(Lx=L, Ly=L, Lz=L, xy=0, xz=0, yz=0, is2D=False) - - @classmethod - def square(cls, L=None): - r"""Construct a 2-dimensional (square) box with equal lengths. - - Args: - L (float): The edge length - - Returns: - :class:`freud.box.Box`: The resulting box object. - """ - # classmethods compiled with cython don't appear to support - # named access to positional arguments, so we keep this to - # recover the behavior - if L is None: - raise TypeError("square() missing 1 required " - "positional argument: L") - return cls(Lx=L, Ly=L, Lz=0, xy=0, xz=0, yz=0, is2D=True) - - @classmethod - def from_box_lengths_and_angles( - cls, L1, L2, L3, alpha, beta, gamma, dimensions=None, - ): - r"""Construct a box from lengths and angles (in radians). - - All the angles provided must be between 0 and :math:`\pi`. - - Args: - L1 (float): - The length of the first lattice vector. - L2 (float): - The length of the second lattice vector. - L3 (float): - The length of the third lattice vector. - alpha (float): - The angle between second and third lattice vector in radians (must be - between 0 and :math:`\pi`). - beta (float): - The angle between first and third lattice vector in radians (must be - between 0 and :math:`\pi`). - gamma (float): - The angle between the first and second lattice vector in radians (must - be between 0 and :math:`\pi`). - dimensions (int): - The number of dimensions (Default value = :code:`None`). - - Returns: - :class:`freud.box.Box`: The resulting box object. - """ - if not 0 < alpha < np.pi: - raise ValueError("alpha must be between 0 and pi.") - if not 0 < beta < np.pi: - raise ValueError("beta must be between 0 and pi.") - if not 0 < gamma < np.pi: - raise ValueError("gamma must be between 0 and pi.") - a1 = np.array([L1, 0, 0]) - a2 = np.array([L2 * np.cos(gamma), L2 * np.sin(gamma), 0]) - a3x = np.cos(beta) - a3y = (np.cos(alpha) - np.cos(beta) * np.cos(gamma)) / np.sin(gamma) - under_sqrt = 1 - a3x**2 - a3y**2 - if under_sqrt < 0: - raise ValueError("The provided angles can not form a valid box.") - a3z = np.sqrt(under_sqrt) - a3 = np.array([L3 * a3x, L3 * a3y, L3 * a3z]) - if dimensions is None: - dimensions = 2 if L3 == 0 else 3 - return cls.from_matrix(np.array([a1, a2, a3]).T, dimensions=dimensions) - - -cdef BoxFromCPP(const freud._box.Box & cppbox): - b = Box(cppbox.getLx(), cppbox.getLy(), cppbox.getLz(), - cppbox.getTiltFactorXY(), cppbox.getTiltFactorXZ(), - cppbox.getTiltFactorYZ(), cppbox.is2D()) - b.periodic = [cppbox.getPeriodicX(), - cppbox.getPeriodicY(), - cppbox.getPeriodicZ()] - return b From d8dda25e612d5697885560d5019113081575e331 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Mon, 3 Jun 2024 21:54:39 +0000 Subject: [PATCH 06/25] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- CMakeLists.txt | 1 - cpp/CMakeLists.txt | 14 +-- cpp/box/Box.h | 95 ++++++++------- cpp/box/CMakeLists.txt | 19 ++- cpp/box/module-box.cc | 2 + cpp/util/utils.h | 1 - freud/CMakeLists.txt | 25 ++-- freud/__init__.py | 54 ++++----- freud/box.py | 254 +++++++++++++++++++++++------------------ freud/util.py | 50 +++++--- 10 files changed, 265 insertions(+), 250 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 836f8f711..881689296 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -79,4 +79,3 @@ set(ignoreMe "${SKBUILD}") add_subdirectory(cpp) add_subdirectory(freud) - diff --git a/cpp/CMakeLists.txt b/cpp/CMakeLists.txt index 4c4c2c1c3..8a625289e 100644 --- a/cpp/CMakeLists.txt +++ b/cpp/CMakeLists.txt @@ -25,15 +25,10 @@ endif() add_subdirectory(box) # commented out for now, uncomment them as the conversion to pybind11 progresses -#add_subdirectory(cluster) -#add_subdirectory(density) -#add_subdirectory(diffraction) -#add_subdirectory(environment) -#add_subdirectory(locality) -#add_subdirectory(order) -#add_subdirectory(parallel) -#add_subdirectory(pmft) -#add_subdirectory(util) +# add_subdirectory(cluster) add_subdirectory(density) +# add_subdirectory(diffraction) add_subdirectory(environment) +# add_subdirectory(locality) add_subdirectory(order) add_subdirectory(parallel) +# add_subdirectory(pmft) add_subdirectory(util) #[[ add_library( @@ -74,4 +69,3 @@ if(CMAKE_EXPORT_COMPILE_COMMANDS) ${PROJECT_SOURCE_DIR}) endif() ]] - diff --git a/cpp/box/Box.h b/cpp/box/Box.h index c87db25fe..27de2a206 100644 --- a/cpp/box/Box.h +++ b/cpp/box/Box.h @@ -6,12 +6,12 @@ #include #include +#include #include #include -#include -#include "utils.h" #include "VectorMath.h" +#include "utils.h" #include @@ -158,7 +158,7 @@ class Box //! Get current stored inverse of L std::array getLinv() const { - return { m_Linv.x, m_Linv.y, m_Linv.z }; + return {m_Linv.x, m_Linv.y, m_Linv.z}; } //! Get tilt factor xy @@ -232,8 +232,8 @@ class Box */ void makeAbsolutePython(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const { - vec3 *vecs_data = (vec3 *)(vecs.data()); - vec3 *out_data = (vec3 *)(out.data()); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { @@ -267,8 +267,8 @@ class Box */ void makeFractionalPython(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const { - vec3 *vecs_data = (vec3 *)(vecs.data()); - vec3 *out_data = (vec3 *)(out.data()); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { @@ -298,11 +298,10 @@ class Box * \param Nvecs Number of vectors \param res Array to save the images */ - void getImages(pybind11_array vecs, unsigned int Nvecs, - pybind11::array_t res) const + void getImages(pybind11_array vecs, unsigned int Nvecs, pybind11::array_t res) const { - vec3 *vecs_data = (vec3 *)(vecs.data()); - vec3 *out_data = (vec3 *)(res.data()); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (res.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { @@ -346,8 +345,8 @@ class Box */ void wrapPython(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const { - vec3 *vecs_data = (vec3 *)(vecs.data()); - vec3 *out_data = (vec3 *)(out.data()); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { @@ -363,11 +362,11 @@ class Box * \param out The array in which to place the wrapped vectors. */ void unwrap(pybind11_array vecs, pybind11_array images, unsigned int Nvecs, - pybind11_array out) const + pybind11_array out) const { - vec3 *vecs_data = (vec3 *)(vecs.data()); - vec3 *images_data = (vec3 *)(images.data()); - vec3 *out_data = (vec3 *)(out.data()); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* images_data = (vec3*) (images.data()); + vec3* out_data = (vec3*) (out.data()); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { @@ -387,7 +386,7 @@ class Box * \param masses Optional array of masses, of length Nvecs * \return Center of mass as a vec3 */ - vec3 centerOfMass(vec3 *vecs, size_t Nvecs, const float *masses) const + vec3 centerOfMass(vec3* vecs, size_t Nvecs, const float* masses) const { // This roughly follows the implementation in // https://en.wikipedia.org/wiki/Center_of_mass#Systems_with_periodic_boundary_conditions @@ -397,8 +396,7 @@ class Box for (size_t i = 0; i < Nvecs; ++i) { vec3 phase(constants::TWO_PI * makeFractional(vecs[i])); - vec3> xi(std::polar(float(1.0), phase.x), - std::polar(float(1.0), phase.y), + vec3> xi(std::polar(float(1.0), phase.x), std::polar(float(1.0), phase.y), std::polar(float(1.0), phase.z)); float mass = masses[i]; total_mass += mass; @@ -406,18 +404,17 @@ class Box } xi_mean /= std::complex(total_mass, 0); - return wrap(makeAbsolute(vec3(std::arg(xi_mean.x), - std::arg(xi_mean.y), - std::arg(xi_mean.z)) + return wrap(makeAbsolute(vec3(std::arg(xi_mean.x), std::arg(xi_mean.y), std::arg(xi_mean.z)) / constants::TWO_PI)); } - std::array centerOfMassPython(pybind11_array vecs, size_t Nvecs, pybind11_array masses) const + std::array centerOfMassPython(pybind11_array vecs, size_t Nvecs, + pybind11_array masses) const { - vec3 *vecs_data = (vec3 *)(vecs.data()); - float *masses_data = (float *)(masses.data()); + vec3* vecs_data = (vec3*) (vecs.data()); + float* masses_data = (float*) (masses.data()); auto com = centerOfMass(vecs_data, Nvecs, masses_data); - return { com.x, com.y, com.z }; + return {com.x, com.y, com.z}; } //! Subtract center of mass from vectors @@ -427,8 +424,8 @@ class Box */ void center(pybind11_array vecs, unsigned int Nvecs, pybind11_array masses) const { - vec3 *vecs_data = (vec3 *)(vecs.data()); - float *masses_data = (float *)(masses.data()); + vec3* vecs_data = (vec3*) (vecs.data()); + float* masses_data = (float*) (masses.data()); vec3 com(centerOfMass(vecs_data, Nvecs, masses_data)); util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { @@ -461,9 +458,9 @@ class Box pybind11_array points, const unsigned int n_points, pybind11_array distances) const { - vec3 *query_points_data = (vec3 *)(query_points.data()); - vec3 *points_data = (vec3 *)(points.data()); - float *distances_data = (float *)(distances.data()); + vec3* query_points_data = (vec3*) (query_points.data()); + vec3* points_data = (vec3*) (points.data()); + float* distances_data = (float*) (distances.data()); if (n_query_points != n_points) { @@ -489,21 +486,21 @@ class Box pybind11_array points, const unsigned int n_points, pybind11_array distances) const { - vec3 *query_points_data = (vec3 *)(query_points.data()); - vec3 *points_data = (vec3 *)(points.data()); - float *distances_data = (float *)(distances.data()); + vec3* query_points_data = (vec3*) (query_points.data()); + vec3* points_data = (vec3*) (points.data()); + float* distances_data = (float*) (distances.data()); - util::forLoopWrapper2D( - 0, n_query_points, 0, n_points, [&](size_t begin_n, size_t end_n, size_t begin_m, size_t end_m) { - for (size_t i = begin_n; i < end_n; ++i) - { - for (size_t j = begin_m; j < end_m; ++j) - { - distances_data[i * n_points + j] = computeDistance( - query_points_data[i], points_data[j]); - } - } - }); + util::forLoopWrapper2D(0, n_query_points, 0, n_points, + [&](size_t begin_n, size_t end_n, size_t begin_m, size_t end_m) { + for (size_t i = begin_n; i < end_n; ++i) + { + for (size_t j = begin_m; j < end_m; ++j) + { + distances_data[i * n_points + j] + = computeDistance(query_points_data[i], points_data[j]); + } + } + }); } //! Get mask of points that fit inside the box. @@ -512,10 +509,10 @@ class Box \param contains_mask Mask of points inside the box. */ void contains(pybind11_array points, const unsigned int n_points, - pybind11_array contains_mask) const + pybind11_array contains_mask) const { - vec3 *points_data = (vec3 *)(points.data()); - bool *contains_mask_data = (bool *)(contains_mask.data()); + vec3* points_data = (vec3*) (points.data()); + bool* contains_mask_data = (bool*) (contains_mask.data()); util::forLoopWrapper(0, n_points, [&](size_t begin, size_t end) { std::transform(&points_data[begin], &points_data[end], &contains_mask_data[begin], diff --git a/cpp/box/CMakeLists.txt b/cpp/box/CMakeLists.txt index 223e6c47e..73a1959ad 100644 --- a/cpp/box/CMakeLists.txt +++ b/cpp/box/CMakeLists.txt @@ -2,23 +2,20 @@ pybind11_add_module(_box SHARED module-box.cc) target_link_libraries(_box PUBLIC TBB::tbb) # this probably isn't needed for box, but may be needed by other modules -#target_compile_definitions( -# box -# # Avoid deprecation warnings for unsupported NumPy API versions. See -# # https://numpy.org/doc/1.19/reference/c-api/deprecations.html -# PRIVATE "NPY_NO_DEPRECATED_API=NPY_1_10_API_VERSION" -# # Default voro++ verbosity is high. -# PRIVATE "VOROPP_VERBOSE=1") +# target_compile_definitions( box # Avoid deprecation warnings for unsupported +# NumPy API versions. See # +# https://numpy.org/doc/1.19/reference/c-api/deprecations.html PRIVATE +# "NPY_NO_DEPRECATED_API=NPY_1_10_API_VERSION" # Default voro++ verbosity is +# high. PRIVATE "VOROPP_VERBOSE=1") if(APPLE) - set_target_properties(_box PROPERTIES INSTALL_RPATH "@loader_path") + set_target_properties(_box PROPERTIES INSTALL_RPATH "@loader_path") else() - set_target_properties(_box PROPERTIES INSTALL_RPATH "\$ORIGIN") + set_target_properties(_box PROPERTIES INSTALL_RPATH "\$ORIGIN") endif() if(_using_conda OR DEFINED ENV{CIBUILDWHEEL}) - set_target_properties(_box - PROPERTIES INSTALL_RPATH_USE_LINK_PATH True) + set_target_properties(_box PROPERTIES INSTALL_RPATH_USE_LINK_PATH True) endif() # install diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc index 37441af84..297759c38 100644 --- a/cpp/box/module-box.cc +++ b/cpp/box/module-box.cc @@ -1,3 +1,5 @@ +// Copyright (c) 2010-2023 The Regents of the University of Michigan +// This file is from the freud project, released under the BSD 3-Clause License. #include #include diff --git a/cpp/util/utils.h b/cpp/util/utils.h index d52cf4b89..3bee6c86b 100644 --- a/cpp/util/utils.h +++ b/cpp/util/utils.h @@ -10,7 +10,6 @@ #include #include - namespace freud { namespace util { //! Clip v if it is outside the range [lo, hi]. diff --git a/freud/CMakeLists.txt b/freud/CMakeLists.txt index f8fc35192..16c8dc557 100644 --- a/freud/CMakeLists.txt +++ b/freud/CMakeLists.txt @@ -20,23 +20,14 @@ if(${PYTHON_VERSION_MAJOR} EQUAL 3 add_compile_options("-Wno-deprecated-declarations") endif() -set(python_files - box.py - util.py) - #cluster.py - #density.py - #diffraction.py - #environment.py - #locality.py - #order.py - #parallel.py - #pmft.py) +set(python_files box.py util.py) +# cluster.py density.py diffraction.py environment.py locality.py order.py +# parallel.py pmft.py) install(FILES ${python_files} DESTINATION freud) -# THIS WILL NEED TO BE MOVED TO THE BUILD OF THE ORDER MODULE -# The SolidLiquid class has an instance of cluster::Cluster as a member, so -# including the header requires the Cluster.h header. Would prefer to inherit -# this information from the _order library, but that's not possible since we're -# linking to libfreud. -#target_include_directories(order PUBLIC ${PROJECT_SOURCE_DIR}/cpp/cluster) +# THIS WILL NEED TO BE MOVED TO THE BUILD OF THE ORDER MODULE The SolidLiquid +# class has an instance of cluster::Cluster as a member, so including the header +# requires the Cluster.h header. Would prefer to inherit this information from +# the _order library, but that's not possible since we're linking to libfreud. +# target_include_directories(order PUBLIC ${PROJECT_SOURCE_DIR}/cpp/cluster) diff --git a/freud/__init__.py b/freud/__init__.py index aadf17ec1..f0da804ef 100644 --- a/freud/__init__.py +++ b/freud/__init__.py @@ -1,51 +1,41 @@ # Copyright (c) 2010-2023 The Regents of the University of Michigan # This file is from the freud project, released under the BSD 3-Clause License. -from . import ( +from . import ( # cluster,; data,; density,; diffraction,; environment,; interface,; locality,; msd,; order,; parallel,; pmft, box, - #cluster, - #data, - #density, - #diffraction, - #environment, - #interface, - #locality, - #msd, - #order, - #parallel, - #pmft, ) from .box import Box -#from .locality import AABBQuery, LinkCell, NeighborList -#from .parallel import NumThreads, get_num_threads, set_num_threads + +# from .locality import AABBQuery, LinkCell, NeighborList +# from .parallel import NumThreads, get_num_threads, set_num_threads # Override TBB's default autoselection. This is necessary because once the # automatic selection runs, the user cannot change it. -#set_num_threads(0) +# set_num_threads(0) __version__ = "3.0.0" __all__ = [ "__version__", "box", - #"cluster", - #"data", - #"density", - #"diffraction", - #"environment", - #"interface", - #"locality", - #"msd", - #"order", - #"parallel", - #"pmft", + # "cluster", + # "data", + # "density", + # "diffraction", + # "environment", + # "interface", + # "locality", + # "msd", + # "order", + # "parallel", + # "pmft", "Box", - #"AABBQuery", - #"LinkCell", - #"NeighborList", - #"get_num_threads", - #"set_num_threads", - #"NumThreads", + # "AABBQuery", + # "LinkCell", + # "NeighborList", + # "get_num_threads", + # "set_num_threads", + # "NumThreads", ] __citation__ = """@article{freud2020, diff --git a/freud/box.py b/freud/box.py index acec10eac..0cf8b0ed6 100644 --- a/freud/box.py +++ b/freud/box.py @@ -12,9 +12,8 @@ import numpy as np -import freud.util - import freud._box +import freud.util logger = logging.getLogger(__name__) @@ -51,42 +50,42 @@ class Box: def __init__(self, Lx, Ly, Lz=0, xy=0, xz=0, yz=0, is2D=None): if is2D is None: - is2D = (Lz == 0) + is2D = Lz == 0 if is2D: if not (Lx and Ly): raise ValueError("Lx and Ly must be nonzero for 2D boxes.") elif Lz != 0 or xz != 0 or yz != 0: warnings.warn( - "Specifying z-dimensions in a 2-dimensional box " - "has no effect!") + "Specifying z-dimensions in a 2-dimensional box " "has no effect!" + ) else: if not (Lx and Ly and Lz): - raise ValueError( - "Lx, Ly, and Lz must be nonzero for 3D boxes.") + raise ValueError("Lx, Ly, and Lz must be nonzero for 3D boxes.") self._cpp_obj = freud._box.Box(Lx, Ly, Lz, xy, xz, yz, is2D) @property def L(self): r""":math:`\left(3, \right)` :class:`numpy.ndarray`: Get or set the box lengths along x, y, and z.""" - return np.asarray([self._cpp_obj.getLx(), - self._cpp_obj.getLy(), - self._cpp_obj.getLz()]) + return np.asarray( + [self._cpp_obj.getLx(), self._cpp_obj.getLy(), self._cpp_obj.getLz()] + ) @L.setter def L(self, value): try: if len(value) != 3: - raise ValueError('setL must be called with a scalar or a list ' - 'of length 3.') + raise ValueError( + "setL must be called with a scalar or a list " "of length 3." + ) except TypeError: # Will fail if object has no length value = (value, value, value) if self.is2D and value[2] != 0: warnings.warn( - "Specifying z-dimensions in a 2-dimensional box " - "has no effect!") + "Specifying z-dimensions in a 2-dimensional box " "has no effect!" + ) self._cpp_obj.setL(value[0], value[1], value[2]) @property @@ -146,10 +145,18 @@ def yz(self, value): def __eq__(self, other): if type(other) != freud.box.Box: return False - return self.Lx == other.Lx and self.Ly == other.Ly and self.Lz == other.Lz \ - and self.xy == other.xy and self.xz == other.xz and self.yz == other.yz \ - and self.is2D == other.is2D and self.periodic_x == other.periodic_x and \ - self.periodic_y == other.periodic_y and self.periodic_z == other.periodic_z + return ( + self.Lx == other.Lx + and self.Ly == other.Ly + and self.Lz == other.Lz + and self.xy == other.xy + and self.xz == other.xz + and self.yz == other.yz + and self.is2D == other.is2D + and self.periodic_x == other.periodic_x + and self.periodic_y == other.periodic_y + and self.periodic_z == other.periodic_z + ) @property def dimensions(self): @@ -199,8 +206,7 @@ def make_absolute(self, fractional_coordinates, out=None): flatten = fractions.ndim == 1 fractions = np.atleast_2d(fractions) fractions = freud.util._convert_array(fractions, shape=(None, 3)) - out = freud.util._convert_array( - out, shape=fractions.shape, allow_copy=False) + out = freud.util._convert_array(out, shape=fractions.shape, allow_copy=False) Np = fractions.shape[0] @@ -228,8 +234,7 @@ def make_fractional(self, absolute_coordinates, out=None): flatten = vecs.ndim == 1 vecs = np.atleast_2d(vecs) vecs = freud.util._convert_array(vecs, shape=(None, 3)) - out = freud.util._convert_array( - out, shape=vecs.shape, allow_copy=False) + out = freud.util._convert_array(out, shape=vecs.shape, allow_copy=False) Np = vecs.shape[0] @@ -315,8 +320,7 @@ def wrap(self, vecs, out=None): flatten = vecs.ndim == 1 vecs = np.atleast_2d(vecs) vecs = freud.util._convert_array(vecs, shape=(None, 3)) - out = freud.util._convert_array( - out, shape=vecs.shape, allow_copy=False) + out = freud.util._convert_array(out, shape=vecs.shape, allow_copy=False) Np = vecs.shape[0] self._cpp_obj.wrap(vecs, Np, out) @@ -351,10 +355,8 @@ def unwrap(self, vecs, imgs, out=None): # Broadcasts (1, 3) to (N, 3) for both arrays vecs, imgs = np.broadcast_arrays(vecs, imgs) vecs = freud.util._convert_array(vecs, shape=(None, 3)).copy() - imgs = freud.util._convert_array( - imgs, shape=vecs.shape, dtype=np.int32) - out = freud.util._convert_array( - out, shape=vecs.shape, allow_copy=False) + imgs = freud.util._convert_array(imgs, shape=vecs.shape, dtype=np.int32) + out = freud.util._convert_array(out, shape=vecs.shape, allow_copy=False) Np = vecs.shape[0] self._cpp_obj.unwrap(vecs, imgs, Np, out) @@ -394,7 +396,7 @@ def center_of_mass(self, vecs, masses=None): Np = vecs.shape[0] if masses is not None: - masses = freud.util._convert_array(masses, shape=(len(vecs), )) + masses = freud.util._convert_array(masses, shape=(len(vecs),)) else: masses = np.ones(Np) @@ -433,7 +435,7 @@ def center(self, vecs, masses=None): Np = vecs.shape[0] if masses is not None: - masses = freud.util._convert_array(masses, shape=(len(vecs), )) + masses = freud.util._convert_array(masses, shape=(len(vecs),)) else: masses = np.ones(Np) @@ -455,19 +457,20 @@ def compute_distances(self, query_points, points): Returns: :math:`\left(N, \right)` :class:`numpy.ndarray`: Array of distances between query points and points. - """ # noqa: E501 + """ # noqa: E501 query_points = freud.util._convert_array( - np.atleast_2d(query_points), shape=(None, 3)) - points = freud.util._convert_array( - np.atleast_2d(points), shape=(None, 3)) + np.atleast_2d(query_points), shape=(None, 3) + ) + points = freud.util._convert_array(np.atleast_2d(points), shape=(None, 3)) n_query_points = query_points.shape[0] n_points = points.shape[0] distances = np.empty(n_query_points, dtype=np.float32) - self._cpp_obj.computeDistances(query_points, n_query_points, points, - n_points, distances) + self._cpp_obj.computeDistances( + query_points, n_query_points, points, n_points, distances + ) return np.asarray(distances) def compute_all_distances(self, query_points, points): @@ -487,16 +490,17 @@ def compute_all_distances(self, query_points, points): Array of distances between query points and points. """ # noqa: E501 query_points = freud.util._convert_array( - np.atleast_2d(query_points), shape=(None, 3)) - points = freud.util._convert_array( - np.atleast_2d(points), shape=(None, 3)) + np.atleast_2d(query_points), shape=(None, 3) + ) + points = freud.util._convert_array(np.atleast_2d(points), shape=(None, 3)) n_query_points = query_points.shape[0] n_points = points.shape[0] distances = np.empty([n_query_points, n_points], dtype=np.float32) - self._cpp_obj.computeAllDistances(query_points, n_query_points, - points, n_points, distances) + self._cpp_obj.computeAllDistances( + query_points, n_query_points, points, n_points, distances + ) return np.asarray(distances) @@ -532,13 +536,11 @@ def contains(self, points): the box, and ``False`` corresponds to points outside the box. """ # noqa: E501 - points = freud.util._convert_array( - np.atleast_2d(points), shape=(None, 3)) + points = freud.util._convert_array(np.atleast_2d(points), shape=(None, 3)) n_points = points.shape[0] - contains_mask = freud.util._convert_array( - np.ones(n_points), dtype=bool) + contains_mask = freud.util._convert_array(np.ones(n_points), dtype=bool) self._cpp_obj.contains(points, n_points, contains_mask) @@ -562,9 +564,13 @@ def cubic(self): def periodic(self): r""":math:`\left(3, \right)` :class:`numpy.ndarray`: Get or set the periodicity of the box in each dimension.""" - return np.asarray([self._cpp_obj.getPeriodicX(), - self._cpp_obj.getPeriodicY(), - self._cpp_obj.getPeriodicZ()]) + return np.asarray( + [ + self._cpp_obj.getPeriodicX(), + self._cpp_obj.getPeriodicY(), + self._cpp_obj.getPeriodicZ(), + ] + ) @periodic.setter def periodic(self, periodic): @@ -616,13 +622,14 @@ def to_dict(self): dict: Box parameters """ return { - 'Lx': self.Lx, - 'Ly': self.Ly, - 'Lz': self.Lz, - 'xy': self.xy, - 'xz': self.xz, - 'yz': self.yz, - 'dimensions': self.dimensions} + "Lx": self.Lx, + "Ly": self.Ly, + "Lz": self.Lz, + "xy": self.xy, + "xz": self.xz, + "yz": self.yz, + "dimensions": self.dimensions, + } def to_matrix(self): r"""Returns the box matrix (3x3). @@ -638,9 +645,13 @@ def to_matrix(self): Returns: :math:`\left(3, 3\right)` :class:`numpy.ndarray`: Box matrix """ - return np.asarray([[self.Lx, self.xy * self.Ly, self.xz * self.Lz], - [0, self.Ly, self.yz * self.Lz], - [0, 0, self.Lz]]) + return np.asarray( + [ + [self.Lx, self.xy * self.Ly, self.xz * self.Lz], + [0, self.Ly, self.yz * self.Lz], + [0, 0, self.Lz], + ] + ) def to_box_lengths_and_angles(self): r"""Return the box lengths and angles. @@ -654,37 +665,45 @@ def to_box_lengths_and_angles(self): (self.xy * self.xz + self.yz) / (np.sqrt(1 + self.xy**2) * np.sqrt(1 + self.xz**2 + self.yz**2)) ) - beta = np.arccos(self.xz/np.sqrt(1+self.xz**2+self.yz**2)) - gamma = np.arccos(self.xy/np.sqrt(1+self.xy**2)) + beta = np.arccos(self.xz / np.sqrt(1 + self.xz**2 + self.yz**2)) + gamma = np.arccos(self.xy / np.sqrt(1 + self.xy**2)) L1 = self.Lx - a2 = [self.Ly*self.xy, self.Ly, 0] - a3 = [self.Lz*self.xz, self.Lz*self.yz, self.Lz] + a2 = [self.Ly * self.xy, self.Ly, 0] + a3 = [self.Lz * self.xz, self.Lz * self.yz, self.Lz] L2 = np.linalg.norm(a2) L3 = np.linalg.norm(a3) return (L1, L2, L3, alpha, beta, gamma) def __repr__(self): - return ("freud.box.{cls}(Lx={Lx}, Ly={Ly}, Lz={Lz}, " - "xy={xy}, xz={xz}, yz={yz}, " - "is2D={is2D})").format(cls=type(self).__name__, - Lx=self.Lx, - Ly=self.Ly, - Lz=self.Lz, - xy=self.xy, - xz=self.xz, - yz=self.yz, - is2D=self.is2D) + return ( + "freud.box.{cls}(Lx={Lx}, Ly={Ly}, Lz={Lz}, " + "xy={xy}, xz={xz}, yz={yz}, " + "is2D={is2D})" + ).format( + cls=type(self).__name__, + Lx=self.Lx, + Ly=self.Ly, + Lz=self.Lz, + xy=self.xy, + xz=self.xz, + yz=self.yz, + is2D=self.is2D, + ) def __str__(self): return repr(self) def __mul__(self, scale): if scale > 0: - return self.__class__(Lx=self.Lx*scale, - Ly=self.Ly*scale, - Lz=self.Lz*scale, - xy=self.xy, xz=self.xz, yz=self.yz, - is2D=self.is2D) + return self.__class__( + Lx=self.Lx * scale, + Ly=self.Ly * scale, + Lz=self.Lz * scale, + xy=self.xy, + xz=self.xz, + yz=self.yz, + is2D=self.is2D, + ) else: raise ValueError("Box can only be multiplied by positive values.") @@ -712,8 +731,10 @@ def plot(self, title=None, ax=None, image=[0, 0, 0], *args, **kwargs): :meth:`matplotlib.axes.Axes.plot`. """ import freud.plot - return freud.plot.box_plot(self, title=title, ax=ax, image=image, - *args, **kwargs) + + return freud.plot.box_plot( + self, title=title, ax=ax, image=image, *args, **kwargs + ) @classmethod def from_box(cls, box, dimensions=None): @@ -758,52 +779,54 @@ def from_box(cls, box, dimensions=None): # Handles freud.box.Box and objects with attributes Lx = box.Lx Ly = box.Ly - Lz = getattr(box, 'Lz', 0) - xy = getattr(box, 'xy', 0) - xz = getattr(box, 'xz', 0) - yz = getattr(box, 'yz', 0) + Lz = getattr(box, "Lz", 0) + xy = getattr(box, "xy", 0) + xz = getattr(box, "xz", 0) + yz = getattr(box, "yz", 0) if dimensions is None: - dimensions = getattr(box, 'dimensions', None) - elif dimensions != getattr(box, 'dimensions', dimensions): + dimensions = getattr(box, "dimensions", None) + elif dimensions != getattr(box, "dimensions", dimensions): raise ValueError( "The provided dimensions argument conflicts with the " - "dimensions attribute of the provided box object.") + "dimensions attribute of the provided box object." + ) except AttributeError: try: # Handle dictionary-like - Lx = box['Lx'] - Ly = box['Ly'] - Lz = box.get('Lz', 0) - xy = box.get('xy', 0) - xz = box.get('xz', 0) - yz = box.get('yz', 0) + Lx = box["Lx"] + Ly = box["Ly"] + Lz = box.get("Lz", 0) + xy = box.get("xy", 0) + xz = box.get("xz", 0) + yz = box.get("yz", 0) if dimensions is None: - dimensions = box.get('dimensions', None) + dimensions = box.get("dimensions", None) else: - if dimensions != box.get('dimensions', dimensions): + if dimensions != box.get("dimensions", dimensions): raise ValueError( "The provided dimensions argument conflicts with " "the dimensions attribute of the provided box " - "object.") + "object." + ) except (IndexError, KeyError, TypeError): if not len(box) in [2, 3, 6]: raise ValueError( "List-like objects must have length 2, 3, or 6 to be " - "converted to freud.box.Box.") + "converted to freud.box.Box." + ) # Handle list-like Lx = box[0] Ly = box[1] Lz = box[2] if len(box) > 2 else 0 xy, xz, yz = box[3:6] if len(box) == 6 else (0, 0, 0) except: # noqa - logger.debug('Supplied box cannot be converted to type ' - 'freud.box.Box.') + logger.debug("Supplied box cannot be converted to type " "freud.box.Box.") raise # Infer dimensions if not provided. if dimensions is None: dimensions = 2 if Lz == 0 else 3 - is2D = (dimensions == 2) + is2D = dimensions == 2 return cls(Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, is2D=is2D) @classmethod @@ -842,9 +865,8 @@ def from_matrix(cls, box_matrix, dimensions=None): xz = yz = 0 if dimensions is None: dimensions = 2 if Lz == 0 else 3 - is2D = (dimensions == 2) - return cls(Lx=Lx, Ly=Ly, Lz=Lz, - xy=xy, xz=xz, yz=yz, is2D=is2D) + is2D = dimensions == 2 + return cls(Lx=Lx, Ly=Ly, Lz=Lz, xy=xy, xz=xz, yz=yz, is2D=is2D) @classmethod def cube(cls, L=None): @@ -877,13 +899,19 @@ def square(cls, L=None): # named access to positional arguments, so we keep this to # recover the behavior if L is None: - raise TypeError("square() missing 1 required " - "positional argument: L") + raise TypeError("square() missing 1 required " "positional argument: L") return cls(Lx=L, Ly=L, Lz=0, xy=0, xz=0, yz=0, is2D=True) @classmethod def from_box_lengths_and_angles( - cls, L1, L2, L3, alpha, beta, gamma, dimensions=None, + cls, + L1, + L2, + L3, + alpha, + beta, + gamma, + dimensions=None, ): r"""Construct a box from lengths and angles (in radians). @@ -932,10 +960,14 @@ def from_box_lengths_and_angles( def BoxFromCPP(cppbox): - b = Box(cppbox.getLx(), cppbox.getLy(), cppbox.getLz(), - cppbox.getTiltFactorXY(), cppbox.getTiltFactorXZ(), - cppbox.getTiltFactorYZ(), cppbox.is2D()) - b.periodic = [cppbox.getPeriodicX(), - cppbox.getPeriodicY(), - cppbox.getPeriodicZ()] + b = Box( + cppbox.getLx(), + cppbox.getLy(), + cppbox.getLz(), + cppbox.getTiltFactorXY(), + cppbox.getTiltFactorXZ(), + cppbox.getTiltFactorYZ(), + cppbox.is2D(), + ) + b.periodic = [cppbox.getPeriodicX(), cppbox.getPeriodicY(), cppbox.getPeriodicZ()] return b diff --git a/freud/util.py b/freud/util.py index 11a1d2866..79e269c16 100644 --- a/freud/util.py +++ b/freud/util.py @@ -8,7 +8,7 @@ import freud.box -class _Compute(object): +class _Compute: r"""Parent class for all compute classes in freud. The primary purpose of this class is to prevent access of uncomputed @@ -40,7 +40,7 @@ def __getattribute__(self, attr): """Compute methods set a flag to indicate that quantities have been computed. Compute must be called before plotting.""" attribute = object.__getattribute__(self, attr) - if attr == 'compute': + if attr == "compute": # Set the attribute *after* computing. This enables # self._called_compute to be used in the compute method itself. compute = attribute @@ -50,11 +50,13 @@ def compute_wrapper(*args, **kwargs): return_value = compute(*args, **kwargs) self._called_compute = True return return_value + return compute_wrapper - elif attr == 'plot': + elif attr == "plot": if not self._called_compute: raise AttributeError( - "The compute method must be called before calling plot.") + "The compute method must be called before calling plot." + ) return attribute @staticmethod @@ -72,17 +74,18 @@ def _computed_property(prop): @wraps(prop) def wrapper(self, *args, **kwargs): if not self._called_compute: - raise AttributeError( - "Property not computed. Call compute first.") + raise AttributeError("Property not computed. Call compute first.") return prop(self, *args, **kwargs) + return wrapper def __str__(self): return repr(self) -def _convert_array(array, shape=None, dtype=np.float32, requirements=("C", ), - allow_copy=True): +def _convert_array( + array, shape=None, dtype=np.float32, requirements=("C",), allow_copy=True +): """Function which takes a given array, checks the dimensions and shape, and converts to a supplied dtype. @@ -113,21 +116,32 @@ def _convert_array(array, shape=None, dtype=np.float32, requirements=("C", ), return_arr = np.require(array, dtype=dtype, requirements=requirements) if not allow_copy and return_arr is not array: - raise ValueError("The provided output array must have dtype " - f"{dtype}, and have the following array flags: " - f"{', '.join(requirements)}.") + raise ValueError( + "The provided output array must have dtype " + f"{dtype}, and have the following array flags: " + f"{', '.join(requirements)}." + ) if shape is not None: if return_arr.ndim != len(shape): - raise ValueError("array.ndim = {}; expected ndim = {}".format( - return_arr.ndim, len(shape))) + raise ValueError( + "array.ndim = {}; expected ndim = {}".format( + return_arr.ndim, len(shape) + ) + ) for i, s in enumerate(shape): if s is not None and return_arr.shape[i] != s: - shape_str = "(" + ", ".join(str(i) if i is not None - else "..." for i in shape) + ")" - raise ValueError('array.shape= {}; expected shape = {}'.format( - return_arr.shape, shape_str)) + shape_str = ( + "(" + + ", ".join(str(i) if i is not None else "..." for i in shape) + + ")" + ) + raise ValueError( + "array.shape= {}; expected shape = {}".format( + return_arr.shape, shape_str + ) + ) return return_arr @@ -153,6 +167,6 @@ def _convert_box(box, dimensions=None): raise if dimensions is not None and box.dimensions != dimensions: - raise ValueError("The box must be {}-dimensional.".format(dimensions)) + raise ValueError(f"The box must be {dimensions}-dimensional.") return box From a9134cb8dcfb8bc4f51f04e9a17a9382e7330ed2 Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Wed, 5 Jun 2024 13:48:00 -0400 Subject: [PATCH 07/25] changes for nanobind --- CMakeLists.txt | 19 +++++++++++++++-- cpp/box/Box.h | 48 +++++++++++++++++++++++------------------- cpp/box/CMakeLists.txt | 2 +- cpp/box/module-box.cc | 10 ++++----- freud/CMakeLists.txt | 12 +++++------ freud/box.py | 2 +- 6 files changed, 56 insertions(+), 37 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 881689296..79b6c4183 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -27,6 +27,16 @@ set(CMAKE_CXX_STANDARD_REQUIRED ON) # https://stackoverflow.com/questions/50600708/combining-cmake-object-libraries-with-shared-libraries set(CMAKE_POSITION_INDEPENDENT_CODE ON) add_subdirectory(CMake) + +# find python +if (CMAKE_VERSION VERSION_LESS 3.18) + set(DEV_MODULE Development) +else() + set(DEV_MODULE Development.Module) +endif() + +find_package(Python 3.8 COMPONENTS Interpreter ${DEV_MODULE} REQUIRED) + find_package_config_first(TBB) if(TBB_FOUND) @@ -36,8 +46,13 @@ if(TBB_FOUND) "[${TBB_LIBRARY}][${TBB_INCLUDE_DIR}]") endif() -# go find pybind11 -find_package(pybind11 CONFIG REQUIRED) +# go find nanobind +execute_process( + COMMAND ${Python_EXECUTABLE} "-m" "nanobind" "--cmake_dir" + OUTPUT_STRIP_TRAILING_WHITESPACE OUTPUT_VARIABLE NB_DIR) +list(APPEND CMAKE_PREFIX_PATH "${NB_DIR}") +message(STATUS ${CMAKE_PREFIX_PATH}) +find_package(nanobind CONFIG REQUIRED) # Fail fast if users have not cloned submodules. if(NOT WIN32) diff --git a/cpp/box/Box.h b/cpp/box/Box.h index 27de2a206..0e2b1e706 100644 --- a/cpp/box/Box.h +++ b/cpp/box/Box.h @@ -6,14 +6,14 @@ #include #include -#include #include #include #include "VectorMath.h" #include "utils.h" -#include +#include +namespace nb = nanobind; /*! \file Box.h \brief Represents simulation boxes and contains helpful wrapping functions. @@ -26,8 +26,8 @@ constexpr float TWO_PI = 2.0 * M_PI; namespace freud { namespace box { -template -using pybind11_array = pybind11::array_t; +template> +using nb_array = nb::ndarray; //! Stores box dimensions and provides common routines for wrapping vectors back into the box /*! Box stores a standard HOOMD simulation box that goes from -L/2 to L/2 in each dimension, allowing Lx, Ly, @@ -156,7 +156,7 @@ class Box } //! Get current stored inverse of L - std::array getLinv() const + std::vector getLinv() const { return {m_Linv.x, m_Linv.y, m_Linv.z}; } @@ -230,7 +230,8 @@ class Box * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void makeAbsolutePython(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const + void makeAbsolutePython(nb_array> vecs, unsigned int Nvecs, + nb_array> out) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); @@ -265,7 +266,8 @@ class Box * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void makeFractionalPython(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const + void makeFractionalPython(nb_array> vecs, unsigned int Nvecs, + nb_array> out) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); @@ -298,7 +300,8 @@ class Box * \param Nvecs Number of vectors \param res Array to save the images */ - void getImages(pybind11_array vecs, unsigned int Nvecs, pybind11::array_t res) const + void getImages(nb_array> vecs, unsigned int Nvecs, + nb_array> res) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (res.data()); @@ -343,7 +346,8 @@ class Box * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void wrapPython(pybind11_array vecs, unsigned int Nvecs, pybind11_array out) const + void wrapPython(nb_array> vecs, unsigned int Nvecs, + nb_array> out) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); @@ -361,8 +365,8 @@ class Box \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void unwrap(pybind11_array vecs, pybind11_array images, unsigned int Nvecs, - pybind11_array out) const + void unwrap(nb_array vecs, nb_array images, unsigned int Nvecs, + nb_array out) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* images_data = (vec3*) (images.data()); @@ -408,8 +412,8 @@ class Box / constants::TWO_PI)); } - std::array centerOfMassPython(pybind11_array vecs, size_t Nvecs, - pybind11_array masses) const + std::vector centerOfMassPython(nb_array vecs, size_t Nvecs, + nb_array> masses) const { vec3* vecs_data = (vec3*) (vecs.data()); float* masses_data = (float*) (masses.data()); @@ -422,7 +426,7 @@ class Box * \param Nvecs Number of vectors * \param masses Optional array of masses, of length Nvecs */ - void center(pybind11_array vecs, unsigned int Nvecs, pybind11_array masses) const + void center(nb_array vecs, unsigned int Nvecs, nb_array> masses) const { vec3* vecs_data = (vec3*) (vecs.data()); float* masses_data = (float*) (masses.data()); @@ -454,9 +458,9 @@ class Box \param distances Pointer to array of length n_query_points containing distances between each point and query_point (overwritten in place). */ - void computeDistances(pybind11_array query_points, const unsigned int n_query_points, - pybind11_array points, const unsigned int n_points, - pybind11_array distances) const + void computeDistances(nb_array query_points, const unsigned int n_query_points, + nb_array points, const unsigned int n_points, + nb_array> distances) const { vec3* query_points_data = (vec3*) (query_points.data()); vec3* points_data = (vec3*) (points.data()); @@ -482,9 +486,9 @@ class Box \param distances Pointer to array of length n_query_points*n_points containing distances between points and query_points (overwritten in place). */ - void computeAllDistances(pybind11_array query_points, const unsigned int n_query_points, - pybind11_array points, const unsigned int n_points, - pybind11_array distances) const + void computeAllDistances(nb_array query_points, const unsigned int n_query_points, + nb_array points, const unsigned int n_points, + nb_array> distances) const { vec3* query_points_data = (vec3*) (query_points.data()); vec3* points_data = (vec3*) (points.data()); @@ -508,8 +512,8 @@ class Box \param n_points The number of points. \param contains_mask Mask of points inside the box. */ - void contains(pybind11_array points, const unsigned int n_points, - pybind11_array contains_mask) const + void contains(nb_array points, const unsigned int n_points, + nb_array> contains_mask) const { vec3* points_data = (vec3*) (points.data()); bool* contains_mask_data = (bool*) (contains_mask.data()); diff --git a/cpp/box/CMakeLists.txt b/cpp/box/CMakeLists.txt index 73a1959ad..94126008a 100644 --- a/cpp/box/CMakeLists.txt +++ b/cpp/box/CMakeLists.txt @@ -1,5 +1,5 @@ # create the target and add the needed properties -pybind11_add_module(_box SHARED module-box.cc) +nanobind_add_module(_box SHARED module-box.cc) target_link_libraries(_box PUBLIC TBB::tbb) # this probably isn't needed for box, but may be needed by other modules # target_compile_definitions( box # Avoid deprecation warnings for unsupported diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc index 297759c38..7f19799e3 100644 --- a/cpp/box/module-box.cc +++ b/cpp/box/module-box.cc @@ -1,18 +1,18 @@ // Copyright (c) 2010-2023 The Regents of the University of Michigan // This file is from the freud project, released under the BSD 3-Clause License. -#include -#include +#include +#include #include "Box.h" using namespace freud::box; -PYBIND11_MODULE(_box, m) +NB_MODULE(_box, m) { - pybind11::class_>(m, "Box") + nanobind::class_(m, "Box") // constructors - .def(pybind11::init()) + .def(nanobind::init()) // getters and setters .def("getLx", &Box::getLx) .def("getLy", &Box::getLy) diff --git a/freud/CMakeLists.txt b/freud/CMakeLists.txt index 16c8dc557..1aecefedf 100644 --- a/freud/CMakeLists.txt +++ b/freud/CMakeLists.txt @@ -13,12 +13,12 @@ include_directories(${NumPy_INCLUDE_DIRS}) # versions concurrently so it shouldn't hide any real issues. For backwards # compatibility with older CMake, I'm using PythonInterp; when we drop support # for CMake < 3.12, we should switch to find_package(Python). -find_package(PythonInterp REQUIRED) -if(${PYTHON_VERSION_MAJOR} EQUAL 3 - AND ${PYTHON_VERSION_MINOR} EQUAL 8 - AND NOT WIN32) - add_compile_options("-Wno-deprecated-declarations") -endif() +#find_package(PythonInterp REQUIRED) +#if(${PYTHON_VERSION_MAJOR} EQUAL 3 +# AND ${PYTHON_VERSION_MINOR} EQUAL 8 +# AND NOT WIN32) +# add_compile_options("-Wno-deprecated-declarations") +#endif() set(python_files box.py util.py) # cluster.py density.py diffraction.py environment.py locality.py order.py diff --git a/freud/box.py b/freud/box.py index 0cf8b0ed6..0811470a1 100644 --- a/freud/box.py +++ b/freud/box.py @@ -398,7 +398,7 @@ def center_of_mass(self, vecs, masses=None): if masses is not None: masses = freud.util._convert_array(masses, shape=(len(vecs),)) else: - masses = np.ones(Np) + masses = np.ones(Np, dtype=np.float32) result = self._cpp_obj.centerOfMass(vecs, Np, masses) return np.asarray(result) From 3888b3f270f9e0f88b1b49308deb657b24856c62 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 5 Jun 2024 17:48:42 +0000 Subject: [PATCH 08/25] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- CMakeLists.txt | 10 +++++++--- cpp/box/Box.h | 13 ++++++------- freud/CMakeLists.txt | 9 +++------ 3 files changed, 16 insertions(+), 16 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 79b6c4183..30adf23b3 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -29,13 +29,16 @@ set(CMAKE_POSITION_INDEPENDENT_CODE ON) add_subdirectory(CMake) # find python -if (CMAKE_VERSION VERSION_LESS 3.18) +if(CMAKE_VERSION VERSION_LESS 3.18) set(DEV_MODULE Development) else() set(DEV_MODULE Development.Module) endif() -find_package(Python 3.8 COMPONENTS Interpreter ${DEV_MODULE} REQUIRED) +find_package( + Python 3.8 + COMPONENTS Interpreter ${DEV_MODULE} + REQUIRED) find_package_config_first(TBB) @@ -49,7 +52,8 @@ endif() # go find nanobind execute_process( COMMAND ${Python_EXECUTABLE} "-m" "nanobind" "--cmake_dir" - OUTPUT_STRIP_TRAILING_WHITESPACE OUTPUT_VARIABLE NB_DIR) + OUTPUT_STRIP_TRAILING_WHITESPACE + OUTPUT_VARIABLE NB_DIR) list(APPEND CMAKE_PREFIX_PATH "${NB_DIR}") message(STATUS ${CMAKE_PREFIX_PATH}) find_package(nanobind CONFIG REQUIRED) diff --git a/cpp/box/Box.h b/cpp/box/Box.h index 0e2b1e706..c8001e2f6 100644 --- a/cpp/box/Box.h +++ b/cpp/box/Box.h @@ -231,7 +231,7 @@ class Box * \param out The array in which to place the wrapped vectors. */ void makeAbsolutePython(nb_array> vecs, unsigned int Nvecs, - nb_array> out) const + nb_array> out) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); @@ -267,7 +267,7 @@ class Box * \param out The array in which to place the wrapped vectors. */ void makeFractionalPython(nb_array> vecs, unsigned int Nvecs, - nb_array> out) const + nb_array> out) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); @@ -301,7 +301,7 @@ class Box \param res Array to save the images */ void getImages(nb_array> vecs, unsigned int Nvecs, - nb_array> res) const + nb_array> res) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (res.data()); @@ -347,7 +347,7 @@ class Box * \param out The array in which to place the wrapped vectors. */ void wrapPython(nb_array> vecs, unsigned int Nvecs, - nb_array> out) const + nb_array> out) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); @@ -365,8 +365,7 @@ class Box \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void unwrap(nb_array vecs, nb_array images, unsigned int Nvecs, - nb_array out) const + void unwrap(nb_array vecs, nb_array images, unsigned int Nvecs, nb_array out) const { vec3* vecs_data = (vec3*) (vecs.data()); vec3* images_data = (vec3*) (images.data()); @@ -413,7 +412,7 @@ class Box } std::vector centerOfMassPython(nb_array vecs, size_t Nvecs, - nb_array> masses) const + nb_array> masses) const { vec3* vecs_data = (vec3*) (vecs.data()); float* masses_data = (float*) (masses.data()); diff --git a/freud/CMakeLists.txt b/freud/CMakeLists.txt index 1aecefedf..6678afa12 100644 --- a/freud/CMakeLists.txt +++ b/freud/CMakeLists.txt @@ -13,12 +13,9 @@ include_directories(${NumPy_INCLUDE_DIRS}) # versions concurrently so it shouldn't hide any real issues. For backwards # compatibility with older CMake, I'm using PythonInterp; when we drop support # for CMake < 3.12, we should switch to find_package(Python). -#find_package(PythonInterp REQUIRED) -#if(${PYTHON_VERSION_MAJOR} EQUAL 3 -# AND ${PYTHON_VERSION_MINOR} EQUAL 8 -# AND NOT WIN32) -# add_compile_options("-Wno-deprecated-declarations") -#endif() +# find_package(PythonInterp REQUIRED) if(${PYTHON_VERSION_MAJOR} EQUAL 3 AND +# ${PYTHON_VERSION_MINOR} EQUAL 8 AND NOT WIN32) +# add_compile_options("-Wno-deprecated-declarations") endif() set(python_files box.py util.py) # cluster.py density.py diffraction.py environment.py locality.py order.py From ee033d3fcc2d74b476ce84655702c77905da6ac3 Mon Sep 17 00:00:00 2001 From: Tommy Waltmann <53307607+tommy-waltmann@users.noreply.github.com> Date: Thu, 6 Jun 2024 11:07:27 -0400 Subject: [PATCH 09/25] Update cpp/box/Box.h Co-authored-by: Jen Bradley <55467578+janbridley@users.noreply.github.com> --- cpp/box/Box.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/cpp/box/Box.h b/cpp/box/Box.h index c8001e2f6..61d4e5c09 100644 --- a/cpp/box/Box.h +++ b/cpp/box/Box.h @@ -228,7 +228,7 @@ class Box /*! \param vecs Vectors of fractional coordinates between 0 and 1 within * parallelepipedal box * \param Nvecs Number of vectors - * \param out The array in which to place the wrapped vectors. + * \param out_data The array in which to place the wrapped vectors. */ void makeAbsolutePython(nb_array> vecs, unsigned int Nvecs, nb_array> out) const From a55c56c7ce34f407d087faa4c262fcdb1cbe053c Mon Sep 17 00:00:00 2001 From: Tommy Waltmann <53307607+tommy-waltmann@users.noreply.github.com> Date: Tue, 11 Jun 2024 15:56:15 -0400 Subject: [PATCH 10/25] Update CMakeLists.txt Co-authored-by: Joshua A. Anderson --- CMakeLists.txt | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 30adf23b3..2471e65fa 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -53,9 +53,7 @@ endif() execute_process( COMMAND ${Python_EXECUTABLE} "-m" "nanobind" "--cmake_dir" OUTPUT_STRIP_TRAILING_WHITESPACE - OUTPUT_VARIABLE NB_DIR) -list(APPEND CMAKE_PREFIX_PATH "${NB_DIR}") -message(STATUS ${CMAKE_PREFIX_PATH}) + OUTPUT_VARIABLE nanobind_ROOT) find_package(nanobind CONFIG REQUIRED) # Fail fast if users have not cloned submodules. From 5fa295c2e3b48615bc3c573b021622f0e2c83d2a Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Tue, 11 Jun 2024 16:02:21 -0400 Subject: [PATCH 11/25] streamline cmake code in python directory --- freud/CMakeLists.txt | 24 ------------------------ 1 file changed, 24 deletions(-) diff --git a/freud/CMakeLists.txt b/freud/CMakeLists.txt index 6678afa12..cdad03428 100644 --- a/freud/CMakeLists.txt +++ b/freud/CMakeLists.txt @@ -1,30 +1,6 @@ -# Need to figure out coverage before CYTHON_FLAGS are set. -option(COVERAGE "Enable coverage" OFF) - -find_package(PythonLibs) -find_package(PythonExtensions REQUIRED) -find_package(NumPy REQUIRED) - -include_directories(${NumPy_INCLUDE_DIRS}) - -# Avoid Cython/Python3.8 minor incompatibility warnings, see -# https://github.com/cython/cython/issues/3474. Note that this option is a bit -# expansive, but it's a temporary fix and we'll be testing on other Python -# versions concurrently so it shouldn't hide any real issues. For backwards -# compatibility with older CMake, I'm using PythonInterp; when we drop support -# for CMake < 3.12, we should switch to find_package(Python). -# find_package(PythonInterp REQUIRED) if(${PYTHON_VERSION_MAJOR} EQUAL 3 AND -# ${PYTHON_VERSION_MINOR} EQUAL 8 AND NOT WIN32) -# add_compile_options("-Wno-deprecated-declarations") endif() - set(python_files box.py util.py) # cluster.py density.py diffraction.py environment.py locality.py order.py # parallel.py pmft.py) install(FILES ${python_files} DESTINATION freud) -# THIS WILL NEED TO BE MOVED TO THE BUILD OF THE ORDER MODULE The SolidLiquid -# class has an instance of cluster::Cluster as a member, so including the header -# requires the Cluster.h header. Would prefer to inherit this information from -# the _order library, but that's not possible since we're linking to libfreud. -# target_include_directories(order PUBLIC ${PROJECT_SOURCE_DIR}/cpp/cluster) From fd6729e6c1ae1c487209ce068ad9171cfa5ae92f Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Tue, 11 Jun 2024 16:02:36 -0400 Subject: [PATCH 12/25] change cmake version requirements --- CMakeLists.txt | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 2471e65fa..d80f37d33 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,6 +1,4 @@ -# CMake 3.12.0 is the oldest version supported by the way freud links TBB to -# object libraries like _cluster. This is also the oldest version tested in CI. -cmake_minimum_required(VERSION 3.12.0) +cmake_minimum_required(VERSION 3.15...3.30) project(freud) From 024b5f45b3a58b303002a22d883e9a3aa8799add Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 11 Jun 2024 20:02:54 +0000 Subject: [PATCH 13/25] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- freud/CMakeLists.txt | 1 - 1 file changed, 1 deletion(-) diff --git a/freud/CMakeLists.txt b/freud/CMakeLists.txt index cdad03428..29ea88c59 100644 --- a/freud/CMakeLists.txt +++ b/freud/CMakeLists.txt @@ -3,4 +3,3 @@ set(python_files box.py util.py) # parallel.py pmft.py) install(FILES ${python_files} DESTINATION freud) - From 73a05b9531cc802c8a67a6e755786558caf2a1f5 Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Tue, 11 Jun 2024 17:07:17 -0400 Subject: [PATCH 14/25] remove unnecessary passing of array sizes and clarify function structure --- cpp/box/Box.h | 219 ++++++++++++++++++++++++++---------------- cpp/box/module-box.cc | 12 +-- freud/box.py | 28 ++---- 3 files changed, 152 insertions(+), 107 deletions(-) diff --git a/cpp/box/Box.h b/cpp/box/Box.h index 61d4e5c09..acc69e310 100644 --- a/cpp/box/Box.h +++ b/cpp/box/Box.h @@ -224,23 +224,28 @@ class Box return v; } + void makeAbsolute(const vec3* vecs, const unsigned int Nvecs, vec3* out) const + { + util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { + for (size_t i = begin; i < end; ++i) + { + out[i] = makeAbsolute(vecs[i]); + } + }); + } + //! Convert fractional coordinates into absolute coordinates in place /*! \param vecs Vectors of fractional coordinates between 0 and 1 within * parallelepipedal box - * \param Nvecs Number of vectors * \param out_data The array in which to place the wrapped vectors. */ - void makeAbsolutePython(nb_array> vecs, unsigned int Nvecs, + void makeAbsolutePython(nb_array> vecs, nb_array> out) const { + unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); - util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { - for (size_t i = begin; i < end; ++i) - { - out_data[i] = makeAbsolute(vecs_data[i]); - } - }); + makeAbsolute(vecs_data, Nvecs, out_data); } //! Convert a point's coordinate from absolute to fractional box coordinates. @@ -261,22 +266,28 @@ class Box return delta; } + void makeFractional(const vec3* vecs, unsigned int Nvecs, vec3* out) const + { + util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { + for (size_t i = begin; i < end; ++i) + { + out[i] = makeFractional(vecs[i]); + } + }); + } + //! Convert point coordinates from absolute to fractional box coordinates. /*! \param vecs Vectors to convert * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void makeFractionalPython(nb_array> vecs, unsigned int Nvecs, + void makeFractionalPython(nb_array> vecs, nb_array> out) const { + unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); - util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { - for (size_t i = begin; i < end; ++i) - { - out_data[i] = makeFractional(vecs_data[i]); - } - }); + makeFractional(vecs_data, Nvecs, out_data); } //! Get periodic image of a vector. @@ -295,22 +306,28 @@ class Box image.z = (int) ((f.z >= float(0.0)) ? f.z + float(0.5) : f.z - float(0.5)); } + void getImages(const vec3* vecs, unsigned int Nvecs, vec3* images) const + { + util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { + for (size_t i = begin; i < end; ++i) + { + getImage(vecs[i], images[i]); + } + }); + } + //! Get the periodic image vectors belongs to /*! \param vecs The vectors to check * \param Nvecs Number of vectors \param res Array to save the images */ - void getImages(nb_array> vecs, unsigned int Nvecs, - nb_array> res) const + void getImagesPython(nb_array> vecs, + nb_array> images) const { + const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); - vec3* out_data = (vec3*) (res.data()); - util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { - for (size_t i = begin; i < end; ++i) - { - getImage(vecs_data[i], out_data[i]); - } - }); + vec3* images_data = (vec3*) (images.data()); + getImages(vecs_data, Nvecs, images_data); } //! Wrap a vector back into the box @@ -341,20 +358,41 @@ class Box return makeAbsolute(v_frac); } + void wrap(const vec3* vecs, unsigned int Nvecs, vec3* out) const + { + util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { + for (size_t i = begin; i < end; ++i) + { + out[i] = wrap(vecs[i]); + } + }); + } + //! Wrap vectors back into the box in place /*! \param vecs Vectors to wrap, updated to the minimum image obeying the periodic settings * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void wrapPython(nb_array> vecs, unsigned int Nvecs, + void wrapPython(nb_array> vecs, nb_array> out) const { + const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); + wrap(vecs_data, Nvecs, out_data); + } + + void unwrap(const vec3* vecs, const vec3* images, unsigned int Nvecs, vec3* out) const + { util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) { - out_data[i] = wrap(vecs_data[i]); + out[i] = vecs[i] + getLatticeVector(0) * float(images[i].x) + + getLatticeVector(1) * float(images[i].y); + if (!m_2d) + { + out[i] += getLatticeVector(2) * float(images[i].z); + } } }); } @@ -365,22 +403,13 @@ class Box \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void unwrap(nb_array vecs, nb_array images, unsigned int Nvecs, nb_array out) const + void unwrapPython(nb_array vecs, nb_array images, nb_array out) const { + const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); vec3* images_data = (vec3*) (images.data()); vec3* out_data = (vec3*) (out.data()); - util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { - for (size_t i = begin; i < end; ++i) - { - out_data[i] = vecs_data[i] + getLatticeVector(0) * float(images_data[i].x) - + getLatticeVector(1) * float(images_data[i].y); - if (!m_2d) - { - out_data[i] += getLatticeVector(2) * float(images_data[i].z); - } - } - }); + unwrap(vecs_data, images_data, Nvecs, out_data); } //! Compute center of mass for vectors @@ -411,32 +440,38 @@ class Box / constants::TWO_PI)); } - std::vector centerOfMassPython(nb_array vecs, size_t Nvecs, + std::vector centerOfMassPython(nb_array vecs, nb_array> masses) const { + const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); float* masses_data = (float*) (masses.data()); auto com = centerOfMass(vecs_data, Nvecs, masses_data); return {com.x, com.y, com.z}; } + void center(vec3* vecs, unsigned int Nvecs, const float* masses) const + { + vec3 com(centerOfMass(vecs, Nvecs, masses)); + util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { + for (size_t i = begin; i < end; ++i) + { + vecs[i] = wrap(vecs[i] - com); + } + }); + } + //! Subtract center of mass from vectors /*! \param vecs Vectors to center * \param Nvecs Number of vectors * \param masses Optional array of masses, of length Nvecs */ - void center(nb_array vecs, unsigned int Nvecs, nb_array> masses) const + void centerPython(nb_array vecs, nb_array> masses) const { + const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); float* masses_data = (float*) (masses.data()); - - vec3 com(centerOfMass(vecs_data, Nvecs, masses_data)); - util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { - for (size_t i = begin; i < end; ++i) - { - vecs_data[i] = wrap(vecs_data[i] - com); - } - }); + center(vecs_data, Nvecs, masses_data); } //! Calculate distance between two points using boundary conditions @@ -449,6 +484,18 @@ class Box return std::sqrt(dot(r_ij, r_ij)); } + void computeDistances(const vec3* query_points, const unsigned int n_query_points, + const vec3* points, const unsigned int n_points, + float* distances) const + { + util::forLoopWrapper(0, n_query_points, [&](size_t begin, size_t end) { + for (size_t i = begin; i < end; ++i) + { + distances[i] = computeDistance(query_points[i], points[i]); + } + }); + } + //! Calculate distances between a set of query points and points. /*! \param query_points Query point positions. \param n_query_points The number of query points. @@ -457,68 +504,61 @@ class Box \param distances Pointer to array of length n_query_points containing distances between each point and query_point (overwritten in place). */ - void computeDistances(nb_array query_points, const unsigned int n_query_points, - nb_array points, const unsigned int n_points, - nb_array> distances) const + void computeDistancesPython(nb_array query_points, nb_array points, + nb_array> distances) const { + const unsigned int n_query_points = query_points.shape(0); vec3* query_points_data = (vec3*) (query_points.data()); + const unsigned int n_points = points.shape(0); vec3* points_data = (vec3*) (points.data()); float* distances_data = (float*) (distances.data()); - if (n_query_points != n_points) { throw std::invalid_argument("The number of query points and points must match."); } - util::forLoopWrapper(0, n_query_points, [&](size_t begin, size_t end) { - for (size_t i = begin; i < end; ++i) - { - distances_data[i] = computeDistance(query_points_data[i], points_data[i]); - } - }); + computeDistances(query_points_data, n_query_points, points_data, n_points, distances_data); } - //! Calculate all pairwise distances between a set of query points and points. - /*! \param query_points Query point positions. - \param n_query_points The number of query points. - \param points Point positions. - \param n_points The number of points. - \param distances Pointer to array of length n_query_points*n_points containing distances between - points and query_points (overwritten in place). - */ - void computeAllDistances(nb_array query_points, const unsigned int n_query_points, - nb_array points, const unsigned int n_points, - nb_array> distances) const + void computeAllDistances(const vec3* query_points, const unsigned int n_query_points, + const vec3* points, const unsigned int n_points, + float* distances) const { - vec3* query_points_data = (vec3*) (query_points.data()); - vec3* points_data = (vec3*) (points.data()); - float* distances_data = (float*) (distances.data()); - util::forLoopWrapper2D(0, n_query_points, 0, n_points, [&](size_t begin_n, size_t end_n, size_t begin_m, size_t end_m) { for (size_t i = begin_n; i < end_n; ++i) { for (size_t j = begin_m; j < end_m; ++j) { - distances_data[i * n_points + j] - = computeDistance(query_points_data[i], points_data[j]); + distances[i * n_points + j] + = computeDistance(query_points[i], points[j]); } } }); } - - //! Get mask of points that fit inside the box. - /*! \param points Point positions. + //! Calculate all pairwise distances between a set of query points and points. + /*! \param query_points Query point positions. + \param n_query_points The number of query points. + \param points Point positions. \param n_points The number of points. - \param contains_mask Mask of points inside the box. + \param distances Pointer to array of length n_query_points*n_points containing distances between + points and query_points (overwritten in place). */ - void contains(nb_array points, const unsigned int n_points, - nb_array> contains_mask) const + void computeAllDistancesPython(nb_array query_points, nb_array points, + nb_array> distances) const { + const unsigned int n_query_points = query_points.shape(0); + vec3* query_points_data = (vec3*) (query_points.data()); + const unsigned int n_points = points.shape(0); vec3* points_data = (vec3*) (points.data()); - bool* contains_mask_data = (bool*) (contains_mask.data()); + float* distances_data = (float*) (distances.data()); + computeAllDistances(query_points_data, n_query_points, points_data, n_points, distances_data); + } + void contains(vec3* points, const unsigned int n_points, + bool* contains_mask) const + { util::forLoopWrapper(0, n_points, [&](size_t begin, size_t end) { - std::transform(&points_data[begin], &points_data[end], &contains_mask_data[begin], + std::transform(&points[begin], &points[end], &contains_mask[begin], [this](const vec3& point) -> bool { vec3 image(0, 0, 0); getImage(point, image); @@ -527,6 +567,19 @@ class Box }); } + //! Get mask of points that fit inside the box. + /*! \param points Point positions. + \param n_points The number of points. + \param contains_mask Mask of points inside the box. + */ + void containsPython(nb_array points, nb_array> contains_mask) const + { + const unsigned int n_points = points.shape(0); + vec3* points_data = (vec3*) (points.data()); + bool* contains_mask_data = (bool*) (contains_mask.data()); + contains(points_data, n_points, contains_mask_data); + } + //! Get the shortest distance between opposite boundary planes of the box /*! The distance between two planes of the lattice is 2 Pi/|b_i|, where * b_1 is the reciprocal lattice vector of the Bravais lattice normal to diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc index 7f19799e3..0f23f1c06 100644 --- a/cpp/box/module-box.cc +++ b/cpp/box/module-box.cc @@ -35,15 +35,15 @@ NB_MODULE(_box, m) .def("is2D", &Box::is2D) .def("set2D", &Box::set2D) .def("getVolume", &Box::getVolume) - .def("center", &Box::center) + .def("center", &Box::centerPython) .def("centerOfMass", &Box::centerOfMassPython) // other stuff .def("makeAbsolute", &Box::makeAbsolutePython) .def("makeFractional", &Box::makeFractionalPython) .def("wrap", &Box::wrapPython) - .def("unwrap", &Box::unwrap) - .def("getImages", &Box::getImages) - .def("computeDistances", &Box::computeDistances) - .def("computeAllDistances", &Box::computeAllDistances) - .def("contains", &Box::contains); + .def("unwrap", &Box::unwrapPython) + .def("getImages", &Box::getImagesPython) + .def("computeDistances", &Box::computeDistancesPython) + .def("computeAllDistances", &Box::computeAllDistancesPython) + .def("contains", &Box::containsPython); } diff --git a/freud/box.py b/freud/box.py index 0811470a1..9b310e23e 100644 --- a/freud/box.py +++ b/freud/box.py @@ -208,9 +208,7 @@ def make_absolute(self, fractional_coordinates, out=None): fractions = freud.util._convert_array(fractions, shape=(None, 3)) out = freud.util._convert_array(out, shape=fractions.shape, allow_copy=False) - Np = fractions.shape[0] - - self._cpp_obj.makeAbsolute(fractions, Np, out) + self._cpp_obj.makeAbsolute(fractions, out) return np.squeeze(out) if flatten else out @@ -236,9 +234,7 @@ def make_fractional(self, absolute_coordinates, out=None): vecs = freud.util._convert_array(vecs, shape=(None, 3)) out = freud.util._convert_array(out, shape=vecs.shape, allow_copy=False) - Np = vecs.shape[0] - - self._cpp_obj.makeFractional(vecs, Np, out) + self._cpp_obj.makeFractional(vecs, out) return np.squeeze(out) if flatten else out @@ -259,8 +255,7 @@ def get_images(self, vecs): vecs = freud.util._convert_array(vecs, shape=(None, 3)) images = np.zeros(vecs.shape, dtype=np.int32) - Np = vecs.shape[0] - self._cpp_obj.getImages(vecs, Np, images) + self._cpp_obj.getImages(vecs, images) return np.squeeze(images) if flatten else images @@ -322,8 +317,7 @@ def wrap(self, vecs, out=None): vecs = freud.util._convert_array(vecs, shape=(None, 3)) out = freud.util._convert_array(out, shape=vecs.shape, allow_copy=False) - Np = vecs.shape[0] - self._cpp_obj.wrap(vecs, Np, out) + self._cpp_obj.wrap(vecs, out) return np.squeeze(out) if flatten else out @@ -358,8 +352,7 @@ def unwrap(self, vecs, imgs, out=None): imgs = freud.util._convert_array(imgs, shape=vecs.shape, dtype=np.int32) out = freud.util._convert_array(out, shape=vecs.shape, allow_copy=False) - Np = vecs.shape[0] - self._cpp_obj.unwrap(vecs, imgs, Np, out) + self._cpp_obj.unwrap(vecs, imgs, out) return np.squeeze(out) if flatten else out @@ -400,7 +393,7 @@ def center_of_mass(self, vecs, masses=None): else: masses = np.ones(Np, dtype=np.float32) - result = self._cpp_obj.centerOfMass(vecs, Np, masses) + result = self._cpp_obj.centerOfMass(vecs, masses) return np.asarray(result) def center(self, vecs, masses=None): @@ -439,7 +432,7 @@ def center(self, vecs, masses=None): else: masses = np.ones(Np) - self._cpp_obj.center(vecs, Np, masses) + self._cpp_obj.center(vecs, masses) return vecs def compute_distances(self, query_points, points): @@ -465,11 +458,10 @@ def compute_distances(self, query_points, points): points = freud.util._convert_array(np.atleast_2d(points), shape=(None, 3)) n_query_points = query_points.shape[0] - n_points = points.shape[0] distances = np.empty(n_query_points, dtype=np.float32) self._cpp_obj.computeDistances( - query_points, n_query_points, points, n_points, distances + query_points, points, distances ) return np.asarray(distances) @@ -499,7 +491,7 @@ def compute_all_distances(self, query_points, points): distances = np.empty([n_query_points, n_points], dtype=np.float32) self._cpp_obj.computeAllDistances( - query_points, n_query_points, points, n_points, distances + query_points, points, distances ) return np.asarray(distances) @@ -542,7 +534,7 @@ def contains(self, points): contains_mask = freud.util._convert_array(np.ones(n_points), dtype=bool) - self._cpp_obj.contains(points, n_points, contains_mask) + self._cpp_obj.contains(points, contains_mask) return np.array(contains_mask).astype(bool) From 72f42d2b747dbacd2b11712ec302549fc4cbb74a Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 11 Jun 2024 21:07:33 +0000 Subject: [PATCH 15/25] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- cpp/box/Box.h | 45 +++++++++++++++++++-------------------------- freud/box.py | 8 ++------ 2 files changed, 21 insertions(+), 32 deletions(-) diff --git a/cpp/box/Box.h b/cpp/box/Box.h index acc69e310..1bdfc1e75 100644 --- a/cpp/box/Box.h +++ b/cpp/box/Box.h @@ -321,8 +321,7 @@ class Box * \param Nvecs Number of vectors \param res Array to save the images */ - void getImagesPython(nb_array> vecs, - nb_array> images) const + void getImagesPython(nb_array> vecs, nb_array> images) const { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); @@ -373,8 +372,7 @@ class Box * \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. */ - void wrapPython(nb_array> vecs, - nb_array> out) const + void wrapPython(nb_array> vecs, nb_array> out) const { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); @@ -440,8 +438,7 @@ class Box / constants::TWO_PI)); } - std::vector centerOfMassPython(nb_array vecs, - nb_array> masses) const + std::vector centerOfMassPython(nb_array vecs, nb_array> masses) const { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); @@ -485,8 +482,7 @@ class Box } void computeDistances(const vec3* query_points, const unsigned int n_query_points, - const vec3* points, const unsigned int n_points, - float* distances) const + const vec3* points, const unsigned int n_points, float* distances) const { util::forLoopWrapper(0, n_query_points, [&](size_t begin, size_t end) { for (size_t i = begin; i < end; ++i) @@ -505,7 +501,7 @@ class Box query_point (overwritten in place). */ void computeDistancesPython(nb_array query_points, nb_array points, - nb_array> distances) const + nb_array> distances) const { const unsigned int n_query_points = query_points.shape(0); vec3* query_points_data = (vec3*) (query_points.data()); @@ -520,20 +516,18 @@ class Box } void computeAllDistances(const vec3* query_points, const unsigned int n_query_points, - const vec3* points, const unsigned int n_points, - float* distances) const - { - util::forLoopWrapper2D(0, n_query_points, 0, n_points, - [&](size_t begin_n, size_t end_n, size_t begin_m, size_t end_m) { - for (size_t i = begin_n; i < end_n; ++i) - { - for (size_t j = begin_m; j < end_m; ++j) - { - distances[i * n_points + j] - = computeDistance(query_points[i], points[j]); - } - } - }); + const vec3* points, const unsigned int n_points, float* distances) const + { + util::forLoopWrapper2D( + 0, n_query_points, 0, n_points, [&](size_t begin_n, size_t end_n, size_t begin_m, size_t end_m) { + for (size_t i = begin_n; i < end_n; ++i) + { + for (size_t j = begin_m; j < end_m; ++j) + { + distances[i * n_points + j] = computeDistance(query_points[i], points[j]); + } + } + }); } //! Calculate all pairwise distances between a set of query points and points. /*! \param query_points Query point positions. @@ -544,7 +538,7 @@ class Box points and query_points (overwritten in place). */ void computeAllDistancesPython(nb_array query_points, nb_array points, - nb_array> distances) const + nb_array> distances) const { const unsigned int n_query_points = query_points.shape(0); vec3* query_points_data = (vec3*) (query_points.data()); @@ -554,8 +548,7 @@ class Box computeAllDistances(query_points_data, n_query_points, points_data, n_points, distances_data); } - void contains(vec3* points, const unsigned int n_points, - bool* contains_mask) const + void contains(vec3* points, const unsigned int n_points, bool* contains_mask) const { util::forLoopWrapper(0, n_points, [&](size_t begin, size_t end) { std::transform(&points[begin], &points[end], &contains_mask[begin], diff --git a/freud/box.py b/freud/box.py index 9b310e23e..5042f12a1 100644 --- a/freud/box.py +++ b/freud/box.py @@ -460,9 +460,7 @@ def compute_distances(self, query_points, points): n_query_points = query_points.shape[0] distances = np.empty(n_query_points, dtype=np.float32) - self._cpp_obj.computeDistances( - query_points, points, distances - ) + self._cpp_obj.computeDistances(query_points, points, distances) return np.asarray(distances) def compute_all_distances(self, query_points, points): @@ -490,9 +488,7 @@ def compute_all_distances(self, query_points, points): n_points = points.shape[0] distances = np.empty([n_query_points, n_points], dtype=np.float32) - self._cpp_obj.computeAllDistances( - query_points, points, distances - ) + self._cpp_obj.computeAllDistances(query_points, points, distances) return np.asarray(distances) From 825ba64509b104688ebb06829952e8d7d438fc98 Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Wed, 12 Jun 2024 14:40:51 -0400 Subject: [PATCH 16/25] remove redundant stuff from python --- freud/box.py | 26 ++++++++++---------------- 1 file changed, 10 insertions(+), 16 deletions(-) diff --git a/freud/box.py b/freud/box.py index 5042f12a1..aab1dfd6b 100644 --- a/freud/box.py +++ b/freud/box.py @@ -202,7 +202,7 @@ def make_absolute(self, fractional_coordinates, out=None): Absolute coordinate vector(s). If ``out`` is provided, a reference to it is returned. """ # noqa: E501 - fractions = np.asarray(fractional_coordinates).copy() + fractions = np.asarray(fractional_coordinates) flatten = fractions.ndim == 1 fractions = np.atleast_2d(fractions) fractions = freud.util._convert_array(fractions, shape=(None, 3)) @@ -228,7 +228,7 @@ def make_fractional(self, absolute_coordinates, out=None): Fractional coordinate vector(s). If ``out`` is provided, a reference to it is returned. """ # noqa: E501 - vecs = np.asarray(absolute_coordinates).copy() + vecs = np.asarray(absolute_coordinates) flatten = vecs.ndim == 1 vecs = np.atleast_2d(vecs) vecs = freud.util._convert_array(vecs, shape=(None, 3)) @@ -253,8 +253,8 @@ def get_images(self, vecs): flatten = vecs.ndim == 1 vecs = np.atleast_2d(vecs) vecs = freud.util._convert_array(vecs, shape=(None, 3)) - images = np.zeros(vecs.shape, dtype=np.int32) + self._cpp_obj.getImages(vecs, images) return np.squeeze(images) if flatten else images @@ -386,12 +386,11 @@ def center_of_mass(self, vecs, masses=None): Center of mass. """ # noqa: E501 vecs = freud.util._convert_array(vecs, shape=(None, 3)) - Np = vecs.shape[0] if masses is not None: masses = freud.util._convert_array(masses, shape=(len(vecs),)) else: - masses = np.ones(Np, dtype=np.float32) + masses = np.ones(vecs.shape[0], dtype=np.float32) result = self._cpp_obj.centerOfMass(vecs, masses) return np.asarray(result) @@ -425,12 +424,11 @@ def center(self, vecs, masses=None): Vectors with center of mass subtracted. """ # noqa: E501 vecs = freud.util._convert_array(vecs, shape=(None, 3)).copy() - Np = vecs.shape[0] if masses is not None: masses = freud.util._convert_array(masses, shape=(len(vecs),)) else: - masses = np.ones(Np) + masses = np.ones(vecs.shape[0], dtype=np.float32) self._cpp_obj.center(vecs, masses) return vecs @@ -457,11 +455,10 @@ def compute_distances(self, query_points, points): ) points = freud.util._convert_array(np.atleast_2d(points), shape=(None, 3)) - n_query_points = query_points.shape[0] - distances = np.empty(n_query_points, dtype=np.float32) + distances = np.empty(query_points.shape[0], dtype=np.float32) self._cpp_obj.computeDistances(query_points, points, distances) - return np.asarray(distances) + return distances def compute_all_distances(self, query_points, points): r"""Calculate distances between all pairs of query points and points, using periodic boundaries. @@ -490,7 +487,7 @@ def compute_all_distances(self, query_points, points): self._cpp_obj.computeAllDistances(query_points, points, distances) - return np.asarray(distances) + return distances def contains(self, points): r"""Compute a boolean array (mask) corresponding to point membership in a box. @@ -525,14 +522,11 @@ def contains(self, points): """ # noqa: E501 points = freud.util._convert_array(np.atleast_2d(points), shape=(None, 3)) - - n_points = points.shape[0] - - contains_mask = freud.util._convert_array(np.ones(n_points), dtype=bool) + contains_mask = freud.util._convert_array(np.ones(points.shape[0]), dtype=bool) self._cpp_obj.contains(points, contains_mask) - return np.array(contains_mask).astype(bool) + return contains_mask @property def cubic(self): From 91e32ef23b0492541b458e8ba753b2356df869ff Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Wed, 12 Jun 2024 16:51:23 -0400 Subject: [PATCH 17/25] move all binding code to its own header --- cpp/box/Box.h | 189 ++++++++++------------------------------- cpp/box/CMakeLists.txt | 2 +- cpp/box/export-box.h | 123 +++++++++++++++++++++++++++ cpp/box/module-box.cc | 21 ++--- 4 files changed, 180 insertions(+), 155 deletions(-) create mode 100644 cpp/box/export-box.h diff --git a/cpp/box/Box.h b/cpp/box/Box.h index 1bdfc1e75..823cabae1 100644 --- a/cpp/box/Box.h +++ b/cpp/box/Box.h @@ -12,9 +12,6 @@ #include "VectorMath.h" #include "utils.h" -#include -namespace nb = nanobind; - /*! \file Box.h \brief Represents simulation boxes and contains helpful wrapping functions. */ @@ -26,9 +23,6 @@ constexpr float TWO_PI = 2.0 * M_PI; namespace freud { namespace box { -template> -using nb_array = nb::ndarray; - //! Stores box dimensions and provides common routines for wrapping vectors back into the box /*! Box stores a standard HOOMD simulation box that goes from -L/2 to L/2 in each dimension, allowing Lx, Ly, Lz, and triclinic tilt factors xy, xz, and yz to be specified independently. @@ -224,6 +218,11 @@ class Box return v; } + //! Convert fractional coordinates into absolute coordinates in place + /*! \param vecs Vectors of fractional coordinates between 0 and 1 within + * parallelepipedal box + * \param out_data The array in which to place the wrapped vectors. + */ void makeAbsolute(const vec3* vecs, const unsigned int Nvecs, vec3* out) const { util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { @@ -234,20 +233,6 @@ class Box }); } - //! Convert fractional coordinates into absolute coordinates in place - /*! \param vecs Vectors of fractional coordinates between 0 and 1 within - * parallelepipedal box - * \param out_data The array in which to place the wrapped vectors. - */ - void makeAbsolutePython(nb_array> vecs, - nb_array> out) const - { - unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* out_data = (vec3*) (out.data()); - makeAbsolute(vecs_data, Nvecs, out_data); - } - //! Convert a point's coordinate from absolute to fractional box coordinates. /*! \param v The vector of the point in absolute coordinates. * \returns The vector of the point in fractional coordinates. @@ -266,6 +251,11 @@ class Box return delta; } + //! Convert point coordinates from absolute to fractional box coordinates. + /*! \param vecs Vectors to convert + * \param Nvecs Number of vectors + * \param out The array in which to place the wrapped vectors. + */ void makeFractional(const vec3* vecs, unsigned int Nvecs, vec3* out) const { util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { @@ -276,20 +266,6 @@ class Box }); } - //! Convert point coordinates from absolute to fractional box coordinates. - /*! \param vecs Vectors to convert - * \param Nvecs Number of vectors - * \param out The array in which to place the wrapped vectors. - */ - void makeFractionalPython(nb_array> vecs, - nb_array> out) const - { - unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* out_data = (vec3*) (out.data()); - makeFractional(vecs_data, Nvecs, out_data); - } - //! Get periodic image of a vector. /*! \param v The vector to check. * \param image The image of a given point. @@ -306,6 +282,11 @@ class Box image.z = (int) ((f.z >= float(0.0)) ? f.z + float(0.5) : f.z - float(0.5)); } + //! Get the periodic image vectors belongs to + /*! \param vecs The vectors to check + * \param Nvecs Number of vectors + \param res Array to save the images + */ void getImages(const vec3* vecs, unsigned int Nvecs, vec3* images) const { util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { @@ -316,19 +297,6 @@ class Box }); } - //! Get the periodic image vectors belongs to - /*! \param vecs The vectors to check - * \param Nvecs Number of vectors - \param res Array to save the images - */ - void getImagesPython(nb_array> vecs, nb_array> images) const - { - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* images_data = (vec3*) (images.data()); - getImages(vecs_data, Nvecs, images_data); - } - //! Wrap a vector back into the box /*! \param v Vector to wrap, updated to the minimum image obeying the periodic settings * \returns Wrapped vector @@ -357,6 +325,11 @@ class Box return makeAbsolute(v_frac); } + //! Wrap vectors back into the box in place + /*! \param vecs Vectors to wrap, updated to the minimum image obeying the periodic settings + * \param Nvecs Number of vectors + * \param out The array in which to place the wrapped vectors. + */ void wrap(const vec3* vecs, unsigned int Nvecs, vec3* out) const { util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { @@ -367,19 +340,12 @@ class Box }); } - //! Wrap vectors back into the box in place - /*! \param vecs Vectors to wrap, updated to the minimum image obeying the periodic settings - * \param Nvecs Number of vectors + //! Unwrap given positions to their absolute location in place + /*! \param vecs Vectors of coordinates to unwrap + * \param images images flags for this point + \param Nvecs Number of vectors * \param out The array in which to place the wrapped vectors. - */ - void wrapPython(nb_array> vecs, nb_array> out) const - { - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* out_data = (vec3*) (out.data()); - wrap(vecs_data, Nvecs, out_data); - } - + */ void unwrap(const vec3* vecs, const vec3* images, unsigned int Nvecs, vec3* out) const { util::forLoopWrapper(0, Nvecs, [&](size_t begin, size_t end) { @@ -395,21 +361,6 @@ class Box }); } - //! Unwrap given positions to their absolute location in place - /*! \param vecs Vectors of coordinates to unwrap - * \param images images flags for this point - \param Nvecs Number of vectors - * \param out The array in which to place the wrapped vectors. - */ - void unwrapPython(nb_array vecs, nb_array images, nb_array out) const - { - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* images_data = (vec3*) (images.data()); - vec3* out_data = (vec3*) (out.data()); - unwrap(vecs_data, images_data, Nvecs, out_data); - } - //! Compute center of mass for vectors /*! \param vecs Vectors to compute center of mass * \param Nvecs Number of vectors @@ -438,15 +389,11 @@ class Box / constants::TWO_PI)); } - std::vector centerOfMassPython(nb_array vecs, nb_array> masses) const - { - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - float* masses_data = (float*) (masses.data()); - auto com = centerOfMass(vecs_data, Nvecs, masses_data); - return {com.x, com.y, com.z}; - } - + //! Subtract center of mass from vectors + /*! \param vecs Vectors to center + * \param Nvecs Number of vectors + * \param masses Optional array of masses, of length Nvecs + */ void center(vec3* vecs, unsigned int Nvecs, const float* masses) const { vec3 com(centerOfMass(vecs, Nvecs, masses)); @@ -458,19 +405,6 @@ class Box }); } - //! Subtract center of mass from vectors - /*! \param vecs Vectors to center - * \param Nvecs Number of vectors - * \param masses Optional array of masses, of length Nvecs - */ - void centerPython(nb_array vecs, nb_array> masses) const - { - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - float* masses_data = (float*) (masses.data()); - center(vecs_data, Nvecs, masses_data); - } - //! Calculate distance between two points using boundary conditions /*! \param r_i Position of first point \param r_j Position of second point @@ -481,6 +415,14 @@ class Box return std::sqrt(dot(r_ij, r_ij)); } + //! Calculate distances between a set of query points and points. + /*! \param query_points Query point positions. + \param n_query_points The number of query points. + \param points Point positions. + \param n_points The number of points. + \param distances Pointer to array of length n_query_points containing distances between each point and + query_point (overwritten in place). + */ void computeDistances(const vec3* query_points, const unsigned int n_query_points, const vec3* points, const unsigned int n_points, float* distances) const { @@ -492,29 +434,14 @@ class Box }); } - //! Calculate distances between a set of query points and points. + //! Calculate all pairwise distances between a set of query points and points. /*! \param query_points Query point positions. \param n_query_points The number of query points. \param points Point positions. \param n_points The number of points. - \param distances Pointer to array of length n_query_points containing distances between each point and - query_point (overwritten in place). + \param distances Pointer to array of length n_query_points*n_points containing distances between + points and query_points (overwritten in place). */ - void computeDistancesPython(nb_array query_points, nb_array points, - nb_array> distances) const - { - const unsigned int n_query_points = query_points.shape(0); - vec3* query_points_data = (vec3*) (query_points.data()); - const unsigned int n_points = points.shape(0); - vec3* points_data = (vec3*) (points.data()); - float* distances_data = (float*) (distances.data()); - if (n_query_points != n_points) - { - throw std::invalid_argument("The number of query points and points must match."); - } - computeDistances(query_points_data, n_query_points, points_data, n_points, distances_data); - } - void computeAllDistances(const vec3* query_points, const unsigned int n_query_points, const vec3* points, const unsigned int n_points, float* distances) const { @@ -529,25 +456,12 @@ class Box } }); } - //! Calculate all pairwise distances between a set of query points and points. - /*! \param query_points Query point positions. - \param n_query_points The number of query points. - \param points Point positions. + + //! Get mask of points that fit inside the box. + /*! \param points Point positions. \param n_points The number of points. - \param distances Pointer to array of length n_query_points*n_points containing distances between - points and query_points (overwritten in place). + \param contains_mask Mask of points inside the box. */ - void computeAllDistancesPython(nb_array query_points, nb_array points, - nb_array> distances) const - { - const unsigned int n_query_points = query_points.shape(0); - vec3* query_points_data = (vec3*) (query_points.data()); - const unsigned int n_points = points.shape(0); - vec3* points_data = (vec3*) (points.data()); - float* distances_data = (float*) (distances.data()); - computeAllDistances(query_points_data, n_query_points, points_data, n_points, distances_data); - } - void contains(vec3* points, const unsigned int n_points, bool* contains_mask) const { util::forLoopWrapper(0, n_points, [&](size_t begin, size_t end) { @@ -560,19 +474,6 @@ class Box }); } - //! Get mask of points that fit inside the box. - /*! \param points Point positions. - \param n_points The number of points. - \param contains_mask Mask of points inside the box. - */ - void containsPython(nb_array points, nb_array> contains_mask) const - { - const unsigned int n_points = points.shape(0); - vec3* points_data = (vec3*) (points.data()); - bool* contains_mask_data = (bool*) (contains_mask.data()); - contains(points_data, n_points, contains_mask_data); - } - //! Get the shortest distance between opposite boundary planes of the box /*! The distance between two planes of the lattice is 2 Pi/|b_i|, where * b_1 is the reciprocal lattice vector of the Bravais lattice normal to diff --git a/cpp/box/CMakeLists.txt b/cpp/box/CMakeLists.txt index 94126008a..7a67bd5bb 100644 --- a/cpp/box/CMakeLists.txt +++ b/cpp/box/CMakeLists.txt @@ -1,5 +1,5 @@ # create the target and add the needed properties -nanobind_add_module(_box SHARED module-box.cc) +nanobind_add_module(_box SHARED module-box.cc export-box.h) target_link_libraries(_box PUBLIC TBB::tbb) # this probably isn't needed for box, but may be needed by other modules # target_compile_definitions( box # Avoid deprecation warnings for unsupported diff --git a/cpp/box/export-box.h b/cpp/box/export-box.h new file mode 100644 index 000000000..b75442dee --- /dev/null +++ b/cpp/box/export-box.h @@ -0,0 +1,123 @@ +#ifndef EXPORT_BOX_H +#define EXPORT_BOX_H + +#include "Box.h" + +#include +namespace nb = nanobind; + + +namespace freud { namespace box { + +template> +using nb_array = nb::ndarray; + + +void makeAbsolutePython(Box box, nb_array> vecs, + nb_array> out) +{ + unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (out.data()); + box.makeAbsolute(vecs_data, Nvecs, out_data); +} + + +void makeFractionalPython(Box box, nb_array> vecs, + nb_array> out) +{ + unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (out.data()); + box.makeFractional(vecs_data, Nvecs, out_data); +} + + +void getImagesPython(Box box, nb_array> vecs, nb_array> images) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* images_data = (vec3*) (images.data()); + box.getImages(vecs_data, Nvecs, images_data); +} + + +void wrapPython(Box box, nb_array> vecs, nb_array> out) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (out.data()); + box.wrap(vecs_data, Nvecs, out_data); +} + + +void unwrapPython(Box box, nb_array vecs, nb_array images, nb_array out) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* images_data = (vec3*) (images.data()); + vec3* out_data = (vec3*) (out.data()); + box.unwrap(vecs_data, images_data, Nvecs, out_data); +} + + +std::vector centerOfMassPython(Box box, nb_array vecs, nb_array> masses) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + float* masses_data = (float*) (masses.data()); + auto com = box.centerOfMass(vecs_data, Nvecs, masses_data); + return {com.x, com.y, com.z}; +} + + +void centerPython(Box box, nb_array vecs, nb_array> masses) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + float* masses_data = (float*) (masses.data()); + box.center(vecs_data, Nvecs, masses_data); +} + + +void computeDistancesPython(Box box, nb_array query_points, nb_array points, + nb_array> distances) +{ + const unsigned int n_query_points = query_points.shape(0); + vec3* query_points_data = (vec3*) (query_points.data()); + const unsigned int n_points = points.shape(0); + vec3* points_data = (vec3*) (points.data()); + float* distances_data = (float*) (distances.data()); + if (n_query_points != n_points) + { + throw std::invalid_argument("The number of query points and points must match."); + } + box.computeDistances(query_points_data, n_query_points, points_data, n_points, distances_data); +} + + +void computeAllDistancesPython(Box box, nb_array query_points, nb_array points, + nb_array> distances) +{ + const unsigned int n_query_points = query_points.shape(0); + vec3* query_points_data = (vec3*) (query_points.data()); + const unsigned int n_points = points.shape(0); + vec3* points_data = (vec3*) (points.data()); + float* distances_data = (float*) (distances.data()); + box.computeAllDistances(query_points_data, n_query_points, points_data, n_points, + distances_data); +} + + +void containsPython(Box box, nb_array points, nb_array> contains_mask) +{ + const unsigned int n_points = points.shape(0); + vec3* points_data = (vec3*) (points.data()); + bool* contains_mask_data = (bool*) (contains_mask.data()); + box.contains(points_data, n_points, contains_mask_data); +} + +}; }; // end namespace freud::box + +#endif + diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc index 0f23f1c06..b9443d5ef 100644 --- a/cpp/box/module-box.cc +++ b/cpp/box/module-box.cc @@ -5,6 +5,7 @@ #include #include "Box.h" +#include "export-box.h" using namespace freud::box; @@ -35,15 +36,15 @@ NB_MODULE(_box, m) .def("is2D", &Box::is2D) .def("set2D", &Box::set2D) .def("getVolume", &Box::getVolume) - .def("center", &Box::centerPython) - .def("centerOfMass", &Box::centerOfMassPython) + .def("center", ¢erPython) + .def("centerOfMass", ¢erOfMassPython) // other stuff - .def("makeAbsolute", &Box::makeAbsolutePython) - .def("makeFractional", &Box::makeFractionalPython) - .def("wrap", &Box::wrapPython) - .def("unwrap", &Box::unwrapPython) - .def("getImages", &Box::getImagesPython) - .def("computeDistances", &Box::computeDistancesPython) - .def("computeAllDistances", &Box::computeAllDistancesPython) - .def("contains", &Box::containsPython); + .def("makeAbsolute", &makeAbsolutePython) + .def("makeFractional", &makeFractionalPython) + .def("wrap", &wrapPython) + .def("unwrap", &unwrapPython) + .def("getImages", &getImagesPython) + .def("computeDistances", &computeDistancesPython) + .def("computeAllDistances", &computeAllDistancesPython) + .def("contains", &containsPython); } From e3e8c061e90bba3cd5478da2edf1e659c94194f2 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 12 Jun 2024 20:52:51 +0000 Subject: [PATCH 18/25] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- cpp/box/export-box.h | 20 +++++--------------- 1 file changed, 5 insertions(+), 15 deletions(-) diff --git a/cpp/box/export-box.h b/cpp/box/export-box.h index b75442dee..968f3600d 100644 --- a/cpp/box/export-box.h +++ b/cpp/box/export-box.h @@ -1,3 +1,6 @@ +// Copyright (c) 2010-2023 The Regents of the University of Michigan +// This file is from the freud project, released under the BSD 3-Clause License. + #ifndef EXPORT_BOX_H #define EXPORT_BOX_H @@ -6,13 +9,11 @@ #include namespace nb = nanobind; - namespace freud { namespace box { template> using nb_array = nb::ndarray; - void makeAbsolutePython(Box box, nb_array> vecs, nb_array> out) { @@ -22,7 +23,6 @@ void makeAbsolutePython(Box box, nb_array> vecs, box.makeAbsolute(vecs_data, Nvecs, out_data); } - void makeFractionalPython(Box box, nb_array> vecs, nb_array> out) { @@ -32,7 +32,6 @@ void makeFractionalPython(Box box, nb_array> vecs, box.makeFractional(vecs_data, Nvecs, out_data); } - void getImagesPython(Box box, nb_array> vecs, nb_array> images) { const unsigned int Nvecs = vecs.shape(0); @@ -41,7 +40,6 @@ void getImagesPython(Box box, nb_array> vecs, nb_array> vecs, nb_array> out) { const unsigned int Nvecs = vecs.shape(0); @@ -50,7 +48,6 @@ void wrapPython(Box box, nb_array> vecs, nb_array vecs, nb_array images, nb_array out) { const unsigned int Nvecs = vecs.shape(0); @@ -60,7 +57,6 @@ void unwrapPython(Box box, nb_array vecs, nb_array images, nb_array< box.unwrap(vecs_data, images_data, Nvecs, out_data); } - std::vector centerOfMassPython(Box box, nb_array vecs, nb_array> masses) { const unsigned int Nvecs = vecs.shape(0); @@ -70,7 +66,6 @@ std::vector centerOfMassPython(Box box, nb_array vecs, nb_array vecs, nb_array> masses) { const unsigned int Nvecs = vecs.shape(0); @@ -79,7 +74,6 @@ void centerPython(Box box, nb_array vecs, nb_array> ma box.center(vecs_data, Nvecs, masses_data); } - void computeDistancesPython(Box box, nb_array query_points, nb_array points, nb_array> distances) { @@ -95,7 +89,6 @@ void computeDistancesPython(Box box, nb_array query_points, nb_array query_points, nb_array points, nb_array> distances) { @@ -104,11 +97,9 @@ void computeAllDistancesPython(Box box, nb_array query_points, nb_array* points_data = (vec3*) (points.data()); float* distances_data = (float*) (distances.data()); - box.computeAllDistances(query_points_data, n_query_points, points_data, n_points, - distances_data); + box.computeAllDistances(query_points_data, n_query_points, points_data, n_points, distances_data); } - void containsPython(Box box, nb_array points, nb_array> contains_mask) { const unsigned int n_points = points.shape(0); @@ -117,7 +108,6 @@ void containsPython(Box box, nb_array points, nb_array> box.contains(points_data, n_points, contains_mask_data); } -}; }; // end namespace freud::box +}; }; // end namespace freud::box #endif - From d00138133f94075e71ac99f8bafc3bcdaafc094f Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Thu, 13 Jun 2024 15:15:14 -0400 Subject: [PATCH 19/25] add namespace and shared pointers --- cpp/box/CMakeLists.txt | 2 +- cpp/box/{export-box.h => export_Box.h} | 44 +++++++++++++------------- cpp/box/module-box.cc | 23 +++++++------- 3 files changed, 35 insertions(+), 34 deletions(-) rename cpp/box/{export-box.h => export_Box.h} (61%) diff --git a/cpp/box/CMakeLists.txt b/cpp/box/CMakeLists.txt index 7a67bd5bb..fbbb0d03e 100644 --- a/cpp/box/CMakeLists.txt +++ b/cpp/box/CMakeLists.txt @@ -1,5 +1,5 @@ # create the target and add the needed properties -nanobind_add_module(_box SHARED module-box.cc export-box.h) +nanobind_add_module(_box SHARED module-box.cc export_Box.h) target_link_libraries(_box PUBLIC TBB::tbb) # this probably isn't needed for box, but may be needed by other modules # target_compile_definitions( box # Avoid deprecation warnings for unsupported diff --git a/cpp/box/export-box.h b/cpp/box/export_Box.h similarity index 61% rename from cpp/box/export-box.h rename to cpp/box/export_Box.h index 968f3600d..4d600d057 100644 --- a/cpp/box/export-box.h +++ b/cpp/box/export_Box.h @@ -9,72 +9,72 @@ #include namespace nb = nanobind; -namespace freud { namespace box { +namespace freud { namespace box { namespace wrap { template> using nb_array = nb::ndarray; -void makeAbsolutePython(Box box, nb_array> vecs, +void makeAbsolute(std::shared_ptr box, nb_array> vecs, nb_array> out) { unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); - box.makeAbsolute(vecs_data, Nvecs, out_data); + box->makeAbsolute(vecs_data, Nvecs, out_data); } -void makeFractionalPython(Box box, nb_array> vecs, +void makeFractional(std::shared_ptr box, nb_array> vecs, nb_array> out) { unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); - box.makeFractional(vecs_data, Nvecs, out_data); + box->makeFractional(vecs_data, Nvecs, out_data); } -void getImagesPython(Box box, nb_array> vecs, nb_array> images) +void getImages(std::shared_ptr box, nb_array> vecs, nb_array> images) { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); vec3* images_data = (vec3*) (images.data()); - box.getImages(vecs_data, Nvecs, images_data); + box->getImages(vecs_data, Nvecs, images_data); } -void wrapPython(Box box, nb_array> vecs, nb_array> out) +void wrap(std::shared_ptr box, nb_array> vecs, nb_array> out) { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); vec3* out_data = (vec3*) (out.data()); - box.wrap(vecs_data, Nvecs, out_data); + box->wrap(vecs_data, Nvecs, out_data); } -void unwrapPython(Box box, nb_array vecs, nb_array images, nb_array out) +void unwrap(std::shared_ptr box, nb_array vecs, nb_array images, nb_array out) { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); vec3* images_data = (vec3*) (images.data()); vec3* out_data = (vec3*) (out.data()); - box.unwrap(vecs_data, images_data, Nvecs, out_data); + box->unwrap(vecs_data, images_data, Nvecs, out_data); } -std::vector centerOfMassPython(Box box, nb_array vecs, nb_array> masses) +std::vector centerOfMass(std::shared_ptr box, nb_array vecs, nb_array> masses) { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); float* masses_data = (float*) (masses.data()); - auto com = box.centerOfMass(vecs_data, Nvecs, masses_data); + auto com = box->centerOfMass(vecs_data, Nvecs, masses_data); return {com.x, com.y, com.z}; } -void centerPython(Box box, nb_array vecs, nb_array> masses) +void center(std::shared_ptr box, nb_array vecs, nb_array> masses) { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); float* masses_data = (float*) (masses.data()); - box.center(vecs_data, Nvecs, masses_data); + box->center(vecs_data, Nvecs, masses_data); } -void computeDistancesPython(Box box, nb_array query_points, nb_array points, +void computeDistances(std::shared_ptr box, nb_array query_points, nb_array points, nb_array> distances) { const unsigned int n_query_points = query_points.shape(0); @@ -86,10 +86,10 @@ void computeDistancesPython(Box box, nb_array query_points, nb_arraycomputeDistances(query_points_data, n_query_points, points_data, n_points, distances_data); } -void computeAllDistancesPython(Box box, nb_array query_points, nb_array points, +void computeAllDistances(std::shared_ptr box, nb_array query_points, nb_array points, nb_array> distances) { const unsigned int n_query_points = query_points.shape(0); @@ -97,17 +97,17 @@ void computeAllDistancesPython(Box box, nb_array query_points, nb_array* points_data = (vec3*) (points.data()); float* distances_data = (float*) (distances.data()); - box.computeAllDistances(query_points_data, n_query_points, points_data, n_points, distances_data); + box->computeAllDistances(query_points_data, n_query_points, points_data, n_points, distances_data); } -void containsPython(Box box, nb_array points, nb_array> contains_mask) +void contains(std::shared_ptr box, nb_array points, nb_array> contains_mask) { const unsigned int n_points = points.shape(0); vec3* points_data = (vec3*) (points.data()); bool* contains_mask_data = (bool*) (contains_mask.data()); - box.contains(points_data, n_points, contains_mask_data); + box->contains(points_data, n_points, contains_mask_data); } -}; }; // end namespace freud::box +}; }; }; // end namespace freud::box #endif diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc index b9443d5ef..5c769742f 100644 --- a/cpp/box/module-box.cc +++ b/cpp/box/module-box.cc @@ -3,9 +3,10 @@ #include #include +#include #include "Box.h" -#include "export-box.h" +#include "export_Box.h" using namespace freud::box; @@ -36,15 +37,15 @@ NB_MODULE(_box, m) .def("is2D", &Box::is2D) .def("set2D", &Box::set2D) .def("getVolume", &Box::getVolume) - .def("center", ¢erPython) - .def("centerOfMass", ¢erOfMassPython) + .def("center", &wrap::center) + .def("centerOfMass", &wrap::centerOfMass) // other stuff - .def("makeAbsolute", &makeAbsolutePython) - .def("makeFractional", &makeFractionalPython) - .def("wrap", &wrapPython) - .def("unwrap", &unwrapPython) - .def("getImages", &getImagesPython) - .def("computeDistances", &computeDistancesPython) - .def("computeAllDistances", &computeAllDistancesPython) - .def("contains", &containsPython); + .def("makeAbsolute", &wrap::makeAbsolute) + .def("makeFractional", &wrap::makeFractional) + .def("wrap", &wrap::wrap) + .def("unwrap", &wrap::unwrap) + .def("getImages", &wrap::getImages) + .def("computeDistances", &wrap::computeDistances) + .def("computeAllDistances", &wrap::computeAllDistances) + .def("contains", &wrap::contains); } From 7af2ccd49cbbd30424a700d5eb1eb7951c8c82e0 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 13 Jun 2024 19:15:32 +0000 Subject: [PATCH 20/25] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- cpp/box/export_Box.h | 19 +++++++++++-------- cpp/box/module-box.cc | 2 +- 2 files changed, 12 insertions(+), 9 deletions(-) diff --git a/cpp/box/export_Box.h b/cpp/box/export_Box.h index 4d600d057..2a082a16b 100644 --- a/cpp/box/export_Box.h +++ b/cpp/box/export_Box.h @@ -15,7 +15,7 @@ template> using nb_array = nb::ndarray; void makeAbsolute(std::shared_ptr box, nb_array> vecs, - nb_array> out) + nb_array> out) { unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); @@ -24,7 +24,7 @@ void makeAbsolute(std::shared_ptr box, nb_array> ve } void makeFractional(std::shared_ptr box, nb_array> vecs, - nb_array> out) + nb_array> out) { unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); @@ -32,7 +32,8 @@ void makeFractional(std::shared_ptr box, nb_array> box->makeFractional(vecs_data, Nvecs, out_data); } -void getImages(std::shared_ptr box, nb_array> vecs, nb_array> images) +void getImages(std::shared_ptr box, nb_array> vecs, + nb_array> images) { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); @@ -40,7 +41,8 @@ void getImages(std::shared_ptr box, nb_array> vecs, box->getImages(vecs_data, Nvecs, images_data); } -void wrap(std::shared_ptr box, nb_array> vecs, nb_array> out) +void wrap(std::shared_ptr box, nb_array> vecs, + nb_array> out) { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); @@ -57,7 +59,8 @@ void unwrap(std::shared_ptr box, nb_array vecs, nb_array images box->unwrap(vecs_data, images_data, Nvecs, out_data); } -std::vector centerOfMass(std::shared_ptr box, nb_array vecs, nb_array> masses) +std::vector centerOfMass(std::shared_ptr box, nb_array vecs, + nb_array> masses) { const unsigned int Nvecs = vecs.shape(0); vec3* vecs_data = (vec3*) (vecs.data()); @@ -75,7 +78,7 @@ void center(std::shared_ptr box, nb_array vecs, nb_array box, nb_array query_points, nb_array points, - nb_array> distances) + nb_array> distances) { const unsigned int n_query_points = query_points.shape(0); vec3* query_points_data = (vec3*) (query_points.data()); @@ -90,7 +93,7 @@ void computeDistances(std::shared_ptr box, nb_array query_points, nb } void computeAllDistances(std::shared_ptr box, nb_array query_points, nb_array points, - nb_array> distances) + nb_array> distances) { const unsigned int n_query_points = query_points.shape(0); vec3* query_points_data = (vec3*) (query_points.data()); @@ -108,6 +111,6 @@ void contains(std::shared_ptr box, nb_array points, nb_arraycontains(points_data, n_points, contains_mask_data); } -}; }; }; // end namespace freud::box +}; }; }; // namespace freud::box::wrap #endif diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc index 5c769742f..5168b850a 100644 --- a/cpp/box/module-box.cc +++ b/cpp/box/module-box.cc @@ -2,8 +2,8 @@ // This file is from the freud project, released under the BSD 3-Clause License. #include -#include #include +#include #include "Box.h" #include "export_Box.h" From e5848c2e522931b3fd4754c56013e51edbe46777 Mon Sep 17 00:00:00 2001 From: Tommy Waltmann <53307607+tommy-waltmann@users.noreply.github.com> Date: Fri, 14 Jun 2024 08:07:42 -0400 Subject: [PATCH 21/25] Update CMakeLists.txt Co-authored-by: Joshua A. Anderson --- CMakeLists.txt | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index d80f37d33..3074eda24 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -53,7 +53,9 @@ execute_process( OUTPUT_STRIP_TRAILING_WHITESPACE OUTPUT_VARIABLE nanobind_ROOT) find_package(nanobind CONFIG REQUIRED) - +if (nanobind_FOUND) + find_package_message(nanobind "Found nanobind: ${nanobind_DIR} ${nanobind_VERSION}" "[${nanobind_DIR},${nanobind_VERSION}]") +endif() # Fail fast if users have not cloned submodules. if(NOT WIN32) string(ASCII 27 Esc) From 111ca13359cdbc18302ac4012d1bf646b9542a61 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 14 Jun 2024 12:07:55 +0000 Subject: [PATCH 22/25] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- CMakeLists.txt | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 3074eda24..95362a546 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -53,8 +53,10 @@ execute_process( OUTPUT_STRIP_TRAILING_WHITESPACE OUTPUT_VARIABLE nanobind_ROOT) find_package(nanobind CONFIG REQUIRED) -if (nanobind_FOUND) - find_package_message(nanobind "Found nanobind: ${nanobind_DIR} ${nanobind_VERSION}" "[${nanobind_DIR},${nanobind_VERSION}]") +if(nanobind_FOUND) + find_package_message( + nanobind "Found nanobind: ${nanobind_DIR} ${nanobind_VERSION}" + "[${nanobind_DIR},${nanobind_VERSION}]") endif() # Fail fast if users have not cloned submodules. if(NOT WIN32) From 2939bf7fc2a048b7a04f2c826c874ca97e1e0751 Mon Sep 17 00:00:00 2001 From: Tommy Waltmann Date: Fri, 14 Jun 2024 08:26:26 -0400 Subject: [PATCH 23/25] move wrapper code to its own translation unit --- CMakeLists.txt | 1 + cpp/box/CMakeLists.txt | 2 +- cpp/box/export_Box.cc | 108 +++++++++++++++++++++++++++++++++++++ cpp/box/export_Box.h | 118 +++++++++-------------------------------- 4 files changed, 134 insertions(+), 95 deletions(-) create mode 100644 cpp/box/export_Box.cc diff --git a/CMakeLists.txt b/CMakeLists.txt index 95362a546..d1bca9668 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -58,6 +58,7 @@ if(nanobind_FOUND) nanobind "Found nanobind: ${nanobind_DIR} ${nanobind_VERSION}" "[${nanobind_DIR},${nanobind_VERSION}]") endif() + # Fail fast if users have not cloned submodules. if(NOT WIN32) string(ASCII 27 Esc) diff --git a/cpp/box/CMakeLists.txt b/cpp/box/CMakeLists.txt index fbbb0d03e..85078e5b2 100644 --- a/cpp/box/CMakeLists.txt +++ b/cpp/box/CMakeLists.txt @@ -1,5 +1,5 @@ # create the target and add the needed properties -nanobind_add_module(_box SHARED module-box.cc export_Box.h) +nanobind_add_module(_box SHARED module-box.cc export_Box.cc) target_link_libraries(_box PUBLIC TBB::tbb) # this probably isn't needed for box, but may be needed by other modules # target_compile_definitions( box # Avoid deprecation warnings for unsupported diff --git a/cpp/box/export_Box.cc b/cpp/box/export_Box.cc new file mode 100644 index 000000000..29068a4d6 --- /dev/null +++ b/cpp/box/export_Box.cc @@ -0,0 +1,108 @@ +// Copyright (c) 2010-2023 The Regents of the University of Michigan +// This file is from the freud project, released under the BSD 3-Clause License. + +#include "export_Box.h" + +namespace nb = nanobind; + +namespace freud { namespace box { namespace wrap { + +void makeAbsolute(std::shared_ptr box, nb_array> vecs, + nb_array> out) +{ + unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (out.data()); + box->makeAbsolute(vecs_data, Nvecs, out_data); +} + +void makeFractional(std::shared_ptr box, nb_array> vecs, + nb_array> out) +{ + unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (out.data()); + box->makeFractional(vecs_data, Nvecs, out_data); +} + +void getImages(std::shared_ptr box, nb_array> vecs, + nb_array> images) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* images_data = (vec3*) (images.data()); + box->getImages(vecs_data, Nvecs, images_data); +} + +void wrap(std::shared_ptr box, nb_array> vecs, + nb_array> out) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* out_data = (vec3*) (out.data()); + box->wrap(vecs_data, Nvecs, out_data); +} + +void unwrap(std::shared_ptr box, nb_array vecs, nb_array images, nb_array out) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + vec3* images_data = (vec3*) (images.data()); + vec3* out_data = (vec3*) (out.data()); + box->unwrap(vecs_data, images_data, Nvecs, out_data); +} + +std::vector centerOfMass(std::shared_ptr box, nb_array vecs, + nb_array> masses) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + float* masses_data = (float*) (masses.data()); + auto com = box->centerOfMass(vecs_data, Nvecs, masses_data); + return {com.x, com.y, com.z}; +} + +void center(std::shared_ptr box, nb_array vecs, nb_array> masses) +{ + const unsigned int Nvecs = vecs.shape(0); + vec3* vecs_data = (vec3*) (vecs.data()); + float* masses_data = (float*) (masses.data()); + box->center(vecs_data, Nvecs, masses_data); +} + +void computeDistances(std::shared_ptr box, nb_array query_points, nb_array points, + nb_array> distances) +{ + const unsigned int n_query_points = query_points.shape(0); + vec3* query_points_data = (vec3*) (query_points.data()); + const unsigned int n_points = points.shape(0); + vec3* points_data = (vec3*) (points.data()); + float* distances_data = (float*) (distances.data()); + if (n_query_points != n_points) + { + throw std::invalid_argument("The number of query points and points must match."); + } + box->computeDistances(query_points_data, n_query_points, points_data, n_points, distances_data); +} + +void computeAllDistances(std::shared_ptr box, nb_array query_points, nb_array points, + nb_array> distances) +{ + const unsigned int n_query_points = query_points.shape(0); + vec3* query_points_data = (vec3*) (query_points.data()); + const unsigned int n_points = points.shape(0); + vec3* points_data = (vec3*) (points.data()); + float* distances_data = (float*) (distances.data()); + box->computeAllDistances(query_points_data, n_query_points, points_data, n_points, distances_data); +} + +void contains(std::shared_ptr box, nb_array points, nb_array> contains_mask) +{ + const unsigned int n_points = points.shape(0); + vec3* points_data = (vec3*) (points.data()); + bool* contains_mask_data = (bool*) (contains_mask.data()); + box->contains(points_data, n_points, contains_mask_data); +} + +}; }; }; // namespace freud::box::wrap + diff --git a/cpp/box/export_Box.h b/cpp/box/export_Box.h index 2a082a16b..7836ebbb6 100644 --- a/cpp/box/export_Box.h +++ b/cpp/box/export_Box.h @@ -7,109 +7,39 @@ #include "Box.h" #include -namespace nb = nanobind; namespace freud { namespace box { namespace wrap { -template> -using nb_array = nb::ndarray; - -void makeAbsolute(std::shared_ptr box, nb_array> vecs, - nb_array> out) -{ - unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* out_data = (vec3*) (out.data()); - box->makeAbsolute(vecs_data, Nvecs, out_data); -} - -void makeFractional(std::shared_ptr box, nb_array> vecs, - nb_array> out) -{ - unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* out_data = (vec3*) (out.data()); - box->makeFractional(vecs_data, Nvecs, out_data); -} - -void getImages(std::shared_ptr box, nb_array> vecs, - nb_array> images) -{ - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* images_data = (vec3*) (images.data()); - box->getImages(vecs_data, Nvecs, images_data); -} - -void wrap(std::shared_ptr box, nb_array> vecs, - nb_array> out) -{ - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* out_data = (vec3*) (out.data()); - box->wrap(vecs_data, Nvecs, out_data); -} - -void unwrap(std::shared_ptr box, nb_array vecs, nb_array images, nb_array out) -{ - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - vec3* images_data = (vec3*) (images.data()); - vec3* out_data = (vec3*) (out.data()); - box->unwrap(vecs_data, images_data, Nvecs, out_data); -} +template> +using nb_array = nanobind::ndarray; + +void makeAbsolute(std::shared_ptr box, nb_array> vecs, + nb_array> out); + +void makeFractional(std::shared_ptr box, nb_array> vecs, + nb_array> out); + +void getImages(std::shared_ptr box, nb_array> vecs, + nb_array> images); + +void wrap(std::shared_ptr box, nb_array> vecs, + nb_array> out); + +void unwrap(std::shared_ptr box, nb_array vecs, nb_array images, nb_array out); std::vector centerOfMass(std::shared_ptr box, nb_array vecs, - nb_array> masses) -{ - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - float* masses_data = (float*) (masses.data()); - auto com = box->centerOfMass(vecs_data, Nvecs, masses_data); - return {com.x, com.y, com.z}; -} - -void center(std::shared_ptr box, nb_array vecs, nb_array> masses) -{ - const unsigned int Nvecs = vecs.shape(0); - vec3* vecs_data = (vec3*) (vecs.data()); - float* masses_data = (float*) (masses.data()); - box->center(vecs_data, Nvecs, masses_data); -} + nb_array> masses); + +void center(std::shared_ptr box, nb_array vecs, nb_array> masses); void computeDistances(std::shared_ptr box, nb_array query_points, nb_array points, - nb_array> distances) -{ - const unsigned int n_query_points = query_points.shape(0); - vec3* query_points_data = (vec3*) (query_points.data()); - const unsigned int n_points = points.shape(0); - vec3* points_data = (vec3*) (points.data()); - float* distances_data = (float*) (distances.data()); - if (n_query_points != n_points) - { - throw std::invalid_argument("The number of query points and points must match."); - } - box->computeDistances(query_points_data, n_query_points, points_data, n_points, distances_data); -} + nb_array> distances); void computeAllDistances(std::shared_ptr box, nb_array query_points, nb_array points, - nb_array> distances) -{ - const unsigned int n_query_points = query_points.shape(0); - vec3* query_points_data = (vec3*) (query_points.data()); - const unsigned int n_points = points.shape(0); - vec3* points_data = (vec3*) (points.data()); - float* distances_data = (float*) (distances.data()); - box->computeAllDistances(query_points_data, n_query_points, points_data, n_points, distances_data); -} - -void contains(std::shared_ptr box, nb_array points, nb_array> contains_mask) -{ - const unsigned int n_points = points.shape(0); - vec3* points_data = (vec3*) (points.data()); - bool* contains_mask_data = (bool*) (contains_mask.data()); - box->contains(points_data, n_points, contains_mask_data); -} + nb_array> distances); + +void contains(std::shared_ptr box, nb_array points, + nb_array> contains_mask); }; }; }; // namespace freud::box::wrap From df19cf91cf3dddc70ffa2bb01430916839163d8a Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 14 Jun 2024 12:26:41 +0000 Subject: [PATCH 24/25] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- cpp/box/export_Box.cc | 1 - cpp/box/export_Box.h | 2 +- 2 files changed, 1 insertion(+), 2 deletions(-) diff --git a/cpp/box/export_Box.cc b/cpp/box/export_Box.cc index 29068a4d6..1410b2e71 100644 --- a/cpp/box/export_Box.cc +++ b/cpp/box/export_Box.cc @@ -105,4 +105,3 @@ void contains(std::shared_ptr box, nb_array points, nb_array box, nb_array query_points, nb_array> distances); void contains(std::shared_ptr box, nb_array points, - nb_array> contains_mask); + nb_array> contains_mask); }; }; }; // namespace freud::box::wrap From 33326d1d596999d33d72fac79467d16bfc9a8b1a Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 19 Jun 2024 13:28:38 +0000 Subject: [PATCH 25/25] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- cpp/box/export_Box.cc | 2 +- cpp/box/export_Box.h | 2 +- cpp/box/module-box.cc | 2 +- freud/util.py | 2 +- 4 files changed, 4 insertions(+), 4 deletions(-) diff --git a/cpp/box/export_Box.cc b/cpp/box/export_Box.cc index 1410b2e71..562796f38 100644 --- a/cpp/box/export_Box.cc +++ b/cpp/box/export_Box.cc @@ -1,4 +1,4 @@ -// Copyright (c) 2010-2023 The Regents of the University of Michigan +// Copyright (c) 2010-2024 The Regents of the University of Michigan // This file is from the freud project, released under the BSD 3-Clause License. #include "export_Box.h" diff --git a/cpp/box/export_Box.h b/cpp/box/export_Box.h index f59f4d33a..fcf319691 100644 --- a/cpp/box/export_Box.h +++ b/cpp/box/export_Box.h @@ -1,4 +1,4 @@ -// Copyright (c) 2010-2023 The Regents of the University of Michigan +// Copyright (c) 2010-2024 The Regents of the University of Michigan // This file is from the freud project, released under the BSD 3-Clause License. #ifndef EXPORT_BOX_H diff --git a/cpp/box/module-box.cc b/cpp/box/module-box.cc index 5168b850a..9fc54f1e5 100644 --- a/cpp/box/module-box.cc +++ b/cpp/box/module-box.cc @@ -1,4 +1,4 @@ -// Copyright (c) 2010-2023 The Regents of the University of Michigan +// Copyright (c) 2010-2024 The Regents of the University of Michigan // This file is from the freud project, released under the BSD 3-Clause License. #include diff --git a/freud/util.py b/freud/util.py index 79e269c16..47b648c93 100644 --- a/freud/util.py +++ b/freud/util.py @@ -1,4 +1,4 @@ -# Copyright (c) 2010-2023 The Regents of the University of Michigan +# Copyright (c) 2010-2024 The Regents of the University of Michigan # This file is from the freud project, released under the BSD 3-Clause License. from functools import wraps