diff --git a/src/wmtk/Mesh.cpp b/src/wmtk/Mesh.cpp index e0cf77b20e..ecc2544a24 100644 --- a/src/wmtk/Mesh.cpp +++ b/src/wmtk/Mesh.cpp @@ -38,6 +38,11 @@ PrimitiveType Mesh::top_simplex_type() const } std::vector Mesh::get_all(PrimitiveType type) const +{ + return get_all(type, false); +} + +std::vector Mesh::get_all(PrimitiveType type, const bool include_deleted) const { ConstAccessor flag_accessor = get_flag_accessor(type); const attribute::CachingAccessor& flag_accessor_indices = flag_accessor.index_access(); @@ -45,9 +50,10 @@ std::vector Mesh::get_all(PrimitiveType type) const long cap = capacity(type); ret.reserve(cap); for (size_t index = 0; index < cap; ++index) { - if ((flag_accessor_indices.const_scalar_attribute(index) & 1)) { + if (flag_accessor_indices.const_scalar_attribute(index) & 1) ret.emplace_back(tuple_from_id(type, index)); - } + else if (include_deleted) + ret.emplace_back(); } return ret; } diff --git a/src/wmtk/Mesh.hpp b/src/wmtk/Mesh.hpp index 9e10b61856..c5f3a741c1 100644 --- a/src/wmtk/Mesh.hpp +++ b/src/wmtk/Mesh.hpp @@ -414,6 +414,15 @@ class Mesh : public std::enable_shared_from_this // hashes for top level simplices (i.e cells) to identify whether tuples // are invalid or not MeshAttributeHandle m_cell_hash_handle; + + + /** + * Generate a vector of Tuples from global vertex/edge/triangle/tetrahedron index + * @param type the type of tuple, can be vertex/edge/triangle/tetrahedron + * @param include_deleted if true returns also the deleted tuples (default false) + * @return vector of Tuples referring to each type + */ + std::vector get_all(PrimitiveType type, const bool include_deleted) const; }; diff --git a/src/wmtk/io/ParaviewWriter.cpp b/src/wmtk/io/ParaviewWriter.cpp index 48d6ecd874..196d4dd9ba 100644 --- a/src/wmtk/io/ParaviewWriter.cpp +++ b/src/wmtk/io/ParaviewWriter.cpp @@ -84,30 +84,35 @@ ParaviewWriter::ParaviewWriter( for (int i = 0; i < 4; ++i) { auto pt = PrimitiveType(i); if (m_enabled[i]) { - const auto tuples = mesh.get_all(pt); + // include deleted tuples so that attributes are aligned + const auto tuples = mesh.get_all(pt, true); cells[i].resize(tuples.size(), i + 1); for (size_t j = 0; j < tuples.size(); ++j) { const auto& t = tuples[j]; - long vid = mesh.id(t, PrimitiveType::Vertex); - cells[i](j, 0) = vid; - if (i > 0) { - auto t1 = mesh.switch_tuple(t, PrimitiveType::Vertex); - - cells[i](j, 1) = mesh.id(t1, PrimitiveType::Vertex); - } - if (i > 1) { - auto t1 = mesh.switch_tuple(t, PrimitiveType::Edge); - auto t2 = mesh.switch_tuple(t1, PrimitiveType::Vertex); - - cells[i](j, 2) = mesh.id(t2, PrimitiveType::Vertex); - } - if (i > 2) { - auto t1 = mesh.switch_tuple(t, PrimitiveType::Face); - auto t2 = mesh.switch_tuple(t1, PrimitiveType::Edge); - auto t3 = mesh.switch_tuple(t2, PrimitiveType::Vertex); - - cells[i](j, 3) = mesh.id(t3, PrimitiveType::Vertex); + if (t.is_null()) { + for (int d = 0; d < i; ++d) cells[i](j, d) = 0; + } else { + long vid = mesh.id(t, PrimitiveType::Vertex); + cells[i](j, 0) = vid; + if (i > 0) { + auto t1 = mesh.switch_tuple(t, PrimitiveType::Vertex); + + cells[i](j, 1) = mesh.id(t1, PrimitiveType::Vertex); + } + if (i > 1) { + auto t1 = mesh.switch_tuple(t, PrimitiveType::Edge); + auto t2 = mesh.switch_tuple(t1, PrimitiveType::Vertex); + + cells[i](j, 2) = mesh.id(t2, PrimitiveType::Vertex); + } + if (i > 2) { + auto t1 = mesh.switch_tuple(t, PrimitiveType::Face); + auto t2 = mesh.switch_tuple(t1, PrimitiveType::Edge); + auto t3 = mesh.switch_tuple(t2, PrimitiveType::Vertex); + + cells[i](j, 3) = mesh.id(t3, PrimitiveType::Vertex); + } } } } diff --git a/tests/test_io.cpp b/tests/test_io.cpp index 5b03b34fb6..0b59fad1f4 100644 --- a/tests/test_io.cpp +++ b/tests/test_io.cpp @@ -6,11 +6,22 @@ #include #include +#include +#include + +#include "tools/DEBUG_TriMesh.hpp" +#include "tools/TriMesh_examples.hpp" + #include #include using namespace wmtk; +using namespace wmtk::tests; + + +constexpr PrimitiveType PV = PrimitiveType::Vertex; +constexpr PrimitiveType PE = PrimitiveType::Edge; TEST_CASE("hdf5_2d", "[io]") { @@ -86,6 +97,69 @@ TEST_CASE("paraview_3d", "[io]") mesh.serialize(writer); } +TEST_CASE("attribute_after_split", "[io]") +{ + DEBUG_TriMesh m = single_triangle_with_position(); + wmtk::MeshAttributeHandle attribute_handle = + m.register_attribute(std::string("test_attribute"), PE, 1); + wmtk::MeshAttributeHandle pos_handle = + m.get_attribute_handle(std::string("position"), PV); + + { + Accessor acc_attribute = m.create_accessor(attribute_handle); + Accessor acc_pos = m.create_accessor(pos_handle); + + const Tuple edge = m.edge_tuple_between_v1_v2(0, 1, 0); + + Eigen::Vector3d p_mid; + { + const Eigen::Vector3d p0 = acc_pos.vector_attribute(edge); + const Eigen::Vector3d p1 = acc_pos.vector_attribute(m.switch_vertex(edge)); + p_mid = 0.5 * (p0 + p1); + } + + { + // set edge(0,1)'s tag as 1 + acc_attribute.scalar_attribute(edge) = 1; + + // all edges hold 0 besides "edge" + for (const Tuple& t : m.get_all(PE)) { + if (m.simplices_are_equal(Simplex::edge(edge), Simplex::edge(t))) { + CHECK(acc_attribute.scalar_attribute(t) == 1); + } else { + CHECK(acc_attribute.scalar_attribute(t) == 0); + } + } + + wmtk::operations::OperationSettings op_settings; + op_settings.split_boundary_edges = true; + op_settings.initialize_invariants(m); + + operations::tri_mesh::EdgeSplit op(m, edge, op_settings); + REQUIRE(op()); + + // set new vertex position + acc_pos.vector_attribute(op.return_tuple()) = p_mid; + } + + // since the default value is 0, there should be no other value in this triangle + for (const Tuple& t : m.get_all(PE)) { + CHECK(acc_attribute.scalar_attribute(t) == 0); + } + } // end of scope for the accessors + + { + Accessor acc_attribute = m.create_accessor(attribute_handle); + for (const Tuple& t : m.get_all(PE)) { + CHECK(acc_attribute.scalar_attribute(t) == 0); + } + } + + // attribute_after_split_edges.hdf contains a 1 in the "test_attribute" + ParaviewWriter writer("attribute_after_split", "position", m, true, true, true, false); + m.serialize(writer); +} + // TEST_CASE("io", "[io][mshio]") // { // using namespace wmtk;