Skip to content

Commit

Permalink
Add parallel versions of minmax functions
Browse files Browse the repository at this point in the history
  • Loading branch information
MaelRL committed Jan 14, 2025
1 parent 306ccad commit 20ecde1
Show file tree
Hide file tree
Showing 3 changed files with 111 additions and 69 deletions.
169 changes: 103 additions & 66 deletions Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/measure.h
Original file line number Diff line number Diff line change
Expand Up @@ -15,20 +15,18 @@

#include <CGAL/license/Polygon_mesh_processing/measure.h>

#include <CGAL/disable_warnings.h>
#include <CGAL/Polygon_mesh_processing/border.h>
#include <CGAL/Polygon_mesh_processing/shape_predicates.h>

#include <CGAL/assertions.h>
#include <CGAL/boost/graph/iterator.h>
#include <CGAL/boost/graph/helpers.h>
#include <CGAL/boost/graph/named_params_helper.h>
#include <CGAL/boost/graph/properties.h>
#include <CGAL/for_each.h>
#include <CGAL/Named_function_parameters.h>
#include <CGAL/boost/graph/named_params_helper.h>

#include <CGAL/Polygon_mesh_processing/border.h>
#include <CGAL/Polygon_mesh_processing/shape_predicates.h>
#include <CGAL/utils_classes.h>

#include <CGAL/Lazy.h> // needed for CGAL::exact(FT)/CGAL::exact(Lazy_exact_nt<T>)
#include <CGAL/utils_classes.h>

#include <boost/container/small_vector.hpp>
#include <boost/graph/graph_traits.hpp>
Expand Down Expand Up @@ -251,6 +249,8 @@ average_edge_length(const PolygonMesh& pmesh,
* returns the shortest and longest edges of a range of edges of a given polygon mesh,
* as well as their respective edge lengths.
*
* @tparam ConcurrencyTag enables sequential versus parallel algorithm. Possible values are
* `Sequential_tag`, `Parallel_tag`, and `Parallel_if_available_tag`.
* @tparam EdgeRange a model of `Range` whose iterator type is `InputIterator` with value type
* `boost::graph_traits<PolygonMesh>::edge_descriptor`.
* @tparam PolygonMesh a model of `HalfedgeGraph`
Expand Down Expand Up @@ -286,7 +286,8 @@ average_edge_length(const PolygonMesh& pmesh,
* @sa `edge_length()`
* @sa `squared_edge_length()`
*/
template<typename EdgeRange,
template<typename ConcurrencyTag = CGAL::Sequential_tag,
typename EdgeRange,
typename PolygonMesh,
typename CGAL_NP_TEMPLATE_PARAMETERS>
#ifdef DOXYGEN_RUNNING
Expand All @@ -303,39 +304,54 @@ minmax_edge_length(const EdgeRange& edge_range,
using parameters::get_parameter;

using edge_descriptor = typename boost::graph_traits<PolygonMesh>::edge_descriptor;
using edge_iterator = typename boost::graph_traits<PolygonMesh>::edge_iterator;

using Geom_traits = typename GetGeomTraits<PolygonMesh, CGAL_NP_CLASS>::type;
using FT = typename Geom_traits::FT;

edge_iterator first = std::cbegin(edge_range), beyond = std::cend(edge_range);
if(first == beyond)
if(std::begin(edge_range) == std::end(edge_range))
{
return std::make_pair(std::make_pair(edge_descriptor(), FT(0)),
std::make_pair(edge_descriptor(), FT(0)));
}

edge_iterator low = first, high = first, eit = first;
FT sq_lo, sq_hi;
sq_lo = sq_hi = squared_edge_length(*eit++, pmesh, np);

for(; eit!=beyond; ++eit)
struct Extremas
{
const FT sq_l = squared_edge_length(*eit, pmesh, np);

if(sq_l < sq_lo) {
low = eit;
sq_lo = sq_l;
}
if(sq_l > sq_hi) {
high = eit;
sq_hi = sq_l;
}
}

CGAL_assertion(low != beyond && high != beyond);
return std::make_pair(std::make_pair(*low, CGAL::approximate_sqrt(sq_lo)),
std::make_pair(*high, CGAL::approximate_sqrt(sq_hi)));
edge_descriptor low{};
edge_descriptor high{};
FT sq_lo{-1};
FT sq_hi{-1};
#ifdef CGAL_LINKED_WITH_TBB
std::mutex mutex;
#endif
};
Extremas extremas;

// Initialize with first edge
auto first = *std::begin(edge_range);
extremas.low = extremas.high = first;
extremas.sq_lo = extremas.sq_hi = squared_edge_length(first, pmesh, np);

CGAL::for_each<ConcurrencyTag>
(edge_range,
[&](const edge_descriptor& e) -> bool
{
const FT sq_l = squared_edge_length(e, pmesh, np);

#ifdef CGAL_LINKED_WITH_TBB
std::lock_guard<std::mutex> lock(extremas.mutex);
#endif
if(extremas.sq_lo < 0 || sq_l < extremas.sq_lo) {
extremas.low = e;
extremas.sq_lo = sq_l;
}
if(extremas.sq_hi < 0 || sq_l > extremas.sq_hi) {
extremas.high = e;
extremas.sq_hi = sq_l;
}
return true;
});

return std::make_pair(std::make_pair(extremas.low, CGAL::approximate_sqrt(extremas.sq_lo)),
std::make_pair(extremas.high, CGAL::approximate_sqrt(extremas.sq_hi)));
}

/*!
Expand All @@ -344,7 +360,8 @@ minmax_edge_length(const EdgeRange& edge_range,
* as well as their respective edge lengths.
* Equivalent to `minmax_edge_length(edges(pmesh), pmesh, np)`
*/
template<typename PolygonMesh,
template<typename ConcurrencyTag = CGAL::Sequential_tag,
typename PolygonMesh,
typename CGAL_NP_TEMPLATE_PARAMETERS>
#ifdef DOXYGEN_RUNNING
std::pair<std::pair<boost::graph_traits<PolygonMesh>::edge_descriptor`, FT>,
Expand All @@ -355,7 +372,7 @@ auto
minmax_edge_length(const PolygonMesh& pmesh,
const CGAL_NP_CLASS& np = parameters::default_values())
{
return minmax_edge_length(edges(pmesh), pmesh, np);
return minmax_edge_length<ConcurrencyTag>(edges(pmesh), pmesh, np);
}

/**
Expand Down Expand Up @@ -1240,6 +1257,8 @@ void match_faces(const PolygonMesh1& m1,
*
* \brief computes the minimum and maximum dihedral angles of a range of edges of a given polygon mesh.
*
* \tparam ConcurrencyTag enables sequential versus parallel algorithm. Possible values are
* `Sequential_tag`, `Parallel_tag`, and `Parallel_if_available_tag`.
* \tparam EdgeRange a model of `Range` whhose iterator type is `InputIterator` with value type
* `boost::graph_traits<PolygonMesh>::edge_descriptor`.
* \tparam TriangleMesh a model of `HalfedgeListGraph`
Expand Down Expand Up @@ -1273,8 +1292,10 @@ void match_faces(const PolygonMesh1& m1,
* \pre There are no degenerate faces in `tmesh`.
*
* \see `detect_sharp_edges()`
* \see `sharp_edges_segmentation()`
*/
template<typename EdgeRange,
template<typename ConcurrencyTag = CGAL::Sequential_tag,
typename EdgeRange,
typename TriangleMesh,
typename CGAL_NP_TEMPLATE_PARAMETERS>
#ifdef DOXYGEN_RUNNING
Expand Down Expand Up @@ -1305,46 +1326,64 @@ minmax_dihedral_angle(const EdgeRange& edge_range,

CGAL_precondition(is_triangle_mesh(tmesh));

edge_descriptor low, high;
FT lo(200), hi(-200);
struct Extremas
{
edge_descriptor low{};
edge_descriptor high{};
FT lo{200};
FT hi{-200};
#ifdef CGAL_LINKED_WITH_TBB
std::mutex mutex;
#endif
};
Extremas extremas;

typename Geom_traits::Compute_approximate_dihedral_angle_3 approx_dh =
gt.compute_approximate_dihedral_angle_3_object();

for(edge_descriptor e : edge_range)
{
CGAL_assertion(is_valid_edge_descriptor(e, tmesh));
if(is_border(e, tmesh))
continue;

const halfedge_descriptor h = halfedge(e, tmesh);
CGAL_assertion(!is_degenerate_triangle_face(h, tmesh));

const vertex_descriptor p = source(h,tmesh);
const vertex_descriptor q = target(h,tmesh);
const vertex_descriptor r = target(next(h,tmesh),tmesh);
const vertex_descriptor s = target(next(opposite(h,tmesh),tmesh),tmesh);

const FT da = approx_dh(get(vpm,p), get(vpm,q), get(vpm,r), get(vpm,s));
if(da < lo) {
low = e;
lo = da;
}
if(da > hi) {
high = e;
hi = da;
}
}
CGAL::for_each<ConcurrencyTag>
(edge_range,
[&](const edge_descriptor& e) -> bool
{
if(is_border(e, tmesh))
return true;

const halfedge_descriptor h = halfedge(e, tmesh);
CGAL_assertion(!is_degenerate_triangle_face(h, tmesh));

const vertex_descriptor p = source(h,tmesh);
const vertex_descriptor q = target(h,tmesh);
const vertex_descriptor r = target(next(h,tmesh),tmesh);
const vertex_descriptor s = target(next(opposite(h,tmesh),tmesh),tmesh);

const FT da = approx_dh(get(vpm,p), get(vpm,q), get(vpm,r), get(vpm,s));

return std::make_pair(std::make_pair(low, lo), std::make_pair(high, hi));
#ifdef CGAL_LINKED_WITH_TBB
std::lock_guard<std::mutex> lock(extremas.mutex);
#endif
if(da < extremas.lo) {
extremas.low = e;
extremas.lo = da;
}
if(da > extremas.hi) {
extremas.high = e;
extremas.hi = da;
}

return true;
});

return std::make_pair(std::make_pair(extremas.low, extremas.lo),
std::make_pair(extremas.high, extremas.hi));
}

/*!
* \ingroup PMP_measure_grp
* computes the minimum and maximum dihedral angles of a given triangle mesh.
* Equivalent to `minmax_dihedral_angle(edges(tmesh), tmesh, np)`
*/
template<typename TriangleMesh,
template<typename ConcurrencyTag = CGAL::Sequential_tag,
typename TriangleMesh,
typename CGAL_NP_TEMPLATE_PARAMETERS>
#ifdef DOXYGEN_RUNNING
std::pair<std::pair<boost::graph_traits<PolygonMesh>::edge_descriptor`, FT>,
Expand All @@ -1355,12 +1394,10 @@ auto
minmax_dihedral_angle(const TriangleMesh& tmesh,
const CGAL_NP_CLASS& np = parameters::default_values())
{
return minmax_dihedral_angle(edges(tmesh), tmesh, np);
return minmax_dihedral_angle<ConcurrencyTag>(edges(tmesh), tmesh, np);
}

} // namespace Polygon_mesh_processing
} // namespace CGAL

#include <CGAL/enable_warnings.h>

#endif // CGAL_POLYGON_MESH_PROCESSING_MEASURE_H
Original file line number Diff line number Diff line change
Expand Up @@ -109,6 +109,7 @@ if(TARGET CGAL::TBB_support)
target_link_libraries(orient_polygon_soup_test PRIVATE CGAL::TBB_support)
target_link_libraries(self_intersection_surface_mesh_test PRIVATE CGAL::TBB_support)
target_link_libraries(test_autorefinement PRIVATE CGAL::TBB_support)
target_link_libraries(measures_test PRIVATE CGAL::TBB_support)
else()
message(STATUS "NOTICE: Intel TBB was not found. Tests will use sequential code.")
endif()
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -109,7 +109,8 @@ void test_pmesh(const Mesh& pmesh)
edge_descriptor longest_edge;
FT longest_edge_length;

auto [shortest_edge_pair, longest_edge_pair] = PMP::minmax_edge_length(pmesh);
auto [shortest_edge_pair, longest_edge_pair] =
PMP::minmax_edge_length<CGAL::Parallel_if_available_tag>(pmesh);
std::tie(shortest_edge, shortest_edge_length) = shortest_edge_pair;
std::tie(longest_edge, longest_edge_length) = longest_edge_pair;

Expand All @@ -124,9 +125,11 @@ void test_pmesh(const Mesh& pmesh)
<< ", length = " << longest_edge_length << std::endl;
for(edge_descriptor e : edges(pmesh))
{
FT edge_length = PMP::edge_length(e, pmesh);
std::cout << "edge length: " << edge_length << std::endl;

// this could be wrong due to numerical errors, but the tested meshes are gentle enough
// such that we don't need a predicate... right?
FT edge_length = PMP::edge_length(e, pmesh);
assert(shortest_edge_length <= edge_length);
assert(longest_edge_length >= edge_length);
}
Expand Down Expand Up @@ -402,7 +405,8 @@ void test_dihedral_angles(const std::string filename,
edge_descriptor edge_with_largest_da;
FT largest_da;

auto [smallest_da_pair, largest_da_pair] = PMP::minmax_dihedral_angle(pmesh);
auto [smallest_da_pair, largest_da_pair] =
PMP::minmax_dihedral_angle<CGAL::Parallel_if_available_tag>(pmesh);
std::tie(edge_with_smallest_da, smallest_da) = smallest_da_pair;
std::tie(edge_with_largest_da, largest_da) = largest_da_pair;

Expand Down

0 comments on commit 20ecde1

Please sign in to comment.