25#include "communication.hpp"
32#include "system/System.hpp"
33#include "thermostat.hpp"
37#include <boost/mpi/collectives.hpp>
52 m_system.lock()->nonbonded_ias->make_particle_type_exist(
type);
55IA_parameters const &ShapeBasedConstraint::get_ia_param(
int type)
const {
56 return m_system.lock()->nonbonded_ias->get_ia_param(
type, part_rep.
type());
60 return all_reduce(
comm_cart, m_local_force, std::plus<>());
64 return all_reduce(
comm_cart, m_outer_normal_force, std::plus<double>());
69 auto global_mindist = std::numeric_limits<double>::infinity();
71 auto const local_mindist = std::accumulate(
72 particles.begin(), particles.end(),
73 std::numeric_limits<double>::infinity(),
74 [
this, &box_geo](
double min,
Particle const &p) {
75 auto const &ia_params = get_ia_param(p.type());
76 if (is_active(ia_params)) {
79 m_shape->calculate_dist(box_geo.folded_position(p.pos()), dist, vec);
80 return std::min(min, dist);
84 boost::mpi::reduce(
comm_cart, local_mindist, global_mindist,
85 boost::mpi::minimum<double>(), 0);
86 return global_mindist;
93 auto const &ia_params = get_ia_param(p.
type());
98 m_shape->calculate_dist(folded_pos, dist, dist_vec);
99 auto &system = *m_system.lock();
100 auto const coulomb_kernel = system.coulomb.pair_force_kernel();
108 outer_normal_vec = -dist_vec / dist;
112 *system.bonded_ias,
get_ptr(coulomb_kernel)) +
117 if (system.thermostat->thermo_switch &
THERMO_DPD) {
119 *system.box_geo, ia_params, dist_vec, dist,
122 system.thermostat->dpd->rng_increment();
125 }
else if (m_penetrable && (dist <= 0)) {
126 if ((!m_only_positive) && (dist < 0)) {
130 *system.bonded_ias,
get_ptr(coulomb_kernel)) +
135 if (system.thermostat->thermo_switch &
THERMO_DPD) {
137 *system.box_geo, ia_params, dist_vec, dist,
140 system.thermostat->dpd->rng_increment();
155 m_local_force -= pf.f;
156 m_outer_normal_force -= outer_normal_vec * pf.f;
161void ShapeBasedConstraint::add_energy(
const Particle &p,
165 auto const &ia_params = get_ia_param(p.
type());
168 auto &system = *m_system.lock();
169 auto const coulomb_kernel = system.coulomb.pair_energy_kernel();
172 m_shape->calculate_dist(folded_pos, dist, vec);
177 }
else if (dist <= 0. and m_penetrable) {
178 if (!m_only_positive and dist < 0.) {
188 if (part_rep.type() >= 0) {
190 p.
type(), part_rep.type(), p.
mol_id(), part_rep.mol_id(), energy);
Vector implementation and trait types for boost qvm interoperability.
Data structures for bonded interactions.
double total_normal_force() const
double min_dist(BoxGeometry const &box_geo, ParticleRange const &particles) const
Calculate the minimum distance between all particle pairs.
Utils::Vector3d total_force() const
Observable for the pressure and energy.
void add_non_bonded_contribution(int type1, int type2, int molid1, int molid2, std::span< const double > data)
boost::mpi::communicator comm_cart
The communicator.
This file contains the defaults for ESPResSo.
const T * get_ptr(std::optional< T > const &opt)
Utils::Vector3d dpd_pair_force(DPDParameters const ¶ms, Utils::Vector3d const &v, double dist, Utils::Vector3d const &noise)
Routines to use DPD as thermostat or pair force .
double calc_non_bonded_pair_energy(Particle const &p1, Particle const &p2, IA_parameters const &ia_params, Utils::Vector3d const &d, double const dist, BondedInteractionsMap const &bonded_ias, Coulomb::ShortRangeEnergyKernel::kernel_type const *coulomb_kernel)
Calculate non-bonded energies between a pair of particles.
This file contains the errorhandling code for severe errors, like a broken bond or illegal parameter ...
#define runtimeErrorMsg()
ParticleForce calc_central_radial_force(IA_parameters const &ia_params, Utils::Vector3d const &d, double const dist)
ParticleForce calc_opposing_force(ParticleForce const &pf, Utils::Vector3d const &d)
ParticleForce calc_non_central_force(Particle const &p1, Particle const &p2, IA_parameters const &ia_params, Utils::Vector3d const &d, double const dist)
static bool is_active(IA_parameters const &data)
Check if a non-bonded interaction is defined.
Various procedures concerning interactions between particles.
constexpr double INACTIVE_CUTOFF
Cutoff for deactivated interactions.
Parameters for non-bonded interactions.
double max_cut
maximal cutoff for this pair of particle types.
Force information on a particle.
Utils::Vector3d torque
torque.
Struct holding all information for one particle.
auto const & type() const
auto const & mol_id() const
Utils::Vector3d thole_pair_force(Particle const &p1, Particle const &p2, IA_parameters const &ia_params, Utils::Vector3d const &d, double dist, BondedInteractionsMap const &bonded_ias, Coulomb::ShortRangeForceKernel::kernel_type const *kernel)
Calculate Thole force.