diff --git a/Heat_method_3/include/CGAL/Heat_method_3/Heat_method_3.h b/Heat_method_3/include/CGAL/Heat_method_3/Heat_method_3.h
index 7f319d40090..4b937237b61 100644
--- a/Heat_method_3/include/CGAL/Heat_method_3/Heat_method_3.h
+++ b/Heat_method_3/include/CGAL/Heat_method_3/Heat_method_3.h
@@ -33,7 +33,7 @@
 #include <CGAL/squared_distance_3.h>
 #include <CGAL/Polygon_mesh_processing/measure.h>
 #include <CGAL/number_utils.h>
-
+#include <CGAL/Heat_method_3/Intrinsic_Delaunay_Triangulation_3.h>
 #include <Eigen/Cholesky>
 #include <Eigen/Sparse>
 
@@ -69,7 +69,6 @@ namespace Heat_method_3 {
             typename Traits,
             typename VertexDistanceMap,
             typename VertexPointMap = typename boost::property_map< TriangleMesh, vertex_point_t>::const_type,
-            typename FaceIndexMap = typename boost::property_map< TriangleMesh, face_index_t>::const_type,
             typename LA = Heat_method_Eigen_traits_3>
   class Heat_method_3
   {
@@ -111,13 +110,13 @@ namespace Heat_method_3 {
   public:
 
     Heat_method_3(const TriangleMesh& tm, VertexDistanceMap vdm)
-      : vertex_id_map(get(Vertex_property_tag(),const_cast<TriangleMesh&>(tm))), face_id_map(get(Face_property_tag(),tm)), tm(tm), vdm(vdm), vpm(get(vertex_point,tm)), fpm(get(boost::face_index,tm))
+      : vertex_id_map(get(Vertex_property_tag(),const_cast<TriangleMesh&>(tm))), face_id_map(get(Face_property_tag(),tm)), v2v(tm),tm(tm), vdm(vdm), vpm(get(vertex_point,tm))
     {
       build();
     }
 
-    Heat_method_3(const TriangleMesh& tm, VertexDistanceMap vdm, VertexPointMap vpm, FaceIndexMap fpm)
-      : tm(tm), vdm(vdm), vpm(vpm), fpm(fpm)
+    Heat_method_3(const TriangleMesh& tm, VertexDistanceMap vdm, VertexPointMap vpm)
+      : v2v(tm), tm(tm), vdm(vdm), vpm(vpm)
     {
       build();
     }
@@ -151,10 +150,13 @@ namespace Heat_method_3 {
        return vertex_id_map;
      }
 
-    bool add_source(vertex_descriptor vd)
+    template <typename VD>
+    bool add_source(VD vd)
     {
       source_change_flag = true;
-      return sources.insert(vd).second;
+      v2v(vd);
+        return true;
+      //return sources.insert(v2v(vd)).second;
     }
 
     /**
@@ -163,7 +165,7 @@ namespace Heat_method_3 {
     bool remove_source(vertex_descriptor vd)
     {
       source_change_flag = true;
-      return (sources.erase(vd) == 1);
+      return (sources.erase(v2v(vd)) == 1);
     }
 
     /**
@@ -541,10 +543,10 @@ namespace Heat_method_3 {
     }
     
     int dimension;
+    V2V<TriangleMesh> v2v;
     const TriangleMesh& tm;
     VertexDistanceMap vdm;
     VertexPointMap vpm;
-    FaceIndexMap fpm;
     std::set<vertex_descriptor> sources;
     double m_time_step;
     Matrix kronecker;
diff --git a/Heat_method_3/include/CGAL/Heat_method_3/Intrinsic_Delaunay_Triangulation_3.h b/Heat_method_3/include/CGAL/Heat_method_3/Intrinsic_Delaunay_Triangulation_3.h
index 84f7c85da48..2a3be8ad7e1 100644
--- a/Heat_method_3/include/CGAL/Heat_method_3/Intrinsic_Delaunay_Triangulation_3.h
+++ b/Heat_method_3/include/CGAL/Heat_method_3/Intrinsic_Delaunay_Triangulation_3.h
@@ -27,7 +27,6 @@
 
 #include <CGAL/disable_warnings.h>
 
-#include <CGAL/Surface_mesh.h>
 #include <CGAL/property_map.h>
 #include <CGAL/double.h>
 #include <CGAL/boost/graph/properties.h>
@@ -359,6 +358,7 @@ namespace Intrinsic_Delaunay_Triangulation_3 {
            copy_face_graph(tmref,tm, std::back_inserter(pairs));
            for(int i=0; i < pairs.size(); i++){
              v2v[pairs[i].second] = pairs[i].first;
+             vtov[pairs[i].first] = pairs[i].second;
            }
            
            vertex_id_map = get(Vertex_property_tag(),const_cast<TriangleMesh&>(tm));
@@ -378,7 +378,6 @@ namespace Intrinsic_Delaunay_Triangulation_3 {
            edge_id_map = get(Edge_property_tag(), const_cast<TriangleMesh&>(tm));
            Index edge_i = 0;
            BOOST_FOREACH(edge_descriptor ed, edges(tm)){
-             std::cout << ed << std::endl;
              mark_edges(edge_i,0)=1;
              edge_lengths(edge_i,0) = Polygon_mesh_processing::edge_length(halfedge(ed,tm),tm);
              put(edge_id_map, ed, edge_i++);
@@ -438,10 +437,52 @@ namespace Intrinsic_Delaunay_Triangulation_3 {
          Eigen::VectorXd edge_lengths;
          Eigen::VectorXd mark_edges;
      public:
-         std::map<vertex_descriptor,vertex_descriptor> v2v;
+       std::map<vertex_descriptor,vertex_descriptor> v2v, vtov;
       };
 
 } // namespace Intrinsic_Delaunay_Triangulation_3
+
+namespace Heat_method_3 {
+
+  
+  template <typename TM>
+  struct V2V {
+
+    V2V(const TM&)
+    {}
+    
+    template <typename T>
+    const T& operator()(const T& t) const
+    {
+      return t;
+    }
+  };
+
+template <typename TM,
+          typename T,
+          typename VDM,
+          typename VPM,
+          typename FIM,
+          typename EIM,
+          typename LA>
+struct V2V<CGAL::Intrinsic_Delaunay_Triangulation_3::Intrinsic_Delaunay_Triangulation_3<TM,T,VDM,VPM,FIM,EIM,LA> >
+{
+  typedef CGAL::Intrinsic_Delaunay_Triangulation_3::Intrinsic_Delaunay_Triangulation_3<TM,T,VDM,VPM,FIM,EIM,LA> Idt;
+  const Idt& idt;
+  typedef typename boost::graph_traits<Idt>::vertex_descriptor Idt_vertex_descriptor;
+  typedef typename boost::graph_traits<TM>::vertex_descriptor TM_vertex_descriptor;
+  
+  V2V(const Idt& idt)
+    : idt(idt)
+  {}
+    
+  Idt_vertex_descriptor operator()(const TM_vertex_descriptor& vd) const
+  {
+    return Idt_vertex_descriptor(idt.vtov.at(vd),idt.triangle_mesh());
+  }
+};
+  } // namespace Heat_method_3
+  
 } // namespace CGAL
 
 namespace boost {
diff --git a/Heat_method_3/test/Heat_method_3/heat_method_intrinsic_surface_mesh_test.cpp b/Heat_method_3/test/Heat_method_3/heat_method_intrinsic_surface_mesh_test.cpp
index fb9244cbcdf..8dec2ef3cfc 100644
--- a/Heat_method_3/test/Heat_method_3/heat_method_intrinsic_surface_mesh_test.cpp
+++ b/Heat_method_3/test/Heat_method_3/heat_method_intrinsic_surface_mesh_test.cpp
@@ -1,5 +1,6 @@
 #include <CGAL/Simple_cartesian.h>
 #include <CGAL/Surface_mesh.h>
+#include <CGAL/Polyhedron_3.h>
 #include <CGAL/Dynamic_property_map.h>
 #include <CGAL/Heat_method_3/Heat_method_3.h>
 #include <CGAL/Heat_method_3/Intrinsic_Delaunay_Triangulation_3.h>
@@ -14,7 +15,8 @@
 typedef CGAL::Simple_cartesian<double>                       Kernel;
 typedef Kernel::Point_3                                      Point;
 typedef Kernel::Point_2                                      Point_2;
-typedef CGAL::Surface_mesh<Point>                            Surface_mesh;
+//typedef CGAL::Surface_mesh<Point>                            Surface_mesh;
+typedef CGAL::Polyhedron_3<Kernel> Surface_mesh;
 
 typedef CGAL::dynamic_vertex_property_t<double> Vertex_distance_tag;
 typedef boost::property_map<Surface_mesh, Vertex_distance_tag >::type Vertex_distance_map;
@@ -123,9 +125,9 @@ int main(int argc, char*argv[])
  
   //source set tests
   Heat_method hm(idt, idt.vertex_distance_map());
-  // hm.add_source(* vertices(idt).first);
 
-  hm.add_source(boost::graph_traits<Idt>::vertex_descriptor(boost::graph_traits<Surface_mesh>::vertex_descriptor(0),sm));
+  hm.add_source(* vertices(sm).first);
+
   hm.update();
 
   BOOST_FOREACH(boost::graph_traits<Surface_mesh>::vertex_descriptor vd, vertices(sm)){
diff --git a/Heat_method_3/test/Heat_method_3/heat_method_surface_mesh_test.cpp b/Heat_method_3/test/Heat_method_3/heat_method_surface_mesh_test.cpp
index bf17b8d3f92..19d7fa7f675 100644
--- a/Heat_method_3/test/Heat_method_3/heat_method_surface_mesh_test.cpp
+++ b/Heat_method_3/test/Heat_method_3/heat_method_surface_mesh_test.cpp
@@ -106,7 +106,6 @@ int main()
 {
   Mesh sm;
   Vertex_distance_map vertex_distance_map = get(Vertex_distance_tag(),sm);
-  bool idf = false;
 
   std::ifstream in("data/pyramid0.off");
   in >> sm;
@@ -115,7 +114,7 @@ int main()
     return 1;
   }
   //source set tests
-  Heat_method hm(sm, vertex_distance_map, idf);
+  Heat_method hm(sm, vertex_distance_map);
   source_set_tests(hm,sm);
   //cotan matrix tests
   const SparseMatrix& M = hm.mass_matrix();
@@ -155,7 +154,7 @@ int main()
     std::cerr << "Problem loading the input data" << std::endl;
     return 1;
   }
-  Heat_method hm2(sm2, vertex_distance_map2, idf);
+  Heat_method hm2(sm2, vertex_distance_map2);
   //Eigen::VectorXd solved_dist_sphere = hm2.distances();
   const SparseMatrix& M2 = hm2.mass_matrix();
   const SparseMatrix& c2 = hm2.cotan_matrix();
@@ -186,7 +185,7 @@ int main()
     std::cerr << "Problem loading the input data" << std::endl;
     return 1;
   }
-  Heat_method hm3(sm3, vertex_distance_map3,idf);
+  Heat_method hm3(sm3, vertex_distance_map3);
   //Eigen::VectorXd solved_dist_sphere = hm2.distances();
   const SparseMatrix& M3 = hm3.mass_matrix();
   const SparseMatrix& c3 = hm3.cotan_matrix();