35#include <boost/mpi/collectives/all_reduce.hpp>
36#include <boost/range/algorithm.hpp>
46 auto save_pos = [](
Particle &p) { p.pos_last_time_step() = p.pos(); };
48 boost::for_each(particles, save_pos);
49 boost::for_each(ghost_particles, save_pos);
60 auto reset_force = [](
Particle &p) { p.rattle_params().correction.fill(0); };
62 boost::for_each(particles, reset_force);
63 boost::for_each(ghost_particles, reset_force);
79 auto const r_ij2 = r_ij.norm2();
81 if (std::abs(1.0 - r_ij2 / ia_params.
d2) > ia_params.
p_tol) {
84 auto const r_ij_dot = r_ij_t * r_ij;
86 0.50 * (ia_params.
d2 - r_ij2) / r_ij_dot / (p1.
mass() + p2.
mass());
88 auto const pos_corr = G * r_ij_t;
106template <
typename Kernel>
110 bool correction =
false;
112 [&correction, &kernel, &box_geo](
Particle &p1,
int bond_id,
116 if (
auto const *bond = boost::get<RigidBond>(&iaparams)) {
117 auto const corrected = kernel(*bond, box_geo, p1, *partners[0]);
135 boost::for_each(particles, [](
Particle &p) {
154 boost::mpi::all_reduce(
comm_cart, repeat_, std::logical_or<bool>());
170 auto const resort_level =
190 auto const v_ij = p1.
v() - p2.
v();
193 auto const v_proj = v_ij * r_ij;
194 if (std::abs(v_proj) > ia_params.
v_tol) {
195 auto const K = v_proj / ia_params.
d2 / (p1.
mass() + p2.
mass());
197 auto const vel_corr = K * r_ij;
214 boost::for_each(particles,
230 boost::mpi::all_reduce(
comm_cart, repeat_, std::logical_or<bool>());
BondedInteractionsMap bonded_ia_params
Field containing the parameters of the bonded ia types.
Data structures for bonded interactions.
mapped_type at(key_type const &key) const
Utils::Vector< T, 3 > get_mi_vector(const Utils::Vector< T, 3 > &a, const Utils::Vector< T, 3 > &b) const
Get the minimum-image vector between two coordinates.
A stripped-down version of std::span from C++17.
boost::mpi::communicator comm_cart
The communicator.
#define SHAKE_MAX_ITERATIONS
Maximal number of iterations in the RATTLE algorithm before it bails out.
This file contains the errorhandling code for severe errors, like a broken bond or illegal parameter ...
#define runtimeErrorMsg()
@ DATA_PART_MOMENTUM
Particle::m.
@ DATA_PART_PROPERTIES
Particle::p.
@ DATA_PART_POSITION
Particle::r.
void save_old_position(const ParticleRange &particles, const ParticleRange &ghost_particles)
copy current position
static bool compute_correction_vector(CellStructure &cs, BoxGeometry const &box_geo, Kernel kernel)
Compute the correction vectors using given kernel.
void correct_velocity_shake(CellStructure &cs, BoxGeometry const &box_geo)
Correction of current velocities using RATTLE algorithm.
static void init_correction_vector(const ParticleRange &particles, const ParticleRange &ghost_particles)
reset correction vectors to zero
void correct_position_shake(CellStructure &cs, BoxGeometry const &box_geo)
Propagate velocity and position while using SHAKE algorithm for bond constraint.
static void apply_positional_correction(const ParticleRange &particles)
Apply positional corrections.
static bool calculate_velocity_correction(RigidBond const &ia_params, BoxGeometry const &box_geo, Particle &p1, Particle &p2)
Calculate the velocity correction for the particles.
static void apply_velocity_correction(ParticleRange const &particles)
Apply velocity corrections.
static bool calculate_positional_correction(RigidBond const &ia_params, BoxGeometry const &box_geo, Particle &p1, Particle &p2)
Calculate the positional correction for the particles.
Definition of the rigid bond data type for the Rattle algorithm.
Describes a cell structure / cell system.
ParticleRange ghost_particles() const
void update_ghosts_and_resort_particle(unsigned data_parts)
Update ghost particles, with particle resort if needed.
void ghosts_update(unsigned data_parts)
Update ghost particles.
bool check_resort_required(Utils::Vector3d const &additional_offset={}) const
Check whether a particle has moved further than half the skin since the last Verlet list update,...
void bond_loop(BondKernel const &bond_kernel)
Bonded pair loop.
void set_resort_particles(Cells::Resort level)
Increase the local resort level at least to level.
ParticleRange local_particles() const
void ghosts_reduce_rattle_correction()
Add rattle corrections from ghost particles to real particles.
Struct holding all information for one particle.
auto const & rattle_params() const
auto const & mass() const
auto const & pos_last_time_step() const
Parameters for the rigid_bond/SHAKE/RATTLE ALGORITHM.
double d2
Square of the length of Constrained Bond.
double v_tol
Velocity Tolerance/Accuracy for termination of RATTLE/SHAKE iterations during velocity corrections.
double p_tol
Positional Tolerance/Accuracy value for termination of RATTLE/SHAKE iterations during position correc...