MeshNavigation.h

/*
 * Copyright 2020 Adobe. All rights reserved.
 * This file is licensed to you under the Apache License, Version 2.0 (the "License");
 * you may not use this file except in compliance with the License. You may obtain a copy
 * of the License at http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software distributed under
 * the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR REPRESENTATIONS
 * OF ANY KIND, either express or implied. See the License for the specific language
 * governing permissions and limitations under the License.
 */
#pragma once
 
#include <lagrange/chain_corners_around_edges.h>
#include <lagrange/chain_corners_around_vertices.h>
#include <lagrange/corner_to_edge_mapping.h>
 
#include <Eigen/Core>
 
namespace lagrange {

template <typename MeshType>
class MeshNavigation
{
public:
    using Scalar = typename MeshType::Scalar;
    using Index = typename MeshType::Index;
    using FacetArray = typename MeshType::FacetArray;
    using IndexArray = Eigen::Matrix<Index, Eigen::Dynamic, 1>;
 
public:
    MeshNavigation(const MeshType& mesh) { initialize(mesh); }
 
public:
    Index get_num_edges() const { return safe_cast<Index>(m_e2c.size()); }

    Index get_edge(Index f, Index lv) const { return m_c2e(f * m_vertex_per_facet + lv); }

    Index get_edge_from_corner(Index c) const { return m_c2e(c); }

    Index get_first_corner_around_edge(Index e) const { return m_e2c[e]; }

    Index get_next_corner_around_edge(Index c) const { return m_next_corner_around_edge[c]; }

    Index get_first_corner_around_vertex(Index v) const { return m_v2c[v]; }

    Index get_next_corner_around_vertex(Index c) const { return m_next_corner_around_vertex[c]; }

    std::array<Index, 2> get_edge_vertices(const FacetArray& facets, Index e) const
    {
        Index c = m_e2c[e];
        if (c == invalid<Index>()) {
            logger().error("Invalid corner id for edge {}", e);
            throw std::runtime_error("Invalid edge id");
        }
        Index f = c / m_vertex_per_facet;
        Index lv = c % m_vertex_per_facet;
        return {{facets(f, lv), facets(f, (lv + 1) % m_vertex_per_facet)}};
    }

    Index get_vertex_opposite_edge(const FacetArray& facets, Index e) const
    {
        la_runtime_assert(m_vertex_per_facet == 3, "This method is only for triangle meshes.");
        Index c = m_e2c[e];
        if (c == invalid<Index>()) {
            logger().error("Invalid corner id for edge {}", e);
            throw std::runtime_error("Invalid edge id");
        }
        Index f = c / m_vertex_per_facet;
        Index lv = c % m_vertex_per_facet;
        return facets(f, (lv + 2) % m_vertex_per_facet);
    }

    Index get_num_facets_around_vertex(Index v) const
    {
        // Count number of incident facets
        Index num_incident_facets = 0;
        foreach_facets_around_vertex(v, [&](Index /*f*/) { ++num_incident_facets; });
        return num_incident_facets;
    }

    Index get_num_facets_around_edge(Index e) const
    {
        // Count number of incident facets
        Index num_incident_facets = 0;
        foreach_facets_around_edge(e, [&](Index /*f*/) { ++num_incident_facets; });
        return num_incident_facets;
    }

    Index get_one_facet_around_edge(Index e) const
    {
        Index c = m_e2c[e];
        if (c != invalid<Index>()) {
            return c / m_vertex_per_facet;
        } else {
            return invalid<Index>();
        }
    }

    Index get_one_corner_around_edge(Index e) const { return m_e2c[e]; }

    Index get_one_corner_around_vertex(Index v) const { return m_v2c[v]; }

    bool is_boundary_edge(Index e) const
    {
        Index c = m_e2c[e];
        assert(c != invalid<Index>());
        return (m_next_corner_around_edge[c] == invalid<Index>());
    }

    bool is_boundary_vertex(Index v) const { return m_is_boundary_vertex[v]; }
 
public:
    template <typename Func>
    void foreach_facets_around_vertex(Index v, Func func) const
    {
        // Loop over incident facets
        for (Index c = m_v2c[v]; c != invalid<Index>(); c = m_next_corner_around_vertex[c]) {
            Index f = c / m_vertex_per_facet;
            func(f);
        }
    }

    template <typename Func>
    void foreach_facets_around_edge(Index e, Func func) const
    {
        // Loop over incident facets
        for (Index c = m_e2c[e]; c != invalid<Index>(); c = m_next_corner_around_edge[c]) {
            Index f = c / m_vertex_per_facet;
            func(f);
        }
    }

    template <typename Func>
    void foreach_corners_around_vertex(Index v, Func func) const
    {
        // Loop over incident facets
        for (Index c = m_v2c[v]; c != invalid<Index>(); c = m_next_corner_around_vertex[c]) {
            func(c);
        }
    }

    template <typename Func>
    void foreach_corners_around_edge(Index e, Func func) const
    {
        // Loop over incident facets
        for (Index c = m_e2c[e]; c != invalid<Index>(); c = m_next_corner_around_edge[c]) {
            func(c);
        }
    }
 
protected:
    void initialize(const MeshType& mesh)
    {
        // Assumed to be constant
        m_vertex_per_facet = mesh.get_vertex_per_facet();
 
        // Compute unique edge ids
        corner_to_edge_mapping(mesh.get_facets(), m_c2e);
 
        // Chain corners around vertices and edges
        chain_corners_around_edges(mesh.get_facets(), m_c2e, m_e2c, m_next_corner_around_edge);
        chain_corners_around_vertices(
            mesh.get_num_vertices(),
            mesh.get_facets(),
            m_v2c,
            m_next_corner_around_vertex);
 
        // Tag boundary vertices
        m_is_boundary_vertex.assign(mesh.get_num_vertices(), false);
        for (Index e = 0; e < get_num_edges(); ++e) {
            if (is_boundary_edge(e)) {
                auto v = get_edge_vertices(mesh.get_facets(), e);
                m_is_boundary_vertex[v[0]] = true;
                m_is_boundary_vertex[v[1]] = true;
            }
        }
    }
 
protected:
    Index m_vertex_per_facet; 
 
    IndexArray m_c2e; 
    IndexArray m_e2c; 
    IndexArray m_next_corner_around_edge; 
    IndexArray m_v2c; 
    IndexArray m_next_corner_around_vertex; 
    std::vector<bool> m_is_boundary_vertex;
};
 
} // namespace lagrange