Merge branch '5.1.x-branch' into 5.2.x-branch

# Conflicts:
#	Tetrahedral_remeshing/examples/Tetrahedral_remeshing/CMakeLists.txt
This commit is contained in:
Laurent Rineau 2020-12-09 15:35:47 +01:00
commit 3fc0ba4435
9 changed files with 200 additions and 31 deletions

View File

@ -3420,6 +3420,7 @@ get_least_square_surface_plane(const Vertex_handle& v,
Bare_point& reference_point,
Surface_patch_index patch_index) const
{
typedef typename C3T3::Triangulation::Triangle Triangle;
typename Gt::Construct_point_3 cp = tr_.geom_traits().construct_point_3_object();
// Get incident facets
@ -3440,45 +3441,47 @@ get_least_square_surface_plane(const Vertex_handle& v,
const Weighted_point& position = tr_.point(v);
// Get adjacent surface points
std::vector<Bare_point> surface_point_vector;
typename Facet_vector::iterator fit = facets.begin(), fend = facets.end();
for ( ; fit != fend; ++fit )
std::vector<Triangle> triangles;
typename C3T3::Facet ref_facet;
for (typename C3T3::Facet f : facets)
{
if ( c3t3_.is_in_complex(*fit) &&
if ( c3t3_.is_in_complex(f) &&
(patch_index == Surface_patch_index() ||
c3t3_.surface_patch_index(*fit) == patch_index) )
c3t3_.surface_patch_index(f) == patch_index) )
{
const Cell_handle cell = fit->first;
const int i = fit->second;
ref_facet = f;
// In the case of a periodic triangulation, the incident facets of a point
// do not necessarily have the same offsets. Worse, the surface centers
// might not have the same offset as their facet. Thus, no solution except
// calling a function 'get_closest_point(p, q)' that simply returns q
// calling a function 'get_closest_triangle(p, t)' that simply returns t
// for a non-periodic triangulation, and checks all possible offsets for
// periodic triangulations
const Bare_point& bp = tr_.get_closest_point(cp(position),
cell->get_facet_surface_center(i));
surface_point_vector.push_back(bp);
Triangle t = c3t3_.triangulation().triangle(f);
Triangle ct = tr_.get_closest_triangle(cp(position), t);
triangles.push_back(ct);
}
}
// In some cases point is not a real surface point
if ( surface_point_vector.empty() )
if ( triangles.empty() )
return boost::none;
// Compute least square fitting plane
Plane_3 plane;
Bare_point point;
CGAL::linear_least_squares_fitting_3(surface_point_vector.begin(),
surface_point_vector.end(),
CGAL::linear_least_squares_fitting_3(triangles.begin(),
triangles.end(),
plane,
point,
Dimension_tag<0>(),
Dimension_tag<2>(),
tr_.geom_traits(),
Default_diagonalize_traits<FT, 3>());
reference_point = surface_point_vector.front();
reference_point = ref_facet.first->get_facet_surface_center(ref_facet.second);
return plane;
}

View File

@ -551,7 +551,6 @@ private:
Cell_circulator current_cell = tr.incident_cells(edge);
Cell_circulator done = current_cell;
CGAL_assertion(c3t3.is_in_complex(current_cell));
// a & b are fixed points
const Weighted_point& wa = tr.point(v);
@ -561,7 +560,6 @@ private:
const Weighted_point& a_b = tr.point(current_cell, current_cell->index(v));
Vector_3 ba = Vector_3(cp(a_b), b);
++current_cell;
CGAL_assertion(c3t3.is_in_complex(current_cell));
CGAL_assertion(current_cell != done);
// c & d are moving points
@ -573,7 +571,6 @@ private:
while ( current_cell != done )
{
CGAL_assertion(c3t3.is_in_complex(current_cell));
Bare_point d = tr.dual(current_cell);
const Weighted_point& a_d = tr.point(current_cell, current_cell->index(v));
Vector_3 da = Vector_3(cp(a_d), d);

View File

@ -97,7 +97,7 @@ public:
{
CGAL_precondition(facet>=0 && facet<4);
char current_bits = bits_;
while (bits_.compare_and_swap(current_bits | (1 << facet), current_bits) != current_bits)
while (bits_.compare_exchange_weak(current_bits, current_bits | (1 << facet)) )
{
current_bits = bits_;
}
@ -108,7 +108,7 @@ public:
{
CGAL_precondition(facet>=0 && facet<4);
char current_bits = bits_;
while (bits_.compare_and_swap(current_bits & (15 & ~(1 << facet)), current_bits) != current_bits)
while (bits_.compare_exchange_weak(current_bits, current_bits & (15 & ~(1 << facet))))
{
current_bits = bits_;
}

View File

@ -59,6 +59,7 @@ public:
typedef typename Geom_traits::FT FT;
typedef typename Base::Bare_point Bare_point;
typedef typename Base::Weighted_point Weighted_point;
typedef typename Base::Triangle Triangle;
typedef typename Base::Vertex_handle Vertex_handle;
typedef typename Base::Cell_handle Cell_handle;
@ -77,11 +78,16 @@ public:
// possibilities). To allow Periodic_3_mesh_3 to use Mesh_3's files,
// each mesh triangulation implements its own version.
Bare_point get_closest_point(const Bare_point& /*p*/, const Bare_point& q) const
const Bare_point& get_closest_point(const Bare_point& /*p*/, const Bare_point& q) const
{
return q;
}
const Triangle& get_closest_triangle(const Bare_point& /*p*/, const Triangle& t) const
{
return t;
}
void set_point(const Vertex_handle v,
const Vector& /*move*/,
const Weighted_point& new_position)

View File

@ -62,6 +62,7 @@ if ( CGAL_FOUND )
create_single_source_cgal_program( "test_mesh_3_issue_1554.cpp" )
create_single_source_cgal_program( "test_mesh_polyhedral_domain_with_features_deprecated.cpp" )
create_single_source_cgal_program( "test_meshing_with_one_step.cpp" )
create_single_source_cgal_program( "test_mesh_cell_base_3.cpp")
foreach(target
test_boost_has_xxx
@ -93,6 +94,7 @@ if ( CGAL_FOUND )
test_c3t3_extract_subdomains_boundaries
test_mesh_3_issue_1554
test_mesh_polyhedral_domain_with_features_deprecated
test_mesh_cell_base_3
test_meshing_with_one_step.cpp)
if(TARGET ${target})
target_link_libraries(${target} PUBLIC CGAL::Eigen_support)
@ -113,6 +115,7 @@ if ( CGAL_FOUND )
test_mesh_capsule_var_distance_bound
test_mesh_3_issue_1554
test_mesh_polyhedral_domain_with_features_deprecated
test_mesh_cell_base_3
)
if(TARGET ${target})
target_link_libraries(${target} PUBLIC CGAL::TBB_support)

View File

@ -0,0 +1,100 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/IO/facets_in_complex_3_to_triangle_mesh.h>
#include <CGAL/Mesh_triangulation_3.h>
#include <CGAL/Mesh_complex_3_in_triangulation_3.h>
#include <CGAL/Mesh_3/tet_soup_to_c3t3.h>
#include <CGAL/Mesh_3/Robust_intersection_traits_3.h>
#include <CGAL/Polyhedral_mesh_domain_with_features_3.h>
#include <CGAL/Mesh_cell_base_3.h>
#include <CGAL/tags.h>
#include <iostream>
#include <fstream>
int main (int argc, char** argv){
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Polyhedral_mesh_domain_with_features_3<K> Polyhedral_mesh_domain;
// Traits classes
typedef CGAL::Kernel_traits<Polyhedral_mesh_domain>::Kernel Robust_intersections_traits;
typedef CGAL::details::Mesh_geom_traits_generator<
Robust_intersections_traits>::type Robust_K;
#ifdef CGAL_LINKED_WITH_TBB
typedef CGAL::Parallel_tag Concurrency_tag;
#else
typedef CGAL::Sequential_tag Concurrency_tag;
#endif
// Triangulation
typedef CGAL::Mesh_cell_base_3<Robust_K, Polyhedral_mesh_domain> Cell_base;
typedef CGAL::Triangulation_cell_base_with_info_3<int, Robust_K, Cell_base> Cell_base_with_info;
typedef CGAL::Mesh_triangulation_3<Polyhedral_mesh_domain,
Robust_intersections_traits,
Concurrency_tag,
CGAL::Default,
Cell_base_with_info>::type Tr;
typedef CGAL::Mesh_complex_3_in_triangulation_3<Tr> C3t3;
// Open file
std::ifstream in (argc > 1 ? argv[1] : "data/elephant.mesh",
std::ios_base::in);
if(!in) {
std::cerr << "Error! Cannot open file " << argv[1] << std::endl;
return 1;
}
C3t3 c3t3;
if(CGAL::build_triangulation_from_file<C3t3::Triangulation, true>(in, c3t3.triangulation()))
{
for( C3t3::Triangulation::Finite_cells_iterator
cit = c3t3.triangulation().finite_cells_begin();
cit != c3t3.triangulation().finite_cells_end();
++cit)
{
CGAL_assertion(cit->subdomain_index() >= 0);
c3t3.add_to_complex(cit, cit->subdomain_index());
for(int i=0; i < 4; ++i)
{
if(cit->surface_patch_index(i)>0)
{
c3t3.add_to_complex(cit, i, cit->surface_patch_index(i));
}
}
}
//if there is no facet in the complex, we add the border facets.
if(c3t3.number_of_facets_in_complex() == 0)
{
for( C3t3::Triangulation::All_cells_iterator
cit = c3t3.triangulation().all_cells_begin();
cit != c3t3.triangulation().all_cells_end();
++cit)
{
if(c3t3.triangulation().is_infinite(cit))
{
for(int i=0; i<4; ++i)
{
if(!c3t3.triangulation().is_infinite(cit, i))
c3t3.add_to_complex(cit, i, 1);
}
}
}
}
}
CGAL::Polyhedron_3<K> poly;
CGAL::facets_in_complex_3_to_triangle_mesh(c3t3, poly);
std::cout << "Graph has " << num_faces(poly) << " faces" << std::endl;
std::cout << "Graph has " << num_vertices(poly) << " vertices" << std::endl;
std::ofstream out("graph.off");
out << poly;
CGAL_assertion(is_valid(poly));
return EXIT_SUCCESS;
}

View File

@ -444,10 +444,10 @@ public:
const Bare_point cq = canonicalize_point(q);
FT min_sq_dist = std::numeric_limits<FT>::infinity();
for(int i = 0; i < 3; ++i) {
for(int j = 0; j < 3; ++j) {
for(int k = 0; k < 3; ++k) {
const Bare_point tcq = construct_point(std::make_pair(cq, Offset(i-1, j-1, k-1)));
for(int i = -1; i < 2; ++i) {
for(int j = -1; j < 2; ++j) {
for(int k = -1; k < 2; ++k) {
const Bare_point tcq = construct_point(std::make_pair(cq, Offset(i, j, k)));
FT sq_dist = geom_traits().compute_squared_distance_3_object()(p, tcq);
if(sq_dist < min_sq_dist)
@ -471,6 +471,48 @@ public:
return cwp(get_closest_point(cp(wp), cp(wq)), cw(wq));
}
Triangle get_closest_triangle(const Bare_point& p, const Triangle& t) const
{
typename Geom_traits::Construct_vector_3 cv = geom_traits().construct_vector_3_object();
typename Geom_traits::Construct_translated_point_3 tr = geom_traits().construct_translated_point_3_object();
typename Geom_traits::Compute_squared_distance_3 csd = geom_traits().compute_squared_distance_3_object();
// It doesn't matter which point we use to canonicalize the triangle as P3M3 is necessarily
// in one cover and we have to look at all the neighboring copies anyway since we do not
// have control of 'p'.
Bare_point canon_p0 = canonicalize_point(t[0]);
Vector_3 move_to_canonical = cv(t[0], canon_p0);
const std::array<Bare_point, 3> ct = { canon_p0,
tr(t[1], move_to_canonical),
tr(t[2], move_to_canonical) };
FT min_sq_dist = std::numeric_limits<FT>::infinity();
Triangle rt;
for(int i = -1; i < 2; ++i) {
for(int j = -1; j < 2; ++j) {
for(int k = -1; k < 2; ++k) {
const Triangle tt(
construct_point(std::make_pair(ct[0], Offset(i, j, k))),
construct_point(std::make_pair(ct[1], Offset(i, j, k))),
construct_point(std::make_pair(ct[2], Offset(i, j, k))));
const FT sq_dist = csd(p, tt);
if(sq_dist == FT(0))
return rt;
if(sq_dist < min_sq_dist) {
rt = tt;
min_sq_dist = sq_dist;
}
}
}
}
return rt;
}
// Warning: This is a periodic version that computes the smallest possible
// distances between p and q, and between p and r FOR ALL POSSIBLE OFFSETS
// before comparing these distances.

View File

@ -19,6 +19,16 @@ endif()
# Use Eigen for Mesh_3
find_package(Eigen3 3.1.0 REQUIRED) #(3.1.0 or greater)
include(CGAL_Eigen_support)
find_package(TBB QUIET)
include(CGAL_TBB_support)
# Concurrent Mesh_3
option(CGAL_ACTIVATE_CONCURRENT_MESH_3 "Activate parallelism in Mesh_3" OFF)
if(CGAL_ACTIVATE_CONCURRENT_MESH_3 OR ENV{CGAL_ACTIVATE_CONCURRENT_MESH_3})
add_definitions(-DCGAL_CONCURRENT_MESH_3)
find_package(TBB REQUIRED)
include(CGAL_TBB_support)
endif()
# Creating entries for all C++ files with "main" routine
# ##########################################################
@ -26,7 +36,9 @@ create_single_source_cgal_program( "tetrahedral_remeshing_example.cpp" )
create_single_source_cgal_program( "tetrahedral_remeshing_with_features.cpp")
create_single_source_cgal_program( "tetrahedral_remeshing_of_one_subdomain.cpp")
create_single_source_cgal_program( "tetrahedral_remeshing_from_mesh.cpp")
create_single_source_cgal_program(
"mesh_and_remesh_polyhedral_domain_with_features.cpp")
target_link_libraries(mesh_and_remesh_polyhedral_domain_with_features
PUBLIC CGAL::Eigen_support)
create_single_source_cgal_program( "mesh_and_remesh_polyhedral_domain_with_features.cpp" )
target_link_libraries(mesh_and_remesh_polyhedral_domain_with_features PUBLIC CGAL::Eigen_support)
if(CGAL_ACTIVATE_CONCURRENT_MESH_3 AND TARGET CGAL::TBB_support)
target_link_libraries(mesh_and_remesh_polyhedral_domain_with_features PRIVATE CGAL::TBB_support)
endif()

View File

@ -14,8 +14,14 @@ typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Mesh_polyhedron_3<K>::type Polyhedron;
typedef CGAL::Polyhedral_mesh_domain_with_features_3<K> Mesh_domain;
#ifdef CGAL_CONCURRENT_MESH_3
typedef CGAL::Parallel_tag Concurrency_tag;
#else
typedef CGAL::Sequential_tag Concurrency_tag;
#endif
// Triangulation for Meshing
typedef CGAL::Mesh_triangulation_3<Mesh_domain, CGAL::Default, CGAL::Default>::type Tr;
typedef CGAL::Mesh_triangulation_3<Mesh_domain, CGAL::Default, Concurrency_tag>::type Tr;
typedef CGAL::Mesh_complex_3_in_triangulation_3<
Tr, Mesh_domain::Corner_index, Mesh_domain::Curve_index> C3t3;