2
0
Эх сурвалжийг харах

int->Index; depreacted pinv; trivial unit test

Alec Jacobson 5 жил өмнө
parent
commit
ad54f2dfe2

+ 7 - 6
include/igl/adjacency_matrix.cpp

@@ -77,18 +77,19 @@ IGL_INLINE void igl::adjacency_matrix(
   typedef Triplet<T> IJV;
   vector<IJV > ijv;
   ijv.reserve(C(C.size()-1)*2);
-  const typename DerivedI::Scalar n = I.maxCoeff()+1;
+  typedef typename DerivedI::Scalar Index;
+  const Index n = I.maxCoeff()+1;
   {
     // loop over polygons
-    for(int p = 0;p<C.size()-1;p++)
+    for(Index p = 0;p<C.size()-1;p++)
     {
       // number of edges
-      const int np = C(p+1)-C(p);
+      const Index np = C(p+1)-C(p);
       // loop over edges
-      for(int c = 0;c<np;c++)
+      for(Index c = 0;c<np;c++)
       {
-        const int i = I(C(p)+c);
-        const int j = I(C(p)+((c+1)%np));
+        const Index i = I(C(p)+c);
+        const Index j = I(C(p)+((c+1)%np));
         ijv.emplace_back(i,j,1);
         ijv.emplace_back(j,i,1);
       }

+ 22 - 21
include/igl/cotmatrix.cpp

@@ -101,39 +101,40 @@ IGL_INLINE void igl::cotmatrix(
   typedef Eigen::Matrix<Scalar,1,3> RowVector3S;
   typedef Eigen::Matrix<Scalar,Eigen::Dynamic,Eigen::Dynamic> MatrixXS;
   typedef Eigen::Matrix<Scalar,Eigen::Dynamic,1> VectorXS;
+  typedef Eigen::Index Index;
   // number of vertices
-  const int n = V.rows();
+  const Index n = V.rows();
   // number of polyfaces
-  const int m = C.size()-1;
+  const Index m = C.size()-1;
   assert(V.cols() == 2 || V.cols() == 3);
   std::vector<Eigen::Triplet<Scalar> > Lfijv;
   std::vector<Eigen::Triplet<Scalar> > Mfijv;
   std::vector<Eigen::Triplet<Scalar> > Pijv;
   // loop over vertices; set identity for original vertices
-  for(int i = 0;i<V.rows();i++) { Pijv.emplace_back(i,i,1); }
+  for(Index i = 0;i<V.rows();i++) { Pijv.emplace_back(i,i,1); }
   // loop over faces
-  for(int p = 0;p<C.size()-1;p++)
+  for(Index p = 0;p<C.size()-1;p++)
   {
     // number of faces/vertices in this simple polygon
-    const int np = C(p+1)-C(p);
+    const Index np = C(p+1)-C(p);
     // Working "local" list of vertices; last vertex is new one
     // this needs to have 3 columns so Eigen doesn't complain about cross
     // products below.
     Eigen::Matrix<Scalar,Eigen::Dynamic,3> X = decltype(X)::Zero(np+1,3);
-    for(int i = 0;i<np;i++){ X.row(i).head(V.cols()) = V.row(I(C(p)+i)); };
+    for(Index i = 0;i<np;i++){ X.row(i).head(V.cols()) = V.row(I(C(p)+i)); };
     // determine weights definig position of inserted vertex
     {
       MatrixXS A = decltype(A)::Zero(np+1,np);
       // My equation (38) would be A w = b.
       VectorXS b = decltype(b)::Zero(np+1);
-      for(int k = 0;k<np;k++)
+      for(Index k = 0;k<np;k++)
       { 
         const RowVector3S Xkp1mk = X.row((k+1)%np)-X.row(k);
         const RowVector3S Xkp1mkck = Xkp1mk.cross(X.row(k));
-        for(int i = 0;i<np;i++)
+        for(Index i = 0;i<np;i++)
         { 
           b(i) -= 2.*(X.row(i).cross(Xkp1mk)).dot(Xkp1mkck);
-          for(int j = 0;j<np;j++)
+          for(Index j = 0;j<np;j++)
           { 
             A(i,j) += 2.*(X.row(j).cross(Xkp1mk)).dot(X.row(i).cross(Xkp1mk));
           }
@@ -145,12 +146,12 @@ IGL_INLINE void igl::cotmatrix(
         Eigen::CompleteOrthogonalDecomposition<Eigen::MatrixXd>(A).solve(b);
       X.row(np) = w.transpose()*X.topRows(np);
       // scatter w into new row of P
-      for(int i = 0;i<np;i++) { Pijv.emplace_back(n+p,I(C(p)+i),w(i)); }
+      for(Index i = 0;i<np;i++) { Pijv.emplace_back(n+p,I(C(p)+i),w(i)); }
     }
     // "local" fan of faces. These could be statically cached, but this will
     // not be the bottleneck.
     Eigen::MatrixXi F(np,3);
-    for(int i = 0;i<np;i++)
+    for(Index i = 0;i<np;i++)
     { 
       F(i,0) = i; 
       F(i,1) = (i+1)%np; 
@@ -167,16 +168,16 @@ IGL_INLINE void igl::cotmatrix(
       Mp = M.diagonal();
     }
     // Scatter into fine Laplacian and mass matrices
-    const auto J = [&n,&np,&p,&I,&C](int i)->int{return i==np?n+p:I(C(p)+i);};
+    const auto J = [&n,&np,&p,&I,&C](Index i)->Index{return i==np?n+p:I(C(p)+i);};
     // Should just build Mf as a vector...
-    for(int i = 0;i<np+1;i++) { Mfijv.emplace_back(J(i),J(i),Mp(i)); }
+    for(Index i = 0;i<np+1;i++) { Mfijv.emplace_back(J(i),J(i),Mp(i)); }
     // loop over faces
-    for(int f = 0;f<np;f++)
+    for(Index f = 0;f<np;f++)
     {
-      for(int c = 0;c<3;c++)
+      for(Index c = 0;c<3;c++)
       {
-        const int i = F(f,(c+1)%3);
-        const int j = F(f,(c+2)%3);
+        const Index i = F(f,(c+1)%3);
+        const Index j = F(f,(c+2)%3);
         // symmetric off-diagonal
         Lfijv.emplace_back(J(i),J(j),K(f,c));
         Lfijv.emplace_back(J(j),J(i),K(f,c));
@@ -202,11 +203,11 @@ IGL_INLINE void igl::cotmatrix(
   MatrixXS Vf = P*V;
   Eigen::MatrixXi Ff(I.size(),3);
   {
-    int f = 0;
-    for(int p = 0;p<C.size()-1;p++)
+    Index f = 0;
+    for(Index p = 0;p<C.size()-1;p++)
     {
-      const int np = C(p+1)-C(p);
-      for(int c = 0;c<np;c++)
+      const Index np = C(p+1)-C(p);
+      for(Index c = 0;c<np;c++)
       {
         Ff(f,0) = I(C(p)+c);
         Ff(f,1) = I(C(p)+(c+1)%np);

+ 3 - 2
include/igl/pinv.h

@@ -1,6 +1,7 @@
 #ifndef IGL_PINV_H
 #define IGL_PINV_H
 #include "igl_inline.h"
+#include "deprecated.h"
 #include <Eigen/Core>
 namespace igl
 {
@@ -16,13 +17,13 @@ namespace igl
   // Obsolete: Use Eigen::CompleteOrthogonalDecomposition<Eigen::MatrixXd>
   // .solve() or .pseudoinverse() instead.
   template <typename DerivedA, typename DerivedX>
-  void pinv(
+  IGL_DEPRECATED void pinv(
     const Eigen::MatrixBase<DerivedA> & A,
     typename DerivedA::Scalar tol,
     Eigen::PlainObjectBase<DerivedX> & X);
   // Wrapper using default tol
   template <typename DerivedA, typename DerivedX>
-  void pinv(
+  IGL_DEPRECATED void pinv(
     const Eigen::MatrixBase<DerivedA> & A,
     Eigen::PlainObjectBase<DerivedX> & X);
 }

+ 7 - 3
include/igl/polygon_corners.cpp

@@ -16,13 +16,17 @@ IGL_INLINE void igl::polygon_corners(
   Eigen::PlainObjectBase<DerivedI> & I,
   Eigen::PlainObjectBase<DerivedC> & C)
 {
-  std::vector<int> vI;vI.reserve(P.size()*4);
+  typedef typename DerivedI::Scalar IType;
+  // JD: Honestly you could do a first loop over P, compute C, and then fill the
+  // entries of I directly. No need for guesses and push_back(), or the extra
+  // copy at the end. That would be more efficient.
+  std::vector<IType> vI;vI.reserve(P.size()*4);
   C.resize(P.size()+1);
   C(0) = 0;
-  for(int p = 0;p<P.size();p++)
+  for(size_t p = 0;p<P.size();p++)
   {
     C(p+1) = C(p)+P[p].size();
-    for(int c = 0;c<P[p].size();c++)
+    for(size_t c = 0;c<P[p].size();c++)
     {
       vI.push_back(P[p][c]);
     }

+ 23 - 22
include/igl/predicates/polygons_to_triangles.cpp

@@ -17,14 +17,15 @@ IGL_INLINE void igl::predicates::polygons_to_triangles(
   Eigen::PlainObjectBase<DerivedF> & F,
   Eigen::PlainObjectBase<DerivedJ> & J)
 {
+  typedef Eigen::Index Index;
   // Each polygon results in #sides-2 triangles. So ∑#sides-2
   F.resize(C(C.size()-1) - (C.size()-1)*2,3);
   J.resize(F.rows());
   {
-    int f = 0;
-    for(int p = 0;p<C.size()-1;p++)
+    Index f = 0;
+    for(Index p = 0;p<C.size()-1;p++)
     {
-      const int np = C(p+1)-C(p);
+      const Index np = C(p+1)-C(p);
       Eigen::MatrixXi pF;
       if(np == 3)
       {
@@ -33,7 +34,7 @@ IGL_INLINE void igl::predicates::polygons_to_triangles(
       {
         // Make little copy of this polygon with an initial fan
         DerivedV pV(np,V.cols());
-        for(int c = 0;c<np;c++)
+        for(Index c = 0;c<np;c++)
         {
           pV.row(c) = V.row(I(C(p)+c));
         }
@@ -69,11 +70,11 @@ IGL_INLINE void igl::predicates::polygons_to_triangles(
         // compute signed area
         {
           double area = 0;
-          for(int c = 0;c<np;c++)
+          for(Index c = 0;c<np;c++)
           {
             area += S((c+0)%np,0)*S((c+1)%np,1) - S((c+1)%np,0)*S((c+0)%np,1);
           }
-          //printf("area: %g\n",area);
+          //prIndexf("area: %g\n",area);
           if(area<0)
           {
             S.col(0) *= -1;
@@ -97,7 +98,7 @@ IGL_INLINE void igl::predicates::polygons_to_triangles(
           //std::cout<<std::endl;
 
           pF.resize(np-2,3);
-          for(int c = 0;c<np;c++)
+          for(Index c = 0;c<np;c++)
           {
             if(c>0 && c<np-1)
             {
@@ -117,13 +118,13 @@ IGL_INLINE void igl::predicates::polygons_to_triangles(
         //  Eigen::MatrixXd pl;
         //  igl::edge_lengths(pV,pF,pl);
 
-        //  typedef Eigen::Matrix<int,Eigen::Dynamic,2> MatrixX2I;
-        //  typedef Eigen::Matrix<int,Eigen::Dynamic,1> VectorXI;
+        //  typedef Eigen::Matrix<Index,Eigen::Dynamic,2> MatrixX2I;
+        //  typedef Eigen::Matrix<Index,Eigen::Dynamic,1> VectorXI;
         //  MatrixX2I E,uE;
         //  VectorXI EMAP;
-        //  std::vector<std::vector<int> > uE2E;
+        //  std::vector<std::vector<Index> > uE2E;
         //  igl::unique_edge_map(pF, E, uE, EMAP, uE2E);
-        //  typedef int Index;
+        //  typedef Index Index;
         //  typedef double Scalar;
         //  const Index num_faces = pF.rows();
         //  std::vector<Index> Q;
@@ -139,11 +140,11 @@ IGL_INLINE void igl::predicates::polygons_to_triangles(
         //    if (uE2E[uei].size() == 2) 
         //    {
         //      double w;
-        //      igl::is_intrinsic_delaunay(pl,uE2E,num_faces,uei,w);
-        //      printf("%d : %0.17f\n",uei,w);
+        //      igl::is_Indexrinsic_delaunay(pl,uE2E,num_faces,uei,w);
+        //      prIndexf("%d : %0.17f\n",uei,w);
         //      if(w<-1e-7) 
         //      {
-        //        printf("  flippin'\n");
+        //        prIndexf("  flippin'\n");
         //        //
         //        //          v1                 v1
         //        //          /|\                / \
@@ -187,7 +188,7 @@ IGL_INLINE void igl::predicates::polygons_to_triangles(
         //          pl(f2,1) = d;
         //          pl(f2,2) = a;
         //        }
-        //        printf("%d,%d %d,%d -> %d,%d\n",uE(uei,0),uE(uei,1),v1,v2,v3,v4);
+        //        prIndexf("%d,%d %d,%d -> %d,%d\n",uE(uei,0),uE(uei,1),v1,v2,v3,v4);
         //        igl::flip_edge(pF, E, uE, EMAP, uE2E, uei);
         //        std::cout<<"  "<<pl.row(f1)<<std::endl;
         //        std::cout<<"  "<<pl.row(f2)<<std::endl;
@@ -195,9 +196,9 @@ IGL_INLINE void igl::predicates::polygons_to_triangles(
         //        //igl::edge_lengths(pV,pF,pl);
         //        // recompute edge lengths of two faces. (extra work on untouched
         //        // edges)
-        //        for(int f : {f1,f2})
+        //        for(Index f : {f1,f2})
         //        {
-        //          for(int c=0;c<3;c++)
+        //          for(Index c=0;c<3;c++)
         //          {
         //            pl(f,c) = 
         //              (pV.row(pF(f,(c+1)%3))-pV.row(pF(f,(c+2)%3))).norm();
@@ -217,19 +218,19 @@ IGL_INLINE void igl::predicates::polygons_to_triangles(
 
 
         //  // check for self-loops (I claim these cannot happen)
-        //  for(int f = 0;f<pF.rows();f++)
+        //  for(Index f = 0;f<pF.rows();f++)
         //  {
-        //    for(int c =0;c<3;c++)
+        //    for(Index c =0;c<3;c++)
         //    {
         //      assert(pF(f,c) != pF(f,(c+1)%3) && "self loops should not exist");
         //    }
         //  }
         //}
       }
-      // Copy into global list
-      for(int i = 0;i<pF.rows();i++)
+      // Copy Indexo global list
+      for(Index i = 0;i<pF.rows();i++)
       {
-        for(int c =0;c<3;c++)
+        for(Index c =0;c<3;c++)
         {
           F(f,c) = I(C(p)+pF(i,c));
         }

+ 26 - 0
tests/include/igl/cotmatrix.cpp

@@ -1,6 +1,32 @@
 #include <test_common.h>
 #include <igl/PI.h>
 #include <igl/cotmatrix.h>
+#include <igl/matrix_to_list.h>
+#include <igl/polygon_corners.h>
+
+TEST_CASE("cotmatrix: poly", "[igl]" )
+{
+  const auto test_case = [](const std::string &param)
+  {
+    Eigen::MatrixXd V;
+    Eigen::MatrixXi F;
+    // Load example mesh: GetParam() will be name of mesh file
+    igl::read_triangle_mesh(test_common::data_path(param), V, F);
+    Eigen::SparseMatrix<double> tL,pL,pM,pP;
+    igl::cotmatrix(V,F,tL);
+    std::vector<std::vector<int> > vF;
+    igl::matrix_to_list(F,vF);
+    // trivial polygon mesh
+    Eigen::VectorXi I,C;
+    igl::polygon_corners(vF,I,C);
+    igl::cotmatrix(V,I,C,pL,pM,pP);
+    REQUIRE (tL.cols() == pL.cols());
+    REQUIRE (tL.rows() == pL.rows());
+    REQUIRE ( tL.isApprox(pL,1e-7) );
+  };
+
+  test_common::run_test_cases(test_common::all_meshes(), test_case);
+}
 
 TEST_CASE("cotmatrix: constant_in_null_space", "[igl]" "[slow]")
 {