fix known issue caused by chain cycle

network_process/interface/fwd_types.hpp

@@ -100,12 +100,12 @@ struct vertex_header_t {
 struct iso_vertex_t {
     vertex_header_t         header{};
     std::array<uint32_t, 4> simplex_vertex_indices{}; /// ? index into a list of tet vertices
-    std::array<uint16_t, 3> subface_indices = {
+    // std::array<uint16_t, 3> subface_indices = {
-        short_invalid_index,
+    //     short_invalid_index,
-        short_invalid_index,
+    //     short_invalid_index,
-        short_invalid_index}; /// list of subfaces whose isosurfaces pass IsoVert (indexed into a global list of
+    //     short_invalid_index}; /// list of subfaces whose isosurfaces pass IsoVert (indexed into a global list of
-    /// subfaces)
+    // /// subfaces)
-    /// so current solution only supports at most 3 iso-surfaces passing through a vertex?
+    // /// so current solution only supports at most 3 iso-surfaces passing through a vertex?
 };
 
 // given the list of vertex indices of a face, return the unique key of the face: (the smallest vert Id,

network_process/src/post_topo/map_chains.cpp

@@ -27,12 +27,13 @@ void map_chain_to_parametric_plane(const baked_blobtree_t&
         const auto subface_index = subface_index_;
         const auto subface_ptr   = subfaces[subface_index];
         const auto subface_type  = subface_types[subface_index];
-        auto  mapping_func = std::bind(internal::get_map_point_to_param_ptr(subface_type), subface_ptr, std::placeholders::_1);
+        auto mapping_func = std::bind(internal::get_map_point_to_param_ptr(subface_type), subface_ptr, std::placeholders::_1);
-        auto& parametric_plane            = parametric_planes[subface_index];
+        const auto cycle  = internal::get_cycle(subface_type);
-        auto& chain_vertices              = parametric_plane.chain_vertices;
+        auto&      parametric_plane            = parametric_planes[subface_index];
-        auto& chain_group_indices         = parametric_plane.chain_group_indices;
+        auto&      chain_vertices              = parametric_plane.chain_vertices;
-        auto& chain_vertex_flags          = parametric_plane.vertex_special_flags;
+        auto&      chain_group_indices         = parametric_plane.chain_group_indices;
-        auto& chain_other_subface_indices = parametric_plane.chain_other_subface_indices;
+        auto&      chain_vertex_flags          = parametric_plane.vertex_special_flags;
+        auto&      chain_other_subface_indices = parametric_plane.chain_other_subface_indices;
         chain_group_indices.reserve(patch_indices.size());
 
         unique_chain_indices.clear();
@@ -42,16 +43,19 @@ void map_chain_to_parametric_plane(const baked_blobtree_t&
             unique_chain_indices.insert(unique_chain_indices.end(), chain_indices.begin(), chain_indices.end());
         }
         std::sort(unique_chain_indices.begin(), unique_chain_indices.end());
-        unique_end_iter = std::unique(unique_chain_indices.begin(), unique_chain_indices.end());
+        unique_end_iter               = std::unique(unique_chain_indices.begin(), unique_chain_indices.end());
-        for (auto iter = unique_chain_indices.begin(); iter != unique_end_iter; ++iter)
-            old_chain_index_to_unique_index[*iter] = std::distance(unique_chain_indices.begin(), iter);
-
         const auto unique_chain_count = std::distance(unique_chain_indices.begin(), unique_end_iter);
+        for (uint32_t i = 0; i < unique_chain_count; ++i) {
+            const auto& chain_index                      = unique_chain_indices[i];
+            old_chain_index_to_unique_index[chain_index] = i;
+        }
+
         chain_vertices.reserve(unique_chain_count);
         chain_vertex_flags.reserve(unique_chain_count);
         chain_other_subface_indices.reserve(unique_chain_count);
-        for (const auto& chain_index : unique_chain_indices) {
+        for (auto iter = unique_chain_indices.begin(); iter != unique_end_iter; ++iter) {
-            const auto& chain = chains[chain_index];
+            const auto  chain_index = *iter;
+            const auto& chain       = chains[chain_index];
 
             auto& chain_vertices_ = chain_vertices.emplace_back();
             chain_vertices_.resize(chain.size());
@@ -72,6 +76,33 @@ void map_chain_to_parametric_plane(const baked_blobtree_t&
             chain_other_subface_indices.emplace_back(other_subface_index);
         }
 
+        // handle chain vertices that cross cycle boundary
+        // it is assumed that cycle is in U direction, and happens at U = 0
+        if (cycle < std::numeric_limits<double>::infinity()) {
+            int32_t     cycle_offset{};
+            const auto& world_to_local          = subface_ptr->world_to_local;
+            auto        unique_chain_index_iter = unique_chain_indices.begin();
+            auto        chain_vertex_iter       = chain_vertices.begin();
+            for (; unique_chain_index_iter != unique_end_iter; ++unique_chain_index_iter, ++chain_vertex_iter) {
+                const auto& chain_index = *unique_chain_index_iter;
+                const auto& chain       = chains[chain_index];
+                auto&       chain_verts = *chain_vertex_iter;
+                for (uint32_t i = 1; i < chain_verts.size(); ++i) {
+                    auto world_position1 = vertices[vertex_old_index_to_unique_index.at(chain[i - 1])];
+                    auto world_position2 = vertices[vertex_old_index_to_unique_index.at(chain[i])];
+                    auto local_position1 = world_to_local.linear() * world_position1 + world_to_local.translation();
+                    auto local_position2 = world_to_local.linear() * world_position2 + world_to_local.translation();
+                    if (local_position1.x() >= 0 && local_position2.x() >= 0) {
+                        if (local_position1.y() <= 0 && local_position2.y() > 0)
+                            ++cycle_offset;
+                        else if (local_position1.y() > 0 && local_position2.y() <= 0)
+                            --cycle_offset;
+                    }
+                    chain_verts[i].x() += cycle * cycle_offset;
+                }
+            }
+        }
+
         for (const auto& patch_index : patch_indices) {
             const auto& chain_indices = chain_of_patch.at(patch_index);
             auto&       chain_group   = chain_group_indices.emplace_back();
@@ -91,12 +122,12 @@ void remap_parametric_plane_vertices(flat_hash_map_mp<uint32_t, parametric_plane
             for (auto& vertex : chain_vertices) uv_bounds = uv_bounds.extend(vertex);
         }
         assert(uv_bounds.sizes().minCoeff() > epsilon);
-        const auto uv_bounds_rpc = uv_bounds.sizes().cwiseInverse().array();
+        // const auto uv_bounds_rpc = uv_bounds.sizes().cwiseInverse().array();
-        for (auto& chain_vertices : parametric_plane.chain_vertices) {
+        // for (auto& chain_vertices : parametric_plane.chain_vertices) {
-            std::transform(chain_vertices.begin(),
+        //     std::transform(chain_vertices.begin(),
-                           chain_vertices.end(),
+        //                    chain_vertices.end(),
-                           chain_vertices.begin(),
+        //                    chain_vertices.begin(),
-                           [&](const Eigen::Vector2d& uv) { return (uv - uv_bounds.min()).array() * uv_bounds_rpc; });
+        //                    [&](const Eigen::Vector2d& uv) { return (uv - uv_bounds.min()).array() * uv_bounds_rpc; });
-        }
+        // }
     }
 }

primitive_process/interface/base/subface.hpp

@@ -95,4 +95,7 @@ def_export_get_function_ptr(eval_dv_constraint);
 def_export_get_function_ptr(eval_intersection);
 
 #undef def_export_get_function_ptr
+
+// for getting cycle for each kind of subface
+PE_API double get_cycle(surface_type);
 } // namespace internal

primitive_process/src/base/subface.cpp

@@ -1,4 +1,5 @@
 #include <base/subface.hpp>
+#include "math/math_defs.hpp"
 
 namespace internal
 {
@@ -34,4 +35,14 @@ def_export_get_function_ptr_impl(eval_dv_constraint);
 def_export_get_function_ptr_impl(eval_intersection);
 
 #undef def_export_get_function_ptr_impl
+
+PE_API double get_cycle(surface_type type)
+{
+    switch (type) {
+        case surface_type::sphere:   return two_pi;
+        case surface_type::cylinder: return two_pi;
+        case surface_type::cone:     return sqrt_2 * pi;
+        default:                     return std::numeric_limits<double>::infinity();
+    }
+}
 } // namespace internal

shared_module/math/math_defs.hpp

@@ -9,6 +9,8 @@ static constexpr auto two_pi      = pi * 2;
 static constexpr auto pi_div_2    = pi / 2;
 static constexpr auto inv_pi      = 1. / pi;
 static constexpr auto inv_two_pi  = 1. / two_pi;
+static constexpr auto sqrt_2      = 1.4142135623730951;
+static constexpr auto sqrt_3      = 1.7320508075688772;
 static const auto     x_direction = Eigen::Vector3d{1.0, 0.0, 0.0};
 
 inline void vec3d_conversion(const vector3d& src, Eigen::Ref<Eigen::Vector3d> dst)