Refactoring of swap and minor c hanges to other operations

components/wildmeshing/wmtk/components/wildmeshing/internal/wildmeshing2d.cpp

@@ -85,7 +85,10 @@
 std::vector<std::pair<std::shared_ptr<Mesh>, std::string>> wildmeshing2d(
     const WildMeshingOptions& options)
 {
-    auto mesh = options.input_mesh;
+    auto mesh = std::dynamic_pointer_cast<TriMesh>(options.input_mesh);
+    if(!bool(mesh)) {
+        throw std::runtime_error("input mesh of wildmeshing3d must be a trimesh");
+    }
 
     if (!mesh->is_connectivity_valid()) {
         throw std::runtime_error("input mesh for wildmeshing connectivity invalid");

components/wildmeshing/wmtk/components/wildmeshing/internal/wildmeshing3d.cpp

@@ -88,7 +88,10 @@
 std::vector<std::pair<std::shared_ptr<Mesh>, std::string>> wildmeshing3d(
     const WildMeshingOptions& options)
 {
-    auto mesh = options.input_mesh;
+    auto mesh = std::dynamic_pointer_cast<TetMesh>(options.input_mesh);
+    if(!bool(mesh)) {
+        throw std::runtime_error("input mesh of wildmeshing3d must be a tetmesh");
+    }
 
     if (!mesh->is_connectivity_valid()) {
         throw std::runtime_error("input mesh for wildmeshing connectivity invalid");
@@ -1775,4 +1778,4 @@
 }
 
 
-} // namespace wmtk::components::internal
+} // namespace wmtk::components::internal

src/wmtk/EdgeMesh.hpp

@@ -44,6 +44,7 @@ class EdgeMesh : public MeshCRTP<EdgeMesh>
         Eigen::Ref<const RowVectors2l> EE,
         Eigen::Ref<const VectorXl> VE);
 
+    using Mesh::is_valid;
     bool is_valid(const Tuple& tuple) const final override;
 
     bool is_connectivity_valid() const override;

src/wmtk/Mesh.cpp

@@ -38,7 +38,7 @@ simplex::IdSimplex Mesh::get_id_simplex(const simplex::Simplex& s) const
 simplex::Simplex Mesh::get_simplex(const simplex::IdSimplex& s) const
 {
     const Tuple& t = tuple_from_id(s.primitive_type(), s.index());
-    return simplex::Simplex(*this, s.primitive_type(), t);
+    return simplex::Simplex(s.primitive_type(), t);
 }
 
 Tuple Mesh::get_tuple_from_id_simplex(const simplex::IdSimplex& s) const
@@ -112,7 +112,16 @@ bool Mesh::is_boundary(const simplex::Simplex& s) const
 
 bool Mesh::is_valid(const Tuple& tuple) const
 {
-    return !tuple.is_null() && !is_removed(tuple);
+    const bool nullity = tuple.is_null();
+    const bool removed = is_removed(tuple);
+    const bool bad = nullity || removed;
+#if !defined(NDEBUG)
+    if(bad) {
+        logger().trace("Mesh::is_valid failed, got nullity:{} removedness:{}", nullity, removed);
+
+    }
+#endif
+    return !bad;
 }
 
 bool Mesh::is_removed(const Tuple& tuple) const
@@ -143,16 +152,7 @@ bool Mesh::is_removed(int64_t index, PrimitiveType pt) const
 
 bool Mesh::is_valid(const simplex::Simplex& s) const
 {
-#if defined(WMTK_ENABLE_SIMPLEX_ID_CACHING)
-    if (!is_valid(s.tuple())) {
-        return false;
-    } else {
-        const int64_t id_tuple = id(s.tuple(), s.primitive_type());
-        return id_tuple == s.m_index;
-    }
-#else
     return is_valid(s.tuple()) && !is_removed(s.tuple(), s.primitive_type());
-#endif
 }
 
 

src/wmtk/Mesh.hpp

@@ -29,6 +29,7 @@
 
 #include "multimesh/attribute/UseParentScopeRAII.hpp"
 
+
 #include "simplex/IdSimplex.hpp"
 #include "simplex/NavigatableSimplex.hpp"
 #include "simplex/Simplex.hpp"

src/wmtk/TetMesh.hpp

@@ -34,6 +34,7 @@ class TetMesh : public MeshCRTP<TetMesh>
     bool is_boundary_edge(const Tuple& tuple) const;
     bool is_boundary_face(const Tuple& tuple) const;
 
+    using Mesh::is_valid;
     bool is_valid(const Tuple& tuple) const final override;
 
     void initialize(

src/wmtk/TriMesh.hpp

@@ -54,6 +54,7 @@ class TriMesh : public MeshCRTP<TriMesh>
     void initialize(Eigen::Ref<const RowVectors3l> F, bool make_free = false);
     void initialize_free(int64_t count);
 
+    using Mesh::is_valid;
     bool is_valid(const Tuple& tuple) const final override;
 
     bool is_connectivity_valid() const final override;

src/wmtk/operations/EdgeSplit.cpp

@@ -20,7 +20,6 @@
     for (const auto& strat : m_new_attr_strategies) {
         if (strat->invalid_state()) {
             all_configured = false;
-            wmtk::logger().warn("Attribute new {} was not configured", strat->name());
         }
     }
     return all_configured;
@@ -86,12 +85,23 @@
     const attribute::MeshAttributeHandle& attribute,
     const std::shared_ptr<const operations::BaseSplitNewAttributeStrategy>& other)
 {
+    // checks through every attr strategy, keeping track of if it's found an attribute to detect duplicate attribute transfers
+    bool done = false;
     for (size_t i = 0; i < m_new_attr_strategies.size(); ++i) {
         if (m_new_attr_strategies[i]->matches_attribute(attribute)) {
+            if(done) {
+                throw std::runtime_error("Two of one attr strat");
+            }
+            auto old = m_new_attr_strategies[i];
             m_new_attr_strategies[i] = other;
-            return;
+            other->invalid_state();
+            done = true;
+            //return;
         }
     }
+    if(done) {
+        return ;
+    }
 
     throw std::runtime_error("unable to find attribute");
 }

src/wmtk/operations/Operation.cpp

@@ -65,9 +65,6 @@ void Operation::add_transfer_strategy(
 
 std::vector<simplex::Simplex> Operation::operator()(const simplex::Simplex& simplex)
 {
-    if (!mesh().is_valid(simplex)) {
-        return {};
-    }
     if (!before(simplex)) {
         return {};
     }
@@ -97,14 +94,9 @@ std::vector<simplex::Simplex> Operation::operator()(const simplex::Simplex& simp
 bool Operation::before(const simplex::Simplex& simplex) const
 {
     // const attribute::Accessor<int64_t> accessor = hash_accessor();
-
-    // if (!mesh().is_valid(
-    //         simplex.tuple(),
-    //         accessor)) { // TODO: chang to is_removed and resurrect later
-    //     return false;
-    // }
-
-    if (mesh().is_removed(simplex.tuple()) || !mesh().is_valid(simplex)) {
+    // we assume the current MeshType's is_valid calls Mesh::is_valid first, which checks if the
+    // simplex is removed or not
+    if (!mesh().is_valid(simplex)) {
         return false;
     }
 
@@ -152,7 +144,6 @@ void Operation::apply_attribute_transfer(const std::vector<simplex::Simplex>& di
     for (const auto& s : direct_mods) {
         if (!s.tuple().is_null()) {
             assert(m_mesh.is_valid(s));
-            assert(m_mesh.get_const_flag_accessor(s.primitive_type()).is_active(s));
             for (const simplex::IdSimplex& ss : simplex::closed_star_iterable(m_mesh, s)) {
                 all.add(ss);
             }
@@ -165,6 +156,7 @@ void Operation::apply_attribute_transfer(const std::vector<simplex::Simplex>& di
     for (const auto& at_ptr : m_attr_transfer_strategies) {
         if (&m_mesh == &(at_ptr->mesh())) {
             for (const simplex::IdSimplex& s : all.simplex_vector()) {
+                assert(m_mesh.get_const_flag_accessor(s.primitive_type()).is_active(s));
                 if (s.primitive_type() == at_ptr->primitive_type()) {
                     at_ptr->run(m_mesh.get_simplex(s));
                 }

src/wmtk/operations/composite/CMakeLists.txt

@@ -1,4 +1,8 @@
 set(SRC_FILES
+
+    EdgeSwap.cpp
+    EdgeSwap.hpp
+
     TriFaceSplit.hpp
     TriFaceSplit.cpp
 

src/wmtk/operations/composite/EdgeSwap.cpp

@@ -0,0 +1,15 @@
+#include "EdgeSwap.hpp"
+
+namespace wmtk::operations::composite {
+EdgeSwap::EdgeSwap(Mesh& m)
+    : EdgeSwap(m, std::make_shared<EdgeSplit>(m), std::make_shared<EdgeCollapse>(m))
+{}
+EdgeSwap::EdgeSwap(
+    Mesh& m,
+    std::shared_ptr<EdgeSplit> split,
+    std::shared_ptr<EdgeCollapse> collapse)
+    : Operation(m)
+    , m_split(std::move(split))
+    , m_collapse(std::move(collapse))
+{}
+} // namespace wmtk::operations::composite

src/wmtk/operations/composite/EdgeSwap.hpp

@@ -0,0 +1,23 @@
+#pragma once
+#include <wmtk/operations/EdgeCollapse.hpp>
+#include <wmtk/operations/EdgeSplit.hpp>
+#include <wmtk/operations/Operation.hpp>
+
+namespace wmtk::operations::composite {
+class EdgeSwap : public Operation
+{
+public:
+    EdgeSwap(Mesh& m);
+    EdgeSwap(Mesh& m, std::shared_ptr<EdgeSplit> split, std::shared_ptr<EdgeCollapse> collapse);
+
+    EdgeSplit& split() { return *m_split; }
+    EdgeCollapse& collapse() { return *m_collapse; }
+    const EdgeSplit& split() const { return *m_split; }
+    const EdgeCollapse& collapse() const { return *m_collapse; }
+
+protected:
+    std::shared_ptr<EdgeSplit> m_split;
+    std::shared_ptr<EdgeCollapse> m_collapse;
+};
+} // namespace wmtk::operations::composite
+

src/wmtk/operations/composite/ProjectOperation.cpp

@@ -22,66 +22,87 @@
     , m_main_op(main_op)
 {
     for (auto& pair : mesh_constaint_pairs) {
-        int64_t count = 0;
-        int64_t index = 0;
-
-        const std::vector<Tuple>& facest =
-            pair.first.mesh().get_all(pair.first.mesh().top_simplex_type());
-
-        const int64_t dim = int64_t(pair.first.mesh().top_simplex_type()) + 1;
-
-        Eigen::MatrixXd vertices(dim * facest.size(), pair.first.dimension());
-        Eigen::MatrixXi faces(facest.size(), dim);
+        add_constraint(pair.second, pair.first);
+    }
+}
 
-        // hugly copy paste
-        if (pair.first.holds<double>()) {
-            const attribute::Accessor<double> accessor =
-                pair.first.mesh().create_const_accessor(pair.first.as<double>());
-
-            for (const auto& f : facest) {
-                auto tmp = faces_single_dimension_tuples(
-                    pair.first.mesh(),
-                    simplex::Simplex(pair.first.mesh(), pair.first.mesh().top_simplex_type(), f),
-                    PrimitiveType::Vertex);
-
-                assert(tmp.size() == dim);
-                for (int64_t j = 0; j < tmp.size(); ++j) {
-                    auto p = accessor.const_vector_attribute(tmp[j]);
-                    faces(index, j) = count;
-                    vertices.row(dim * index + j) = p;
-
-                    ++count;
-                }
-                ++index;
+void ProjectOperation::add_constraint(
+    const attribute::MeshAttributeHandle& mah,
+    const attribute::MeshAttributeHandle& projection_mah)
+{
+    int64_t count = 0;
+    int64_t index = 0;
+
+    const std::vector<Tuple>& facest =
+        projection_mah.mesh().get_all(projection_mah.mesh().top_simplex_type());
+
+    const int64_t dim = int64_t(projection_mah.mesh().top_simplex_type()) + 1;
+
+    Eigen::MatrixXd vertices(dim * facest.size(), projection_mah.dimension());
+    Eigen::MatrixXi faces(facest.size(), dim);
+
+    // hugly copy paste
+    if (projection_mah.holds<double>()) {
+        const attribute::Accessor<double> accessor =
+            projection_mah.mesh().create_const_accessor(projection_mah.as<double>());
+
+        for (const auto& f : facest) {
+            auto tmp = faces_single_dimension_tuples(
+                projection_mah.mesh(),
+                simplex::Simplex(
+                    projection_mah.mesh(),
+                    projection_mah.mesh().top_simplex_type(),
+                    f),
+                PrimitiveType::Vertex);
+
+            assert(tmp.size() == dim);
+            for (int64_t j = 0; j < tmp.size(); ++j) {
+                auto p = accessor.const_vector_attribute(tmp[j]);
+                faces(index, j) = count;
+                vertices.row(dim * index + j) = p;
+
+                ++count;
             }
-        } else {
-            const attribute::Accessor<Rational> accessor =
-                pair.first.mesh().create_const_accessor(pair.first.as<Rational>());
-
-
-            for (const auto& f : facest) {
-                auto tmp = faces_single_dimension_tuples(
-                    pair.first.mesh(),
-                    simplex::Simplex(pair.first.mesh(), pair.first.mesh().top_simplex_type(), f),
-                    PrimitiveType::Vertex);
-
-                assert(tmp.size() == dim);
-                for (int64_t j = 0; j < tmp.size(); ++j) {
-                    auto p = accessor.const_vector_attribute(tmp[j]).cast<double>();
-                    faces(index, j) = count;
-                    vertices.row(dim * index + j) = p;
-
-                    ++count;
-                }
-                ++index;
+            ++index;
+        }
+    } else {
+        const attribute::Accessor<Rational> accessor =
+            projection_mah.mesh().create_const_accessor(projection_mah.as<Rational>());
+
+
+        for (const auto& f : facest) {
+            auto tmp = faces_single_dimension_tuples(
+                projection_mah.mesh(),
+                simplex::Simplex(
+                    projection_mah.mesh(),
+                    projection_mah.mesh().top_simplex_type(),
+                    f),
+                PrimitiveType::Vertex);
+
+            assert(tmp.size() == dim);
+            for (int64_t j = 0; j < tmp.size(); ++j) {
+                auto p = accessor.const_vector_attribute(tmp[j]).cast<double>();
+                faces(index, j) = count;
+                vertices.row(dim * index + j) = p;
+
+                ++count;
             }
+            ++index;
         }
+    }
 
-        auto bvh = std::make_shared<SimpleBVH::BVH>();
-        bvh->init(vertices, faces, 1e-10);
+    auto bvh = std::make_shared<SimpleBVH::BVH>();
+    bvh->init(vertices, faces, 1e-10);
+    m_bvh.emplace_back(mah, bvh);
+}
 
-        m_bvh.emplace_back(pair.second, bvh);
-    }
+void ProjectOperation::add_constraint(
+    const attribute::MeshAttributeHandle& mah,
+    std::shared_ptr<SimpleBVH::BVH> bvh)
+{
+    assert(mah.is_valid());
+    assert(bool(bvh));
+    m_bvh.emplace_back(mah, bvh);
 }
 
 std::vector<simplex::Simplex> ProjectOperation::execute(const simplex::Simplex& simplex)

src/wmtk/operations/composite/ProjectOperation.hpp

@@ -22,12 +22,20 @@ class ProjectOperation : public AttributesUpdate
 
     ProjectOperation(
         std::shared_ptr<Operation> main_op,
-        const std::vector<MeshConstrainPair>& mesh_constaint_pairs);
+        const std::vector<MeshConstrainPair>& mesh_constaint_pairs = {});
 
     std::vector<simplex::Simplex> execute(const simplex::Simplex& simplex) override;
     PrimitiveType primitive_type() const override { return m_main_op->primitive_type(); }
 
 
+    void add_constraint(
+        const attribute::MeshAttributeHandle& mah,
+        std::shared_ptr<SimpleBVH::BVH> bvh);
+
+    void add_constraint(
+        const attribute::MeshAttributeHandle& mah,
+        const attribute::MeshAttributeHandle& bvh_mesh_attribute);
+
 private:
     using BVHConstrainPair =
         std::pair<attribute::MeshAttributeHandle, std::shared_ptr<SimpleBVH::BVH>>;