Rewrite collision detection interface

[jngrad]

The collision detection offers several protocols for binding particles together. Some of these protocols rely on specific features, such as VIRTUAL_SITES and VIRTUAL_SITES_RELATIVE. These protocols have different parameters, although with some overlap. For historical reasons, these parameters have been packed into a single global struct Collision_parameters, and there is only one class constructor in C++ and in Cython (the system.collision_detection.set_params() method), where all parameters are mixed.

This design decision has several major drawbacks:

• we need to manually curate a python dict or a C++ std::map containing for each protocol the set of parameters that are valid
  • every time we set a new protocol from the Python script, we need to check the input parameters against that dictionary
  • every time we get the list of active parameters from the Python script, the complete parameter struct is retrieved by the Python interface and stored dict, then all parameters that are irrelevant to the currently active protocol are deleted from the dict
• the struct parameters are marked as read/write in the script interface, but the Cython interface had to disable the __setattr__() method to make them read-only, since we cannot guarantee that changing one parameter will leave the struct in a valid state
• to add a new protocol, one has to modify code that affects all existing protocols, because their code have been interleaved in common functions with if/else acting on the enum value that encodes which protocol is currently active
• calling set_params() with invalid arguments overwrites the global state of the feature with invalid values that throw an exception, but leaves the system in an indeterminate state until the next call of set_params() with valid arguments (fixed by fe12f7425d32aa48ad712c6b464f0fb672bb506b)

The collision detection interface could be redesigned by splitting the global struct Collision_parameters into individual structs, one for each protocol. This would make them completely independent from one another, would significantly reduce the complexity of the C++ script interface class, remove the need to curate the list of valid parameters in a std::map or dict, and allow the user to instantiate multiple collision detection objects in the same pypresso script (only one of them would be activate at a time). The bond breakage feature can serve as a blueprint for the new collision detection interface.

[jngrad]

The collision detection feature is explained in more details in the user guide section Creating bonds when particles collide. We would like to refactor the feature to make it easier to extend with new methods. @christophlohrmann has a research project where particles, upon collision with a wall, remain bound to the wall by a newly created particle inside the wall. This has applications in bacterial biofilm formation, but could also be used to model soot formation on surfaces, or binding of organic molecules to metallic catalysts. Offline discussions with @RudolfWeeber also lead to the decision to remove the angle bond formation method.

The task can be broken down as follows:

1. remove the bind_three_particles collision method
   1. use textual search to find and erase all occurrences of this method in the following files:
      • src/core/collision.hpp
      • src/core/collision.cpp
      • src/script_interface/collision_detection/CollisionDetection.hpp
      • src/python/espressomd/collision_detection.py
      • testsuite/python/collision_detection.py
      • testsuite/python/collision_detection_interface.py
      • doc/sphinx/advanced_methods.rst
   2. use compiler warnings about unused variables to erase all parameters of that method that are not used in other methods
   3. use cmake . -D ESPRESSO_CTEST_ARGS=-j$(nproc) once to set the project to use all CPU cores during testing
   4. use make check_python_skip_long to run the testsuite and make sure the method is no longer used
   5. optionally, use code coverage analysis to verify there is no leftover code
2. write one parameter struct per collision method
   1. split Collision_parameters into multiple structs CollisionOff, CollisionBindCenters, CollisionBindVS, CollisionGlueToSurf, containing only the subset of parameters that are needed (the compiler will notify you if a parameter is missing), and split the Collision_parameters::initialize() method accordingly
   2. change the type of the global variable collision_params from Collision_parameters to a std::variant<CollisionOff,CollisionBindCenters,CollisionBindVS,CollisionGlueToSurf>
   3. use std::visit(visitor, collision_params) to execute a collision method
      • here visitor is an instance of a struct or class that carries out an operation using the underlying object in the collision_params (it's passed as first argument to the overloaded functions of the visitor)
      • for details please refer to the Boost tutorial or look for examples in the ESPResSo C++ core
      • we use std::variant and std::visit rather than the corresponding Boost functionality, because variants are now part of the C++ standard library
   4. to adapt the script interface, use src/script_interface/lees_edwards/LeesEdwards.hpp, src/script_interface/lees_edwards/initialize.cpp and src/python/espressomd/lees_edwards.py as templates

[RudolfWeeber]

Results of offline discussion

• the detection of collisions is in the hot-loop, i.e., this part is performance critical
• Therefore, there, only particle types and distances should be checked
• The actual action should be handled when processing the queue

Rough sketch of final design

I think, the final interface should look as follows

Python

center_bond = ...
vs_bond = ...

system.bonded_inter.add(c_bond)
system.bonded_inter.add(vs_bond)
bind_at_point_of_collision = espressomd.collision_detection.BindAtPointOfCollision(
    vs_placement=0.5, center_bond=center_bond, vs_bond=vs_bond)
change_types = espressomd.collision_detection.ChangeTypes(from=(0,0), to=(1,1))
bind_centers = espressomd.collision_detection.BindCenters(bond=center_bond)

system.collision_detection.add(
  p_types=(0,0), distance=1.1*sigma,
  actions=[bind_at_point_of_collision,change_particle_types])
)

or p_types=(-1,-1) (or None in Python to match all particle types).

C++

I.e., on the C++ side, this is a

struct CollisionRule {
   std::pair<size_t,size_t> particle_types;
   double distance;
   std::vector<CollisionAction> actions;
}

std::vector<CollisionRule> collision_rules;

struct Collision {
   int p_id1;
   int p_id2;
   std::vector<CollisionAction> actions;
}

std::vector<Collision> collision_queue;

The individual actions should support (on the C++ side)

• a function bool are_bound(Particle const&, Particle const&) which determines if two particle are already bound. This would typically be done by checking the presence of the center bond
• a function bool run(CellStructure&, Particle&, Particle&) which performs the action

the detect_collision(CellStructure&, Particle&, Particle& ) called in the hot loop would

• compare the particle types and distance between the particles against the parameters in collision_rules
• if the particle pair matches one or more rules, call the are_bound() function on the pair. If that returns false, queue the pair

Helper function

bool types_match(std::pair<size_t,size_t> target_types, size_t type1, size_t type2)

• returns true if the types in type1 and type2 match the types in target_types (which come from the collision_rule)
• Note that particles can come in in any order
• Special handling of -1, i.e., ignore type
• This function can be written and tested independently from the rest

Implementation in steps

So intermediate steps could be:

1. delete the angle bond method
2. keep the cutoff and one CollisionAction (bind_centers, bind_at_point_of_collision, etc.) as global parameter
3. but split the parameters for bind_centers (center bond), bind_at_point_of_collision (center bond, vs bond, vs placement), glue to surface (bonds, vs placement, new particle types) into separate structs
4. write the particle type filtering function
5. implement the collision_rules vector and use it in detect_collision()

This can be merged in steps:

• step 1.
• step 2.+3.
• step 4.
• step 5.

[jngrad]

Here are the updated instructions:

1. extract methods CollisionDetection::initialize() and CollisionDetection::handle_collisions() and the associated data members into a self-contained struct called MergedProtocol, of which a shared pointer will be stored in CollisionDetection as std::shared_ptr<MergedProtocol> m_protocol
   • the source code is in src/core/collision.hpp and src/core/collision.cpp
2. split the struct MergedProtocol into multiple structs called Off, BindCenters, BindVS, GlueToSurf, and store them inside a namespace Collision in a new file src/core/collision/protocols.hpp; store the CollisionDetection class inside a namespace Collision into a new file src/core/collision/collision.hpp
   1. these structs should contain only the subset of parameters they need (the compiler will notify you if a data member is missing or is unused)
   2. split the initialize() and handle_collisions() methods accordingly
   3. define using ActiveProtocol = std::variant<Off,BindCenters,BindVS,GlueToSurf>
   4. change the type of the shared pointer CollisionDetection::m_protocol from std::shared_ptr<MergedProtocol> to std::shared_ptr<ActiveProtocol>
   5. use std::visit(visitor, m_protocol) to execute a collision protocol
      • here visitor is an instance of a struct or class that carries out an operation using the underlying object in the m_protocol (it's passed as first argument to the overloaded functions of the visitor)
      • for details please refer to the Boost tutorial or look for examples in the ESPResSo C++ core
      • we use std::variant and std::visit rather than the corresponding Boost functionality, because variants are now part of the C++ standard library
      • use src/core/lees_edwards/lees_edwards.hpp and src/core/lees_edwards/protocols.hpp as a guide
3. adapt the script interface to the new C++ API
   • use src/script_interface/lees_edwards/LeesEdwards.hpp, src/script_interface/lees_edwards/initialize.cpp and src/python/espressomd/lees_edwards.py as a guide