Merge pull request #39 from tovrstra/random-cell

Random cell

cellcutoff/cell.cpp

@@ -23,6 +23,7 @@
 
 #include <algorithm>
 #include <cmath>
+#include <random>
 #include <stdexcept>
 #include <vector>
 
@@ -344,6 +345,43 @@ void Cell::bars_cutoff_low(SphereSlice* slice, int ivec, std::vector<int>* bars)
 }
 
 
+// Free functions
+
+
+Cell* create_random_cell(unsigned int seed, const int nvec,
+    const double scale, const double ratio, const bool cuboid) {
+  // Range check
+  if ((nvec < 0) || (nvec > 3))
+    throw std::domain_error("A random cell must be 0D, 1D, 2D or 2D periodic.");
+  if (nvec == 0)
+    return new Cell();
+  // Randomly construct a cell till a decent one (sufficient volume) is found.
+  double vecs[9];
+  std::minstd_rand gen(seed);
+  while (true) {
+    // Fill the array with random data for given seed.
+    std::uniform_real_distribution<double> dis(-scale, scale);
+    for (int i=0; i < 9; ++i)
+      vecs[i] = dis(gen);
+    if (cuboid) {
+      vecs[1] = 0.0;
+      vecs[2] = 0.0;
+      if (nvec > 1) {
+        vecs[3] = 0.0;
+        vecs[5] = 0.0;
+      }
+      if (nvec > 2) {
+        vecs[6] = 0.0;
+        vecs[7] = 0.0;
+      }
+    }
+    if (fabs(vec3::triple(vecs, vecs + 3, vecs + 6)) > pow(ratio*scale, 3)) {
+      return new Cell(vecs, nvec);
+    }
+  }
+}
+
+
 }  // namespace cellcutoff
 
 // vim: textwidth=90 et ts=2 sw=2

cellcutoff/cell.h

@@ -30,6 +30,7 @@
 #define CELLCUTOFF_CELL_H_
 
 #include <exception>
+#include <memory>
 #include <string>
 #include <vector>
 
@@ -353,6 +354,35 @@ class Cell {
 };
 
 
+// Free functions
+
+
+/** @brief
+        Return a randomized but sensible unit cell.
+
+    @param seed
+        The random seed used to initialize the cell vectors.
+
+    @param nvec
+        The dimensionality of the cell.
+
+    @param scale
+        Determines the overall size of the random cell. The components of the cell vectors
+        are sampled from a uniform random distribution over the interval [-scale, scale].
+
+    @param ratio
+        Controls the minimal volume of the cell. Random cell vectors are tried until the
+        volume is above (ratio*scale)**nvec. When a suitable volume is found, these
+        cell vectors are used. The closer to 2, the more cubic the cell. High values
+        of this parameter may make this function slow because many trials will be needed.
+
+    @param cuboid
+        When true, an orthorhombic cell is constructed.
+  */
+Cell* create_random_cell(const unsigned int seed, const int nvec, const double scale = 1.0,
+  const double ratio = 0.1, const bool cuboid = false);
+
+
 }  // namespace cellcutoff
 
 

cellcutoff/tests/common.cpp

@@ -99,36 +99,6 @@ unsigned int fill_random_permutation(const unsigned int seed, int* array,
 }
 
 
-std::unique_ptr<cl::Cell> create_random_cell_nvec(unsigned int seed, const int nvec,
-    const double scale, const double ratio, const bool cuboid) {
-  // Range check
-  if ((nvec < 0) || (nvec > 3))
-    throw std::domain_error("A random cell must be 0D, 1D, 2D or 2D periodic.");
-  if (nvec == 0)
-    return std::unique_ptr<cl::Cell>(new cl::Cell());
-  // Randomly construct a cell till a decent one (sufficient volume) is found.
-  double vecs[9];
-  while (true) {
-    seed = fill_random_double(seed, vecs, 9, -scale, scale);
-    if (cuboid) {
-      vecs[1] = 0.0;
-      vecs[2] = 0.0;
-      if (nvec > 1) {
-        vecs[3] = 0.0;
-        vecs[5] = 0.0;
-      }
-      if (nvec > 2) {
-        vecs[6] = 0.0;
-        vecs[7] = 0.0;
-      }
-    }
-    if (fabs(vec3::triple(vecs, vecs + 3, vecs + 6)) > pow(ratio*scale, 3)) {
-      return std::unique_ptr<cl::Cell>(new cl::Cell(vecs, nvec));
-    }
-  }
-}
-
-
 unsigned int random_point(unsigned int seed,  const double* center,
     const double cutoff, double* point, double* norm) {
   seed = fill_random_double(seed, point, 3, -cutoff, cutoff);
@@ -146,8 +116,6 @@ TEST(CommonTest, domain) {
   EXPECT_THROW(fill_random_int(0, nullptr, 1, 1, 0), std::domain_error);
   EXPECT_THROW(fill_random_permutation(0, nullptr, 0), std::domain_error);
   EXPECT_THROW(fill_random_permutation(0, nullptr, -1), std::domain_error);
-  EXPECT_THROW(create_random_cell_nvec(1277, -1, 1, false), std::domain_error);
-  EXPECT_THROW(create_random_cell_nvec(1274, 4, 1, false), std::domain_error);
 }
 
 // vim: textwidth=90 et ts=2 sw=2

cellcutoff/tests/common.h

@@ -22,8 +22,6 @@
 #ifndef CELLCUTOFF_TESTS_COMMON_H_
 #define CELLCUTOFF_TESTS_COMMON_H_
 
-#include <memory>
-
 #include "cellcutoff/cell.h"
 
 
@@ -57,30 +55,6 @@ unsigned int fill_random_int(const unsigned int seed, int* array, const int size
 //! Fills and array of int with a random permutation
 unsigned int fill_random_permutation(const unsigned int seed, int* array, const int size);
 
-/** @brief
-        Return a randomized but sensible unit cell.
-
-    @param seed
-        The random seed used to initialize the cell vectors.
-
-    @param nvec
-        The dimensionality of the cell.
-
-    @param scale
-        Determines the overall size of the random cell. The components of the cell vectors
-        are sampled from a uniform random distribution over the interval [-scale, scale].
-
-    @param ratio
-        Controls the minimal volume of the cell. Random cell vectors are tried until the
-        volume is above (ratio*scale)**nvec. When a suitable volume is found, these
-        cell vectors are used.
-
-    @param cuboid
-        When true, an orthorhombic cell is constructed.
-  */
-std::unique_ptr<cl::Cell> create_random_cell_nvec(const unsigned int seed, const int nvec,
-    const double scale = 1.0, const double ratio = 0.1, const bool cuboid = false);
-
 //! Compute a random point in a cubic box centered around center. Also computes distance.
 unsigned int random_point(const unsigned int seed,  const double* center,
     const double cutoff, double* point, double* norm);

cellcutoff/tests/test_cell.cpp

@@ -72,7 +72,7 @@ class CellTest : public ::testing::Test {
 
   std::unique_ptr<cl::Cell> create_random_cell(const unsigned int seed,
       const double scale = 1.0, const double ratio = 0.1, const bool cuboid = false) {
-    return create_random_cell_nvec(seed, nvec, scale, ratio, cuboid);
+    return std::unique_ptr<cl::Cell>(cl::create_random_cell(seed, nvec, scale, ratio, cuboid));
   }
 
   int nvec;
@@ -1282,6 +1282,18 @@ TEST_F(CellTest3, bars_cutoff_corners) {
   EXPECT_LE(NREP*14, nbar_total);
 }
 
+
+// Free functions: create_random_cell
+// ----------------------------------
+
+
+
+TEST(RandomCell, domain) {
+  EXPECT_THROW(cl::create_random_cell(1277, -1, 1, false), std::domain_error);
+  EXPECT_THROW(cl::create_random_cell(1274, 4, 1, false), std::domain_error);
+}
+
+
 // Instantiation of parameterized tests
 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 

cellcutoff/tests/test_decomposition.cpp

@@ -24,6 +24,7 @@
 
 #include <gtest/gtest.h>
 
+#include <cellcutoff/cell.h>
 #include <cellcutoff/decomposition.h>
 
 #include "cellcutoff/tests/common.h"
@@ -173,7 +174,7 @@ TEST(DecompositionTest, assign_icell_example_shape) {
 TEST(DecompositionTest, assign_icell_random_wrap) {
   for (int irep = 0; irep < NREP; ++irep) {
     // Get a random 3D cell
-    std::unique_ptr<cl::Cell> cell(create_random_cell_nvec(irep*NPOINT, 3, 2));
+    std::unique_ptr<cl::Cell> cell(cl::create_random_cell(irep*NPOINT, 3, 2));
     // Get a subcell
     double threshold = 0.2;
     int shape[3] = {-1, -1, -1};
@@ -272,7 +273,7 @@ TEST(DecompositionTest, random_cell_map) {
       fill_random_double(ipoint+3157, cart, 3, -5.0, 5.0);
       points.push_back(cl::Point(cart));
     }
-    std::unique_ptr<cl::Cell> subcell(create_random_cell_nvec(irep*NPOINT, 3));
+    std::unique_ptr<cl::Cell> subcell(cl::create_random_cell(irep*NPOINT, 3));
     cl::assign_icell(*subcell, points.data(), points.size(), sizeof(cl::Point));
     cl::sort_by_icell(points.data(), points.size(), sizeof(cl::Point));
     std::unique_ptr<cl::CellMap> cell_map(

cellcutoff/tests/test_iterators.cpp

@@ -43,7 +43,7 @@ class BarIteratorTestP : public ::testing::TestWithParam<int> {
 
   std::unique_ptr<cl::Cell> create_random_cell(const unsigned int seed,
       const double scale = 1.0, const double ratio = 0.1, const bool cuboid = false) {
-    return create_random_cell_nvec(seed, nvec, scale, ratio, cuboid);
+    return std::unique_ptr<cl::Cell>(cl::create_random_cell(seed, nvec, scale, ratio, cuboid));
   }
 
   int nvec;

cellcutoff/tests/test_usage.cpp

@@ -52,7 +52,7 @@ class UsageTestP : public ::testing::TestWithParam<int> {
 
   std::unique_ptr<cl::Cell> create_random_cell(const unsigned int seed,
       const double scale = 1.0, const double ratio = 0.1, const bool cuboid = false) {
-    return create_random_cell_nvec(seed, nvec, scale, ratio, cuboid);
+    return std::unique_ptr<cl::Cell>(cl::create_random_cell(seed, nvec, scale, ratio, cuboid));
   }
 
   int nvec;

python-cellcutoff/cellcutoff/__init__.py

@@ -21,4 +21,4 @@
 
 
 from .version import __version__
-from .ext import Cell
+from .ext import Cell, create_random_cell

python-cellcutoff/cellcutoff/cell.pxd

@@ -49,3 +49,6 @@ cdef extern from "cellcutoff/cell.h" namespace "cellcutoff":
 
         size_t ranges_cutoff(const double* center, double cutoff, int* ranges_begin,
             int* ranges_end) const;
+
+    Cell* create_random_cell(const unsigned int seed, const int nvec,
+        const double scale, const double ratio, const bool cuboid)

python-cellcutoff/cellcutoff/ext.pyx

@@ -26,11 +26,12 @@ cimport numpy as np
 np.import_array()
 
 from libc.string cimport memcpy
+from libcpp cimport bool
 
-cimport cellcutoff.cell
+cimport cellcutoff.cell as cell
 
 
-__all__ = ['Cell']
+__all__ = ['Cell', 'create_random_cell']
 
 
 def check_array_arg(name, arg, expected_shape):
@@ -267,3 +268,38 @@ cdef class Cell(object):
         cdef np.ndarray[int, ndim=1] ranges_end = np.zeros(3, np.intc)
         self._this.ranges_cutoff(&center[0], cutoff, &ranges_begin[0], &ranges_end[0])
         return ranges_begin, ranges_end
+
+
+def create_random_cell(int seed, int nvec, double scale=10.0,
+                       double ratio=0.1, bool cuboid=False) -> Cell:
+    """Return a sensible random cell.
+
+    Parameters
+    ----------
+    seed
+        A random seed, for reproducile testing.
+    nvec
+        The dimensionality of the cell.
+    scale
+        Determines the overall size of the random cell. The components of the
+        cell vectors are sampled from a uniform random distribution over the
+        interval [-scale, scale].
+    ratio
+        Controls the minimal volume of the cell. Random cell vectors are tried
+        until the volume is above (ratio*scale)**nvec. When a suitable volume is
+        found, these cell vectors are used. The closer to 2, the more cubic the
+        cell. High values of this parameter may make this function slow because
+        many trials will be needed.
+    cuboid
+        When True, a random orthorhombic cell is generated.
+
+    Returns
+    -------
+    cell
+        The random cell object.
+
+    """
+    cdef cell.Cell* cpp_cell = cell.create_random_cell(seed, nvec, scale, ratio, cuboid)
+    cdef Cell result = Cell.__new__(Cell, None, initvoid=True)
+    result._this = cpp_cell
+    return result