Skip to content

Software design

Jump to bottom
Jean-Noël Grad edited this page Dec 6, 2022 · 16 revisions

This page documents the way ESPResSo is designed and how it interfaces with libraries.

Table of Contents

Core

MPI parallel code

ESPResSo uses a MPI callback mechanism based on a visitor pattern. At the start of the program, global variables from the head node are communicated to the worker nodes. The script interface sends and receives data from the core on the head node. When particles are created, the head node distributes them to the corresponding worker nodes based on the currently active domain decomposition.

Functions in the script interface that change the global state of the system in the core or accesses particles are implemented using callbacks. There are also "event" callbacks that get triggered when the global state of the system changes, e.g. the box dimensions in NpT simulations, to synchronize the global state on every MPI rank.

Below is one possible implementation of a callback to calculate the total linear momentum of the system. The momentum calculation is written in the callback calc_linear_momentum_local(). It is registered as a callback and tagged for a reduction operation using the macro REGISTER_CALLBACK_REDUCTION(). The calc_linear_momentum() function is responsible for calling the callbacks on all ranks. This function is exported in the Python interface, where it can be called by the user on the head node. 

#include "Particle.hpp"
#include "cells.hpp"
#include "communication.hpp"
#include <utils/Vector.hpp>

/** Callback to calculate the linear momentum of the cell on this node. */
static Utils::Vector3d calc_linear_momentum_local() {
  auto const particles = cell_structure.local_particles();
  auto const momentum =
      std::accumulate(particles.begin(), particles.end(), Utils::Vector3d{},
                      [](Utils::Vector3d &m, Particle const &p) {
                        return m + p.p.mass * p.m.v;
                      });
  return momentum;
}

REGISTER_CALLBACK_REDUCTION(calc_linear_momentum_local,
                            std::plus<Utils::Vector3d>())

/** Calculate the total linear momentum of the system.
 *  Must be called on the head node.
 */
Utils::Vector3d calc_linear_momentum() {
  Utils::Vector3d const result = mpi_call(::Communication::Result::reduction,
                                          std::plus<Utils::Vector3d>(),
                                          calc_linear_momentum_local);
  return result;
}

There are specialized macros and mpi_call() overloads to accommodate for a wide range of collective operations. See the doxygen documentation of header files MpiCallbacks.hpp and communication.hpp for more details.

Exception mechanism

Non-parallel functions

Functions called exclusively on the head node should throw an exception. The corresponding import in the Python interface should decorate the symbol with except +, so that the C++ exception gets transformed into a Python exception. 

# sd.pxd
cdef extern from "sd.hpp":
    void set_sd_viscosity(double eta) except +
/* sd.hpp */
void set_sd_viscosity(double eta);

/* sd.cpp */
#include <stdexcept>
#include <string>

double sd_viscosity;
void set_sd_viscosity(double eta) {
  if (eta < 0.0)
    throw std::runtime_error("invalid viscosity: " + std::to_string(eta));
  sd_viscosity = eta;
}

Parallel functions

Functions called by MPI callbacks cannot throw exceptions. Instead they should use the runtime errors logging mechanism to log the error message. All Python functions that indirectly call such a callback must add a call to handle_errors(), which will halt the flow of the program if the logger contains error messages (by wrapping them in a Python exception). 

/* p3m.hpp */
void p3m_tune();

/* p3m.cpp */
#include "p3m.hpp"
#include "errorhandling.hpp"
#include "communication.hpp"

void p3m_tune_local() {
  int error = sanity_checks();
  if (error) {
    runtimeErrorMsg() << "sanity checks failed";
    return;
  }
  /* implementation not shown */
}

REGISTER_CALLBACK(p3m_tune_local)

void p3m_tune() {
  mpi_call_all(p3m_tune_local);
}
# electrostatics.pxd
cdef extern from "p3m.hpp":
    void p3m_tune()

# electrostatics.pyx
from . cimport electrostatics  # get symbols from the pxd file
from .utils cimport handle_errors

def tune():
    p3m_tune()
    # convert the logged error messages into a Python exception
    handle_errors("P3M tuning failed")

CUDA functions

CUDA functions return an error code upon failure. This code can be converted to a runtime error using throw cuda_runtime_error_cuda(error_code), or alternatively via the convenience macro CUDA_CHECK(). 

/* cuda_init.cu */
#include "cuda_utils.cuh"
#include <cuda.h>

int cuda_get_n_gpus() {
  int deviceCount;
  CUDA_CHECK(cudaGetDeviceCount(&deviceCount))
  return deviceCount;
}

In an MPI environment, the error can be caught via its parent class cuda_runtime_error. 

/* gpu_callbacks.cpp */
#include "cuda_init.hpp"
#include "cuda_utils.hpp"

int n_devices_local_node() {
  int n_devices;
  try {
    n_devices = cuda_get_n_gpus();
  } catch (cuda_runtime_error const &err) {
    n_devices = 0;
  }
  return n_devices;
}

Cython interface

Whenever possible, Cython should only check variable types and the core should check the variable values. This avoids unnecessary code duplication between the core and interface for range checking.

Exceptions need to be checked in the testsuite using with self.assertRaises(RuntimeError). The test case should also check the system isn't left in an undefined state, i.e. the exceptions are raised before values in the core are changed. See for example commit 8977e50.

Submodules

Internal libraries

Functionality that doesn't depend on the state of the system can be extracted from the EspressoCore target and stored in a submodule of src/. This is used for example by the EspressoUtils module, which is a dependency of EspressoCore.

Such submodules can adopt the Pitchfork layout sketched below. With this layout, passing the compiler flag -Isrc/library/include makes the public header file library/feature.hpp visible to consumers of that library while the implementation details in src/library/src/feature_impl.hpp remain hidden. 

src/
 \--- library/
       |--- include/
       |     \--- library/
       |           \--- feature.hpp
       |--- src/
       |     |--- feature_impl.hpp
       |     \--- feature_impl.cpp
       \--- tests/
             \--- feature_test.cpp

External libraries

External libraries can be added to the project after approval from the core team.

CMake

CMake-based projects that provide functionality to ESPResSo are included using the FetchContent method. This approach enables pinning the library version and patching the source code. 

if(ESPRESSO_BUILD_WITH_WALBERLA)
  include(FetchContent)
  FetchContent_Declare(
    walberla
    GIT_REPOSITORY https://i10git.cs.fau.de/walberla/walberla
    GIT_TAG        64a4d68
    PATCH_COMMAND  patch -p0 --ignore-whitespace -i
                   ${PROJECT_SOURCE_DIR}/cmake/waLBerlaFunctions.patch
  )
  FetchContent_Populate(walberla)
  set(WALBERLA_BUILD_TESTS off CACHE BOOL "")
  set(CMAKE_POSITION_INDEPENDENT_CODE on CACHE BOOL "")
  add_subdirectory("${walberla_SOURCE_DIR}" "${walberla_BINARY_DIR}")
  set(WALBERLA_LIBS walberla::core walberla::blockforest walberla::boundary walberla::field walberla::lbm)
endif()

Another method consists in loading the library from the system via environment variables, such as $PKG_CONFIG_PATH or $LD_LIBRARY_PATH. 

find_package(NumPy REQUIRED)
find_program(IPYTHON_EXECUTABLE NAMES jupyter ipython3 ipython)
pkg_check_modules(SCAFACOS scafacos)

C++ bridge

For better separation of concerns, the external library should not expose its public interface to the EspressoCore directly, but instead rely on a private submodule. With this approach, the external library is a private dependency of a submodule, and that submodule is a private dependency of EspressoCore. The submodule should use a C++ design pattern to hide the details of the external library.

Using the Pitchfork layout outlined in section Internal libraries as an example, one could declare an abstract base class in feature.hpp and derive a concrete class in feature_impl.hpp. An instance of the class would be stored in EspressoCore, however EspressoCore wouldn't have access to the header files of the external library. 

/* /usr/include/external_library/ExternalFeature.hpp */
class ExternalFeature {
public:
  ExternalFeature(int argument): m_argument(argument) {}
  double energy(int position) const;
private:
  int m_argument;
};

/* src/library/include/library/feature.hpp */
class FeatureBase {
public:
  virtual double get_energy(int position) const = 0;
};
FeatureBase *new_feature(int argument);

/* src/library/src/feature_impl.hpp */
#include <library/feature.hpp>
#include <ExternalFeature.hpp>

class Feature: public FeatureBase {
public:
  Feature(int argument): m_feature(ExternalFeature(argument)) {}
  double get_energy(int position) const override;
private:
  ExternalFeature m_feature;
};

/* src/library/src/feature_impl.cpp */
#include <library/feature.hpp>
#include "feature_impl.hpp"

FeatureBase *new_feature(int argument) {
  FeatureBase *instance = new Feature(argument);
  return instance;
}

double Feature::get_energy(int position) const {
  return m_feature.energy(position);
}

/* src/core/feature_instance.hpp */
#include <library/feature.hpp>

FeatureBase *feature();

/* src/core/feature_instance.cpp */
#include "communication.hpp"
#include "feature_instance.hpp"
#include <library/feature.hpp>
#include <stdexcept>

namespace {
FeatureBase *feature_instance = nullptr;
}

FeatureBase *feature() {
  if (!feature_instance)
    throw std::runtime_error("feature is uninitialized.");
  return feature_instance;
}

void init_feature_local(int argument) {
  feature_instance = new_feature(argument);
}

REGISTER_CALLBACK(init_feature_local)

void init_feature(int argument) {
  Communication::mpiCallbacks().call_all(init_feature_local, argument);
}

/* src/core/analysis.hpp */
double feature_energy(int position);

/* src/core/analysis.cpp */
#include "analysis.hpp"
#include "feature_instance.hpp"

double feature_energy(int position) {
  return feature()->get_energy(position);
}

Script interface

Exporting symbols from the core

C++ functions and classes from the ESPResSo core are exported in the Python interface using the src/python/espressomd/*.pxd files. These exported symbols can then be used within cdef functions in src/python/espressomd/*.pyx files. 

# interface.pxd
cdef extern from "cuda_init.hpp":
    int cuda_get_device()

# interface.pyx
from . cimport interface  # get symbols from the pxd file

cdef device(self):
    """Get device id."""
    dev = cuda_get_device()
    if dev == -1:
        raise Exception("cuda device get error")
    return dev

Likewise, struct and class declarations can be copy-pasted in .pxd files to allow returning a copy of a core object in a Cython function.

Using handles to core and script interface objects

Classes declared in the core can be instantiated in the core and managed by the script interface, or be directly instantiated in the script interface. This is achieved using "glue" classes in src/script_interface, which derive from AutoParameters. These classes set up getters and setters to members of the managed object (via add_parameters in the constructor) and provide a mechanism to forward function calls (via do_call_method()). Values are automatically converted from the python types (int, float, string, etc.) to the corresponding C++ types (int, double, std::string, etc.) and passed between the Python interface and script interface via a boost::variant.

Python classes derived from ScriptInterfaceHelper provide the Python interface. These classes name the AutoParameters class to instantiate (field _so_name), which will in turn instantiate the managed object. The Python class must also explicitly list which methods of the managed object are visible to the user (field _so_bind_methods). Python classes also control on which MPI nodes the AutoParameters object is instantiated (field _so_creation_policy): either on all nodes (GLOBAL creation policy, default) or only on the head node (LOCAL creation policy).

Below is an example with observable classes. The Observable parent class has a shape() method that can be called directly from an instance of the observable. The AutoParameters class has an additional calculate() method that is not meant to be used directly by the user, but can be called from within the class with call_method("calculate"). This approach is used when the result of a function needs to be post-processed, in this case to reshape the flat array into an N-dimensional array. The script interface class instantiates the managed object and provides a method observable() that returns a shared pointer to it. This is used by Accumulator classes to get a handle to the observable they are accumulating data from, and optionally forward it to the core global variable auto_update_accumulators (a vector of shared pointers to script interface observables) when the user calls system.auto_update_accumulators.add(my_accumulator). 

# observables.py
from .script_interface import ScriptInterfaceHelper, script_interface_register

@script_interface_register
class Observable(ScriptInterfaceHelper):
    """
    Base class for all observables.

    Methods
    -------
    shape()
        Return the shape of the observable.
    """
    _so_name = "Observables::Observable"
    _so_bind_methods = ("shape",)
    _so_creation_policy = "LOCAL"
    def calculate(self):
        return np.array(self.call_method("calculate")).reshape(self.shape())

@script_interface_register
class ComVelocity(Observable):
    """Calculates the center of mass velocity for particles with given ids."""
    _so_name = "Observables::ComVelocity"
/* src/script_interface/observables/initialize.hpp */
#include "script_interface/ObjectHandle.hpp"
#include <utils/Factory.hpp>

namespace ScriptInterface {
namespace Observables {
void initialize(Utils::Factory<ObjectHandle> *om);
}
}

/* src/script_interface/observables/initialize.cpp */
#include "script_interface/ScriptInterface.hpp"
#include "core/observables/Observable.hpp"
#include "core/observables/ComVelocity.hpp"
#include "initialize.hpp"
#include <utils/Factory.hpp>
#include <memory>

namespace ScriptInterface {
namespace Observables {

/** Base class for script interfaces to core Observables classes. */
class Observable : public ObjectHandle {
public:
  virtual std::shared_ptr<::Observables::Observable> observable() const = 0;
  Variant do_call_method(std::string const &method,
                         VariantMap const &parameters) override {
    if (method == "calculate")
      return observable()->operator()();
    if (method == "shape")
      return observable()->shape();
    return {};
  }
};

/** Base class for script interfaces to particle-based observables.
 *  @tparam CorePidObs Any particle-based observable core class.
 */
template <typename CorePidObs>
class PidObservable : public AutoParameters<PidObservable<CorePidObs>, Observable> {
public:
  PidObservable() {
    this->add_parameters({
      // AutoParameter object for CorePidObs::ids in the form {name, setter lambda, getter lambda}
      {"ids",
       AutoParameter::read_only,
       [this]() { return m_observable->ids(); }
      }
    });
  }
  void do_construct(VariantMap const &params) override {
    // instantiate the managed object in the script interface
    m_obs = make_shared_from_args<CorePidObs, std::vector<int>>(params, "ids");
  }
  std::shared_ptr<::Observables::Observable> observable() const override {
    return m_obs;
  }
private:
  // pointer to the managed object
  std::shared_ptr<CorePidObs> m_obs;
};

/** Register script interface bridge to core observable. */
void initialize(Utils::Factory<ObjectHandle> *om) {
  om->register_new<PidObservable<::Observables::ComVelocity>>("Observables::ComVelocity");
}

} /* namespace Observables */
} /* namespace ScriptInterface */

The AutoParameters framework

The AutoParameters framework provides tools to register members of the wrapped core class and expose them to the python script interface. It is possible to configure access by passing pointers to the core class setter/getter functions or by passing custom lambdas. Values can also be marked as read-only. This is achieved by constructing a list of AutoParameter objects in the add_parameters() function. The code sample below illustrates some of the AutoParameter constructor overloads (see the full list). 

/* file src/core/MyClass.hpp */
class MyClass {
private:
  int m_id;
  double m_radius;
  double m_height;
  double m_phi;

public:
  int &get_id_as_reference() { return m_id; }
  void set_radius(double const &radius) { m_radius = radius; }
  double get_radius() const { return m_radius; }
  void set_phi(double phi) { m_phi = phi; }
  double get_phi() const { return m_phi; }
  double const &get_height_const_ref() const { return m_height; }
};

/* file src/script_interface/MyClass.hpp */
#include "core/MyClass.hpp"
#include "script_interface/auto_parameters/AutoParameters.hpp"
#include <cmath>
#include <memory>

namespace ScriptInterface {
class MyClass : public AutoParameters<MyClass> {
private:
  std::shared_ptr<::MyClass> m_obj;
  constexpr static double deg2rad = M_PI / 180.;

public:
  MyClass() {
    add_parameters({
        // class member 'id' has a getter that returns a reference
        {
           "id",
           m_obj->get_id_as_reference()
        },
        // class member 'radius' has a getter and setter
        {
           "radius",
            m_obj,
            &::MyClass::set_radius,
            &::MyClass::get_radius
        },
        // class member access by lambda (e.g. for unit conversion)
        {
           "phi", 
           [this](Variant const &v) { m_obj->set_phi(get_value<double>(v) * deg2rad); },
           [this]() { return m_obj->get_phi() / deg2rad; }
        },
        // immutable class member
        {
           "height",
           AutoParameter::read_only,
           [this]() { return m_obj->get_height_const_ref(); }
        }
     });
  }
};
} // namespace ScriptInterface

Checkpointing

To enable checkpointing, Python classes need to define special methods to serialize and deserialize themselves. This is usually achieved by implementing adequate __getstate__(self) and __setstate__(self, params) methods. The result of the __getstate__ can be any data type, usually a list or dict, and will be written to a checkpoint file by the pickle module. When loading the system from a checkpoint file, pickle will read the value and pass it to the class as the params argument to the __setstate__ function to instantiate an object. During this instantiation, the __init__(self) method is not called.

Clone this wiki locally