CompositeEnvelopes.hpp

#pragma once
 
#include <wmtk/SurfaceTagAttributes.h>
#include <wmtk/envelope/Envelope.hpp>
#include <wmtk/utils/Logger.hpp>
 
#include <cstdint>
#include <map>
#include <memory>
#include <mutex>
#include <shared_mutex>
#include <utility>
#include <vector>
 
namespace wmtk {
 
class IntersectionEnvelope final : public Envelope
{
public:
    explicit IntersectionEnvelope(std::vector<std::shared_ptr<Envelope>> members)
        : Envelope(/*exact=*/false)
        , m_members(std::move(members))
    {}
 
    bool is_outside(const std::array<Eigen::Vector3d, 3>& tris) const override
    {
        for (const auto& e : m_members) {
            if (e->is_outside(tris)) return true;
        }
        return false;
    }
 
    bool is_outside(const Eigen::Vector3d& pts) const override
    {
        for (const auto& e : m_members) {
            if (e->is_outside(pts)) return true;
        }
        return false;
    }
 
    bool is_outside(const std::array<Eigen::Vector2d, 2>& edge) const override
    {
        for (const auto& e : m_members) {
            if (e->is_outside(edge)) return true;
        }
        return false;
    }
 
    bool is_outside(const Eigen::Vector2d& pts) const override
    {
        for (const auto& e : m_members) {
            if (e->is_outside(pts)) return true;
        }
        return false;
    }
 
    const std::vector<std::shared_ptr<Envelope>>& members() const { return m_members; }
 
private:
    std::vector<std::shared_ptr<Envelope>> m_members;
};
 
class UnionEnvelope final : public Envelope
{
public:
    explicit UnionEnvelope(std::vector<std::shared_ptr<Envelope>> members)
        : Envelope(/*exact=*/false)
        , m_members(std::move(members))
    {}
 
    bool is_outside(const std::array<Eigen::Vector3d, 3>& tris) const override
    {
        for (const auto& e : m_members) {
            if (!e->is_outside(tris)) return false;
        }
        return true;
    }
 
    bool is_outside(const Eigen::Vector3d& pts) const override
    {
        for (const auto& e : m_members) {
            if (!e->is_outside(pts)) return false;
        }
        return true;
    }
 
    bool is_outside(const std::array<Eigen::Vector2d, 2>& edge) const override
    {
        for (const auto& e : m_members) {
            if (!e->is_outside(edge)) return false;
        }
        return true;
    }
 
    bool is_outside(const Eigen::Vector2d& pts) const override
    {
        for (const auto& e : m_members) {
            if (!e->is_outside(pts)) return false;
        }
        return true;
    }
 
private:
    std::vector<std::shared_ptr<Envelope>> m_members;
};
 
class InputEnvelopes
{
public:
    void set(std::vector<std::shared_ptr<Envelope>> members)
    {
        if (members.size() > size_t(SurfaceTagAttributes::max_inputs)) {
            log_and_throw_error(
                "{} inputs; a tracked simplex records its inputs in {} bits, so at most {} are "
                "supported",
                members.size(),
                SurfaceTagAttributes::max_inputs,
                SurfaceTagAttributes::max_inputs);
        }
        m_members = std::move(members);
        m_intersections.clear();
    }
 
    size_t size() const { return m_members.size(); }
    bool empty() const { return m_members.empty(); }
    const std::shared_ptr<Envelope>& member(size_t i) const { return m_members[i]; }
 
    std::shared_ptr<Envelope> containment(const InputSet set) const
    {
        check(set);
        if ((set & (set - 1)) == 0) return m_members[lowest_bit(set)];
        // Read-mostly: a set's intersection is built once, then every check only looks it up,
        // concurrently. An exclusive lock on every lookup serialized the operations' threads.
        {
            std::shared_lock<std::shared_mutex> lock(m_mutex);
            const auto it = m_intersections.find(set);
            if (it != m_intersections.end()) return it->second;
        }
        std::unique_lock<std::shared_mutex> lock(m_mutex);
        std::shared_ptr<Envelope>& e = m_intersections[set];
        if (!e) {
            std::vector<std::shared_ptr<Envelope>> members;
            for (size_t i = 0; i < m_members.size(); ++i) {
                if (set & (InputSet(1) << i)) members.push_back(m_members[i]);
            }
            e = std::make_shared<IntersectionEnvelope>(std::move(members));
        }
        return e;
    }
 
    template <typename Point>
    std::shared_ptr<Envelope> pull(const InputSet set, const Point& p) const
    {
        check(set);
        if ((set & (set - 1)) == 0) return m_members[lowest_bit(set)];
        std::shared_ptr<Envelope> best;
        double worst = -1.;
        for (size_t i = 0; i < m_members.size(); ++i) {
            if (!(set & (InputSet(1) << i))) continue;
            const double d2 = m_members[i]->squared_distance(p);
            if (d2 > worst) {
                worst = d2;
                best = m_members[i];
            }
        }
        return best;
    }
 
private:
    void check(const InputSet set) const
    {
        if (set == 0) {
            log_and_throw_error(
                "a tracked simplex has an empty input set: the insertion or an operation did "
                "not carry it");
        }
        if (m_members.size() < size_t(SurfaceTagAttributes::max_inputs) &&
            (set >> m_members.size())) {
            log_and_throw_error(
                "input set {:#x} names an input beyond the {} that have envelopes",
                set,
                m_members.size());
        }
    }
    static size_t lowest_bit(const InputSet set)
    {
        size_t i = 0;
        while (!(set & (InputSet(1) << i))) ++i;
        return i;
    }
 
    std::vector<std::shared_ptr<Envelope>> m_members;
    mutable std::map<InputSet, std::shared_ptr<Envelope>> m_intersections;
    mutable std::shared_mutex m_mutex;
};
 
} // namespace wmtk