From 40a3d27e9b1a692ccee280cfebe9017392c070cf Mon Sep 17 00:00:00 2001
From: tomvanmele <tom.v.mele@gmail.com>
Date: Sun, 19 Nov 2023 20:07:06 +0100
Subject: [PATCH] loop and sqrt3

---
 src/compas_cgal/subdivision.py | 30 ++++++++++++++++--
 src/subdivision.cpp            | 56 ++++++++++++++++++++++++++++++----
 src/subdivision.h              | 16 ++++++++--
 3 files changed, 90 insertions(+), 12 deletions(-)

diff --git a/src/compas_cgal/subdivision.py b/src/compas_cgal/subdivision.py
index 0790e12..d2ce7ba 100644
--- a/src/compas_cgal/subdivision.py
+++ b/src/compas_cgal/subdivision.py
@@ -6,7 +6,7 @@
 
 
 def mesh_subdivide_catmull_clark(mesh: VerticesFaces, k=1) -> VerticesFacesNumpy:
-    """Subdivide a mesh withe the Catmull Clark scheme.
+    """Subdivide a mesh with the Catmull Clark scheme.
 
     Parameters
     ----------
@@ -38,7 +38,7 @@ def mesh_subdivide_catmull_clark(mesh: VerticesFaces, k=1) -> VerticesFacesNumpy
 
 
 def mesh_subdivide_loop(mesh: VerticesFaces, k=1) -> VerticesFacesNumpy:
-    """Subdivide a mesh withe the Loop scheme.
+    """Subdivide a mesh with the Loop scheme.
 
     Parameters
     ----------
@@ -52,4 +52,28 @@ def mesh_subdivide_loop(mesh: VerticesFaces, k=1) -> VerticesFacesNumpy:
     :attr:`compas_cgal.types.VerticesFacesNumpy`
 
     """
-    raise NotImplementedError
+    V, F = mesh
+    V = np.asarray(V, dtype=np.float64)
+    F = np.asarray(F, dtype=np.int32)
+    return subdivision.subd_loop(V, F, k)
+
+
+def mesh_subdivide_sqrt3(mesh: VerticesFaces, k=1) -> VerticesFacesNumpy:
+    """Subdivide a mesh with the Sqrt3 scheme.
+
+    Parameters
+    ----------
+    mesh : :attr:`compas_cgal.types.VerticesFaces`
+        The mesh to remesh.
+    k : int, optional
+        The number of subdivision steps.
+
+    Returns
+    -------
+    :attr:`compas_cgal.types.VerticesFacesNumpy`
+
+    """
+    V, F = mesh
+    V = np.asarray(V, dtype=np.float64)
+    F = np.asarray(F, dtype=np.int32)
+    return subdivision.subd_sqrt3(V, F, k)
diff --git a/src/subdivision.cpp b/src/subdivision.cpp
index b2ea038..984937c 100644
--- a/src/subdivision.cpp
+++ b/src/subdivision.cpp
@@ -1,11 +1,10 @@
 #include "subdivision.h"
 #include <CGAL/subdivision_method_3.h>
 
-
 std::tuple<compas::RowMatrixXd, compas::RowMatrixXi>
 subd_catmullclark(
-    compas::RowMatrixXd & V,
-    std::vector< std::vector<int> > & faces,
+    compas::RowMatrixXd &V,
+    std::vector<std::vector<int>> &faces,
     unsigned int k)
 {
     Mesh mesh = compas::ngon_from_vertices_and_faces(V, faces);
@@ -18,8 +17,40 @@ subd_catmullclark(
     return R;
 };
 
+std::tuple<compas::RowMatrixXd, compas::RowMatrixXi>
+subd_loop(
+    compas::RowMatrixXd &V,
+    compas::RowMatrixXi &F,
+    unsigned int k)
+{
+    Mesh mesh = compas::mesh_from_vertices_and_faces(V, F);
+
+    CGAL::Subdivision_method_3::Loop_subdivision(mesh, CGAL::parameters::number_of_iterations(k));
+
+    mesh.collect_garbage();
+
+    std::tuple<compas::RowMatrixXd, compas::RowMatrixXi> R = compas::mesh_to_vertices_and_faces(mesh);
+    return R;
+};
+
+std::tuple<compas::RowMatrixXd, compas::RowMatrixXi>
+subd_sqrt3(
+    compas::RowMatrixXd &V,
+    compas::RowMatrixXi &F,
+    unsigned int k)
+{
+    Mesh mesh = compas::mesh_from_vertices_and_faces(V, F);
 
-void init_subdivision(pybind11::module & m) {
+    CGAL::Subdivision_method_3::Sqrt3_subdivision(mesh, CGAL::parameters::number_of_iterations(k));
+
+    mesh.collect_garbage();
+
+    std::tuple<compas::RowMatrixXd, compas::RowMatrixXi> R = compas::mesh_to_vertices_and_faces(mesh);
+    return R;
+};
+
+void init_subdivision(pybind11::module &m)
+{
     pybind11::module submodule = m.def_submodule("subdivision");
 
     submodule.def(
@@ -27,6 +58,19 @@ void init_subdivision(pybind11::module & m) {
         &subd_catmullclark,
         pybind11::arg("V").noconvert(),
         pybind11::arg("faces").noconvert(),
-        pybind11::arg("k")
-    );
+        pybind11::arg("k"));
+
+    submodule.def(
+        "subd_loop",
+        &subd_loop,
+        pybind11::arg("V").noconvert(),
+        pybind11::arg("F").noconvert(),
+        pybind11::arg("k"));
+
+    submodule.def(
+        "subd_sqrt3",
+        &subd_sqrt3,
+        pybind11::arg("V").noconvert(),
+        pybind11::arg("F").noconvert(),
+        pybind11::arg("k"));
 };
diff --git a/src/subdivision.h b/src/subdivision.h
index 40eed55..071770f 100644
--- a/src/subdivision.h
+++ b/src/subdivision.h
@@ -3,12 +3,22 @@
 
 #include <compas.h>
 
-
 std::tuple<compas::RowMatrixXd, compas::RowMatrixXi>
 subd_catmullclark(
-    compas::RowMatrixXd & V,
-    std::vector< std::vector<int> > & faces,
+    compas::RowMatrixXd &V,
+    std::vector<std::vector<int>> &faces,
+    unsigned int k);
+
+std::tuple<compas::RowMatrixXd, compas::RowMatrixXi>
+subd_loop(
+    compas::RowMatrixXd &V,
+    compas::RowMatrixXi &F,
     unsigned int k);
 
+std::tuple<compas::RowMatrixXd, compas::RowMatrixXi>
+subd_sqrt3(
+    compas::RowMatrixXd &V,
+    compas::RowMatrixXi &F,
+    unsigned int k);
 
 #endif /* COMPAS_SUBD_H */