Skip to content

Software design

Jump to bottom
Jean-Noël Grad edited this page Jan 27, 2021 · 15 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

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;
}

Functions called by MPI callbacks cannot throw exceptions. Instead they should use the runtime errors logging mechanism to store a string containing 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")

Submodules

Internal libraries

Functionality that doesn't depend on the state of the system can be extracted from the EspressoCore target and store in a submodule of src/. This is used for example by the EspressoUtils module, which is a dependency of EspressoCore. Such submodules typically adopt the Pitchfork layout, i.e. 

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

With this layout, passing the compiler flag -Isrc/library/include will make the public header file library/feature.hpp visible to consumers of that library while the implementation details in src/library/src remain private.

External libraries

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(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(WITH_WALBERLA)

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, 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 objects

When exposing C++ classes in Python, it is sometimes necessary to write a bridge between the core implementation and the script interface. This is required when the Python interface needs to access the original core object through a pointer instead of working on a local copy of the core object. This is achieved by the "glue" classes in src/script_interface, which derive from AutoParameters. Corresponding Python classes derived from ScriptInterfaceHelper can then be written. They name the AutoParameters class to instantiate (field _so_name), which will in turn instantiate the core object and store a shared pointer to it. The Python class must also explicitly list the subset of getters and setters (field _so_bind_methods) that are visible to the user. 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 core via a boost::variant.

A class derived from AutoParameters stores a shared pointer to a core object and sets up getters and optionally setters to members of the core object. The AutoParameters constructor takes a vector of AutoParameter objects. The AutoParameter is simply a wrapper for getters and setters constructed as lambdas that deal with the boost::variant.

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. 

# 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.cpp */
#include <script_interface/ObjectHandle.hpp>
#include <utils/Factory.hpp>

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

/* src/script_interface/observables/initialize.cpp */
#include "script_interface/ScriptInterface.hpp"
#include "core/observables/Observable.hpp"
#include "core/observables/ComVelocity.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 (name, setter, getter) */
      {"ids",
       [this](Variant const &v) { m_obs->ids() = get_value<std::vector<int>>(v); },
       [this]() { return m_obs->ids(); }
      }
    });
  }
  void do_construct(VariantMap const &params) override {
    m_obs = make_shared_from_args<CorePidObs, std::vector<int>>(params, "ids");
  }
  std::shared_ptr<::Observables::Observable> observable() const override {
    return m_obs;
  }
private:
  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 */

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