From cd53dc3fc652be64782016c61898d0ce41cc679d Mon Sep 17 00:00:00 2001 From: Jose Agustin Arguedas Leiva <agustin.arguedas@ds.mpg.de> Date: Tue, 2 Mar 2021 10:44:47 +0100 Subject: [PATCH 1/2] moved variables into protected to use in derived classes --- cpp/particles/p2p/p2p_ghost_collisions.hpp | 15 ++++++++------- 1 file changed, 8 insertions(+), 7 deletions(-) diff --git a/cpp/particles/p2p/p2p_ghost_collisions.hpp b/cpp/particles/p2p/p2p_ghost_collisions.hpp index 166ac5fc..4671b91a 100644 --- a/cpp/particles/p2p/p2p_ghost_collisions.hpp +++ b/cpp/particles/p2p/p2p_ghost_collisions.hpp @@ -45,17 +45,12 @@ class p2p_ghost_collisions const double pi = atan(1.0)*4; const double pi2 = atan(1.0)*8; // depending on the particle shape, we will be deciding whether particles intersect in different ways - enum particle_shape {CYLINDER, SPHERE, DISK}; - particle_shape current_particle_shape; - - // description for cylinders: - double cylinder_width; - double cylinder_length; // description for disks: double disk_width; - bool isActive; protected: + enum particle_shape {CYLINDER, SPHERE, DISK}; + particle_shape current_particle_shape; void add_colliding_pair(partsize_t idx_part1, partsize_t idx_part2) { // store colliding particle ids in order, to be able to identify pairs more easily @@ -70,6 +65,12 @@ class p2p_ghost_collisions this->collision_pairs_local.push_back(idx_part_small); this->collision_pairs_local.push_back(idx_part_big); } + // description for cylinders: + double cylinder_width; + double cylinder_length; + + bool isActive; + bool synchronisation; std::vector <partsize_t> collision_pairs_local; std::vector <partsize_t> collision_pairs_global; -- GitLab From 99a1922ed7b765a4c5e8b9c83411b58401ecf00b Mon Sep 17 00:00:00 2001 From: Jose Agustin Arguedas Leiva <agustin.arguedas@ds.mpg.de> Date: Wed, 3 Mar 2021 13:34:07 +0100 Subject: [PATCH 2/2] passed functionality in p2p cylinders to a function --- cpp/particles/p2p/p2p_ghost_collisions.hpp | 55 ++++++++++++---------- 1 file changed, 30 insertions(+), 25 deletions(-) diff --git a/cpp/particles/p2p/p2p_ghost_collisions.hpp b/cpp/particles/p2p/p2p_ghost_collisions.hpp index 4671b91a..84bb63b0 100644 --- a/cpp/particles/p2p/p2p_ghost_collisions.hpp +++ b/cpp/particles/p2p/p2p_ghost_collisions.hpp @@ -90,6 +90,16 @@ public: this->collision_pairs_local.reserve(src.collision_pairs_local.capacity()); } + double rod_distance(double px, double py, double pz, double qx, double qy, double qz, double t, double s, double x, double y, double z){ + double x_dist, y_dist, z_dist; + double min_distance; + x_dist = this->cylinder_length*0.5*t*px-this->cylinder_length*0.5*s*qx+x; + y_dist = this->cylinder_length*0.5*t*py-this->cylinder_length*0.5*s*qy+y; + z_dist = this->cylinder_length*0.5*t*pz-this->cylinder_length*0.5*s*qz+z; + min_distance = sqrt(x_dist*x_dist+y_dist*y_dist+z_dist*z_dist); + return min_distance; + } + template <int size_particle_rhs> void init_result_array(real_number /*rhs*/[], const partsize_t /*nbParticles*/) const{ } @@ -209,7 +219,8 @@ public: break; case CYLINDER: { - double pq, xp, xq, t, s, x_dist, y_dist, z_dist ,min_distance; + double pq, xp, xq, t, s, min_distance, min_distance_current; + double px,py, pz, qx, qy, qz; const double x = xseparation; const double y = yseparation; @@ -218,12 +229,17 @@ public: if( r <= this->cylinder_length ) { - /* p and q are the orientation vectors of the first and second particles. */ + px = pos_part1[IDXC_X+3]; + py = pos_part1[IDXC_Y+3]; + pz = pos_part1[IDXC_Z+3]; + qx = pos_part2[IDXC_X+3]; + qy = pos_part2[IDXC_Y+3]; + qz = pos_part2[IDXC_Z+3]; /* pq, xp, xq are scalar products of these vectors with x, relative position */ - pq = pos_part1[IDXC_X+3] * pos_part2[IDXC_X+3] + pos_part1[IDXC_Y+3] * pos_part2[IDXC_Y+3] + pos_part1[IDXC_Z+3] * pos_part2[IDXC_Z+3]; - xp = x * pos_part1[IDXC_X+3] + y * pos_part1[IDXC_Y+3] + z * pos_part1[IDXC_Z+3]; - xq = x * pos_part2[IDXC_X+3] + y * pos_part2[IDXC_Y+3] + z * pos_part2[IDXC_Z+3]; + pq = px * qx + py * qy + pz * qz; + xp = x * px + y * py + z * pz; + xq = x * qx + y * qy + z * qz; /* t and s parametrize the two rods. Find min distance: */ assert(this->cylinder_length > 0); t = 2.0/(this->cylinder_length*(pq*pq-1.0))*(xp-pq*xq); @@ -232,10 +248,7 @@ public: if( abs(t)<=1.0 and abs(s)<=1.0 ) { /* Get minimal distance in case of both t and s in {-1,1}. Else: check edges */ - x_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_X+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_X+3]+x; - y_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Y+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Y+3]+y; - z_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Z+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Z+3]+z; - min_distance = sqrt(x_dist*x_dist+y_dist*y_dist+z_dist*z_dist); + min_distance = this->rod_distance(px, py, pz, qx, qy, qz, t, s, x, y, z); } else { @@ -244,34 +257,26 @@ public: t = 1.0; s = t*pq+2.0/this->cylinder_length*xq; if( abs(s)>1.0 ) { s = s / abs(s) ;} - x_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_X+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_X+3]+x; - y_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Y+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Y+3]+y; - z_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Z+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Z+3]+z; - min_distance = fmin( sqrt(x_dist*x_dist+y_dist*y_dist+z_dist*z_dist), min_distance ); + min_distance_current = this->rod_distance(px, py, pz, qx, qy, qz, t, s, x, y, z); + min_distance = fmin( min_distance_current, min_distance ); /* t fixed at -1, find min along s */ t = -1.0; s = t*pq+2.0/this->cylinder_length*xq; if( abs(s)>1.0 ) { s = s / abs(s) ;} - x_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_X+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_X+3]+x; - y_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Y+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Y+3]+y; - z_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Z+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Z+3]+z; - min_distance = fmin( sqrt(x_dist*x_dist+y_dist*y_dist+z_dist*z_dist), min_distance ); + min_distance_current = this->rod_distance(px, py, pz, qx, qy, qz, t, s, x, y, z); + min_distance = fmin( min_distance_current, min_distance ); /* s fixed at 1, find min along t */ s = 1.0; t = s*pq-2.0/this->cylinder_length*xp; if( abs(t)>1.0 ) { t = t / abs(t) ;} - x_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_X+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_X+3]+x; - y_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Y+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Y+3]+y; - z_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Z+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Z+3]+z; - min_distance = fmin( sqrt(x_dist*x_dist+y_dist*y_dist+z_dist*z_dist), min_distance ); + min_distance_current = this->rod_distance(px, py, pz, qx, qy, qz, t, s, x, y, z); + min_distance = fmin( min_distance_current, min_distance ); /* s fixed at -1, find min along t */ s = -1.0; t = s*pq-2.0/this->cylinder_length*xp; if( abs(t)>1.0 ) { t = t / abs(t) ;} - x_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_X+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_X+3]+x; - y_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Y+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Y+3]+y; - z_dist = this->cylinder_length*0.5*t*pos_part1[IDXC_Z+3]-this->cylinder_length*0.5*s*pos_part2[IDXC_Z+3]+z; - min_distance = fmin( sqrt(x_dist*x_dist+y_dist*y_dist+z_dist*z_dist), min_distance ); + min_distance_current = this->rod_distance(px, py, pz, qx, qy, qz, t, s, x, y, z); + min_distance = fmin( min_distance_current, min_distance ); } /* If cylinders overlap count it as a collision */ if( min_distance<=this->cylinder_width ) -- GitLab