diff --git a/src/remesher.cpp b/src/remesher.cpp
index 84e5feb1..1253d4b6 100644
--- a/src/remesher.cpp
+++ b/src/remesher.cpp
@@ -6,7 +6,31 @@
 #ifdef WITH_CUDA
 #include <cuda_runtime.h>
 #endif
+#include <CGAL/boost/graph/graph_traits_Surface_mesh.h>
+#include <CGAL/Polygon_mesh_processing/remesh.h>
+#include <CGAL/Polygon_mesh_processing/border.h>
+#include <boost/function_output_iterator.hpp>
+#include <fstream>
+#include <vector>
 #include "remesher.h"
+#include "booleanmesh.h"
+
+typedef boost::graph_traits<CgalMesh>::halfedge_descriptor halfedge_descriptor;
+typedef boost::graph_traits<CgalMesh>::edge_descriptor     edge_descriptor;
+namespace PMP = CGAL::Polygon_mesh_processing;
+
+struct halfedge2edge
+{
+  halfedge2edge(const CgalMesh& m, std::vector<edge_descriptor>& edges)
+    : m_mesh(m), m_edges(edges)
+  {}
+  void operator()(const halfedge_descriptor& h) const
+  {
+    m_edges.push_back(edge(h, m_mesh));
+  }
+  const CgalMesh& m_mesh;
+  std::vector<edge_descriptor>& m_edges;
+};
 
 using namespace qflow;
 
@@ -36,6 +60,31 @@ const std::vector<std::pair<QUuid, QUuid>> &Remesher::getRemeshedVertexSources()
     return m_remeshedVertexSources;
 }
 
+void Remesher::isotropicRemesh(float targetEdgeLength, unsigned int iterationNum)
+{
+    CgalMesh *mesh = buildCgalMesh<CgalKernel>(m_vertices, m_triangles);
+    if (nullptr == mesh)
+        return;
+    
+    std::vector<edge_descriptor> border;
+    PMP::border_halfedges(faces(*mesh),
+        *mesh,
+        boost::make_function_output_iterator(halfedge2edge(*mesh, border)));
+    PMP::split_long_edges(border, targetEdgeLength, *mesh);
+    
+    PMP::isotropic_remeshing(faces(*mesh),
+        targetEdgeLength,
+        *mesh,
+        PMP::parameters::number_of_iterations(iterationNum)
+        .protect_constraints(true));
+    
+    m_vertices.clear();
+    m_triangles.clear();
+    fetchFromCgalMesh<CgalKernel>(mesh, m_vertices, m_triangles);
+    
+    delete mesh;
+}
+
 void Remesher::remesh()
 {
     Parametrizer field;
@@ -44,6 +93,8 @@ void Remesher::remesh()
     cudaFree(0);
 #endif
     
+    isotropicRemesh();
+    
     field.V.resize(3, m_vertices.size());
     field.F.resize(3, m_triangles.size());
     for (decltype(m_vertices.size()) i = 0; i < m_vertices.size(); i++) {
diff --git a/src/remesher.h b/src/remesher.h
index 7aa3a720..e1585567 100644
--- a/src/remesher.h
+++ b/src/remesher.h
@@ -27,6 +27,7 @@ private:
     std::vector<std::pair<QVector3D, float>> m_nodes;
     std::vector<std::pair<QUuid, QUuid>> m_sourceIds;
     void resolveSources();
+    void isotropicRemesh(float targetEdgeLength=0.021, unsigned int iterationNum=3);
 };
 
 #endif