Remove GPU implementation of Stokesian Dynamics

doc/sphinx/installation.rst

@@ -491,13 +491,9 @@ using a CMake flag (see :ref:`Options and Variables`).
 - ``GSL`` Enables features relying on the GNU Scientific Library, e.g.
   :meth:`espressomd.cluster_analysis.Cluster.fractal_dimension`.
 
-- ``STOKESIAN_DYNAMICS`` Enables the Stokesian Dynamics feature for CPU
+- ``STOKESIAN_DYNAMICS`` Enables the Stokesian Dynamics feature
   (see :ref:`Stokesian Dynamics`). Requires BLAS and LAPACK.
 
-- ``STOKESIAN_DYNAMICS_GPU`` Enables the Stokesian Dynamics feature for GPU
-  (see :ref:`Stokesian Dynamics`). Requires thrust/cuBLAS/cuSolver.
-  Requires ``EXPERIMENTAL_FEATURES``.
-
 
 
 .. _Configuring:

samples/dancing.py

@@ -31,41 +31,26 @@
 
 parser = argparse.ArgumentParser(epilog=__doc__)
 group = parser.add_mutually_exclusive_group()
-group.add_argument('--cpu', action='store_true', help='Use CPU implementation')
-group.add_argument('--gpu', action='store_true', help='Use GPU implementation')
-group = parser.add_mutually_exclusive_group()
 group.add_argument('--ft', action='store_true', help='Use FT approximation')
 group.add_argument('--fts', action='store_true', help='Use FTS approximation')
 args = parser.parse_args()
 
-required_features = ["STOKESIAN_DYNAMICS"]
-
-if args.gpu:
-    print("Using GPU implementation")
-    required_features.append("CUDA")
-    required_features.append("STOKESIAN_DYNAMICS_GPU")
-    sd_device = "gpu"
-else:
-    print("Using CPU implementation")
-    sd_device = "cpu"
-
 if args.ft:
     print("Using FT approximation method")
     sd_method = "ft"
 else:
     print("Using FTS approximation method")
     sd_method = "fts"
 
-espressomd.assert_features(required_features)
+espressomd.assert_features("STOKESIAN_DYNAMICS")
 
 system = espressomd.System(box_l=[10, 10, 10])
 system.time_step = 1.5
 system.cell_system.skin = 0.4
 system.periodicity = [False, False, False]
 
 system.integrator.set_stokesian_dynamics(
-    viscosity=1.0, radii={0: 1.0}, device=sd_device,
-    approximation_method=sd_method)
+    viscosity=1.0, radii={0: 1.0}, approximation_method=sd_method)
 
 system.part.add(pos=[-5, 0, 0], rotation=[1, 1, 1])
 system.part.add(pos=[0, 0, 0], rotation=[1, 1, 1])
@@ -94,8 +79,7 @@
     plt.plot(data[:, i], data[:, i + 1], linestyle='dashed')
 
 plt.title("Trajectory of sedimenting spheres\nsolid line: simulation "
-          "({} on {}), dashed line: paper (FTS)"
-          .format(sd_method.upper(), sd_device.upper()))
+          "({}), dashed line: paper (FTS)".format(sd_method.upper()))
 plt.xlabel("x")
 plt.ylabel("y")
 plt.show()

src/config/features.def

@@ -26,10 +26,6 @@ ENGINE                          implies ROTATION, EXTERNAL_FORCES
 PARTICLE_ANISOTROPY             implies ROTATION
 STOKESIAN_DYNAMICS              requires BLAS and LAPACK
 STOKESIAN_DYNAMICS              implies ROTATION
-# Needs cuBLAS, which is expected to be included with CUDA
-STOKESIAN_DYNAMICS_GPU          requires CUDA
-STOKESIAN_DYNAMICS_GPU          implies ROTATION
-STOKESIAN_DYNAMICS_GPU          requires EXPERIMENTAL_FEATURES
 
 /* Rotation */
 ROTATION
@@ -111,4 +107,3 @@ GSL external
 BLAS external
 LAPACK external
 STOKESIAN_DYNAMICS external
-STOKESIAN_DYNAMICS_GPU external

src/core/CMakeLists.txt

@@ -108,12 +108,7 @@ if(WITH_UNIT_TESTS)
   add_subdirectory(unit_tests)
 endif(WITH_UNIT_TESTS)
 
-if(STOKESIAN_DYNAMICS OR STOKESIAN_DYNAMICS_GPU)
+if(STOKESIAN_DYNAMICS)
   add_subdirectory(stokesian_dynamics)
-  if(STOKESIAN_DYNAMICS)
-    target_link_libraries(EspressoCore PRIVATE StokesianDynamics::sd_cpu)
-  endif()
-  if(STOKESIAN_DYNAMICS_GPU)
-    target_link_libraries(EspressoCore PRIVATE StokesianDynamics::sd_gpu)
-  endif()
+  target_link_libraries(EspressoCore PRIVATE StokesianDynamics::sd_cpu)
 endif()

src/core/stokesian_dynamics/sd_interface.cpp

@@ -19,14 +19,8 @@
 
 #include "sd_interface.hpp"
 
-#if defined(STOKESIAN_DYNAMICS) || defined(STOKESIAN_DYNAMICS_GPU)
 #ifdef STOKESIAN_DYNAMICS
 #include "stokesian_dynamics/sd_cpu.hpp"
-#endif
-
-#ifdef STOKESIAN_DYNAMICS_GPU
-#include "stokesian_dynamics/sd_gpu.hpp"
-#endif
 
 #include "thermostat.hpp"
 
@@ -37,6 +31,7 @@
 #include <boost/range/algorithm.hpp>
 #include <boost/serialization/is_bitwise_serializable.hpp>
 
+#include <string>
 #include <vector>
 
 namespace {
@@ -63,8 +58,6 @@ struct SD_particle_data {
 
 double sd_viscosity = -1.0;
 
-enum { CPU, GPU, INVALID } device = INVALID;
-
 std::unordered_map<int, double> radius_dict;
 
 double sd_kT = 0.0;
@@ -118,35 +111,6 @@ void set_sd_viscosity(double eta) {
 
 double get_sd_viscosity() { return sd_viscosity; }
 
-void set_sd_device(std::string const &dev) {
-#ifdef STOKESIAN_DYNAMICS
-  if (dev == "cpu") {
-    device = CPU;
-    return;
-  }
-#endif
-#ifdef STOKESIAN_DYNAMICS_GPU
-  if (dev == "gpu") {
-    device = GPU;
-
-    return;
-  }
-#endif
-
-  throw std::runtime_error("Invalid device " + dev);
-}
-
-std::string get_sd_device() {
-  switch (device) {
-  case CPU:
-    return "cpu";
-  case GPU:
-    return "gpu";
-  default:
-    return "invalid";
-  }
-}
-
 void set_sd_radius_dict(std::unordered_map<int, double> const &x) {
   /* Check that radii are positive */
   for (auto const &kv : x) {
@@ -224,27 +188,10 @@ void propagate_vel_pos_sd(const ParticleRange &particles,
       ++i;
     }
 
-    switch (device) {
-#ifdef STOKESIAN_DYNAMICS
-    case CPU:
-      v_sd = sd_cpu(x_host, f_host, a_host, n_part, sd_viscosity,
-                    std::sqrt(sd_kT / time_step),
-                    static_cast<std::size_t>(stokesian.rng_counter()),
-                    static_cast<std::size_t>(stokesian.rng_seed()), sd_flags);
-      break;
-#endif
-
-#ifdef STOKESIAN_DYNAMICS_GPU
-    case GPU:
-      v_sd = sd_gpu(x_host, f_host, a_host, n_part, sd_viscosity,
-                    std::sqrt(sd_kT / time_step),
-                    static_cast<std::size_t>(stokesian.rng_counter()),
-                    static_cast<std::size_t>(stokesian.rng_seed()), sd_flags);
-      break;
-#endif
-    default:
-      throw std::runtime_error("Stokesian dynamics device not set.");
-    }
+    v_sd = sd_cpu(x_host, f_host, a_host, n_part, sd_viscosity,
+                  std::sqrt(sd_kT / time_step),
+                  static_cast<std::size_t>(stokesian.rng_counter()),
+                  static_cast<std::size_t>(stokesian.rng_seed()), sd_flags);
   } else { // if (this_node == 0)
     v_sd.resize(particles.size() * 6);
   } // if (this_node == 0) {...} else

src/core/stokesian_dynamics/sd_interface.hpp

@@ -32,15 +32,11 @@
 
 #include <boost/mpi/communicator.hpp>
 
-#include <string>
 #include <unordered_map>
 
 void set_sd_viscosity(double eta);
 double get_sd_viscosity();
 
-void set_sd_device(std::string const &dev);
-std::string get_sd_device();
-
 void set_sd_radius_dict(std::unordered_map<int, double> const &x);
 std::unordered_map<int, double> get_sd_radius_dict();
 

src/core/thermostat.cpp

@@ -74,7 +74,7 @@ ThermalizedBondThermostat thermalized_bond = {};
 #ifdef DPD
 DPDThermostat dpd = {};
 #endif
-#if defined(STOKESIAN_DYNAMICS) || defined(STOKESIAN_DYNAMICS_GPU)
+#ifdef STOKESIAN_DYNAMICS
 StokesianThermostat stokesian = {};
 #endif
 
@@ -85,7 +85,7 @@ REGISTER_THERMOSTAT_CALLBACKS(thermalized_bond)
 #ifdef DPD
 REGISTER_THERMOSTAT_CALLBACKS(dpd)
 #endif
-#if defined(STOKESIAN_DYNAMICS) || defined(STOKESIAN_DYNAMICS_GPU)
+#ifdef STOKESIAN_DYNAMICS
 REGISTER_THERMOSTAT_CALLBACKS(stokesian)
 #endif
 
@@ -133,7 +133,7 @@ void philox_counter_increment() {
     dpd.rng_increment();
   }
 #endif
-#if defined(STOKESIAN_DYNAMICS) || defined(STOKESIAN_DYNAMICS_GPU)
+#ifdef STOKESIAN_DYNAMICS
   if (thermo_switch & THERMO_SD) {
     stokesian.rng_increment();
   }

src/core/thermostat.hpp

@@ -316,7 +316,7 @@ struct ThermalizedBondThermostat : public BaseThermostat {};
 struct DPDThermostat : public BaseThermostat {};
 #endif
 
-#if defined(STOKESIAN_DYNAMICS) || defined(STOKESIAN_DYNAMICS_GPU)
+#ifdef STOKESIAN_DYNAMICS
 /** %Thermostat for Stokesian dynamics. */
 struct StokesianThermostat : public BaseThermostat {
   StokesianThermostat() { rng_initialize(0); }
@@ -343,7 +343,7 @@ NEW_THERMOSTAT(thermalized_bond)
 #ifdef DPD
 NEW_THERMOSTAT(dpd)
 #endif
-#if defined(STOKESIAN_DYNAMICS) || defined(STOKESIAN_DYNAMICS_GPU)
+#ifdef STOKESIAN_DYNAMICS
 NEW_THERMOSTAT(stokesian)
 #endif
 
@@ -355,7 +355,7 @@ extern ThermalizedBondThermostat thermalized_bond;
 #ifdef DPD
 extern DPDThermostat dpd;
 #endif
-#if defined(STOKESIAN_DYNAMICS) || defined(STOKESIAN_DYNAMICS_GPU)
+#ifdef STOKESIAN_DYNAMICS
 extern StokesianThermostat stokesian;
 #endif
 

src/python/espressomd/integrate.pxd

@@ -43,20 +43,17 @@ IF NPT:
                                              cbool zdir_rescale, cbool cubic_box)
 
 cdef extern from "stokesian_dynamics/sd_interface.hpp":
-    IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
+    IF STOKESIAN_DYNAMICS:
         void set_sd_viscosity(double eta)
         double get_sd_viscosity()
 
-        void set_sd_device(const string & dev)
-        string get_sd_device()
-
         void set_sd_radius_dict(const unordered_map[int, double] & radius_dict)
         unordered_map[int, double] get_sd_radius_dict()
 
         void set_sd_flags(int flg)
         int get_sd_flags()
 
-IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
+IF STOKESIAN_DYNAMICS:
     cpdef enum flags:
         NONE = 0,
         SELF_MOBILITY = 1 << 0,

src/python/espressomd/integrate.pyx

@@ -406,7 +406,7 @@ cdef class BrownianDynamics(Integrator):
         integrate_set_bd()
 
 
-IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
+IF STOKESIAN_DYNAMICS:
     cdef class StokesianDynamics(Integrator):
         """
         Stokesian Dynamics integrator.
@@ -417,8 +417,6 @@ IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
             Bulk viscosity.
         radii : :obj:`dict`
             Dictionary that maps particle types to radii.
-        device : :obj:`str`, optional, \{'cpu', 'gpu'\}
-            Device to execute on.
         approximation_method : :obj:`str`, optional, \{'ft', 'fts'\}
             Chooses the method of the mobility approximation.
             ``'fts'`` is more accurate. Default is ``'fts'``.
@@ -432,20 +430,15 @@ IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
         """
 
         def default_params(self):
-            IF STOKESIAN_DYNAMICS:
-                sd_device_str = "cpu"
-            ELIF STOKESIAN_DYNAMICS_GPU:
-                sd_device_str = "gpu"
             return {"lubrication": False, "approximation_method": "fts",
-                    "self_mobility": True, "pair_mobility": True,
-                    "device": sd_device_str}
+                    "self_mobility": True, "pair_mobility": True}
 
         def valid_keys(self):
             """All parameters that can be set.
 
             """
-            return {"radii", "viscosity", "device", "lubrication",
-                    "approximation_method", "self_mobility", "pair_mobility"}
+            return {"radii", "viscosity", "lubrication", "approximation_method",
+                    "self_mobility", "pair_mobility"}
 
         def required_keys(self):
             """Parameters that have to be set.
@@ -457,9 +450,6 @@ IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
             check_type_or_throw_except(
                 self._params["viscosity"], 1, float,
                 "viscosity must be a number")
-            check_type_or_throw_except(
-                self._params["device"], 1, str,
-                "device must be a string")
             check_type_or_throw_except(
                 self._params["radii"], 1, dict,
                 "radii must be a dictionary")
@@ -486,7 +476,6 @@ IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
         def _set_params_in_es_core(self):
             integrate_set_sd()
             set_sd_radius_dict(self._params["radii"])
-            set_sd_device(to_char_pointer(self._params["device"].lower()))
             set_sd_viscosity(self._params["viscosity"])
             fl = flags.NONE
             if self._params["lubrication"]:

src/python/espressomd/thermostat.pxd

@@ -61,7 +61,7 @@ cdef extern from "thermostat.hpp":
     IF DPD:
         cdef cppclass DPDThermostat(BaseThermostat):
             pass
-    IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
+    IF STOKESIAN_DYNAMICS:
         cdef cppclass StokesianThermostat(BaseThermostat):
             pass
 
@@ -71,27 +71,27 @@ cdef extern from "thermostat.hpp":
     ThermalizedBondThermostat thermalized_bond
     IF DPD:
         DPDThermostat dpd
-    IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
+    IF STOKESIAN_DYNAMICS:
         StokesianThermostat stokesian
 
     void langevin_set_rng_seed(stdint.uint32_t seed)
     void brownian_set_rng_seed(stdint.uint32_t seed)
     void npt_iso_set_rng_seed(stdint.uint32_t seed)
     IF DPD:
         void dpd_set_rng_seed(stdint.uint32_t seed)
-    IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
+    IF STOKESIAN_DYNAMICS:
         void stokesian_set_rng_seed(stdint.uint32_t seed)
 
     void langevin_set_rng_counter(stdint.uint64_t counter)
     void brownian_set_rng_counter(stdint.uint64_t counter)
     void npt_iso_set_rng_counter(stdint.uint64_t counter)
     IF DPD:
         void dpd_set_rng_counter(stdint.uint64_t counter)
-    IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
+    IF STOKESIAN_DYNAMICS:
         void stokesian_set_rng_counter(stdint.uint64_t counter)
 
 cdef extern from "stokesian_dynamics/sd_interface.hpp":
-    IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU):
+    IF STOKESIAN_DYNAMICS:
         void set_sd_kT(double kT)
         double get_sd_kT()