diff --git a/Installation/changes.html b/Installation/changes.html index 9021b6dc11b..4006638bc59 100644 --- a/Installation/changes.html +++ b/Installation/changes.html @@ -148,6 +148,14 @@ and src/ directories). + +

3D Mesh Generation

+ +

Release 4.4

@@ -158,7 +166,7 @@ and src/ directories).
  • Additional supported platforms: @@ -298,7 +306,7 @@ and src/ directories). @@ -340,7 +348,7 @@ and src/ directories). constructor added to CGAL::Object. However, it is recommended to upgrade your code. The previous behavior can be restored by defining the - macro CGAL_INTERSECTION_VERSION to 1. + macro CGAL_INTERSECTION_VERSION to1.
  • 2D Arrangements

    @@ -663,7 +671,7 @@ and src/ directories).
  • Additional supported platforms:
  • Improved configuration for essential and optional external third party software
  • diff --git a/Mesh_3/applications/CMakeLists.txt b/Mesh_3/applications/CMakeLists.txt index 94144bcbafa..7324e7ef359 100644 --- a/Mesh_3/applications/CMakeLists.txt +++ b/Mesh_3/applications/CMakeLists.txt @@ -38,7 +38,7 @@ if ( CGAL_FOUND ) endif() create_single_source_cgal_program( "mesh_polyhedral_domain.cpp" ) create_single_source_cgal_program( "output_distribution_to_stdout.cpp" ) -# create_single_source_cgal_program( "mesher_tester.cpp" "../examples/Mesh_3/implicit_functions.cpp" ) + create_single_source_cgal_program( "mesher_tester.cpp" "../examples/Mesh_3/implicit_functions.cpp" ) create_single_source_cgal_program( "mesher_tester_image.cpp" ) create_single_source_cgal_program( "mesher_tester_polyhedron.cpp" ) else() diff --git a/Mesh_3/applications/mesh_implicit_domains.cpp b/Mesh_3/applications/mesh_implicit_domains.cpp index 0902c39bded..0926a94ecea 100644 --- a/Mesh_3/applications/mesh_implicit_domains.cpp +++ b/Mesh_3/applications/mesh_implicit_domains.cpp @@ -26,15 +26,15 @@ //****************************************************************************** -#include +#include "debug.h" #include #include #include #include -#include -#include +#include +#include #include #include "../examples/Mesh_3/implicit_functions.h" @@ -51,10 +51,10 @@ using namespace CGAL::parameters; // Domain struct K: public CGAL::Exact_predicates_inexact_constructions_kernel {}; typedef FT_to_point_function_wrapper Function; -typedef CGAL::Mesh_3::Implicit_vector_to_labeled_function_wrapper +typedef CGAL::Implicit_multi_domain_to_labeling_function_wrapper Function_wrapper; typedef Function_wrapper::Function_vector Function_vector; -typedef CGAL::Mesh_3::Labeled_mesh_domain_3 Mesh_domain; +typedef CGAL::Labeled_mesh_domain_3 Mesh_domain; // Triangulation typedef CGAL::Mesh_triangulation_3::type Tr; @@ -95,7 +95,7 @@ void set_function(Function_vector& v, { if ( vm.count(function_name) ) { - v.push_back(&f); + v.push_back(f); std::cout << function_name << " "; } } diff --git a/Mesh_3/applications/mesher_tester.cpp b/Mesh_3/applications/mesher_tester.cpp index 3e89a3932d3..ee09c482d57 100644 --- a/Mesh_3/applications/mesher_tester.cpp +++ b/Mesh_3/applications/mesher_tester.cpp @@ -14,8 +14,8 @@ #include #include // Implicit domain -#include -#include +#include +#include #include #include "../examples/Mesh_3/implicit_functions.h" @@ -63,9 +63,9 @@ typedef CGAL::Polyhedron_3 Polyhedron; typedef CGAL::Polyhedral_mesh_domain_3 Polyhedral_domain; // Implicit domain typedef FT_to_point_function_wrapper I_Function; -typedef CGAL::Mesh_3::Implicit_vector_to_labeled_function_wrapper I_Function_wrapper; +typedef CGAL::Implicit_multi_domain_to_labeling_function_wrapper I_Function_wrapper; typedef I_Function_wrapper::Function_vector I_Function_vector; -typedef CGAL::Mesh_3::Labeled_mesh_domain_3 Implicit_domain; +typedef CGAL::Labeled_mesh_domain_3 Implicit_domain; // 3D Image typedef CGAL::Image_3 Image; typedef CGAL::Labeled_image_mesh_domain_3 Image_domain; @@ -121,7 +121,7 @@ void set_implicit_function(I_Function_vector& v, { if(vm.count(function_name)) { - v.push_back(&f); + v.push_back(f); std::cout << function_name << " "; } } @@ -398,7 +398,7 @@ void mesh(const std::string& data, const po::variables_map& vm) //Load the domain std::cout << "****** [" << filename << "] Create domain..."; std::flush(std::cout); - Domain_builder domain_builder(it->string()); + Domain_builder domain_builder(it->path().string()); std::cout <<"done (" << timer.time() << "s) ******\n\n"; while(std::getline(file_param,line_param)) diff --git a/Mesh_3/applications/mesher_tester.h b/Mesh_3/applications/mesher_tester.h index 9f5d9960603..12dd8c36823 100644 --- a/Mesh_3/applications/mesher_tester.h +++ b/Mesh_3/applications/mesher_tester.h @@ -570,7 +570,7 @@ void mesh(const std::string& data, const std::string& output_dir, const int nb_t //Load the domain std::stringstream cout_loc; cout_loc << "+ [" << filename << "] Create domain..."; - std::tr1::shared_ptr > pdomain_builder(new Domain_builder(it->string())); + std::tr1::shared_ptr > pdomain_builder(new Domain_builder(it->path().string())); cout_loc << "done (" << timer.time() << "s)\n"; std::cout << cout_loc.str(); diff --git a/Mesh_3/demo/Mesh_3/C3t3_type.h b/Mesh_3/demo/Mesh_3/C3t3_type.h index b8338d067d9..f17042b778f 100644 --- a/Mesh_3/demo/Mesh_3/C3t3_type.h +++ b/Mesh_3/demo/Mesh_3/C3t3_type.h @@ -42,7 +42,7 @@ typedef CGAL::Polyhedral_mesh_domain_with_features_3 Image_mesh_domain; typedef Wrapper Function_wrapper; -typedef CGAL::Mesh_3::Labeled_mesh_domain_3 Function_mesh_domain; +typedef CGAL::Labeled_mesh_domain_3 Function_mesh_domain; // Triangulation typedef CGAL::Mesh_triangulation_3::type Tr; diff --git a/Mesh_3/doc/Mesh_3/CGAL/Implicit_to_labeling_function_wrapper.h b/Mesh_3/doc/Mesh_3/CGAL/Implicit_to_labeling_function_wrapper.h new file mode 100644 index 00000000000..d506857fde8 --- /dev/null +++ b/Mesh_3/doc/Mesh_3/CGAL/Implicit_to_labeling_function_wrapper.h @@ -0,0 +1,64 @@ +namespace CGAL { + +/*! +\ingroup PkgMesh_3Domains + +The class `Implicit_multi_domain_to_labeling_function_wrapper` is an helping class to get a function with integer values +labeling the components of a multi-domain. The multidomain is described through a set of function {fi(p), i=1, ...n}. +Each component corresponds to a sign vector [s1, s2, ..., sn] where si is the sign of the function fi(p) at a point p of the component. +This wrapper class can be passed to `Labeled_mesh_domain_3` as first template parameter. + +\par Example +For example, the multidomain described by the three funstions [f1,f2,f3] and the two sign vectors [-,-,+] and [+,-,+] + includes two components.
    +The first one matches the locus of points satisfying f1(p)<0 and f2(p)<0 and f3(p)>0.
    +The second one matches the locus of points satisfying f1(p)>0 and f2(p)<0 and f3(p)>0.
    + +\tparam Function provides the definition of the function. +This parameter stands for a model of the concept ImplicitFunction described in the surface mesh generation package. +The number types Function::FT and BGT::FT are required to match. + +\sa `Implicit_mesh_domain_3`. + +*/ +template +class Implicit_multi_domain_to_labeling_function_wrapper +{ +public: + /// \name Types + /// @{ + //! + typedef std::vector Function_vector; + //! + typedef typename Function::Point Point_3; + /// @} + + /// \name Creation + /// @{ + /*! + * \brief Construction from a vector of implicit functions. + * \param implicit_functions the vector of implicit functions. + * + * Poistions vectors are built automatically so that the union of components equals the union of the functions. + */ + Implicit_multi_domain_to_labeling_function_wrapper (const Function_vector& implicit_functions); + + /*! + * \brief Construction from a vector of implicit functions and a vector of vector of signs. + * \param implicit_functions the vector of implicit functions. + * \param positions_vectors the vector of vector of signs. Each vector of position describes a component. + * \sa `Sign` + */ + Implicit_multi_domain_to_labeling_function_wrapper (const Function_vector& implicit_functions, const std::vector >& positions_vectors); + + /*! + * \brief Construction from a vector of implicit functions and a vector of strings. + * \param implicit_functions the vector of implicit functions. + * \param positions_strings the vector of string. The strings contained in this vector must contain '+' or '-' only. Each string (vector of positions) describes a component. + */ + Implicit_multi_domain_to_labeling_function_wrapper (const Function_vector& implicit_functions, const std::vector& positions_strings); +/// @} + +}; /* end Implicit_multi_domain_to_labeling_function_wrapper */ +/// @} +} /* end namespace CGAL */ diff --git a/Mesh_3/doc/Mesh_3/CGAL/Labeled_mesh_domain_3.h b/Mesh_3/doc/Mesh_3/CGAL/Labeled_mesh_domain_3.h new file mode 100644 index 00000000000..8755edeeb8c --- /dev/null +++ b/Mesh_3/doc/Mesh_3/CGAL/Labeled_mesh_domain_3.h @@ -0,0 +1,96 @@ +namespace CGAL { + +/*! +\ingroup PkgMesh_3Domains + +\brief The class `Labeled_mesh_domain_3` implements indexed domains. + +This class is a model of concept `MeshDomain_3`. + +Any boundary facet is labeled , a, where b!=0. + +This class includes a function that provides, the subdomain index of any +query point. An intersection between a segment and bounding +surfaces is detected when both segment endpoints are associated with different +values of subdomain indices. The intersection is then constructed by bisection. +The bisection stops when the query segment is shorter than an error bound +`e` given by the product of the +length of the diagonal of the bounding box (in world coordinates), or the radius of the bounding sphere, and +a relative error bound passed as argument to the constructor of `Labeled_mesh_domain_3`. + +Implicit_mesh_domain_3 is a heir of Labeled_mesh_domain_3. It uses a satisfying labeling function if there is only one component to mesh. + +\tparam LabelingFunction is the type of the input function.
    +This parameter stands for a model of the concept ImplicitFunction described in the surface mesh generation package.
    +Labeling function f must return 0 if the point isn't located in any subdomain. Usually, the return type of labeling functions are integer.
    +Let p be a Point. +
      +
    • f(p)=0 means that p is outside domain.
    • +
    • f(p)=a, a!=0 means that p is inside subdomain a.
    • +
    +Implicit_multi_domain_to_labeling_function_wrapper is a good candidate for this template parameter +if there are several components to mesh. + +\tparam BGT is a geometric traits class that provides +the basic operations to implement +intersection tests and intersection computations +through a bisection method. This parameter must be instantiated +with a model of the concept `BisectionGeometricTraits_3`. + +\cgalModels MeshDomain_3 + +\sa `Implicit_mesh_domain_3` +\sa `Implicit_multi_domain_to_labeling_function_wrapper` +\sa `CGAL::make_mesh_3()`. + +*/ +template +class Labeled_mesh_domain_3 +{ +public: + +/// \name Creation +/// @{ + +/*! +\brief Construction from a labeling function, a bounding Sphere and a relative error bound. +\param f the labeling function. +\param bounding_sphere the bounding sphere of the meshable space. +\param relative_error_bound is the relative error bound used to compute intersection points between the implicit surface and query segments. The +bisection is stopped when the length of the intersected segment is less than the product of `relative_error_bound` by the radius of +`bounding_sphere`. +*/ +Labeled_mesh_domain_3(const LabelingFunction& f, + const Sphere_3& bounding_sphere, + const FT& relative_error_bound = FT(1e-3)); + +/*! +\brief Construction from a labeling function, a bounding box and a relative error bound. +\param f the labeling function. +\param bbox the bounding box of the meshable space. +\param relative_error_bound is the relative error bound used to compute intersection points between the implicit surface and query segments. The +bisection is stopped when the length of the intersected segment is less than the product of `relative_error_bound` by the diagonal of +`bounding_box`. +*/ +Labeled_mesh_domain_3(const LabelingFunction& f, + const Bbox_3& bbox, + const FT& relative_error_bound = FT(1e-3)); + +/*! +\brief Construction from a function, a bounding Iso_cuboid_3 and a relative error bound. +\param f the function. +\param bbox the bounding box of the meshable space. +\param relative_error_bound is the relative error bound used to compute intersection points between the implicit surface and query segments. The +bisection is stopped when the length of the intersected segment is less than the product of `relative_error_bound` by the diagonal of +`bounding_box`. +*/ +Labeled_mesh_domain_3(const LabelingFunction& f, + const Iso_cuboid_3& bbox, + const FT& relative_error_bound = FT(1e-3)); + +/// @} + +}; /* end Labeled_mesh_domain_3 */ +} /* end namespace CGAL */ diff --git a/Mesh_3/examples/Mesh_3/CMakeLists.txt b/Mesh_3/examples/Mesh_3/CMakeLists.txt index 6dd0cd0837f..e973b4f6bb4 100644 --- a/Mesh_3/examples/Mesh_3/CMakeLists.txt +++ b/Mesh_3/examples/Mesh_3/CMakeLists.txt @@ -35,7 +35,10 @@ if ( CGAL_FOUND ) create_single_source_cgal_program( "mesh_implicit_sphere.cpp" ) create_single_source_cgal_program( "mesh_implicit_sphere_variable_size.cpp" ) create_single_source_cgal_program( "mesh_two_implicit_spheres_with_balls.cpp" ) -# create_single_source_cgal_program( "mesh_implicit_domains.cpp" "implicit_functions.cpp" ) + create_single_source_cgal_program( "mesh_implicit_domains_2.cpp" "implicit_functions.cpp" ) + create_single_source_cgal_program( "mesh_cubes_intersection.cpp" ) + create_single_source_cgal_program( "mesh_cubes_intersection_with_features.cpp" ) + create_single_source_cgal_program( "mesh_implicit_domains.cpp" "implicit_functions.cpp" ) create_single_source_cgal_program( "mesh_polyhedral_domain.cpp" ) create_single_source_cgal_program( "mesh_polyhedral_domain_with_features.cpp" ) if( WITH_CGAL_ImageIO ) diff --git a/Mesh_3/examples/Mesh_3/implicit_functions.h b/Mesh_3/examples/Mesh_3/implicit_functions.h index c8cfd5bc9bd..e9a3ce0d05b 100644 --- a/Mesh_3/examples/Mesh_3/implicit_functions.h +++ b/Mesh_3/examples/Mesh_3/implicit_functions.h @@ -24,12 +24,13 @@ double sphere_function (double x, double y, double z) // (c=(0,0,0), r=Sq_radius -template -class FT_to_point_function_wrapper : public std::unary_function +template +class FT_to_point_function_wrapper : public std::unary_function { typedef FT (*Implicit_function)(FT, FT, FT); Implicit_function function; public: + typedef P Point; FT_to_point_function_wrapper(Implicit_function f) : function(f) {} FT operator()(Point p) const { return function(p.x(), p.y(), p.z()); } }; diff --git a/Mesh_3/examples/Mesh_3/mesh_cubes_intersection.cpp b/Mesh_3/examples/Mesh_3/mesh_cubes_intersection.cpp new file mode 100644 index 00000000000..3b9da29830e --- /dev/null +++ b/Mesh_3/examples/Mesh_3/mesh_cubes_intersection.cpp @@ -0,0 +1,97 @@ + +//****************************************************************************** +// File Description : +// Outputs to out.mesh a mesh of implicit domains. +//****************************************************************************** + + + +#include + +#include +#include +#include + +#include +#include +#include + +// IO +#include + +using namespace CGAL::parameters; + +// Domain +typedef CGAL::Exact_predicates_inexact_constructions_kernel K; +typedef K::Point_3 Point; +typedef K::FT FT; +typedef FT (*Function)(const Point&); +typedef CGAL::Implicit_multi_domain_to_labeling_function_wrapper + Function_wrapper; +typedef Function_wrapper::Function_vector Function_vector; +typedef CGAL::Labeled_mesh_domain_3 Mesh_domain; + +// Triangulation +typedef CGAL::Mesh_triangulation_3::type Tr; +typedef CGAL::Mesh_complex_3_in_triangulation_3 C3t3; + +// Mesh Criteria +typedef CGAL::Mesh_criteria_3 Mesh_criteria; +typedef Mesh_criteria::Facet_criteria Facet_criteria; +typedef Mesh_criteria::Cell_criteria Cell_criteria; + + +double cube_function_1 (const Point& p) +{ + if( p.x() >= 0 && p.x() <= 2 && + p.y() >= 0 && p.y() <= 2 && + p.z() >= 0 && p.z() <= 2 ) + return -1.; + return 1.; +} + +double cube_function_2 (const Point& p) +{ + if( p.x() >= 1 && p.x() <= 3 && + p.y() >= 1 && p.y() <= 3 && + p.z() >= 1 && p.z() <= 3 ) + return -1.; + return 1.; +} + +int main() +{ + // Define functions + Function f1 = cube_function_1; + Function f2 = cube_function_2; + + Function_vector v; + v.push_back(f1); + v.push_back(f2); + + std::vector vps; + vps.push_back("--"); + + // Domain (Warning: Sphere_3 constructor uses square radius !) + Mesh_domain domain(Function_wrapper(v, vps), K::Sphere_3(CGAL::ORIGIN, 5.*5.)); + + // Set mesh criteria + Mesh_criteria criteria(edge_size = 0.15, + facet_angle = 30, facet_size = 0.2, + cell_radius_edge_ratio = 2, cell_size = 0.4); + + // Mesh generation + C3t3 c3t3 = CGAL::make_mesh_3(domain, criteria, no_exude(), no_perturb()); + + // Perturbation (maximum cpu time: 10s, targeted dihedral angle: default) + CGAL::perturb_mesh_3(c3t3, domain, time_limit = 10); + + // Exudation + CGAL::exude_mesh_3(c3t3,12); + + // Output + std::ofstream medit_file("out_cubes_intersection.mesh"); + CGAL::output_to_medit(medit_file, c3t3); + + return 0; +} diff --git a/Mesh_3/examples/Mesh_3/mesh_cubes_intersection_with_features.cpp b/Mesh_3/examples/Mesh_3/mesh_cubes_intersection_with_features.cpp new file mode 100644 index 00000000000..0ab5d0d7ffc --- /dev/null +++ b/Mesh_3/examples/Mesh_3/mesh_cubes_intersection_with_features.cpp @@ -0,0 +1,184 @@ + +//****************************************************************************** +// File Description : +// Outputs to out.mesh a mesh of implicit domains. +//****************************************************************************** + + + +#include + +#include +#include +#include + +#include +#include +#include +#include + +// IO +#include + +using namespace CGAL::parameters; + +// Domain +typedef CGAL::Exact_predicates_inexact_constructions_kernel K; +typedef K::Point_3 Point; +typedef K::FT FT; +typedef FT (*Function)(const Point&); +typedef CGAL::Implicit_multi_domain_to_labeling_function_wrapper + Function_wrapper; +typedef Function_wrapper::Function_vector Function_vector; +typedef CGAL::Labeled_mesh_domain_3 Mesh_domain; +typedef CGAL::Mesh_domain_with_polyline_features_3 Mesh_domain_with_features; + +// Triangulation +typedef CGAL::Mesh_triangulation_3::type Tr; +typedef CGAL::Mesh_complex_3_in_triangulation_3 C3t3; + +// Mesh Criteria +typedef CGAL::Mesh_criteria_3 Mesh_criteria; +typedef Mesh_criteria::Facet_criteria Facet_criteria; +typedef Mesh_criteria::Cell_criteria Cell_criteria; + + +double cube_function_1 (const Point& p) +{ + if( p.x() >= 0 && p.x() <= 2 && + p.y() >= 0 && p.y() <= 2 && + p.z() >= 0 && p.z() <= 2 ) + return -1.; + return 1.; +} + +double cube_function_2 (const Point& p) +{ + if( p.x() >= 1 && p.x() <= 3 && + p.y() >= 1 && p.y() <= 3 && + p.z() >= 1 && p.z() <= 3 ) + return -1.; + return 1.; +} + +// Polyline +typedef std::vector Polyline_3; +typedef std::list Polylines; + +void create_polylines (Polylines& polylines) +{ + { + Polyline_3 polyline; + polyline.push_back(Point(1,1,1)); + polyline.push_back(Point(2,1,1)); + polylines.push_back(polyline); + } + { + Polyline_3 polyline; + polyline.push_back(Point(2,1,1)); + polyline.push_back(Point(2,2,1)); + polylines.push_back(polyline); + } + { + Polyline_3 polyline; + polyline.push_back(Point(2,2,1)); + polyline.push_back(Point(1,2,1)); + polylines.push_back(polyline); + } + { + Polyline_3 polyline; + polyline.push_back(Point(1,2,1)); + polyline.push_back(Point(1,1,1)); + polylines.push_back(polyline); + } +//---------- + { + Polyline_3 polyline; + polyline.push_back(Point(1,1,2)); + polyline.push_back(Point(2,1,2)); + polylines.push_back(polyline); + } + { + Polyline_3 polyline; + polyline.push_back(Point(2,1,2)); + polyline.push_back(Point(2,2,2)); + polylines.push_back(polyline); + } + { + Polyline_3 polyline; + polyline.push_back(Point(2,2,2)); + polyline.push_back(Point(1,2,2)); + polylines.push_back(polyline); + } + { + Polyline_3 polyline; + polyline.push_back(Point(1,2,2)); + polyline.push_back(Point(1,1,2)); + polylines.push_back(polyline); + } + //---------- + { + Polyline_3 polyline; + polyline.push_back(Point(1,1,1)); + polyline.push_back(Point(1,1,2)); + polylines.push_back(polyline); + } + { + Polyline_3 polyline; + polyline.push_back(Point(2,1,1)); + polyline.push_back(Point(2,1,2)); + polylines.push_back(polyline); + } + { + Polyline_3 polyline; + polyline.push_back(Point(2,2,1)); + polyline.push_back(Point(2,2,2)); + polylines.push_back(polyline); + } + { + Polyline_3 polyline; + polyline.push_back(Point(1,2,1)); + polyline.push_back(Point(1,2,2)); + polylines.push_back(polyline); + } +} + +int main() +{ + // Define functions + Function f1 = cube_function_1; + Function f2 = cube_function_2; + + Function_vector v; + v.push_back(f1); + v.push_back(f2); + + std::vector vps; + vps.push_back("--"); + + // Domain (Warning: Sphere_3 constructor uses square radius !) + Mesh_domain_with_features domain(Function_wrapper(v, vps), K::Sphere_3(CGAL::ORIGIN, 5.*5.)); + Polylines polylines; + create_polylines(polylines); + domain.add_features(polylines.begin(),polylines.end()); + + // Set mesh criteria + Mesh_criteria criteria(edge_size = 0.15, + facet_angle = 30, facet_size = 0.2, + cell_radius_edge_ratio = 2, cell_size = 0.4); + + // Mesh generation + C3t3 c3t3 = CGAL::make_mesh_3(domain, criteria, no_exude(), no_perturb()); + + // Perturbation (maximum cpu time: 10s, targeted dihedral angle: default) + CGAL::perturb_mesh_3(c3t3, domain, time_limit = 10); + + // Exudation + CGAL::exude_mesh_3(c3t3,12); + + // Output + std::ofstream medit_file("out_cubes_intersection_with_features.mesh"); + CGAL::output_to_medit(medit_file, c3t3); + + return 0; +} diff --git a/Mesh_3/examples/Mesh_3/mesh_implicit_domains.cpp b/Mesh_3/examples/Mesh_3/mesh_implicit_domains.cpp index 5e62f176986..b3a69f54507 100644 --- a/Mesh_3/examples/Mesh_3/mesh_implicit_domains.cpp +++ b/Mesh_3/examples/Mesh_3/mesh_implicit_domains.cpp @@ -8,15 +8,15 @@ -#include +#include "debug.h" #include #include #include #include -#include -#include +#include +#include #include #include "implicit_functions.h" @@ -28,10 +28,10 @@ using namespace CGAL::parameters; // Domain typedef CGAL::Exact_predicates_inexact_constructions_kernel K; typedef FT_to_point_function_wrapper Function; -typedef CGAL::Mesh_3::Implicit_vector_to_labeled_function_wrapper +typedef CGAL::Implicit_multi_domain_to_labeling_function_wrapper Function_wrapper; typedef Function_wrapper::Function_vector Function_vector; -typedef CGAL::Mesh_3::Labeled_mesh_domain_3 Mesh_domain; +typedef CGAL::Labeled_mesh_domain_3 Mesh_domain; // Triangulation typedef CGAL::Mesh_triangulation_3::type Tr; @@ -51,8 +51,8 @@ int main() Function f3(&sphere_function<2>); Function_vector v; - v.push_back(&f1); - v.push_back(&f2); + v.push_back(f1); + v.push_back(f2); //v.push_back(&f3); // Domain (Warning: Sphere_3 constructor uses square radius !) diff --git a/Mesh_3/examples/Mesh_3/mesh_implicit_domains_2.cpp b/Mesh_3/examples/Mesh_3/mesh_implicit_domains_2.cpp new file mode 100644 index 00000000000..af6efcccdfe --- /dev/null +++ b/Mesh_3/examples/Mesh_3/mesh_implicit_domains_2.cpp @@ -0,0 +1,80 @@ + +//****************************************************************************** +// File Description : +// Outputs to out.mesh a mesh of implicit domains. These domains are defined +// by a vector of functions. Each n-uplet of sign of function values defines a +// subdomain. +//****************************************************************************** + + + +#include + +#include +#include +#include + +#include +#include +#include +#include "implicit_functions.h" + +// IO +#include + +using namespace CGAL::parameters; + +// Domain +typedef CGAL::Exact_predicates_inexact_constructions_kernel K; +typedef FT_to_point_function_wrapper Function; +typedef CGAL::Implicit_multi_domain_to_labeling_function_wrapper + Function_wrapper; +typedef Function_wrapper::Function_vector Function_vector; +typedef CGAL::Labeled_mesh_domain_3 Mesh_domain; + +// Triangulation +typedef CGAL::Mesh_triangulation_3::type Tr; +typedef CGAL::Mesh_complex_3_in_triangulation_3 C3t3; + +// Mesh Criteria +typedef CGAL::Mesh_criteria_3 Mesh_criteria; +typedef Mesh_criteria::Facet_criteria Facet_criteria; +typedef Mesh_criteria::Cell_criteria Cell_criteria; + + +int main() +{ + // Define functions + Function f1(&torus_function); + Function f2(&sphere_function<5>); + + Function_vector v; + v.push_back(f1); + v.push_back(f2); + + std::vector vps; + vps.push_back("--"); + + // Domain (Warning: Sphere_3 constructor uses square radius !) + Mesh_domain domain(Function_wrapper(v, vps), K::Sphere_3(CGAL::ORIGIN, 5.*5.)); + + // Set mesh criteria + Facet_criteria facet_criteria(30, 0.2, 0.02); // angle, size, approximation + Cell_criteria cell_criteria(2., 0.4); // radius-edge ratio, size + Mesh_criteria criteria(facet_criteria, cell_criteria); + + // Mesh generation + C3t3 c3t3 = CGAL::make_mesh_3(domain, criteria, no_exude(), no_perturb()); + + // Perturbation (maximum cpu time: 10s, targeted dihedral angle: default) + CGAL::perturb_mesh_3(c3t3, domain, time_limit = 10); + + // Exudation + CGAL::exude_mesh_3(c3t3,12); + + // Output + std::ofstream medit_file("out.mesh"); + CGAL::output_to_medit(medit_file, c3t3); + + return 0; +} diff --git a/Mesh_3/examples/Mesh_3/mesh_polyhedral_edge_tolerance_region.cpp b/Mesh_3/examples/Mesh_3/mesh_polyhedral_edge_tolerance_region.cpp index 2b6fbacffa0..f1165de5d2e 100644 --- a/Mesh_3/examples/Mesh_3/mesh_polyhedral_edge_tolerance_region.cpp +++ b/Mesh_3/examples/Mesh_3/mesh_polyhedral_edge_tolerance_region.cpp @@ -15,7 +15,7 @@ #include #include -#include +#include #include #include @@ -28,7 +28,7 @@ typedef CGAL::Exact_predicates_inexact_constructions_kernel K; typedef CGAL::Polyhedron_3 Polyhedron; typedef CGAL::Mesh_3::Polyhedral_edge_tolerance_to_labeled_function_wrapper Polyhedral_wrapper; -typedef CGAL::Mesh_3::Labeled_mesh_domain_3 Mesh_domain; +typedef CGAL::Labeled_mesh_domain_3 Mesh_domain; // Triangulation typedef CGAL::Mesh_triangulation_3::type Tr; diff --git a/Mesh_3/examples/Mesh_3/mesh_polyhedral_implicit_function.cpp b/Mesh_3/examples/Mesh_3/mesh_polyhedral_implicit_function.cpp index b260c1bf2d1..e852417fb01 100644 --- a/Mesh_3/examples/Mesh_3/mesh_polyhedral_implicit_function.cpp +++ b/Mesh_3/examples/Mesh_3/mesh_polyhedral_implicit_function.cpp @@ -13,7 +13,7 @@ #include #include -#include +#include #include #include @@ -26,7 +26,7 @@ struct K: public CGAL::Exact_predicates_inexact_constructions_kernel {}; typedef CGAL::Polyhedron_3 Polyhedron; typedef CGAL::Mesh_3::Polyhedral_to_labeled_function_wrapper Polyhedral_wrapper; -typedef CGAL::Mesh_3::Labeled_mesh_domain_3 Mesh_domain; +typedef CGAL::Labeled_mesh_domain_3 Mesh_domain; // Triangulation typedef CGAL::Mesh_triangulation_3::type Tr; diff --git a/Mesh_3/examples/Mesh_3/mesh_polyhedral_surface_tolerance_region.cpp b/Mesh_3/examples/Mesh_3/mesh_polyhedral_surface_tolerance_region.cpp index 826d1432a62..d731dc81c55 100644 --- a/Mesh_3/examples/Mesh_3/mesh_polyhedral_surface_tolerance_region.cpp +++ b/Mesh_3/examples/Mesh_3/mesh_polyhedral_surface_tolerance_region.cpp @@ -14,7 +14,7 @@ #include #include -#include +#include #include #include @@ -27,7 +27,7 @@ typedef CGAL::Exact_predicates_inexact_constructions_kernel K; typedef CGAL::Polyhedron_3 Polyhedron; typedef CGAL::Mesh_3::Polyhedral_tolerance_to_labeled_function_wrapper Polyhedral_wrapper; -typedef CGAL::Mesh_3::Labeled_mesh_domain_3 Mesh_domain; +typedef CGAL::Labeled_mesh_domain_3 Mesh_domain; // Triangulation typedef CGAL::Mesh_triangulation_3::type Tr; diff --git a/Mesh_3/include/CGAL/Implicit_mesh_domain_3.h b/Mesh_3/include/CGAL/Implicit_mesh_domain_3.h index 16238012a35..f8c6cae7283 100644 --- a/Mesh_3/include/CGAL/Implicit_mesh_domain_3.h +++ b/Mesh_3/include/CGAL/Implicit_mesh_domain_3.h @@ -31,8 +31,8 @@ # pragma warning(disable:4180) // qualifier applied to function type has no meaning; ignored #endif -#include -#include +#include +#include namespace CGAL { @@ -45,13 +45,13 @@ namespace CGAL { */ template > + class Wrapper = Implicit_to_labeling_function_wrapper > class Implicit_mesh_domain_3 - : public Mesh_3::Labeled_mesh_domain_3 + : public Labeled_mesh_domain_3 { public: /// Base type - typedef Mesh_3::Labeled_mesh_domain_3 Base; + typedef Labeled_mesh_domain_3 Base; /// Public types typedef typename Base::Sphere_3 Sphere_3; diff --git a/Mesh_3/include/CGAL/Implicit_to_labeling_function_wrapper.h b/Mesh_3/include/CGAL/Implicit_to_labeling_function_wrapper.h new file mode 100644 index 00000000000..cf6b7df1c44 --- /dev/null +++ b/Mesh_3/include/CGAL/Implicit_to_labeling_function_wrapper.h @@ -0,0 +1,286 @@ +// Copyright (c) 2009 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org). +// You can redistribute it and/or modify it under the terms of the GNU +// General Public License as published by the Free Software Foundation, +// either version 3 of the License, or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL$ +// $Id$ +// +// +// Author(s) : Stéphane Tayeb, Aymeric PELLE +// +//****************************************************************************** +// File Description : +// Implicit_to_labeling_function_wrapper and +// Implicit_vector_to_labeling_function_wrapper class declaration +// and implementation. +// +// See classes description to have more information. +//****************************************************************************** + +#ifndef CGAL_IMPLICIT_TO_LABELING_FUNCTION_WRAPPER_H +#define CGAL_IMPLICIT_TO_LABELING_FUNCTION_WRAPPER_H + + +#if defined(BOOST_MSVC) +# pragma warning(push) +# pragma warning(disable:4180) // qualifier applied to function type has no meaning; ignored +#endif + +#include + +namespace CGAL { + +#include + +/** + * @class Implicit_to_labeling_function_wrapper + * + * This class is designed to wrap an implicit function which describes a domain + * by [p is inside if f(p)<0] to a function which takes its values into {0,1}. + * f(p)=0 means that p is outside the domain. + */ +template +class Implicit_to_labeling_function_wrapper +{ +public: + // Types + typedef int return_type; + typedef typename BGT::Point_3 Point_3; + + /// Constructor + Implicit_to_labeling_function_wrapper(const Function_& f) + : r_f_(f) {} + + // Default copy constructor and assignment operator are ok + + /// Destructor + ~Implicit_to_labeling_function_wrapper() {} + + /// Operator () + return_type operator()(const Point_3& p, const bool = true) const + { + return ( (r_f_(p)<0) ? 1 : 0 ); + } + +private: + /// Function to wrap + const Function_& r_f_; + +}; // end class Implicit_to_labeling_function_wrapper + + + +/** + * \deprecated + * + * @class Implicit_vector_to_labeling_function_wrapper + * + * Wraps a set of implicit function [f1,f2,...] to one function F which + * takes its values into N. + * + * Let p be a point. + * F(p) = 0b000000(f2(p)<0)(f1(p)<0) + * + * It can handle at most 8 functions. + */ +template +class Implicit_vector_to_labeling_function_wrapper +{ +public: + // Types + typedef int return_type; + typedef std::vector Function_vector; + typedef typename BGT::Point_3 Point_3; + + /// Constructor + Implicit_vector_to_labeling_function_wrapper(const std::vector& v) + : function_vector_(v) {} + + // Default copy constructor and assignment operator are ok + + /// Destructor + ~Implicit_vector_to_labeling_function_wrapper() {} + + /// Operator () + return_type operator()(const Point_3& p, const bool = true) const + { + int nb_func = function_vector_.size(); + if ( nb_func > 8 ) + { + CGAL_error_msg("We support at most 8 functions !"); + } + + char bits = 0; + for ( int i = 0 ; i < nb_func ; ++i ) + { + // Insert value into bits : we compute fi(p) and insert result at + // bit i of bits + bits |= ( ((*function_vector_[i])(p) < 0) << i ); + } + + return ( static_cast(bits) ); + } + +private: + /// Functions to wrap + const Function_vector function_vector_; + +}; // end class Implicit_to_labeling_function_wrapper + +template +class Implicit_multi_domain_to_labeling_function_wrapper +{ + template + class Implicit_function_traits + { + public: + typedef typename T_::Point Point; + }; + + template + class Implicit_function_traits + { + public: + typedef typename boost::remove_reference< + typename boost::remove_cv< Point_ >::type>::type Point; + }; + +public: + typedef int return_type; + typedef ImplicitFunction Function; + typedef typename Implicit_function_traits::Point Point_3; + typedef std::vector Function_vector; + +private: + std::vector funcs; + std::vector > bmasks; + +public: + Implicit_multi_domain_to_labeling_function_wrapper (const Function_vector& vf, const std::vector >& vps) + : funcs(vf), bmasks(vps.size(), boost::dynamic_bitset<>(funcs.size() * 2, false)) + { + assert(funcs.size()); + + std::size_t mask_index = 0; + for (std::vector >::const_iterator mask_iter = vps.begin(), mask_end_iter = vps.end(); + mask_iter != mask_end_iter; + ++mask_iter) + { + const std::vector& mask = *mask_iter; + assert(funcs.size() == mask.size()); + boost::dynamic_bitset<>& bmask = bmasks[mask_index++]; + + std::size_t bit_index = 0; + for (std::vector::const_iterator iter = mask.begin(), endIter = mask.end(); iter != endIter; ++iter) + { + std::string::value_type character = *iter; + assert(character == POSITIVE || character == NEGATIVE); + + bmask[bit_index] = character == POSITIVE; + ++bit_index; + bmask[bit_index] = character == NEGATIVE; + ++bit_index; + } + } + std::sort(bmasks.begin(), bmasks.end()); + } + + Implicit_multi_domain_to_labeling_function_wrapper (const Function_vector& vf) + : funcs(vf) + { + assert(funcs.size()); + + bmasks.reserve((1 << funcs.size()) - 1); + bmasks.push_back(boost::dynamic_bitset<>(std::string("10"))); + bmasks.push_back(boost::dynamic_bitset<>(std::string("01"))); + + for (std::size_t i = 0; i < funcs.size()-1; ++i) + { + std::size_t c_size = bmasks.size(); + for (std::size_t index = 0; index < c_size; ++index) + { + boost::dynamic_bitset<> aux = bmasks[index]; + aux.push_back(true); + aux.push_back(false); + bmasks.push_back(aux); + bmasks[index].push_back(false); + bmasks[index].push_back(true); + } + } + bmasks.pop_back(); + std::sort(bmasks.begin(), bmasks.end()); + } + + Implicit_multi_domain_to_labeling_function_wrapper (const Function_vector& vf, const std::vector& vps) + : funcs(vf), bmasks(vps.size(), boost::dynamic_bitset<>(funcs.size() * 2, false)) + { + assert(funcs.size()); + + std::size_t mask_index = 0; + for (std::vector::const_iterator mask_iter = vps.begin(), mask_end_iter = vps.end(); + mask_iter != mask_end_iter; + ++mask_iter) + { + const std::string& mask_str = *mask_iter; + assert(funcs.size() == mask_str.length()); + boost::dynamic_bitset<>& bmask = bmasks[mask_index++]; + + std::size_t bit_index = 0; + for (std::string::const_iterator iter = mask_str.begin(), endIter = mask_str.end(); iter != endIter; ++iter) + { + std::string::value_type character = *iter; + assert(character == '+' || character == '-'); + + bmask[bit_index] = character == '+'; + ++bit_index; + bmask[bit_index] = character == '-'; + ++bit_index; + } + } + std::sort(bmasks.begin(), bmasks.end()); + } + + return_type operator() (const Point_3& p, const bool = true) const + { + boost::dynamic_bitset<> bmask(funcs.size() * 2, false); + + std::size_t i = 0; + for (typename std::vector::const_iterator iter = funcs.begin(), endIter = funcs.end(); + iter != endIter; + ++iter) + { + const Function& function = *iter; + + double fres = function(p); + bmask[i] = fres > 0; + ++i; + bmask[i] = fres < 0; + ++i; + } + + std::vector >::const_iterator iter = std::lower_bound(bmasks.begin(), bmasks.end(), bmask); + if (iter != bmasks.end() && *iter == bmask) + return static_cast(1 + (iter - bmasks.begin())); + return 0; + } +}; + +} // end namespace CGAL + + + +#if defined(BOOST_MSVC) +# pragma warning(pop) +#endif + +#endif // CGAL_IMPLICIT_TO_LABELING_FUNCTION_WRAPPER_H diff --git a/Mesh_3/include/CGAL/Labeled_image_mesh_domain_3.h b/Mesh_3/include/CGAL/Labeled_image_mesh_domain_3.h index 74fba59dfd7..b4a2b350949 100644 --- a/Mesh_3/include/CGAL/Labeled_image_mesh_domain_3.h +++ b/Mesh_3/include/CGAL/Labeled_image_mesh_domain_3.h @@ -28,7 +28,7 @@ #define CGAL_LABELED_IMAGE_MESH_DOMAIN_3_H -#include +#include #include @@ -43,10 +43,10 @@ template > class Labeled_image_mesh_domain_3 -: public Mesh_3::Labeled_mesh_domain_3 +: public Labeled_mesh_domain_3 { public: - typedef Mesh_3::Labeled_mesh_domain_3 Base; + typedef Labeled_mesh_domain_3 Base; typedef typename Base::Sphere_3 Sphere_3; typedef typename Base::FT FT; diff --git a/Mesh_3/include/CGAL/Labeled_mesh_domain_3.h b/Mesh_3/include/CGAL/Labeled_mesh_domain_3.h new file mode 100644 index 00000000000..6ac61cc3b35 --- /dev/null +++ b/Mesh_3/include/CGAL/Labeled_mesh_domain_3.h @@ -0,0 +1,605 @@ +// Copyright (c) 2009 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org). +// You can redistribute it and/or modify it under the terms of the GNU +// General Public License as published by the Free Software Foundation, +// either version 3 of the License, or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL$ +// $Id$ +// +// +// Author(s) : Stéphane Tayeb, Aymeric PELLE +// +//****************************************************************************** +// File Description : +// class Labeled_mesh_domain_3. See class description. +//****************************************************************************** + +#ifndef CGAL_LABELED_MESH_DOMAIN_3_H +#define CGAL_LABELED_MESH_DOMAIN_3_H + +#include + +#include + +#include +#include + +#include +#include +#include +#include +#include + +namespace CGAL { + +/** + * \class Labeled_mesh_domain_3 + * + * Function f must take his values into N. + * Let p be a Point. + * - f(p)=0 means that p is outside domain. + * - f(p)=a, a!=0 means that p is inside subdomain a. + * + * Any boundary facet is labelled , a, where b!=0. + */ +template +class Labeled_mesh_domain_3 +{ +public: + /// Geometric object types + typedef typename BGT::Point_3 Point_3; + typedef typename BGT::Segment_3 Segment_3; + typedef typename BGT::Ray_3 Ray_3; + typedef typename BGT::Line_3 Line_3; + typedef typename BGT::Vector_3 Vector_3; + typedef typename BGT::Sphere_3 Sphere_3; + typedef CGAL::Bbox_3 Bbox_3; + +protected: + typedef typename BGT::Iso_cuboid_3 Iso_cuboid_3; + +public: + // Kernel_traits compatibility + typedef BGT R; + + //------------------------------------------------------- + // Index Types + //------------------------------------------------------- + /// Type of indexes for cells of the input complex + typedef typename Function::return_type Subdomain_index; + typedef boost::optional Subdomain; + /// Type of indexes for surface patch of the input complex + typedef std::pair Surface_patch_index; + typedef boost::optional Surface_patch; + /// Type of indexes to characterize the lowest dimensional face of the input + /// complex on which a vertex lie + typedef boost::variant Index; + typedef CGAL::cpp11::tuple Intersection; + + + typedef typename BGT::FT FT; + typedef BGT Geom_traits; + + /** + * @brief Constructor + */ + Labeled_mesh_domain_3(const Function& f, + const Sphere_3& bounding_sphere, + const FT& error_bound = FT(1e-3)); + + Labeled_mesh_domain_3(const Function& f, + const Bbox_3& bbox, + const FT& error_bound = FT(1e-3)); + + Labeled_mesh_domain_3(const Function& f, + const Iso_cuboid_3& bbox, + const FT& error_bound = FT(1e-3)); + + /// Destructor + virtual ~Labeled_mesh_domain_3() {} + + /** + * Constructs a set of \ccc{n} points on the surface, and output them to + * the output iterator \ccc{pts} whose value type is required to be + * \ccc{std::pair}. + */ + struct Construct_initial_points + { + Construct_initial_points(const Labeled_mesh_domain_3& domain) + : r_domain_(domain) {} + + template + OutputIterator operator()(OutputIterator pts, const int n = 12) const; + + private: + const Labeled_mesh_domain_3& r_domain_; + }; + + /// Returns Construct_initial_points object + Construct_initial_points construct_initial_points_object() const + { + return Construct_initial_points(*this); + } + + /** + * Returns true if point~\ccc{p} is in the domain. If \ccc{p} is in the + * domain, the parameter index is set to the index of the subdomain + * including $p$. It is set to the default value otherwise. + */ + struct Is_in_domain + { + Is_in_domain(const Labeled_mesh_domain_3& domain) : r_domain_(domain) {} + + Subdomain operator()(const Point_3& p) const + { + // f(p)==0 means p is outside the domain + Subdomain_index index = (r_domain_.function_)(p); + if ( Subdomain_index() == index ) + return Subdomain(); + else + return Subdomain(index); + } + private: + const Labeled_mesh_domain_3& r_domain_; + }; + + /// Returns Is_in_domain object + Is_in_domain is_in_domain_object() const { return Is_in_domain(*this); } + + /** + * Returns true is the element \ccc{type} intersect properly any of the + * surface patches describing the either the domain boundary or some + * subdomain boundary. + * \ccc{Type} is either \ccc{Segment_3}, \ccc{Ray_3} or \ccc{Line_3}. + * Parameter index is set to the index of the intersected surface patch + * if \ccc{true} is returned and to the default \ccc{Surface_patch_index} + * value otherwise. + */ + struct Do_intersect_surface + { + Do_intersect_surface(const Labeled_mesh_domain_3& domain) + : r_domain_(domain) {} + + Surface_patch operator()(const Segment_3& s) const + { + return this->operator()(s.source(), s.target()); + } + + Surface_patch operator()(const Ray_3& r) const + { + return clip_to_segment(r); + } + + Surface_patch operator()(const Line_3& l) const + { + return clip_to_segment(l); + } + + private: + /// Returns true if points \c a & \c b do not belong to the same subdomain + /// \c index is set to the surface index of subdomains f(a), f(b) + Surface_patch operator()(const Point_3& a, const Point_3& b) const + { + // If f(a) != f(b), then [a,b] intersects some surface. Here we consider + // [a,b] intersects surface_patch labelled (or ). + // It may be false, further rafinement will improve precision + const Subdomain_index value_a = r_domain_.function_(a); + const Subdomain_index value_b = r_domain_.function_(b); + + if ( value_a != value_b ) + return Surface_patch(r_domain_.make_surface_index(value_a, value_b)); + else + return Surface_patch(); + } + + /** + * Clips \c query to a segment \c s, and call operator()(s) + */ + template + Surface_patch clip_to_segment(const Query& query) const + { + typename cpp11::result_of::type + clipped = CGAL::intersection(query, r_domain_.bbox_); + + if(clipped) +#if CGAL_INTERSECTION_VERSION > 1 + if(const Segment_3* s = boost::get(&*clipped)) + return this->operator()(*s); +#else + if(const Segment_3* s = object_cast(&clipped)) + return this->operator()(*s); +#endif + + return Surface_patch(); + } + + private: + const Labeled_mesh_domain_3& r_domain_; + }; + + /// Returns Do_intersect_surface object + Do_intersect_surface do_intersect_surface_object() const + { + return Do_intersect_surface(*this); + } + + /** + * Returns a point in the intersection of the primitive \ccc{type} + * with some boundary surface. + * \ccc{Type1} is either \ccc{Segment_3}, \ccc{Ray_3} or \ccc{Line_3}. + * The integer \ccc{dimension} is set to the dimension of the lowest + * dimensional face in the input complex containing the returned point, and + * \ccc{index} is set to the index to be stored at a mesh vertex lying + * on this face. + */ + struct Construct_intersection + { + Construct_intersection(const Labeled_mesh_domain_3& domain) + : r_domain_(domain) {} + + Intersection operator()(const Segment_3& s) const + { +#ifndef CGAL_MESH_3_NO_LONGER_CALLS_DO_INTERSECT_3 + CGAL_precondition(r_domain_.do_intersect_surface_object()(s)); +#endif // NOT CGAL_MESH_3_NO_LONGER_CALLS_DO_INTERSECT_3 + return this->operator()(s.source(),s.target()); + } + + Intersection operator()(const Ray_3& r) const + { + return clip_to_segment(r); + } + + Intersection operator()(const Line_3& l) const + { + return clip_to_segment(l); + } + + private: + /** + * Returns a point in the intersection of [a,b] with the surface + * \c a must be the source point, and \c b the out point. It's important + * because it drives bisection cuts. + * Indeed, the returned point is the first intersection from \c [a,b] + * with a subdomain surface. + */ + Intersection operator()(const Point_3& a, const Point_3& b) const + { + // Functors + typename BGT::Compute_squared_distance_3 squared_distance = + BGT().compute_squared_distance_3_object(); + typename BGT::Construct_midpoint_3 midpoint = + BGT().construct_midpoint_3_object(); + + // Non const points + Point_3 p1 = a; + Point_3 p2 = b; + Point_3 mid = midpoint(p1, p2); + + // Cannot be const: those values are modified below. + Subdomain_index value_at_p1 = r_domain_.function_(p1); + Subdomain_index value_at_p2 = r_domain_.function_(p2); + Subdomain_index value_at_mid = r_domain_.function_(mid,true); + + // If both extremities are in the same subdomain, + // there is no intersection. + // This should not happen... + if( value_at_p1 == value_at_p2 ) + { + return Intersection(); + } + + // Construct the surface patch index and index from the values at 'a' + // and 'b'. Even if the bissection find out a different pair of + // values, the reported index will be constructed from the initial + // values. + const Surface_patch_index sp_index = + r_domain_.make_surface_index(value_at_p1, value_at_p2); + const Index index = r_domain_.index_from_surface_patch_index(sp_index); + + // Else lets find a point (by bisection) + // Bisection ends when the point is near than error bound from surface + while(true) + { + // If the two points are enough close, then we return midpoint + if ( squared_distance(p1, p2) < r_domain_.squared_error_bound_ ) + { + CGAL_assertion(value_at_p1 != value_at_p2); + return Intersection(mid, index, 2); + } + + // Else we must go on + // Here we consider that p1(a) is the source point. Thus, we keep p1 and + // change p2 if f(p1)!=f(p2). + // That allows us to find the first intersection from a of [a,b] with + // a surface. + if ( value_at_p1 != value_at_mid ) + { + p2 = mid; + value_at_p2 = value_at_mid; + } + else + { + p1 = mid; + value_at_p1 = value_at_mid; + } + + mid = midpoint(p1, p2); + value_at_mid = r_domain_.function_(mid,true); + } + } + + /// Clips \c query to a segment \c s, and call operator()(s) + template + Intersection clip_to_segment(const Query& query) const + { + typename cpp11::result_of::type + clipped = CGAL::intersection(query, r_domain_.bbox_); + + if(clipped) +#if CGAL_INTERSECTION_VERSION > 1 + if(const Segment_3* s = boost::get(&*clipped)) + return this->operator()(*s); +#else + if(const Segment_3* s = object_cast(&clipped)) + return this->operator()(*s); +#endif + + return Intersection(); + } + + private: + const Labeled_mesh_domain_3& r_domain_; + }; + + /// Returns Construct_intersection object + Construct_intersection construct_intersection_object() const + { + return Construct_intersection(*this); + } + + /** + * Returns the index to be stored in a vertex lying on the surface identified + * by \c index. + */ + Index index_from_surface_patch_index(const Surface_patch_index& index) const + { return Index(index); } + + /** + * Returns the index to be stored in a vertex lying in the subdomain + * identified by \c index. + */ + Index index_from_subdomain_index(const Subdomain_index& index) const + { return Index(index); } + + /** + * Returns the \c Surface_patch_index of the surface patch + * where lies a vertex with dimension 2 and index \c index. + */ + Surface_patch_index surface_patch_index(const Index& index) const + { return boost::get(index); } + + /** + * Returns the index of the subdomain containing a vertex + * with dimension 3 and index \c index. + */ + Subdomain_index subdomain_index(const Index& index) const + { return boost::get(index); } + + // ----------------------------------- + // Backward Compatibility + // ----------------------------------- +#ifndef CGAL_MESH_3_NO_DEPRECATED_SURFACE_INDEX + typedef Surface_patch_index Surface_index; + + Index index_from_surface_index(const Surface_index& index) const + { return index_from_surface_patch_index(index); } + + Surface_index surface_index(const Index& index) const + { return surface_patch_index(index); } +#endif // CGAL_MESH_3_NO_DEPRECATED_SURFACE_INDEX + // ----------------------------------- + // End backward Compatibility + // ----------------------------------- + + +private: + /// Returns Surface_patch_index from \c i and \c j + Surface_patch_index make_surface_index(const Subdomain_index i, + const Subdomain_index j) const + { + if ( i < j ) return Surface_patch_index(i,j); + else return Surface_patch_index(j,i); + } + + /// Returns squared error bound from \c bbox and \c error + FT squared_error_bound(const Iso_cuboid_3& bbox, const FT& error) const + { + typename BGT::Compute_squared_distance_3 squared_distance = + BGT().compute_squared_distance_3_object(); + return squared_distance((bbox.min)(), (bbox.max)())*error*error/4; + } + + /// Returns squared error bound from \c sphere and \c error + FT squared_error_bound(const Sphere_3& sphere, const FT& error) const + { + typename BGT::Compute_squared_radius_3 squared_radius = + BGT().compute_squared_radius_3_object(); + return squared_radius(sphere)*error*error; + } + + /// Returns the bounding sphere of an Iso_cuboid_3 + Sphere_3 bounding_sphere(const Iso_cuboid_3& bbox) const + { + typename BGT::Construct_sphere_3 sphere = BGT().construct_sphere_3_object(); + return sphere((bbox.min)(), (bbox.max)()); + } + + /// Returns and Iso_cuboid_3 from a Bbox_3 + Iso_cuboid_3 iso_cuboid(const Bbox_3& bbox) + { + const Point_3 p_min(bbox.xmin(), bbox.ymin(), bbox.zmin()); + const Point_3 p_max(bbox.xmax(), bbox.ymax(), bbox.zmax()); + + return Iso_cuboid_3(p_min,p_max); + } + +protected: + /// Returns bounding box + const Iso_cuboid_3& bounding_box() const { return bbox_; } + +private: + /// The function which answers subdomain queries + const Function function_; + /// The bounding box + const Iso_cuboid_3 bbox_; + /// Error bound relative to sphere radius + FT squared_error_bound_; + + +private: + // Disabled copy constructor & assignment operator + typedef Labeled_mesh_domain_3 Self; + Labeled_mesh_domain_3(const Self& src); + Self& operator=(const Self& src); + +}; // end class Labeled_mesh_domain_3 + + + + +//------------------------------------------------------- +// Method implementation +//------------------------------------------------------- +template +Labeled_mesh_domain_3::Labeled_mesh_domain_3( + const F& f, + const Sphere_3& bounding_sphere, + const FT& error_bound ) +: function_(f) +, bbox_(iso_cuboid(bounding_sphere.bbox())) +, squared_error_bound_(squared_error_bound(bounding_sphere,error_bound)) +{ + // TODO : CGAL_ASSERT(0 < f(bounding_sphere.get_center()) ) ? +} + +template +Labeled_mesh_domain_3::Labeled_mesh_domain_3( + const F& f, + const Bbox_3& bbox, + const FT& error_bound ) +: function_(f) +, bbox_(iso_cuboid(bbox)) +, squared_error_bound_(squared_error_bound(bbox_,error_bound)) +{ + // TODO : CGAL_ASSERT(0 < f(bounding_sphere.get_center()) ) ? +} + +template +Labeled_mesh_domain_3::Labeled_mesh_domain_3( + const F& f, + const Iso_cuboid_3& bbox, + const FT& error_bound ) +: function_(f) +, bbox_(bbox) +, squared_error_bound_(squared_error_bound(bbox_,error_bound)) +{ + // TODO : CGAL_ASSERT(0 < f( bbox.get_center()) ) ? +} + + + +template +template +OutputIterator +Labeled_mesh_domain_3::Construct_initial_points::operator()( + OutputIterator pts, + const int nb_points) const +{ + // Create point_iterator on and in bounding_sphere + typedef Random_points_on_sphere_3 Random_points_on_sphere_3; + typedef Random_points_in_sphere_3 Random_points_in_sphere_3; + + + const FT squared_radius = BGT().compute_squared_radius_3_object()( + r_domain_.bounding_sphere(r_domain_.bbox_)); + + const double radius = std::sqrt(CGAL::to_double(squared_radius)); + + Random_points_on_sphere_3 random_point_on_sphere(radius); + Random_points_in_sphere_3 random_point_in_sphere(radius); + + // Get some functors + typename BGT::Construct_segment_3 segment_3 = + BGT().construct_segment_3_object(); + typename BGT::Construct_vector_3 vector_3 = + BGT().construct_vector_3_object(); + typename BGT::Construct_translated_point_3 translate = + BGT().construct_translated_point_3_object(); + typename BGT::Construct_center_3 center = BGT().construct_center_3_object(); + + // Get translation from origin to sphere center + Point_3 center_pt = center(r_domain_.bounding_sphere(r_domain_.bbox_)); + const Vector_3 sphere_translation = vector_3(CGAL::ORIGIN, center_pt); + + // Create nb_point points + int n = nb_points; +#ifdef CGAL_MESH_3_VERBOSE + std::cerr << "construct initial points:\n"; +#endif + while ( 0 != n ) + { + // Get a random segment + const Point_3 random_point = translate(*random_point_on_sphere, + sphere_translation); + const Segment_3 random_seg = segment_3(center_pt, random_point); + + // Add the intersection to the output if it exists + Surface_patch surface = r_domain_.do_intersect_surface_object()(random_seg); + if ( surface ) + { + const Point_3 intersect_pt = CGAL::cpp11::get<0>( + r_domain_.construct_intersection_object()(random_seg)); + *pts++ = std::make_pair(intersect_pt, + r_domain_.index_from_surface_patch_index(*surface)); + --n; + +#ifdef CGAL_MESH_3_VERBOSE + std::cerr << boost::format("\r \r" + "%1%/%2% initial point(s) found...") + % (nb_points - n) + % nb_points; +#endif + } + else + { + // Get a new random point into sphere as center of object + // It may be necessary if the center of the domain is empty, e.g. torus + // In general case, it is good for input point dispersion + ++random_point_in_sphere; + center_pt = translate(*random_point_in_sphere, sphere_translation); + } + ++random_point_on_sphere; + } + +#ifdef CGAL_MESH_3_VERBOSE + std::cerr << "\n"; +#endif + return pts; +} + + +} // end namespace CGAL + +#endif // LABELLED_MESH_TRAITS_3_H_ diff --git a/Mesh_3/include/CGAL/Mesh_3/Implicit_to_labeled_function_wrapper.h b/Mesh_3/include/CGAL/Mesh_3/Implicit_to_labeled_function_wrapper.h index ef0eab9f6cb..9f16d01426e 100644 --- a/Mesh_3/include/CGAL/Mesh_3/Implicit_to_labeled_function_wrapper.h +++ b/Mesh_3/include/CGAL/Mesh_3/Implicit_to_labeled_function_wrapper.h @@ -20,8 +20,8 @@ // //****************************************************************************** // File Description : -// Implicit_to_labeled_function_wrapper and -// Implicit_vector_to_labeled_function_wrapper class declaration +// Implicit_to_labeling_function_wrapper and +// Implicit_vector_to_labeling_function_wrapper class declaration // and implementation. // // See classes description to have more information. @@ -30,12 +30,17 @@ #ifndef CGAL_MESH_3_IMPLICIT_TO_LABELED_FUNCTION_WRAPPER_H #define CGAL_MESH_3_IMPLICIT_TO_LABELED_FUNCTION_WRAPPER_H - #if defined(BOOST_MSVC) # pragma warning(push) # pragma warning(disable:4180) // qualifier applied to function type has no meaning; ignored #endif +#define CGAL_DEPRECATED_HEADER "" +#define CGAL_REPLACEMENT_HEADER "" +#include + +#include + namespace CGAL { namespace Mesh_3 { @@ -81,6 +86,8 @@ private: /** + * \deprecated + * * @class Implicit_vector_to_labeled_function_wrapper * * Wraps a set of implicit function [f1,f2,...] to one function F which @@ -112,7 +119,7 @@ public: /// Operator () return_type operator()(const Point_3& p, const bool = true) const { - int nb_func = function_vector_.size(); + int nb_func = static_cast(function_vector_.size()); if ( nb_func > 8 ) { CGAL_error_msg("We support at most 8 functions !"); @@ -135,7 +142,6 @@ private: }; // end class Implicit_to_labeled_function_wrapper - } // end namespace Mesh_3 } // end namespace CGAL diff --git a/Mesh_3/include/CGAL/Mesh_3/Labeled_mesh_domain_3.h b/Mesh_3/include/CGAL/Mesh_3/Labeled_mesh_domain_3.h index 92062a7c1c7..5279ed3e92b 100644 --- a/Mesh_3/include/CGAL/Mesh_3/Labeled_mesh_domain_3.h +++ b/Mesh_3/include/CGAL/Mesh_3/Labeled_mesh_domain_3.h @@ -26,6 +26,10 @@ #ifndef CGAL_MESH_3_LABELED_MESH_DOMAIN_3_H #define CGAL_MESH_3_LABELED_MESH_DOMAIN_3_H +#define CGAL_DEPRECATED_HEADER "" +#define CGAL_REPLACEMENT_HEADER "" +#include + #include #include @@ -68,6 +72,9 @@ public: typedef typename BGT::Sphere_3 Sphere_3; typedef CGAL::Bbox_3 Bbox_3; + typedef typename BGT::Iso_cuboid_3 Iso_cuboid_3; + +public: // Kernel_traits compatibility typedef BGT R; @@ -408,9 +415,6 @@ public: // ----------------------------------- -private: - typedef typename BGT::Iso_cuboid_3 Iso_cuboid_3; - private: /// Returns Surface_patch_index from \c i and \c j Surface_patch_index make_surface_index(const Subdomain_index i, @@ -504,7 +508,6 @@ Labeled_mesh_domain_3::Labeled_mesh_domain_3( } - template template OutputIterator diff --git a/Mesh_3/test/Mesh_3/CMakeLists.txt b/Mesh_3/test/Mesh_3/CMakeLists.txt index 547412ea87f..9b791b1cde3 100644 --- a/Mesh_3/test/Mesh_3/CMakeLists.txt +++ b/Mesh_3/test/Mesh_3/CMakeLists.txt @@ -25,10 +25,15 @@ if ( CGAL_FOUND ) create_single_source_cgal_program( "test_backward_compatibility.cpp" ) create_single_source_cgal_program( "test_boost_has_xxx.cpp" ) create_single_source_cgal_program( "test_c3t3.cpp" ) + create_single_source_cgal_program( "test_implicit_multi_domain_to_labeling_function_wrapper.cpp" ) + create_single_source_cgal_program( "test_mesh_3_implicit_vector_to_labeled_function_wrapper.cpp" ) create_single_source_cgal_program( "test_c3t3_io.cpp" ) create_single_source_cgal_program( "test_c3t3_with_features.cpp" ) create_single_source_cgal_program( "test_criteria.cpp" ) create_single_source_cgal_program( "test_domain_with_polyline_features.cpp" ) + create_single_source_cgal_program( "test_mesh_3_labeled_mesh_domain_3.cpp" ) + create_single_source_cgal_program( "test_mesh_implicit_domains.cpp" "implicit_functions.cpp" ) + create_single_source_cgal_program( "test_labeled_mesh_domain_3.cpp" ) create_single_source_cgal_program( "test_mesh_criteria_creation.cpp" ) if(CGAL_ImageIO_USE_ZLIB) create_single_source_cgal_program( "test_meshing_3D_image.cpp" ) diff --git a/Mesh_3/test/Mesh_3/debug.h b/Mesh_3/test/Mesh_3/debug.h new file mode 100644 index 00000000000..690f03b204f --- /dev/null +++ b/Mesh_3/test/Mesh_3/debug.h @@ -0,0 +1,23 @@ +#ifdef BENCH +# define NDEBUG +#endif + +#define CGAL_MESH_3_DEBUG_BEFORE_CONFLICTS 0 +#define CGAL_MESH_3_DEBUG_AFTER_NO_INSERTION 0 +#define CGAL_MESH_3_DEBUG_AFTER_INSERTION 0 +#define CGAL_MESH_3_DEBUG_DOUBLE_MAP 0 + +// #define CGAL_MESHES_NO_OUTPUT + +// #define CGAL_SURFACE_MESHER_EDGES_DEBUG_INTERSECTION +// #define CGAL_SURFACE_MESHER_DEBUG_POLYHEDRAL_SURFACE_CONSTRUCTION +// #define CGAL_SURFACE_MESHER_VERBOSE +// #define CGAL_POLYHEDRAL_SURFACE_VERBOSE_CONSTRUCTION +// #define CGAL_SURFACE_MESHER_EDGES_DEBUG_INSERTIONS +// #define CGAL_MESHES_DEBUG_REFINEMENT_POINTS +// #define CGAL_DEBUG_OCTREE_CONSTRUCTION + +#define CGAL_MESH_3_VERBOSE +#define CGAL_MESH_3_IO_VERBOSE + +// #define OPTIMIZE_INTERSECTION diff --git a/Mesh_3/test/Mesh_3/implicit_functions.cpp b/Mesh_3/test/Mesh_3/implicit_functions.cpp new file mode 100644 index 00000000000..0ede1606d8b --- /dev/null +++ b/Mesh_3/test/Mesh_3/implicit_functions.cpp @@ -0,0 +1,303 @@ +///////////////// Code for functions of famous surfaces ///////////////// + + +// Sphere (r=1) + +double sphere_function (double x, double y, double z) { + double x2=x*x, y2=y*y, z2=z*z; + + return x2+y2+z2-1; +} + + +// Ellipsoid (r=1) + +double ellipsoid_function (double x, double y, double z) { + double x2=x*x, y2=y*y, z2=z*z; + + return x2+2*y2+4*z2-1; +} + + +// Torus (r=2) + +double torus_function (double x, double y, double z) { + double x2=x*x, y2=y*y, z2=z*z; + double x4=x2*x2, y4=y2*y2, z4=z2*z2; + + return x4 + y4 + z4 + 2 *x2* y2 + 2* + x2*z2 + 2*y2* z2 - 5 *x2 + 4* y2 - 5*z2+4; +} + + +// "Chair" (r=6) + +double chair_function (double x, double y, double z) { + double x2=x*x, y2=y*y, z2=z*z; + double x4=x2*x2, y4=y2*y2, z4=z2*z2; + + return x4-1.2*x2*y2+3.6*x2*z2-7.50*x2+y4+3.6*y2*z2-7.50*y2+.2*z4-7.50*z2+64.0625-16.0*z*y2+16.0*x2*z; +} + + +// "Tanglecube" (r=4) + +double tanglecube_function (double x, double y, double z) { + double x2=x*x, y2=y*y, z2=z*z; + double x4=x2*x2, y4=y2*y2, z4=z2*z2; + + return x4 - 5*x2 + y4 - 5*y2 + z4 - 5*z2 + 11.8; +} + +double cube_function (double x, double y, double z){ + if( x < 1 && -x < 1 && + y < 1 && -y < 1 && + z < 1 && -z < 1 ) + return -1.; + return 1.; +} + +// Barth's octic surface (degree 8) + +double octic_function (double x, double y, double z) { // r=2 + double x2=x*x, y2=y*y, z2=z*z; + double x4=x2*x2, y4=y2*y2, z4=z2*z2; + double x6=x4*x2, y6=y4*y2, z6=z4*z2; + double x8=x4*x4, y8=y4*y4, z8=z4*z4; + + return 43.30495169 *x2* y2 + 43.30495169 *x2* z2 + 43.30495169 *y2 * z2 + 44.36067980 *x6* y2 + 44.36067980* x6 * z2 + 66.54101970* x4* y4 + 66.54101970* x4 * z4 + 44.36067980 *x2 * y6 - 11.70820393* x2 - 11.70820393* y2 - 11.70820393* z2 + 37.65247585 *x4 + 37.65247585 *y4 + 37.65247585* z4 + 11.09016995* x8 + 11.09016995 *y8 + 11.09016995* z8 + 133.0820394 *x2* y4* z2 + 133.0820394 *x2 * y2 * z4 + 44.36067980* x2 * z6 + 44.36067980 *y6 * z2 + 66.54101970 *y4 * z4 + 44.36067980 *y2 * z6 + 133.0820394* x4* y2 * z2 - 91.95742756 *x4 * y2 - 91.95742756 *x4 *z2 - 91.95742756* x2 * y4 - 91.95742756 *x2 *z4 - 91.95742756* y4 * z2 - 183.9148551* x2 *y2 *z2 - 30.65247585 *x6 - 30.65247585* y6 - 91.95742756 *y2 * z4 - 30.65247585* z6 + 1.618033988; +} + + +// "Heart" + +double heart_function (double x, double y, double z) { // radius = 2 + + return (2*x*x+y*y+z*z-1)*(2*x*x+y*y+z*z-1)*(2*x*x+y*y+z*z-1) - (0.1*x*x+y*y)*z*z*z; +} + + +// Klein's bottle + +double klein_function (double x, double y, double z) { // radius = 4 + + return (x*x+y*y+z*z+2*y-1)*((x*x+y*y+z*z-2*y-1)*(x*x+y*y+z*z-2*y-1)-8*z*z)+16*x*z*(x*x+y*y+z*z-2*y-1); +} + + +// Ring + +double ring_function (double x, double y, double z) { // radius = ? + double e=0.1; + + double f1 = x*x+y*y+z*z-1; + double f2 = x; + + f1 = f1*f1-e*e; + f2 = f2*f2-e*e; + + if (f1 < 0 && f2 < 0) + return -1.; + else if (f1 > 0 || f2 > 0) + return 1.; + else + return 0.; +} + + +// False knot + +double false_knot_function (double x, double y, double z) { // radius = 1 + double d=1.2, e=0.1; + + double f1 = x*(x*x-z*z)-2*x*z*z-y*y+d*d-x*x-z*z-y*y; + + double m2 = z*(x*x-z*z)+2*x*x*z; + double f2 = 4*y*y*(d*d-x*x-z*z-y*y) - m2*m2; + + f1 = f1*f1-e*e; + f2 = f2*f2-e*e; + + if (f1 < 0 && f2 < 0) + return -1.; + else if (f1 > 0 || f2 > 0) + return 1.; + else + return 0.; +} + + +// Knot 1 + +void puiss(double& x, double& y, int n) { + + double xx = 1, yy = 0; + + while(n>0) { + if (n&1) { + double xxx = xx, yyy = yy; + xx = xxx*x - yyy*y; + yy = xxx*y + yyy*x; + } + + double xxx = x, yyy = y; + x=xxx*xxx-yyy*yyy; + y=2*xxx*yyy; + + n/=2; + } + + x = xx; + y = yy; +} + + + +double knot1_function (double a, double b, double c) { // radius = 4 + double e=0.1; + + double x, y, z, t, den; + den=1+a*a+b*b+c*c; + x=2*a/den; + y=2*b/den; + z=2*c/den; + t=(1-a*a-b*b-c*c)/den; + + double x3=x, y3=y, z2=z, t2=t; + puiss(x3,y3,3); + puiss(z2,t2,2); + + double f1 = z2-x3; + double f2 = t2-y3; + + f1 = f1*f1; + f2 = f2*f2; + e=e*e/(den-1); + + if (f1+f2-e < 0) + return -1.; + else if (f1+f2-e > 0) + return 1.; + else + return 0.; +} + +/* +double knot1_function (double x, double y, double z) { // radius = 4 + double e=0.1; + + double c1, c2, c3, c4, den; + den=1+x*x+y*y+z*z; + c1=2*x/den; + c2=2*y/den; + c3=2*z/den; + c4=(1-x*x-y*y-z*z)/den; + + double f1 = c1*(c1*c1-c2*c2)-2*c1*c2*c2-c3*c3+c4*c4; + double f2 = c2*(c1*c1-c2*c2)+2*c1*c1*c2-2*c3*c4; + + f1 = f1*f1; + f2 = f2*f2; + e=e*e/(den-1); + + if (f1+f2-e < 0) + return -1.; + else if (f1+f2-e > 0) + return 1.; + else + return 0.; +} +*/ + +double knot1_surf1_function (double x, double y, double z) { // radius = 4 + double c1, c2, c3, c4, den; + den=1+x*x+y*y+z*z; + c1=2*x/den; + c2=2*y/den; + c3=2*z/den; + c4=(1-x*x-y*y-z*z)/den; + + return c1*(c1*c1-c2*c2)-2*c1*c2*c2-c3*c3+c4*c4; +} + + +double knot1_surf2_function (double x, double y, double z) { // radius = 4 + double c1, c2, c3, c4, den; + den=1+x*x+y*y+z*z; + c1=2*x/den; + c2=2*y/den; + c3=2*z/den; + c4=(1-x*x-y*y-z*z)/den; + + return c2*(c1*c1-c2*c2)+2*c1*c1*c2-2*c3*c4; +} + + +// Knot 2 + +double knot2_function (double a, double b, double c) { // radius = 4 + double e=0.025; + + double x, y, z, t, den; + den=1+a*a+b*b+c*c; + x=2*a/den; + y=2*b/den; + z=2*c/den; + t=(1-a*a-b*b-c*c)/den; + + double x7=x, y7=y, x13=x, y13=y; + puiss(x7,y7,7); + puiss(x13,y13,13); + + double z3=z, t3=t; + puiss(z3,t3,3); + + double f1t = (z3-x7)*(z3-x7+100*x13) - (t3-y7)*(t3-y7+100*y13); + double f2t = (z3-x7)*(t3-y7+100*y13) + (t3-y7)*(z3-x7+100*x13); + double f1 = f1t*(z3-x7-100*x13) - f2t*(t3-y7-100*y13); + double f2 = f1t*(t3-y7-100*y13) + f2t*(z3-x7-100*x13); + + f1 = f1*f1; + f2 = f2*f2; + e=e*e/(den-1); + + if (f1+f2-e < 0) + return -1.; + else if (f1+f2-e > 0) + return 1.; + else + return 0.; +} + + +// Knot 3 + + +double knot3_function (double a, double b, double c) { // radius = 4 + double e=0.025; + + double x, y, z, t, den; + den=1+a*a+b*b+c*c; + x=2*a/den; + y=2*b/den; + z=2*c/den; + t=(1-a*a-b*b-c*c)/den; + + double x19=x, y19=y, z17=z, t17=t; + puiss(x19,y19,19); + puiss(z17,t17,17); + + double f1 = z17-x19; + double f2 = t17-y19; + + f1 = f1*f1; + f2 = f2*f2; + e=e*e/(den-1); + + if (f1+f2-e < 0) + return -1.; + else if (f1+f2-e > 0) + return 1.; + else + return 0.; +} diff --git a/Mesh_3/test/Mesh_3/implicit_functions.h b/Mesh_3/test/Mesh_3/implicit_functions.h new file mode 100644 index 00000000000..e9a3ce0d05b --- /dev/null +++ b/Mesh_3/test/Mesh_3/implicit_functions.h @@ -0,0 +1,36 @@ +///////////////// Definitions of several famous surfaces ///////////////// +double sphere_function (double, double, double); // (c=(0,0,0), r=1) +double ellipsoid_function (double, double, double); // (c=(0,0,0), r=1) +double torus_function (double, double, double); // (c=(0,0,0), r=2) +double chair_function (double, double, double); // (c=(0,0,0), r=6) +double tanglecube_function (double, double, double); // (c=(0,0,0), r=4) +double octic_function (double, double, double); // (c=(0,0,0), r=2) +double heart_function (double, double, double); // (c=(0,0,0), r=2) +double klein_function (double, double, double); // (c=(0,0,0), r=4) +double ring_function (double, double, double); // (c=(0,0,0), r=?) +double false_knot_function (double, double, double); // (c=(0,0,0), r=1) +double knot1_function (double, double, double); // (c=(0,0,0), r=4) +double knot2_function (double, double, double); // (c=(0,0,0), r=4) +double knot3_function (double, double, double); // (c=(0,0,0), r=4) +double cube_function (double, double, double); // (c=(0,0,0), r=2) + + +template +double sphere_function (double x, double y, double z) // (c=(0,0,0), r=Sq_radius) +{ + double x2=x*x, y2=y*y, z2=z*z; + return (x2+y2+z2)/Sq_radius - 1; +} + + + +template +class FT_to_point_function_wrapper : public std::unary_function +{ + typedef FT (*Implicit_function)(FT, FT, FT); + Implicit_function function; +public: + typedef P Point; + FT_to_point_function_wrapper(Implicit_function f) : function(f) {} + FT operator()(Point p) const { return function(p.x(), p.y(), p.z()); } +}; diff --git a/Mesh_3/test/Mesh_3/test_implicit_multi_domain_to_labeling_function_wrapper.cpp b/Mesh_3/test/Mesh_3/test_implicit_multi_domain_to_labeling_function_wrapper.cpp new file mode 100644 index 00000000000..7607ea50cd0 --- /dev/null +++ b/Mesh_3/test/Mesh_3/test_implicit_multi_domain_to_labeling_function_wrapper.cpp @@ -0,0 +1,141 @@ +#include +#include + + +typedef CGAL::Exact_predicates_inexact_constructions_kernel K_e_i; +typedef K_e_i::Point_3 Point_3; +typedef double (Function) (const Point_3&); + +typedef CGAL::Implicit_multi_domain_to_labeling_function_wrapper Labeling_function; +typedef Labeling_function::Function_vector Function_vector; + + +double cube_function_1 (const Point_3& p) +{ + if( p.x() > 0 && p.x() < 2 && + p.y() > 0 && p.y() < 2 && + p.z() > 0 && p.z() < 2 ) + return -1.; + return 1.; +} + +double cube_function_2 (const Point_3& p) +{ + if( p.x() > 1 && p.x() < 3 && + p.y() > 1 && p.y() < 3 && + p.z() > 1 && p.z() < 3 ) + return -1.; + return 1.; +} + + +void test_constructor_vfunc () +{ + Function_vector vfunc; + vfunc.push_back(cube_function_1); + vfunc.push_back(cube_function_2); + + Labeling_function lfunc(vfunc); + + Point_3 p1(0.5, 0.5, 0.5); + Point_3 p2(2.5, 2.5, 2.5); + Point_3 p3(1.5, 1.5, 1.5); + Point_3 p_out(4., 4., 4.); + + Labeling_function::return_type rp1 = lfunc(p1); + Labeling_function::return_type rp2 = lfunc(p2); + Labeling_function::return_type rp3 = lfunc(p3); + Labeling_function::return_type rp_out = lfunc(p_out); + + assert(rp1 != 0); + assert(rp2 != 0); + assert(rp3 != 0); + assert(rp_out == 0); + + assert(rp1 != rp2); + assert(rp1 != rp3); + assert(rp2 != rp3); +} + +void test_constructor_vfunc_vvpos () +{ + Function_vector vfunc; + vfunc.push_back(cube_function_1); + vfunc.push_back(cube_function_2); + + std::vector > vvpos; + std::vector vpos1; + vpos1.push_back(CGAL::NEGATIVE); + vpos1.push_back(CGAL::POSITIVE); + vvpos.push_back(vpos1); + std::vector vpos2; + vpos2.push_back(CGAL::POSITIVE); + vpos2.push_back(CGAL::NEGATIVE); + vvpos.push_back(vpos2); + std::vector vpos3; + vpos3.push_back(CGAL::NEGATIVE); + vpos3.push_back(CGAL::NEGATIVE); + vvpos.push_back(vpos3); + + Labeling_function lfunc(vfunc, vvpos); + + Point_3 p1(0.5, 0.5, 0.5); + Point_3 p2(2.5, 2.5, 2.5); + Point_3 p3(1.5, 1.5, 1.5); + Point_3 p_out(4., 4., 4.); + + Labeling_function::return_type rp1 = lfunc(p1); + Labeling_function::return_type rp2 = lfunc(p2); + Labeling_function::return_type rp3 = lfunc(p3); + Labeling_function::return_type rp_out = lfunc(p_out); + + assert(rp1 != 0); + assert(rp2 != 0); + assert(rp3 != 0); + assert(rp_out == 0); + + assert(rp1 != rp2); + assert(rp1 != rp3); + assert(rp2 != rp3); +} + +void test_constructor_vfunc_vstr () +{ + Function_vector vfunc; + vfunc.push_back(cube_function_1); + vfunc.push_back(cube_function_2); + + std::vector vstr; + vstr.push_back("-+"); + vstr.push_back("+-"); + vstr.push_back("--"); + + Labeling_function lfunc(vfunc, vstr); + + Point_3 p1(0.5, 0.5, 0.5); + Point_3 p2(2.5, 2.5, 2.5); + Point_3 p3(1.5, 1.5, 1.5); + Point_3 p_out(4., 4., 4.); + + Labeling_function::return_type rp1 = lfunc(p1); + Labeling_function::return_type rp2 = lfunc(p2); + Labeling_function::return_type rp3 = lfunc(p3); + Labeling_function::return_type rp_out = lfunc(p_out); + + assert(rp1 != 0); + assert(rp2 != 0); + assert(rp3 != 0); + assert(rp_out == 0); + + assert(rp1 != rp2); + assert(rp1 != rp3); + assert(rp2 != rp3); +} + +int main () +{ + test_constructor_vfunc(); + test_constructor_vfunc_vvpos(); + test_constructor_vfunc_vstr(); + return 0; +} diff --git a/Mesh_3/test/Mesh_3/test_labeled_mesh_domain_3.cpp b/Mesh_3/test/Mesh_3/test_labeled_mesh_domain_3.cpp new file mode 100644 index 00000000000..62142278690 --- /dev/null +++ b/Mesh_3/test/Mesh_3/test_labeled_mesh_domain_3.cpp @@ -0,0 +1,271 @@ +// Copyright (c) 2009 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org). +// You can redistribute it and/or modify it under the terms of the GNU +// General Public License as published by the Free Software Foundation, +// either version 3 of the License, or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL$ +// $Id$ +// +// +// Author(s) : Stephane Tayeb +// +//****************************************************************************** +// File Description : +//****************************************************************************** + +#include "test_meshing_utilities.h" +#include +#include + + +template +struct LM3_tester +{ + typedef typename K::Point_3 Point_3; + typedef typename K::FT FT; + + struct Function + { + typedef Point_3 Point; + typedef FT (*Func)(const Point_3&); + + Function (Func f) : f_(f) {} + FT operator() (const Point_3& p) const { return f_(p); } + + private: + Func f_; + }; + + typedef CGAL::Implicit_to_labeling_function_wrapper Function_wrapper; + typedef CGAL::Labeled_mesh_domain_3 Mesh_domain; + + static FT shape_function (const Point_3& p) + { + if (p.x() < 0) + return -1; + const FT x2=p.x()*p.x(), y2=p.y()*p.y(), z2=p.z()*p.z(); + return x2 + y2 + z2 - 1; + } + + static FT sphere_function (const Point_3& p) + { + const FT x2=p.x()*p.x(), y2=p.y()*p.y(), z2=p.z()*p.z(); + return x2 + y2 + z2 - 1; + } + + void operator() () const + { + typedef typename K::Sphere_3 Sphere_3; + typedef typename K::Iso_cuboid_3 Iso_cuboid_3; + + test_domain(Sphere_3(CGAL::ORIGIN, 4.)); + test_domain(CGAL::Bbox_3(-2.,-2.,-2., 2.,2.,2.)); + test_domain(Iso_cuboid_3(Point_3(-2.,-2.,-2.), Point_3(2.,2.,2.))); + } + +private: + template + void test_domain (const BoundingShape& bounding_shape) const + { + FT error_bound(1e-3); + + Function f_sphere(&sphere_function); + Function_wrapper wrapper_1(f_sphere); + Mesh_domain domain(wrapper_1, bounding_shape, error_bound); + test_construct_initial_points(domain, error_bound); + + Function f_shape(&shape_function); + Function_wrapper wrapper_2(f_shape); + Mesh_domain domain_2(wrapper_2, bounding_shape, error_bound); + test_is_in_domain(domain_2); + test_do_intersect_surface(domain_2); + test_construct_intersection(domain_2); + } + + void test_construct_initial_points (const Mesh_domain& domain, FT error_bound) const + { + typedef typename Mesh_domain::Construct_initial_points Construct_initial_points; + typedef typename Mesh_domain::Index Index; + typedef typename std::vector >::const_iterator Points_const_iterator; + typedef typename K::Compute_squared_distance_3 Compute_squared_distance_3; + + Compute_squared_distance_3 squared_distance = K().compute_squared_distance_3_object(); + + Construct_initial_points construct_initial_points = domain.construct_initial_points_object(); + std::vector > points; + int point_count = 12; + construct_initial_points(std::back_inserter(points), point_count); + + for (Points_const_iterator iter = points.begin(), end_iter = points.end(); iter != end_iter; ++iter) + { + const Point_3& p = iter->first; + + FT sd = squared_distance(CGAL::ORIGIN, p); + FT diff = sd - 1; + if (diff < FT(0.)) + diff = -diff; + assert(diff <= error_bound); + } + } + + void test_is_in_domain (const Mesh_domain& domain) const + { + typedef typename Mesh_domain::Is_in_domain Is_in_domain; + typedef typename Mesh_domain::Subdomain Subdomain; + typedef typename Mesh_domain::Subdomain_index Subdomain_index; + + Is_in_domain is_in_domain = domain.is_in_domain_object(); + + { + Subdomain subdomain = is_in_domain(Point_3(CGAL::ORIGIN)); + assert(subdomain); + Subdomain_index subdomain_index = *subdomain; + assert(subdomain_index == 1); + } + + { + Subdomain subdomain = is_in_domain(Point_3(1.5, 0., 0.)); + assert(!subdomain); + } + } + + void test_do_intersect_surface (const Mesh_domain& domain) const + { + typedef typename Mesh_domain::Do_intersect_surface Do_intersect_surface; + typedef typename Mesh_domain::Surface_patch Surface_patch; + typedef typename Mesh_domain::Surface_patch_index Surface_patch_index; + typedef typename Mesh_domain::Segment_3 Segment_3; + typedef typename Mesh_domain::Ray_3 Ray_3; + typedef typename Mesh_domain::Line_3 Line_3; + typedef typename Mesh_domain::Vector_3 Vector_3; + + Do_intersect_surface do_intersect_surface = domain.do_intersect_surface_object(); + + { + Segment_3 s(Point_3(CGAL::ORIGIN), Point_3(1.5, 0., 0.)); + Surface_patch p = do_intersect_surface(s); + assert(p); + Surface_patch_index pi = *p; + assert(pi.first == 0 && pi.second == 1); + } + + { + Segment_3 s(Point_3(1.5, 1.5, 0.), Point_3(1.5, 0., 0.)); + Surface_patch p = do_intersect_surface(s); + assert(!p); + } + + { + Ray_3 r(Point_3(CGAL::ORIGIN), Vector_3(1., 0., 0.)); + Surface_patch p = do_intersect_surface(r); + assert(p); + Surface_patch_index pi = *p; + assert(pi.first == 0 && pi.second == 1); + } + + { + Ray_3 r(Point_3(1.5, 0., 0.), Vector_3(0., 1., 0.)); + Surface_patch p = do_intersect_surface(r); + assert(!p); + } + + { + Line_3 l(Point_3(CGAL::ORIGIN), Point_3(1.5, 0., 0.)); + Surface_patch p = do_intersect_surface(l); + assert(p); + Surface_patch_index pi = *p; + assert(pi.first == 0 && pi.second == 1); + } + + { + Line_3 l(Point_3(1.5, 0., 0.), Point_3(1.5, 0.5, 0.)); + Surface_patch p = do_intersect_surface(l); + assert(!p); + } + } + + + void test_construct_intersection (const Mesh_domain& domain) const + { + typedef typename Mesh_domain::Construct_intersection Construct_intersection; + typedef typename Mesh_domain::Intersection Intersection; + typedef typename Mesh_domain::Subdomain_index Subdomain_index; + typedef typename Mesh_domain::Surface_patch_index Surface_patch_index; + typedef typename Mesh_domain::Index Index; + typedef typename Mesh_domain::Segment_3 Segment_3; + typedef typename Mesh_domain::Ray_3 Ray_3; + typedef typename Mesh_domain::Line_3 Line_3; + typedef typename Mesh_domain::Vector_3 Vector_3; + + Construct_intersection construct_intersection = domain.construct_intersection_object(); + + { + Segment_3 s(Point_3(CGAL::ORIGIN), Point_3(1.5, 0., 0.)); + Intersection i = construct_intersection(s); + assert(CGAL::cpp11::get<0>(i) != Point_3(1., 0., 0.)); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 2); + } + + { + Segment_3 s(Point_3(1.5, 1.5, 0.), Point_3(1.5, 0., 0.)); + Intersection i = construct_intersection(s); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 0); + } + + { + Ray_3 r(Point_3(CGAL::ORIGIN), Vector_3(1., 0., 0.)); + Intersection i = construct_intersection(r); + assert(CGAL::cpp11::get<0>(i) != Point_3(1., 0., 0.)); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 2); + } + + { + Ray_3 r(Point_3(1.5, 0., 0.), Vector_3(0., 1., 0.)); + Intersection i = construct_intersection(r); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 0); + } + + { + Line_3 l(Point_3(CGAL::ORIGIN), Point_3(1.5, 0., 0.)); + Intersection i = construct_intersection(l); + assert(CGAL::cpp11::get<0>(i) != Point_3(1., 0., 0.)); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 2); + } + + { + Line_3 l(Point_3(1.5, 0., 0.), Point_3(1.5, 0.5, 0.)); + Intersection i = construct_intersection(l); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 0); + } + } +}; + + +int main () +{ + LM3_tester test_epic; + test_epic(); + + return EXIT_SUCCESS; +} diff --git a/Mesh_3/test/Mesh_3/test_mesh_3_implicit_vector_to_labeled_function_wrapper.cpp b/Mesh_3/test/Mesh_3/test_mesh_3_implicit_vector_to_labeled_function_wrapper.cpp new file mode 100644 index 00000000000..3a45601d56c --- /dev/null +++ b/Mesh_3/test/Mesh_3/test_mesh_3_implicit_vector_to_labeled_function_wrapper.cpp @@ -0,0 +1,66 @@ +#define CGAL_NO_DEPRECATION_WARNINGS 1 + +#include +#include + + +typedef CGAL::Exact_predicates_inexact_constructions_kernel K_e_i; +typedef K_e_i::Point_3 Point_3; +typedef double (Function) (const Point_3&); + +typedef CGAL::Mesh_3::Implicit_vector_to_labeled_function_wrapper Labeling_function; +typedef Labeling_function::Function_vector Function_vector; + + +double cube_function_1 (const Point_3& p) +{ + if( p.x() > 0 && p.x() < 2 && + p.y() > 0 && p.y() < 2 && + p.z() > 0 && p.z() < 2 ) + return -1.; + return 1.; +} + +double cube_function_2 (const Point_3& p) +{ + if( p.x() > 1 && p.x() < 3 && + p.y() > 1 && p.y() < 3 && + p.z() > 1 && p.z() < 3 ) + return -1.; + return 1.; +} + + +void test_constructor_vfunc () +{ + Function_vector vfunc; + vfunc.push_back(cube_function_1); + vfunc.push_back(cube_function_2); + + Labeling_function lfunc(vfunc); + + Point_3 p1(0.5, 0.5, 0.5); + Point_3 p2(2.5, 2.5, 2.5); + Point_3 p3(1.5, 1.5, 1.5); + Point_3 p_out(4., 4., 4.); + + Labeling_function::return_type rp1 = lfunc(p1); + Labeling_function::return_type rp2 = lfunc(p2); + Labeling_function::return_type rp3 = lfunc(p3); + Labeling_function::return_type rp_out = lfunc(p_out); + + assert(rp1 != 0); + assert(rp2 != 0); + assert(rp3 != 0); + assert(rp_out == 0); + + assert(rp1 != rp2); + assert(rp1 != rp3); + assert(rp2 != rp3); +} + +int main () +{ + test_constructor_vfunc(); + return 0; +} diff --git a/Mesh_3/test/Mesh_3/test_mesh_3_labeled_mesh_domain_3.cpp b/Mesh_3/test/Mesh_3/test_mesh_3_labeled_mesh_domain_3.cpp new file mode 100644 index 00000000000..2a104c037f3 --- /dev/null +++ b/Mesh_3/test/Mesh_3/test_mesh_3_labeled_mesh_domain_3.cpp @@ -0,0 +1,270 @@ +// Copyright (c) 2009 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org). +// You can redistribute it and/or modify it under the terms of the GNU +// General Public License as published by the Free Software Foundation, +// either version 3 of the License, or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL$ +// $Id$ +// +// +// Author(s) : Aymeric PELLE +// +//****************************************************************************** +// File Description : +//****************************************************************************** + +#define CGAL_NO_DEPRECATION_WARNINGS 1 + +#include "test_meshing_utilities.h" +#include +#include + +template +struct LM3_tester +{ + typedef typename K::Point_3 Point_3; + typedef typename K::FT FT; + + struct Function + { + typedef Point_3 Point; + typedef FT (*Func)(const Point_3&); + + Function (Func f) : f_(f) {} + FT operator() (const Point_3& p) const { return f_(p); } + + private: + Func f_; + }; + + typedef CGAL::Mesh_3::Implicit_to_labeled_function_wrapper Function_wrapper; + typedef CGAL::Mesh_3::Labeled_mesh_domain_3 Mesh_domain; + + static FT shape_function (const Point_3& p) + { + if (p.x() < 0) + return -1; + const FT x2=p.x()*p.x(), y2=p.y()*p.y(), z2=p.z()*p.z(); + return x2 + y2 + z2 - 1; + } + + static FT sphere_function (const Point_3& p) + { + const FT x2=p.x()*p.x(), y2=p.y()*p.y(), z2=p.z()*p.z(); + return x2 + y2 + z2 - 1; + } + + void operator() () const + { + typedef typename K::Sphere_3 Sphere_3; + + test_domain(Sphere_3(CGAL::ORIGIN, 4.)); + test_domain(CGAL::Bbox_3(-2.,-2.,-2., 2.,2.,2.)); + } + +private: + template + void test_domain (const BoundingShape& bounding_shape) const + { + FT error_bound(1e-3); + + Function f_sphere(&sphere_function); + Function_wrapper wrapper_1(f_sphere); + Mesh_domain domain(wrapper_1, bounding_shape, error_bound); + test_construct_initial_points(domain, error_bound); + + Function f_shape(&shape_function); + Function_wrapper wrapper_2(f_shape); + Mesh_domain domain_2(wrapper_2, bounding_shape, error_bound); + test_is_in_domain(domain_2); + test_do_intersect_surface(domain_2); + test_construct_intersection(domain_2); + } + + void test_construct_initial_points (const Mesh_domain& domain, FT error_bound) const + { + typedef typename Mesh_domain::Construct_initial_points Construct_initial_points; + typedef typename Mesh_domain::Index Index; + typedef typename std::vector >::const_iterator Points_const_iterator; + typedef typename K::Compute_squared_distance_3 Compute_squared_distance_3; + + Compute_squared_distance_3 squared_distance = K().compute_squared_distance_3_object(); + + Construct_initial_points construct_initial_points = domain.construct_initial_points_object(); + std::vector > points; + int point_count = 12; + construct_initial_points(std::back_inserter(points), point_count); + + for (Points_const_iterator iter = points.begin(), end_iter = points.end(); iter != end_iter; ++iter) + { + const Point_3& p = iter->first; + + FT sd = squared_distance(CGAL::ORIGIN, p); + FT diff = sd - 1; + if (diff < FT(0.)) + diff = -diff; + assert(diff <= error_bound); + } + } + + void test_is_in_domain (const Mesh_domain& domain) const + { + typedef typename Mesh_domain::Is_in_domain Is_in_domain; + typedef typename Mesh_domain::Subdomain Subdomain; + typedef typename Mesh_domain::Subdomain_index Subdomain_index; + + Is_in_domain is_in_domain = domain.is_in_domain_object(); + + { + Subdomain subdomain = is_in_domain(Point_3(CGAL::ORIGIN)); + assert(subdomain); + Subdomain_index subdomain_index = *subdomain; + assert(subdomain_index == 1); + } + + { + Subdomain subdomain = is_in_domain(Point_3(1.5, 0., 0.)); + assert(!subdomain); + } + } + + void test_do_intersect_surface (const Mesh_domain& domain) const + { + typedef typename Mesh_domain::Do_intersect_surface Do_intersect_surface; + typedef typename Mesh_domain::Surface_patch Surface_patch; + typedef typename Mesh_domain::Surface_patch_index Surface_patch_index; + typedef typename Mesh_domain::Segment_3 Segment_3; + typedef typename Mesh_domain::Ray_3 Ray_3; + typedef typename Mesh_domain::Line_3 Line_3; + typedef typename Mesh_domain::Vector_3 Vector_3; + + Do_intersect_surface do_intersect_surface = domain.do_intersect_surface_object(); + + { + Segment_3 s(Point_3(CGAL::ORIGIN), Point_3(1.5, 0., 0.)); + Surface_patch p = do_intersect_surface(s); + assert(p); + Surface_patch_index pi = *p; + assert(pi.first == 0 && pi.second == 1); + } + + { + Segment_3 s(Point_3(1.5, 1.5, 0.), Point_3(1.5, 0., 0.)); + Surface_patch p = do_intersect_surface(s); + assert(!p); + } + + { + Ray_3 r(Point_3(CGAL::ORIGIN), Vector_3(1., 0., 0.)); + Surface_patch p = do_intersect_surface(r); + assert(p); + Surface_patch_index pi = *p; + assert(pi.first == 0 && pi.second == 1); + } + + { + Ray_3 r(Point_3(1.5, 0., 0.), Vector_3(0., 1., 0.)); + Surface_patch p = do_intersect_surface(r); + assert(!p); + } + + { + Line_3 l(Point_3(CGAL::ORIGIN), Point_3(1.5, 0., 0.)); + Surface_patch p = do_intersect_surface(l); + assert(p); + Surface_patch_index pi = *p; + assert(pi.first == 0 && pi.second == 1); + } + + { + Line_3 l(Point_3(1.5, 0., 0.), Point_3(1.5, 0.5, 0.)); + Surface_patch p = do_intersect_surface(l); + assert(!p); + } + } + + + void test_construct_intersection (const Mesh_domain& domain) const + { + typedef typename Mesh_domain::Construct_intersection Construct_intersection; + typedef typename Mesh_domain::Intersection Intersection; + typedef typename Mesh_domain::Subdomain_index Subdomain_index; + typedef typename Mesh_domain::Surface_patch_index Surface_patch_index; + typedef typename Mesh_domain::Index Index; + typedef typename Mesh_domain::Segment_3 Segment_3; + typedef typename Mesh_domain::Ray_3 Ray_3; + typedef typename Mesh_domain::Line_3 Line_3; + typedef typename Mesh_domain::Vector_3 Vector_3; + + Construct_intersection construct_intersection = domain.construct_intersection_object(); + + { + Segment_3 s(Point_3(CGAL::ORIGIN), Point_3(1.5, 0., 0.)); + Intersection i = construct_intersection(s); + assert(CGAL::cpp11::get<0>(i) != Point_3(1., 0., 0.)); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 2); + } + + { + Segment_3 s(Point_3(1.5, 1.5, 0.), Point_3(1.5, 0., 0.)); + Intersection i = construct_intersection(s); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 0); + } + + { + Ray_3 r(Point_3(CGAL::ORIGIN), Vector_3(1., 0., 0.)); + Intersection i = construct_intersection(r); + assert(CGAL::cpp11::get<0>(i) != Point_3(1., 0., 0.)); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 2); + } + + { + Ray_3 r(Point_3(1.5, 0., 0.), Vector_3(0., 1., 0.)); + Intersection i = construct_intersection(r); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 0); + } + + { + Line_3 l(Point_3(CGAL::ORIGIN), Point_3(1.5, 0., 0.)); + Intersection i = construct_intersection(l); + assert(CGAL::cpp11::get<0>(i) != Point_3(1., 0., 0.)); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 2); + } + + { + Line_3 l(Point_3(1.5, 0., 0.), Point_3(1.5, 0.5, 0.)); + Intersection i = construct_intersection(l); + Index ii = CGAL::cpp11::get<1>(i); + assert(boost::get(&ii)); + assert(CGAL::cpp11::get<2>(i) == 0); + } + } +}; + + +int main () +{ + LM3_tester test_epic; + test_epic(); + + return EXIT_SUCCESS; +} diff --git a/Mesh_3/test/Mesh_3/test_mesh_implicit_domains.cpp b/Mesh_3/test/Mesh_3/test_mesh_implicit_domains.cpp new file mode 100644 index 00000000000..ee002ca0c61 --- /dev/null +++ b/Mesh_3/test/Mesh_3/test_mesh_implicit_domains.cpp @@ -0,0 +1,81 @@ + +//****************************************************************************** +// File Description : +// Outputs to out.mesh a mesh of implicit domains. These domains are defined +// by a vector of functions. Each n-uplet of sign of function values defines a +// subdomain. +//****************************************************************************** + + +#define CGAL_NO_DEPRECATION_WARNINGS 1 + +#include "debug.h" +#include + +#include +#include +#include + +#include +#include +#include +#include "implicit_functions.h" + +// IO +#include + +using namespace CGAL::parameters; + +// Domain +typedef CGAL::Exact_predicates_inexact_constructions_kernel K; +typedef FT_to_point_function_wrapper Function; +typedef CGAL::Mesh_3::Implicit_vector_to_labeled_function_wrapper + Function_wrapper; +typedef Function_wrapper::Function_vector Function_vector; +typedef CGAL::Mesh_3::Labeled_mesh_domain_3 Mesh_domain; + +// Triangulation +typedef CGAL::Mesh_triangulation_3::type Tr; +typedef CGAL::Mesh_complex_3_in_triangulation_3 C3t3; + +// Mesh Criteria +typedef CGAL::Mesh_criteria_3 Mesh_criteria; +typedef Mesh_criteria::Facet_criteria Facet_criteria; +typedef Mesh_criteria::Cell_criteria Cell_criteria; + + +int main() +{ + // Define functions + Function f1(&torus_function); + Function f2(&sphere_function<5>); + Function f3(&sphere_function<2>); + + Function_vector v; + v.push_back(&f1); + v.push_back(&f2); + //v.push_back(&f3); + + // Domain (Warning: Sphere_3 constructor uses square radius !) + Mesh_domain domain(v, K::Sphere_3(CGAL::ORIGIN, 5.*5.), 1e-6); + + // Set mesh criteria + Facet_criteria facet_criteria(30, 0.2, 0.02); // angle, size, approximation + Cell_criteria cell_criteria(2., 0.4); // radius-edge ratio, size + Mesh_criteria criteria(facet_criteria, cell_criteria); + + // Mesh generation + C3t3 c3t3 = CGAL::make_mesh_3(domain, criteria, no_exude(), no_perturb()); + + // Perturbation (maximum cpu time: 10s, targeted dihedral angle: default) + CGAL::perturb_mesh_3(c3t3, domain, time_limit = 10); + + // Exudation + CGAL::exude_mesh_3(c3t3,12); + + // Output + std::ofstream medit_file("out.mesh"); + CGAL::output_to_medit(medit_file, c3t3); + + return 0; +} diff --git a/Surface_mesher/doc/Surface_mesher/Concepts/ImplicitFunction.h b/Surface_mesher/doc/Surface_mesher/Concepts/ImplicitFunction.h index 1df71a27491..1af848f7d31 100644 --- a/Surface_mesher/doc/Surface_mesher/Concepts/ImplicitFunction.h +++ b/Surface_mesher/doc/Surface_mesher/Concepts/ImplicitFunction.h @@ -7,7 +7,7 @@ The concept `ImplicitFunction` describes a function object whose `operator()` computes the values of a function \f$ f : \mathbb{R}^3 \longrightarrow \mathbb{R}\f$. -\cgalHasModel `CGAL::Gray_level_image_3` +\cgalHasModel CGAL::Gray_level_image_3 \cgalHasModel any pointer to a function of type `FT (*)(Point)` \sa `CGAL::Implicit_surface_3` @@ -20,7 +20,7 @@ public: /// \name Types /// @{ - +///The following types aren't required for any pointer to a function of type `FT (*)(Point)`. /*! Number type */