diff --git a/cpp/particles/p2p/p2p_computer.hpp b/cpp/particles/p2p/p2p_computer.hpp
index 47d193e095403d329cd704e8623567c0f080a39b..ec8b2bcd35d765939a04d178ed65e3b60888a245 100644
--- a/cpp/particles/p2p/p2p_computer.hpp
+++ b/cpp/particles/p2p/p2p_computer.hpp
@@ -71,11 +71,16 @@ public:
template <int size_particle_positions, int size_particle_rhs>
void compute_interaction(const partsize_t /*idx_part1*/,
- const real_number pos_part1[], real_number rhs_part1[],
+ const real_number pos_part1[],
+ real_number rhs_part1[],
const partsize_t /*idx_part2*/,
- const real_number pos_part2[], real_number rhs_part2[],
- const real_number dist_pow2, const real_number cutoff,
- const real_number /*xshift_coef*/, const real_number /*yshift_coef*/, const real_number /*zshift_coef*/) const{
+ const real_number pos_part2[],
+ real_number rhs_part2[],
+ const real_number dist_pow2,
+ const real_number cutoff,
+ const real_number /*xseparation*/,
+ const real_number /*yseparation*/,
+ const real_number /*zseparation*/) const{
static_assert(size_particle_positions == 6, "This kernel works only with 6 values for one particle's position+orientation");
static_assert(size_particle_rhs == 6, "This kernel works only with 6 values per particle's rhs");
diff --git a/cpp/particles/p2p/p2p_computer_empty.hpp b/cpp/particles/p2p/p2p_computer_empty.hpp
index ce4e282c3bdb7ff4bfc479b5a0b16d106f0492c1..03d8cd64995ed6a190f1d444d2846d230985dec8 100644
--- a/cpp/particles/p2p/p2p_computer_empty.hpp
+++ b/cpp/particles/p2p/p2p_computer_empty.hpp
@@ -44,8 +44,11 @@ public:
const real_number /*pos_part1*/[], real_number /*rhs_part1*/[],
const partsize_t /*idx_part2*/,
const real_number /*pos_part2*/[], real_number /*rhs_part2*/[],
- const real_number /*dist_pow2*/, const real_number /*cutoff*/,
- const real_number /*xshift_coef*/, const real_number /*yshift_coef*/, const real_number /*zshift_coef*/) const{
+ const real_number /*dist_pow2*/,
+ const real_number /*cutoff*/,
+ const real_number /*xseparation*/,
+ const real_number /*yseparation*/,
+ const real_number /*zseparation*/) const{
}
void merge(const p2p_computer_empty&){}
diff --git a/cpp/particles/p2p/p2p_ghost_collisions.hpp b/cpp/particles/p2p/p2p_ghost_collisions.hpp
index b89fe8574d20620cfca8320b4f9fbebc0ef85405..13a7fb842b010f29eceb74be803738f35cf27f3b 100644
--- a/cpp/particles/p2p/p2p_ghost_collisions.hpp
+++ b/cpp/particles/p2p/p2p_ghost_collisions.hpp
@@ -216,11 +216,11 @@ public:
const partsize_t idx_part2,
const real_number pos_part2[],
real_number /*rhs_part2*/[],
- const real_number /*dist_pow2*/,
+ const real_number dist_pow2,
const real_number /*cutoff*/,
- const real_number /*xshift_coef*/,
- const real_number /*yshift_coef*/,
- const real_number /*zshift_coef*/){
+ const real_number xseparation,
+ const real_number yseparation,
+ const real_number zseparation){
switch(this->current_particle_shape)
{
case SPHERE:
@@ -233,19 +233,14 @@ public:
break;
case CYLINDER:
{
- double x, y, z, pq, xp, xq, t, s, x_dist, y_dist, z_dist ,min_distance;
+ double pq, xp, xq, t, s, x_dist, y_dist, z_dist ,min_distance;
- /* Relative position vector of the two particles (x,y,z)^T */
- /* complicated usage of fmod, fabs and sgn because C++ fmod does not do what it should. */
- x = std::fmod(pos_part2[0],pi2)-std::fmod(pos_part1[0],pi2);
- y = std::fmod(pos_part2[1],pi2)-std::fmod(pos_part1[1],pi2);
- z = std::fmod(pos_part2[2],pi2)-std::fmod(pos_part1[2],pi2);
+ const double x = xseparation;
+ const double y = yseparation;
+ const double z = zseparation;
+ const double r = sqrt(dist_pow2);
- x = ( std::fmod( std::fabs(x) +pi ,pi2) - pi ) * sgn(x) ;
- y = ( std::fmod( std::fabs(y) +pi ,pi2) - pi ) * sgn(y) ;
- z = ( std::fmod( std::fabs(z) +pi ,pi2) - pi ) * sgn(z) ;
-
- if( sqrt(x*x+y*y+z*z) <= this->cylinder_length )
+ if( r <= this->cylinder_length )
{
/* p and q are the orientation vectors of the first and second particles. */
@@ -384,7 +379,9 @@ public:
D1 = sqrt(x1x1+x0x0+t*t-2.0*x0x1-2.0*x1p0*t);
D2 = sqrt(x2x2+x0x0+t*t-2.0*x0x2-2.0*x2p0*t);
/* Check if a collision is possible. Exit if true. */
- if( D1<=this->disk_width and D2<=this->disk_width ){
+ if ((D1 <= (this->disk_width/2)) and
+ (D2 <= (this->disk_width/2)))
+ {
std::pair <partsize_t, partsize_t> single_collision_pair(idx_part1, idx_part2);
this->collision_pairs_local.insert(single_collision_pair);
assert(idx_part1!=idx_part2);
@@ -397,7 +394,9 @@ public:
D1 = sqrt(x1x1+x0x0+t*t-2.0*x0x1-2.0*x1p0*t);
D2 = sqrt(x2x2+x0x0+t*t-2.0*x0x2-2.0*x2p0*t);
/* Check if a collision is possible. Exit if true. */
- if( D1<=this->disk_width and D2<=this->disk_width ){
+ if ((D1 <= (this->disk_width/2)) and
+ (D2 <= (this->disk_width/2)))
+ {
std::pair <partsize_t, partsize_t> single_collision_pair(idx_part1, idx_part2);
this->collision_pairs_local.insert(single_collision_pair);
assert(idx_part1!=idx_part2);
diff --git a/cpp/particles/p2p/p2p_merge_collisions.hpp b/cpp/particles/p2p/p2p_merge_collisions.hpp
index 96fff9c2d7f6fa9f23d21ff346f1b8b819c9c1d5..cf93f3e8001699b3cb8de2db4c043409bae4861b 100644
--- a/cpp/particles/p2p/p2p_merge_collisions.hpp
+++ b/cpp/particles/p2p/p2p_merge_collisions.hpp
@@ -47,11 +47,16 @@ public:
template <int size_particle_positions, int size_particle_rhs>
void compute_interaction(const partsize_t idx_part1,
- const real_number /*pos_part1*/[], real_number /*rhs_part1*/[],
+ const real_number /*pos_part1*/[],
+ real_number /*rhs_part1*/[],
const partsize_t idx_part2,
- const real_number /*pos_part2*/[], real_number /*rhs_part2*/[],
- const real_number /*dist_pow2*/, const real_number /*cutoff*/,
- const real_number /*xshift_coef*/, const real_number /*yshift_coef*/, const real_number /*zshift_coef*/){
+ const real_number /*pos_part2*/[],
+ real_number /*rhs_part2*/[],
+ const real_number /*dist_pow2*/,
+ const real_number /*cutoff*/,
+ const real_number /*xseparation*/,
+ const real_number /*yseparation*/,
+ const real_number /*zseparation*/){
mergedParticles.emplace_back(std::max(idx_part1,idx_part2));
}