74#include <boost/optional.hpp>
75#include <boost/variant.hpp>
85 auto force_factor = 0.;
95#ifdef LENNARD_JONES_GENERIC
142 pf.f += force_factor * d;
216 if (dist < ia_params.
max_cut) {
235 auto const q1q2 = p1.
q() * p2.
q();
236 if (q1q2 != 0. and coulomb_kernel !=
nullptr) {
237 pf.f += (*coulomb_kernel)(q1q2, d, dist);
240 (*elc_kernel)(p1, p2, q1q2);
261 ia_params, d, dist, dist2);
273 if (dipoles_kernel) {
274 pf += (*dipoles_kernel)(p1, p2, d, dist, dist2);
298 if (
auto const *iap = boost::get<FeneBond>(&iaparams)) {
299 return iap->force(dx);
301 if (
auto const *iap = boost::get<HarmonicBond>(&iaparams)) {
302 return iap->force(dx);
304 if (
auto const *iap = boost::get<QuarticBond>(&iaparams)) {
305 return iap->force(dx);
308 if (
auto const *iap = boost::get<BondedCoulomb>(&iaparams)) {
309 return iap->force(p1.
q() * p2.
q(), dx);
311 if (
auto const *iap = boost::get<BondedCoulombSR>(&iaparams)) {
312 return iap->force(dx, *kernel);
315#ifdef BOND_CONSTRAINT
316 if (boost::get<RigidBond>(&iaparams)) {
321 if (
auto const *iap = boost::get<TabulatedDistanceBond>(&iaparams)) {
322 return iap->force(dx);
325 if (boost::get<VirtualBond>(&iaparams)) {
337 if (
auto const *iap = boost::get<ThermalizedBond>(&iaparams)) {
338 auto result = iap->forces(p1, p2, dx);
340 auto const &forces = result.get();
342 p1.
force() += std::get<0>(forces);
343 p2.
force() += std::get<1>(forces);
350 p1.
force() += result.get();
351 p2.
force() -= result.get();
362inline boost::optional<
363 std::tuple<Utils::Vector3d, Utils::Vector3d, Utils::Vector3d>>
369 if (
auto const *iap = boost::get<AngleHarmonicBond>(&iaparams)) {
370 return iap->forces(vec1, vec2);
372 if (
auto const *iap = boost::get<AngleCosineBond>(&iaparams)) {
373 return iap->forces(vec1, vec2);
375 if (
auto const *iap = boost::get<AngleCossquareBond>(&iaparams)) {
376 return iap->forces(vec1, vec2);
379 if (
auto const *iap = boost::get<TabulatedAngleBond>(&iaparams)) {
380 return iap->forces(vec1, vec2);
383 if (
auto const *iap = boost::get<IBMTriel>(&iaparams)) {
384 return iap->calc_forces(vec1, vec2);
393 if (boost::get<OifGlobalForcesBond>(&iaparams)) {
399 auto const &forces = result.get();
401 p1.
force() += std::get<0>(forces);
402 p2.
force() += std::get<1>(forces);
403 p3.
force() += std::get<2>(forces);
416 if (
auto const *iap = boost::get<OifLocalForcesBond>(&iaparams)) {
417 return iap->calc_forces(box_geo, p1, p2, p3, p4);
419 if (
auto const *iap = boost::get<IBMTribend>(&iaparams)) {
420 return iap->calc_forces(box_geo, p1, p2, p3, p4);
426 if (
auto const *iap = boost::get<DihedralBond>(&iaparams)) {
427 return iap->forces(v12, v23, v34);
430 if (
auto const *iap = boost::get<TabulatedDihedralBond>(&iaparams)) {
431 return iap->forces(v12, v23, v34);
444 auto const &forces = result.get();
446 p1.
force() += std::get<0>(forces);
447 p2.
force() += std::get<1>(forces);
448 p3.
force() += std::get<2>(forces);
449 p4.
force() += std::get<3>(forces);
464 if (partners.
size() == 1) {
467 p1.
id(), {{partners[0]->id(), std::nullopt}}, bond_id, d)) {
471 if (partners.size() == 2) {
473 box_geo.get_mi_vector(partners[0]->
pos(), partners[1]->
pos()).norm();
474 if (bond_breakage.check_and_handle_breakage(
475 p1.id(), {{partners[0]->id(), partners[1]->id()}}, bond_id, d)) {
493 *partners[1], *partners[2]);
Vector implementation and trait types for boost qvm interoperability.
__shared__ int pos[MAXDEPTH *THREADS5/WARPSIZE]
Routines to calculate the Born-Meyer-Huggins-Tosi-Fumi potential between particle pairs.
double BMHTF_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate BMHTF force factor.
BondedInteractionsMap bonded_ia_params
Field containing the parameters of the bonded ia types.
Data structures for bonded interactions.
int number_of_partners(Bonded_IA_Parameters const &iaparams)
Return the number of bonded partners for the specified bond.
boost::variant< NoneBond, FeneBond, HarmonicBond, QuarticBond, BondedCoulomb, BondedCoulombSR, AngleHarmonicBond, AngleCosineBond, AngleCossquareBond, DihedralBond, TabulatedDistanceBond, TabulatedAngleBond, TabulatedDihedralBond, ThermalizedBond, RigidBond, IBMTriel, IBMVolCons, IBMTribend, OifGlobalForcesBond, OifLocalForcesBond, VirtualBond > Bonded_IA_Parameters
Variant in which to store the parameters of an individual bonded interaction.
Routines to calculate the Buckingham potential between particle pairs.
double buck_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate Buckingham force factor.
bool check_and_handle_breakage(int particle_id, BondPartners const &bond_partners, int bond_type, double distance)
Check if the bond between the particles should break, if yes, queue it.
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.
std::shared_ptr< DPDThermostat > dpd
int thermo_switch
Bitmask of currently active thermostats.
A stripped-down version of std::span from C++17.
DEVICE_QUALIFIER constexpr size_type size() const
This file contains the defaults for ESPResSo.
Routines to calculate the Gaussian potential between particle pairs.
double gaussian_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate Gaussian force factor.
__device__ void vector_product(float const *a, float const *b, float *out)
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 .
This file contains the errorhandling code for severe errors, like a broken bond or illegal parameter ...
bool do_nonbonded(Particle const &p1, Particle const &p2)
Determine if the non-bonded interactions between p1 and p2 should be calculated.
void npt_add_virial_force_contribution(const Utils::Vector3d &force, const Utils::Vector3d &d)
Update the NpT virial.
boost::optional< std::tuple< Utils::Vector3d, Utils::Vector3d, Utils::Vector3d > > calc_bonded_three_body_force(Bonded_IA_Parameters const &iaparams, BoxGeometry const &box_geo, Particle const &p1, Particle const &p2, Particle const &p3)
bool add_bonded_force(Particle &p1, int bond_id, Utils::Span< Particle * > partners, BondBreakage::BondBreakage &bond_breakage, BoxGeometry const &box_geo, Coulomb::ShortRangeForceKernel::kernel_type const *kernel)
ParticleForce calc_central_radial_force(IA_parameters const &ia_params, Utils::Vector3d const &d, double const dist)
void add_non_bonded_pair_force(Particle &p1, Particle &p2, Utils::Vector3d const &d, double dist, double dist2, IA_parameters const &ia_params, Thermostat::Thermostat const &thermostat, BoxGeometry const &box_geo, Coulomb::ShortRangeForceKernel::kernel_type const *coulomb_kernel, Dipoles::ShortRangeForceKernel::kernel_type const *dipoles_kernel, Coulomb::ShortRangeForceCorrectionsKernel::kernel_type const *elc_kernel)
Calculate non-bonded forces between a pair of particles and update their forces and torques.
ParticleForce calc_opposing_force(ParticleForce const &pf, Utils::Vector3d const &d)
bool add_bonded_two_body_force(Bonded_IA_Parameters const &iaparams, Particle &p1, Particle &p2, BoxGeometry const &box_geo, Coulomb::ShortRangeForceKernel::kernel_type const *kernel)
boost::optional< std::tuple< Utils::Vector3d, Utils::Vector3d, Utils::Vector3d, Utils::Vector3d > > calc_bonded_four_body_force(Bonded_IA_Parameters const &iaparams, BoxGeometry const &box_geo, Particle const &p1, Particle const &p2, Particle const &p3, Particle const &p4)
bool add_bonded_three_body_force(Bonded_IA_Parameters const &iaparams, BoxGeometry const &box_geo, Particle &p1, Particle &p2, Particle &p3)
boost::optional< Utils::Vector3d > calc_bond_pair_force(Particle const &p1, Particle const &p2, Bonded_IA_Parameters const &iaparams, Utils::Vector3d const &dx, Coulomb::ShortRangeForceKernel::kernel_type const *kernel)
Compute the bonded interaction force between particle pairs.
ParticleForce calc_central_radial_charge_force(Particle const &p1, Particle const &p2, IA_parameters const &ia_params, Utils::Vector3d const &d, double const dist, Coulomb::ShortRangeForceKernel::kernel_type const *coulomb_kernel)
bool add_bonded_four_body_force(Bonded_IA_Parameters const &iaparams, BoxGeometry const &box_geo, Particle &p1, Particle &p2, Particle &p3, Particle &p4)
ParticleForce calc_non_central_force(Particle const &p1, Particle const &p2, IA_parameters const &ia_params, Utils::Vector3d const &d, double const dist)
Routines to calculate the Gay-Berne potential between particle pairs.
ParticleForce gb_pair_force(Utils::Quaternion< double > const &qi, Utils::Quaternion< double > const &qj, IA_parameters const &ia_params, Utils::Vector3d const &d, double dist)
Calculate Gay-Berne force and torques.
Routines to calculate the hat potential between particle pairs.
double hat_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate hat force factor.
Routines to calculate the Hertzian potential between particle pairs.
double hertzian_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate Hertzian force factor.
Routines to calculate the Lennard-Jones potential between particle pairs.
double lj_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate Lennard-Jones force factor.
Routines to calculate the Lennard-Jones with cosine tail potential between particle pairs.
double ljcos2_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate Lennard-Jones cosine squared force factor.
Routines to calculate the Lennard-Jones+cosine potential between particle pairs.
double ljcos_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate Lennard-Jones cosine force factor.
Routines to calculate the generalized Lennard-Jones potential between particle pairs.
double ljgen_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate Lennard-Jones force factor.
Routines to calculate the Morse potential between particle pairs.
double morse_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate Morse force factor.
Various procedures concerning interactions between particles.
Routines to calculate the energy and/or force for particle pairs via interpolation of lookup tables.
double tabulated_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate a non-bonded pair force factor by linear interpolation from a table.
Routines to calculate the OIF global forces energy or/and and force for a particle triple (triangle f...
Routines to calculate the OIF local forces for a particle quadruple (two neighboring triangles with c...
Routines to calculate the smooth step potential between particle pairs.
double SmSt_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate smooth step force factor.
Routines to calculate the soft-sphere potential between particle pairs.
double soft_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate soft-sphere force factor.
Exception indicating that a bond with an unexpected number of partners was encountered.
Exception indicating that a bond type was unknown.
Solver::ShortRangeForceCorrectionsKernel kernel_type
Solver::ShortRangeForceKernel kernel_type
Solver::ShortRangeForceKernel kernel_type
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 & quat() const
auto const & force_and_torque() const
auto const & force() const
Routines to calculate the Thole damping potential between particle pairs.
Utils::Vector3d thole_pair_force(Particle const &p1, Particle const &p2, IA_parameters const &ia_params, Utils::Vector3d const &d, double dist, Coulomb::ShortRangeForceKernel::kernel_type const *kernel)
Calculate Thole force.
Routines to calculate the Weeks-Chandler-Andersen potential between particle pairs.
double wca_pair_force_factor(IA_parameters const &ia_params, double dist)
Calculate WCA force factor.