diff --git a/src/main_logic/main.cpp b/src/main_logic/main.cpp index f1a44a4..15fdbf7 100644 --- a/src/main_logic/main.cpp +++ b/src/main_logic/main.cpp @@ -72,7 +72,7 @@ int main() if(!polyhedronDispatchFailed) { while (true) - { + { get_mapnodes_reflections(simulation_map, mapnodes_reflections); for (RunIterationFunc f : iteration_runners) { diff --git a/src/main_logic/simulation_logic.cu b/src/main_logic/simulation_logic.cu index b230564..6bfdd4c 100644 --- a/src/main_logic/simulation_logic.cu +++ b/src/main_logic/simulation_logic.cu @@ -1,7 +1,5 @@ #ifdef COMPILE_FOR_CPU - #include - #endif //COMPILE_FOR_CPU #include diff --git a/src/simulation_objects/geometric/Face.cu b/src/simulation_objects/geometric/Face.cu index c4c2e10..797b167 100644 --- a/src/simulation_objects/geometric/Face.cu +++ b/src/simulation_objects/geometric/Face.cu @@ -1,5 +1,9 @@ #include +#ifdef COMPILE_FOR_CPU +#include +#endif //COMPILE_FOR_CPU + #include "Face.cuh" #include "../MapNode.cuh" #include "Polyhedron.cuh" @@ -110,6 +114,42 @@ __host__ __device__ bool operator==(const Face &a, const Face &b) return true; } +__host__ __device__ bool Face::contains_point(SpacePoint p) +{ + /* + * Calculate the area of the face as sum of areas of triangles with vertices `vertices[0]`, `vertices[i]`, + * `vertices[i + 1]` for all `i` > 1. + * WARNING: only works for convex polyhedrons. + */ + double face_area_by_triangles = 0; + for(int i = 1; i < n_of_vertices - 1; ++i) + { + SpacePoint a = vertices[i] - vertices[0]; + SpacePoint b = vertices[i + 1] - vertices[0]; + + // Magnitude of the cross product `a % b` is the area of the parallelogram, its half is the area of the triangle + face_area_by_triangles += get_distance(origin, a % b) / 2; + } + + /* + * Calculate the area of the face as sum of areas of triangles with vertices `p` (point to be checked), + * `vertices[i]`, `vertices[i + 1]` for all `i` + */ + double face_area_with_point = 0; + for(int i = 0; i < n_of_vertices; ++i) + { + SpacePoint this_vertex = vertices[i], next_vertex = vertices[(i + 1) % n_of_vertices]; + SpacePoint a = this_vertex - p, b = next_vertex - p; + + // Magnitude of the cross product `a % b` is the area of the parallelogram, its half is the area of the triangle + face_area_with_point += get_distance(origin, a % b) / 2; + } + + // The face contains point if and only if the two areas are equal to each other + return std::abs(face_area_by_triangles - face_area_with_point) < eps; +} + + __host__ __device__ bool operator!=(const Face &a, const Face &b) { return !(a == b); diff --git a/src/simulation_objects/geometric/Face.cuh b/src/simulation_objects/geometric/Face.cuh index ff80a43..470570c 100644 --- a/src/simulation_objects/geometric/Face.cuh +++ b/src/simulation_objects/geometric/Face.cuh @@ -33,7 +33,8 @@ public: * * @param vertices Array of polyhedron vertices that belong to the face. * Must be ordered in a special way (see note below) - * @param n_of_vertices Number of vertices on the face + * @param n_of_vertices Number of elements in the `vertices` array (i.e. number of vertices, INCLUDING the repeated + * one (see below)) * * @note The vertices order matters: looking on a face from outside the polyhedron, some vertex (let's call * it A) must be saved to `vertices[0]`. It's neighbour clockwise - vertex B (B is next to A clockwise) must be @@ -125,6 +126,21 @@ public: */ __host__ __device__ friend bool operator==(const Face &a, const Face &b); + + /** + * Checks whether the given point belongs to the face via area calculation. + * + * Checks if the given point belongs to the face by calculating actual area of the face and area of the face as if + * the point was on it. If the areas match (equal +- epsilon), the point is considered to belong to the face. + * + * @param p Point to be checked + * + * @returns `true` if the point belongs to the face, `false` otherwise + * + * @warning The method only works correctly if the face is a **convex** polygon + */ + __host__ __device__ bool contains_point(SpacePoint p); + private: /// Pointer to some node laying on the face MapNode *node; diff --git a/src/simulation_objects/geometric/Polyhedron.cu b/src/simulation_objects/geometric/Polyhedron.cu index 3eaebbf..0819596 100644 --- a/src/simulation_objects/geometric/Polyhedron.cu +++ b/src/simulation_objects/geometric/Polyhedron.cu @@ -72,10 +72,17 @@ __host__ __device__ Face *Polyhedron::find_face_by_point(SpacePoint point) const for(int i = 0; i < n_of_faces; ++i) { Face *face = &faces[i]; + + /* + // Old version of the check. Probably will be investigated later? SpacePoint normal = (point - face->get_vertices()[0]) % (face->get_vertices()[1] - face->get_vertices()[0]); normal = normal / get_distance(normal, origin); if(normal * face->get_normal() >= 1 - eps) return face; + */ + + if(face->contains_point(point)) + return face; } return &faces[0]; } @@ -96,14 +103,14 @@ __host__ __device__ double Polyhedron::calculate_square_of_surface() double square = 0; for(int i = 0; i < n_of_faces; ++i) { - // Cause first vertex of face repeats again in the end the condition is `j < faces[i].get_n_of_vertices() - 2` + // Start with the vertex `1`, not `0`, because the 0th vertex is repeated in the end of the vertices array for(int j = 1; j < faces[i].get_n_of_vertices() - 1; ++j) { SpacePoint a = faces[i].get_vertices()[j + 1] - faces[i].get_vertices()[0]; SpacePoint b = faces[i].get_vertices()[j] - faces[i].get_vertices()[0]; double sign_of_square = (a % b) * faces[i].get_normal(); - sign_of_square /= abs(sign_of_square); + sign_of_square /= std::abs(sign_of_square); square += sign_of_square * (a ^ b) / 2; } } @@ -191,6 +198,5 @@ __host__ __device__ SpacePoint get_projected_vector_end(SpacePoint vector_start, vector_end, current_face, &intersection_edge_vertex_id); } - // If vector AB does not intersect any edge of face, `vector_end` equals `b` return vector_end; }