Reducing memory consumption

Skin_surface_3/examples/Skin_surface_3/NGHK_skin_surface_subdiv.C

@@ -6,8 +6,7 @@
 #include <CGAL/subdivide_skin_surface_mesh_3.h>
 #include <list>
 
-#include <fstream>
+#include "skin_surface_writer.h"
-#include <CGAL/IO/Polyhedron_iostream.h>
 
 typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
 typedef CGAL::Regular_triangulation_euclidean_traits_3<K>   Traits;
@@ -36,6 +35,7 @@ int main(int argc, char *argv[]) {
   CGAL::subdivide_skin_surface_mesh_3(p, skin_surface);
 
   std::ofstream out("mesh.off");
+  write_polyhedron_with_normals(p, skin_surface, out);
   out << p;
 
   return 0;

Skin_surface_3/examples/Skin_surface_3/makefile

@@ -13,7 +13,7 @@ include $(CGAL_MAKEFILE)
 #                    compiler flags
 #---------------------------------------------------------------------#
 
-CXXFLAGS = \
+CXXFLAGS = -g\
            -I.\
            -Idsrpdb/include\
            -I../../include \
@@ -30,7 +30,7 @@ LIBPATH = \
           $(TESTSUITE_LIBPATH) \
           $(CGAL_LIBPATH)
 
-LDFLAGS = \
+LDFLAGS = -g \
           $(TESTSUITE_LDFLAGS) \
           $(LONG_NAME_PROBLEM_LDFLAGS) \
           $(CGAL_LDFLAGS)
@@ -40,17 +40,12 @@ LDFLAGS = \
 #---------------------------------------------------------------------#
 
 all:            \
-                implicit_mixed_complex$(EXE_EXT) \
                 skin_surface_pdb_reader$(EXE_EXT) \
                 skin_surface_simple$(EXE_EXT) \
                 skin_surface_subdiv$(EXE_EXT) \
                 skin_surface_subdiv_with_normals$(EXE_EXT) \
-                Surface_mesher_skin_surface$(EXE_EXT) \
                 union_of_balls_simple$(EXE_EXT) 
 
-implicit_mixed_complex$(EXE_EXT): implicit_mixed_complex$(OBJ_EXT)
-	$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)implicit_mixed_complex implicit_mixed_complex$(OBJ_EXT) $(LDFLAGS)
-
 NGHK_skin_surface_simple$(EXE_EXT): NGHK_skin_surface_simple$(OBJ_EXT)
 	$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)NGHK_skin_surface_simple NGHK_skin_surface_simple$(OBJ_EXT) $(LDFLAGS)
 
@@ -75,19 +70,14 @@ skin_surface_subdiv$(EXE_EXT): skin_surface_subdiv$(OBJ_EXT)
 skin_surface_subdiv_with_normals$(EXE_EXT): skin_surface_subdiv_with_normals$(OBJ_EXT)
 	$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)skin_surface_subdiv_with_normals skin_surface_subdiv_with_normals$(OBJ_EXT) $(LDFLAGS)
 
-Surface_mesher_skin_surface$(EXE_EXT): Surface_mesher_skin_surface$(OBJ_EXT)
-	$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)Surface_mesher_skin_surface Surface_mesher_skin_surface$(OBJ_EXT) $(LDFLAGS)
-
 union_of_balls_simple$(EXE_EXT): union_of_balls_simple$(OBJ_EXT)
 	$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)union_of_balls_simple union_of_balls_simple$(OBJ_EXT) $(LDFLAGS)
 
 clean: \
-                   implicit_mixed_complex.clean \
                    skin_surface_pdb_reader.clean \
                    skin_surface_simple.clean \
                    skin_surface_subdiv.clean \
                    skin_surface_subdiv_with_normals.clean \
-                   Surface_mesher_skin_surface.clean \
                    union_of_balls_simple.clean 
 
 #---------------------------------------------------------------------#

Skin_surface_3/examples/Skin_surface_3/skin_surface_pdb_reader.C

@@ -10,8 +10,7 @@
 #include <extract_balls_from_pdb.h>
 
 #include <list>
-#include <fstream>
+#include "skin_surface_writer.h"
-#include <CGAL/IO/Polyhedron_iostream.h>
 
 typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
 typedef CGAL::Regular_triangulation_euclidean_traits_3<K>   Traits;
@@ -25,20 +24,26 @@ int main(int argc, char *argv[]) {
     std::cout << "Usage: " << argv[0] << " <pdb-file>" << std::endl;
     return 1;
   }
+
   std::list<Weighted_point> l;
   FT                        shrinkfactor = 0.5;
 
   extract_balls_from_pdb(argv[1], K(), std::back_inserter(l));
+  std::cout << "Read pdb" << std::endl;
 
-  Skin_surface_3 skin_surface(l.begin(), l.end(), shrinkfactor);
+  Skin_surface_3 skin_surface(l.begin(), l.end(), shrinkfactor, true, Traits(), true);
+  std::cout << "Constructed Skin_surface_3" << std::endl;
 
   Polyhedron p;
   CGAL::mesh_skin_surface_3(skin_surface, p);
+  std::cout << "Meshed Skin_surface_3" << std::endl;
 
   //  CGAL::subdivide_skin_surface_mesh_3(p, skin_surface);
+  //  std::cout << "Subdivided Skin_surface_3" << std::endl;
 
   std::ofstream out("mesh.off");
-  out << p;
+  write_polyhedron_with_normals(p, skin_surface, out);
+  //out << p; 
 
   return 0;
 }

Skin_surface_3/include/CGAL/Compute_anchor_3.h

@@ -265,7 +265,6 @@ compute_anchor_del(Facet const &f) {
   CGAL_assertion(!reg.is_infinite(f));
   equiv_anchors.clear();
 
-  Simplex s;
   int i;
   Sign result;
   bool contains_center = true;
@@ -289,7 +288,7 @@ compute_anchor_del(Facet const &f) {
     if (i==3) e.third = ((f.second+2)&3); 
     else e.third = ((f.second+3)&3);
 
-    s = anchor_del(e);
+    Simplex s = anchor_del(e);
     if (s.dimension() == 1) {
       equiv_anchors.clear();
       return s;
@@ -315,7 +314,6 @@ compute_anchor_del(Cell_handle const ch) {
   equiv_anchors.clear();
   
   Simplex s;
-
   bool contains_center = true;
   Sign result;
   for (int i=0; (i<4) && contains_center; i++) {
@@ -428,9 +426,9 @@ compute_anchor_vor (Vertex_handle const v) {
       } else if (s.dimension() == 3) {
         //  s lies on a Voronoi vertex
         Cell_handle ch=s;
+	CGAL_assertion(ch != Cell_handle());
         int index = ch->index(v);
         for (int i=1; (i<4) && (found); i++) {
-          Simplex tmp;
           tmp = anchor_vor(Facet(ch, (index+i)&3));
 	  found = (tmp == s);
         }
@@ -463,6 +461,7 @@ compute_anchor_vor (Vertex_handle const v) {
       s = tmp;
     }
   }
+  CGAL_assertion(false);
   return Simplex();
 }
 
@@ -561,6 +560,7 @@ Compute_anchor_3<RegularTriangulation3>::compute_anchor_vor (Facet const &f) {
   
   Sign side;
   
+  CGAL_assertion(f.first != Cell_handle());
   if (!reg.is_infinite(f.first)) {
     side = test_anchor(f.first, f.second);
     if (side==NEGATIVE) {
@@ -571,6 +571,7 @@ Compute_anchor_3<RegularTriangulation3>::compute_anchor_vor (Facet const &f) {
   }
 
   Cell_handle neighbor = f.first->neighbor(f.second);
+  CGAL_assertion(neighbor != Cell_handle());
   if (!reg.is_infinite(neighbor)) {
     int n_index = neighbor->index(f.first);
     side = test_anchor(neighbor, n_index);

Skin_surface_3/include/CGAL/Marching_tetrahedra_traits_skin_surface_3.h

@@ -61,10 +61,10 @@ private:
   HDS_RT value(const Cell_handle &ch, const HDS_point &p) const {
     return ch->surf->value(p);
   }
-  HDS_RT value(const Cell_handle &ch, const Triang_point &p) const {
+//   HDS_RT value(const Cell_handle &ch, const Triang_point &p) const {
-    // NGHK: Remove the to_double later ...
+//     // NGHK: Remove the to_double later ...
-    return CGAL::to_double(ch->surf->value(converter(p)));
+//     return CGAL::to_double(ch->surf->value(converter(p)));
-  }
+//   }
 
   Converter converter;
   HDS_RT iso_value;

Skin_surface_3/include/CGAL/Skin_surface_3.h

@@ -77,7 +77,7 @@ public:
 //                                           Triangulated_mixed_complex_tds;
 
   // defining the triangulated mixed complex:
-  typedef Exact_predicates_exact_constructions_kernel    TMC_Traits;
+  typedef Exact_predicates_inexact_constructions_kernel    TMC_Traits;
 public:
 
 #ifdef CGAL_SKIN_SURFACE_USE_EXACT_IMPLICIT_SURFACE
@@ -100,13 +100,23 @@ public:
   template < class WP_iterator >
   Skin_surface_3(WP_iterator begin, WP_iterator end, 
 		 FT shrink_factor,
+		 bool grow_balls = true,
 		 Gt gt = Gt(),
 		 bool verbose = false
 		 ) 
-    : regular(), gt(gt), shrink(shrink_factor), verbose(verbose) {
+    : gt(gt), shrink(shrink_factor), verbose(verbose) {
 
+    CGAL_assertion(begin != end);
+
+    Regular regular;
+    if (grow_balls) {
+      for (; begin != end; begin++) {
+	regular.insert(Weighted_point(*begin, begin->weight()/shrink));
+      }
+    } else {
       regular.insert(begin, end);
-    construct_bounding_box();
+    }
+    construct_bounding_box(regular);
   
     if (verbose) {
       std::cerr << "Triangulation ready" << std::endl;
@@ -123,6 +133,12 @@ public:
       std::cerr << "Vertices: " << _tmc.number_of_vertices() << std::endl;
       std::cerr << "Cells:    " << _tmc.number_of_cells() << std::endl;
     }
+//     std::ofstream out("vertices.txt");
+//     for (typename Triangulated_mixed_complex::Finite_vertices_iterator 
+// 	   vit = _tmc.finite_vertices_begin();
+// 	 vit != _tmc.finite_vertices_end(); vit ++) {
+//       out << vit->point().x().exact() << std::endl;
+//     }
   }
   const Triangulated_mixed_complex &triangulated_mixed_complex() const {
     return _tmc;
@@ -132,9 +148,8 @@ public:
     return _tmc.locate(p);
   }
 private:
-  void construct_bounding_box();
+  void construct_bounding_box(Regular &regular);
 
-  Regular regular;
   Gt &gt;
   FT shrink;
   Triangulated_mixed_complex _tmc;
@@ -144,7 +159,7 @@ private:
 template <class SkinSurfaceTraits_3> 
 void 
 Skin_surface_3<SkinSurfaceTraits_3>::
-construct_bounding_box() 
+construct_bounding_box(Regular &regular) 
 {
   typedef typename Regular::Finite_vertices_iterator Finite_vertices_iterator;
   typedef typename Regular::Geom_traits     GT;

Skin_surface_3/include/CGAL/Skin_surface_quadratic_surface_3.h

@@ -65,10 +65,15 @@ public:
     return inv_conv(compute_gradient(xp));
   }
 
-  Vector compute_gradient(Point const &x) {
+  template <class Input_point>
-    std::cout << "NGHK: NOT YET IMPLEMENTED" << std::endl;
+  Vector compute_gradient(Input_point const &x) {
-    // NGHK: TODO:
+    typedef Cartesian_converter<typename Input_point::R, K> Converter;
-    return (x-p);
+    
+    Vector v = Converter()(x) - p;
+
+    return Vector(2*Q[0]*v.x() + Q[1]*v.y() + Q[3]*v.z(),
+		  Q[1]*v.x() + 2*Q[2]*v.y() + Q[4]*v.z(),
+		  Q[3]*v.x() + Q[4]*v.y() + 2*Q[5]*v.z());
   }
 	
   /// Construct the intersection point with the segment (p0,p1)

Skin_surface_3/include/CGAL/marching_tetrahedra_3.h

@@ -120,6 +120,7 @@ public:
 
 
   bool is_inside(T_cell_iterator const ch, int i) {
+    return (traits.sign(ch,i) == POSITIVE);
     T_vertex_map_it it = triang_vertex_signs.find(ch->vertex(i));
     
     if (it == triang_vertex_signs.end()) {

Skin_surface_3/include/CGAL/triangulate_mixed_complex_3.h

@@ -29,8 +29,6 @@
 // NGHK: move this one to SkinSurfaceTraits
 #include <CGAL/Compute_anchor_3.h>
 
-#include <CGAL/Union_find.h>
-
 CGAL_BEGIN_NAMESPACE
 
 
@@ -96,6 +94,9 @@ private:
 
   typedef Compute_anchor_3<Regular>                       Compute_anchor;
   typedef std::pair<Rt_Simplex,Rt_Simplex>                Symb_anchor;
+
+  typedef std::map<Rt_Simplex, Rt_Simplex>       Anchor_map;
+  typedef typename Anchor_map::iterator          Anchor_map_iterator;
 public:
 
   Mixed_complex_triangulator_3(
@@ -151,16 +152,6 @@ private:
 
     triangulation_incr_builder.end_triangulation();
     
-    std::cout << map_del.size() << " vs. " << map_vor.size() << std::endl;
-    std::cout << anchor_del.size() << " vs. " << anchor_vor.size() << std::endl;
-    std::cout << "Union_find: " << anchor_del.number_of_sets () << "vs. " << anchor_del.size () << std::endl;
-    std::cout << "Union_find: " << anchor_vor.number_of_sets () << "vs. " << anchor_vor.size () << std::endl;
-    std::cout << anchors.size() << std::endl;
-
-    anchor_del.clear();
-    anchor_vor.clear();
-    map_del.clear();
-    map_vor.clear();
     anchors.clear();
   }
 
@@ -174,10 +165,19 @@ private:
   void construct_anchor_vor(Rt_Simplex const &sVor);
   void construct_anchors();
   Rt_Simplex get_anchor_del(Rt_Simplex const &sDel) {
-    return *anchor_del.find(map_del[sDel]);
+    return find_anchor(anchor_del2, sDel);
   }
   Rt_Simplex get_anchor_vor(Rt_Simplex const &sVor) {
-    return *anchor_vor.find(map_vor[sVor]);
+    return find_anchor(anchor_vor2, sVor);
+  }  
+  Rt_Simplex find_anchor(Anchor_map &a_map, Rt_Simplex const&s) {
+    Anchor_map_iterator it = a_map.find(s);
+    CGAL_assertion(it != a_map.end());
+    if (it->second == s) return s;
+
+    Rt_Simplex result = find_anchor(a_map,it->second);
+    it->second = result;
+    return result;
   }
   void construct_vertices();
   
@@ -225,18 +225,12 @@ private:
   struct Index_v {
     Unique_hash_map < Rt_Vertex_handle, Tmc_Vertex_handle > V;
   };
-  // Facets on the border of the simplicial complex:
-  // name is given by (dim del,dim vor)
 
   // index to vertex
   Unique_hash_map < Rt_Cell_handle, Index_c4 > index_03;
   
-  typedef Union_find<Rt_Simplex>                 Union_find_anchor;
-  typedef typename Union_find_anchor::handle     Union_find_anchor_handle;
-  typedef typename Union_find_anchor::iterator   Union_find_anchor_iterator;
 
-  Union_find_anchor                              anchor_del, anchor_vor;
+  Anchor_map                                     anchor_del2, anchor_vor2;
-  std::map<Rt_Simplex, Union_find_anchor_handle> map_del, map_vor;
   std::map<Symb_anchor, Tmc_Vertex_handle>        anchors;
 };
 
@@ -261,20 +255,28 @@ Mixed_complex_triangulator_3<
   TriangulatedMixedComplex_3,
   TriangulatedMixedComplexObserver_3>::
 construct_anchor_del(Rt_Simplex const &sDel) {
-  Union_find_anchor_handle handle = anchor_del.make_set(sDel);
-  map_del[sDel] = handle;
-  
   Rt_Simplex s = compute_anchor_obj.anchor_del(sDel);
-  if (sDel != s) {
+  anchor_del2.insert(std::pair<Rt_Simplex,Rt_Simplex>(sDel, s));
-    anchor_del.unify_sets(handle, map_del[s]);
-  }
 
   // degenerate simplices:
   if (compute_anchor_obj.is_degenerate()) {
-    typename Compute_anchor::Simplex_iterator it;
+    s = find_anchor(anchor_del2, sDel);
-    for (it = compute_anchor_obj.equivalent_anchors_begin();
+    typename Compute_anchor::Simplex_iterator degenerate_it;
-	 it != compute_anchor_obj.equivalent_anchors_end(); it++) {
+    for (degenerate_it = compute_anchor_obj.equivalent_anchors_begin();
-      anchor_del.unify_sets(handle, map_del[*it]);
+	 degenerate_it != compute_anchor_obj.equivalent_anchors_end(); 
+	 degenerate_it++) {
+      Anchor_map_iterator it = anchor_del2.find(*degenerate_it);
+      CGAL_assertion(it != anchor_del2.end());
+      Anchor_map_iterator it2 = anchor_del2.find(it->second);
+      CGAL_assertion(it2 != anchor_del2.end());
+      // Merge sets:
+      while (it != it2) {
+	it->second = s;
+	it = it2;
+	it2 = anchor_del2.find(it->second);
+	CGAL_assertion(it2 != anchor_del2.end());
+      }
+      it->second = s;
     }
   }
 }
@@ -289,27 +291,31 @@ Mixed_complex_triangulator_3<
   TriangulatedMixedComplex_3,
   TriangulatedMixedComplexObserver_3>::
 construct_anchor_vor(Rt_Simplex const &sVor) {
-  Union_find_anchor_handle handle = anchor_vor.make_set(sVor);
-  map_vor[sVor] = handle;
-
   Rt_Simplex s = compute_anchor_obj.anchor_vor(sVor);
-  if (sVor != s) {
+  anchor_vor2.insert(std::pair<Rt_Simplex,Rt_Simplex>(sVor, s));
-    anchor_vor.unify_sets(handle, map_vor[s]);
-  }
 
   // degenerate simplices:
   if (compute_anchor_obj.is_degenerate()) {
-    typename Compute_anchor::Simplex_iterator it;
+    s = find_anchor(anchor_vor2, sVor);
-    for (it = compute_anchor_obj.equivalent_anchors_begin();
+    typename Compute_anchor::Simplex_iterator degenerate_it;
-	 it != compute_anchor_obj.equivalent_anchors_end(); it++) {
+    for (degenerate_it = compute_anchor_obj.equivalent_anchors_begin();
-      typename std::map<Rt_Simplex, Union_find_anchor_handle>::iterator h_it;
+	 degenerate_it != compute_anchor_obj.equivalent_anchors_end(); 
-      h_it = map_vor.find(*it);
+	 degenerate_it++) {
       // Possibly not found for 2 Voronoi vertices with the same center,
       // If the first vertex is inserted and the second is already found.
-      if (h_it != map_vor.end()) {
+      // see compute_anchor_obj.anchor_vor(Cell_handle)
-	anchor_vor.unify_sets(handle, (*h_it).second);
+      Anchor_map_iterator it = anchor_vor2.find(*degenerate_it);
-      } else {
+      if (it != anchor_vor2.end()) {
-	CGAL_assertion(s.dimension() == 3);
+	Anchor_map_iterator it2 = anchor_vor2.find(it->second);
+	CGAL_assertion(it2 != anchor_vor2.end());
+	// Merge sets:
+	while (it != it2) {
+	  it->second = s;
+	  it = it2;
+	  it2 = anchor_vor2.find(it->second);
+	  CGAL_assertion(it2 != anchor_vor2.end());
+	}
+	it->second = s;
       }
     }
   }
@@ -334,9 +340,7 @@ construct_anchors() {
   // Compute anchor points:
   for (vit=regular.finite_vertices_begin();
        vit!=regular.finite_vertices_end(); vit++) {
-    s = Rt_Simplex(vit);
+    construct_anchor_del(Rt_Simplex(vit));
-    construct_anchor_del(s);
-    CGAL_assertion(s.dimension() == 0);
   }  
   for (eit=regular.finite_edges_begin();
        eit!=regular.finite_edges_end(); eit++) {
@@ -371,6 +375,7 @@ construct_anchors() {
   }
   for (vit=regular.finite_vertices_begin();
        vit!=regular.finite_vertices_end(); vit++) {
+    CGAL_assertion(vit->cell() != Rt_Cell_handle());
     s = Rt_Simplex(vit);
     construct_anchor_vor(s);
     CGAL_assertion(s.dimension() == 0);
@@ -404,7 +409,7 @@ construct_vertices() {
   if (verbose) std::cout << "construct_anchors" << std::endl;
   construct_anchors();
 
-  if (verbose) std::cout << "1" << std::endl;
+  if (verbose) std::cout << "9 ";
   // anchor dimDel=0, dimVor=3
   for (cit=regular.finite_cells_begin();
        cit!=regular.finite_cells_end(); cit++) {
@@ -419,7 +424,7 @@ construct_vertices() {
     }
   }
 
-  if (verbose) std::cout << "2" << std::endl;
+  if (verbose) std::cout << "8 ";
   // anchor dimDel=1, dimVor=3
   for (cit=regular.finite_cells_begin(); cit!=regular.finite_cells_end(); cit++) {
     sVor = get_anchor_vor(Rt_Simplex(cit));
@@ -435,7 +440,7 @@ construct_vertices() {
     }
   }
 
-  if (verbose) std::cout << "3" << std::endl;
+  if (verbose) std::cout << "7 ";
   // anchor dimDel=2, dimVor=3 and dimDel=0, dimVor=2
   for (fit=regular.finite_facets_begin(); fit!=regular.finite_facets_end(); fit++) {
     // anchor dimDel=2, dimVor=3
@@ -473,7 +478,7 @@ construct_vertices() {
     }
   }
 	
-  if (verbose) std::cout << "4" << std::endl;
+  if (verbose) std::cout << "6 ";
   // anchor dimDel=0, dimVor=1
   for (eit=regular.finite_edges_begin(); eit!=regular.finite_edges_end(); eit++) {
     sVor = get_anchor_vor(*eit);
@@ -495,7 +500,7 @@ construct_vertices() {
     }
   }
 	
-  if (verbose) std::cout << "5" << std::endl;
+  if (verbose) std::cout << "5 ";
   // anchor dimDel=3, dimVor=3
   for (cit=regular.finite_cells_begin(); cit!=regular.finite_cells_end(); cit++) {
     sDel = get_anchor_del(Rt_Simplex(cit));
@@ -508,7 +513,7 @@ construct_vertices() {
   }
 
 
-  if (verbose) std::cout << "6" << std::endl;
+  if (verbose) std::cout << "4 ";
   // anchor dimDel=0, dimVor=0
   for (vit=regular.finite_vertices_begin(); vit!=regular.finite_vertices_end(); vit++) {
     sDel = get_anchor_del(Rt_Simplex(vit));
@@ -520,7 +525,7 @@ construct_vertices() {
     }
   }
 	
-  if (verbose) std::cout << "7" << std::endl;
+  if (verbose) std::cout << "3 ";
   // anchor dimDel=1, dimVor=2
   for (fit=regular.finite_facets_begin(); fit!=regular.finite_facets_end(); fit++) {
     c1 = fit->first;
@@ -542,7 +547,7 @@ construct_vertices() {
     }
   }
 	
-  if (verbose) std::cout << "8" << std::endl;
+  if (verbose) std::cout << "2 ";
   // anchor dimDel=2, dimVor=2
   for (fit=regular.finite_facets_begin(); fit!=regular.finite_facets_end(); fit++) {
     c1 = fit->first;
@@ -557,7 +562,7 @@ construct_vertices() {
     }
   }
 	
-  if (verbose) std::cout << "9" << std::endl;
+  if (verbose) std::cout << "1" << std::endl;
   // anchor dimDel=1, dimVor=1
   for (eit=regular.finite_edges_begin(); eit!=regular.finite_edges_end(); eit++) {
     v1 = eit->first->vertex(eit->second);

Skin_surface_3/test/Skin_surface_3/surface_test.C

@@ -55,7 +55,7 @@ int main(int argc, char *argv[]) {
 	   cell != cells.end(); cell++) {
 	if (!tmc.is_infinite(*cell)) {
 	  Quadratic_surface::RT val2 = (*cell)->surf->value(vit->point());
- 	  CGAL_assertion(val == val2);
+ 	  CGAL_assertion(std::abs(val - val2) < 10e-6);
 	}
       }
     }