PhysicsOperator.cpp

Go to the documentation of this file.

#include "PhysicsOperator.hpp"
 
#include "MayaFlux/Nodes/Network/NodeNetwork.hpp"
 
#include "MayaFlux/Journal/Archivist.hpp"
 
#include "MayaFlux/Kinesis/Tendency/ForceFields.hpp"
#include "MayaFlux/Transitive/Reflect/EnumReflect.hpp"
 
namespace MayaFlux::Nodes::Network {
 
PhysicsOperator::PhysicsOperator()
{
    MF_DEBUG(Journal::Component::Nodes, Journal::Context::NodeProcessing,
        "PhysicsOperator created");
}
 
//-----------------------------------------------------------------------------
// GraphicsOperator Interface (Basic Initialization)
//-----------------------------------------------------------------------------
 
void PhysicsOperator::initialize(const std::vector<PointVertex>& vertices)
{
    if (vertices.empty()) {
        MF_WARN(Journal::Component::Nodes, Journal::Context::NodeProcessing,
            "Cannot initialize PhysicsOperator with zero vertices");
        return;
    }
 
    m_collections.clear();
    add_collection(vertices);
 
    MF_DEBUG(Journal::Component::Nodes, Journal::Context::NodeProcessing,
        "PhysicsOperator initialized with {} points", vertices.size());
}
 
//-----------------------------------------------------------------------------
// Advanced Initialization (Multiple Collections)
//-----------------------------------------------------------------------------
 
void PhysicsOperator::initialize_collections(
    const std::vector<std::vector<PointVertex>>& collections)
{
    for (const auto& collection : collections) {
        add_collection(collection);
    }
 
    MF_DEBUG(Journal::Component::Nodes, Journal::Context::NodeProcessing,
        "PhysicsOperator initialized with {} collections",
        collections.size());
}
 
void PhysicsOperator::add_collection(
    const std::vector<PointVertex>& vertices,
    float mass_multiplier)
{
    if (vertices.empty()) {
        MF_WARN(Journal::Component::Nodes, Journal::Context::NodeProcessing,
            "Cannot add collection with zero vertices");
        return;
    }
 
    CollectionGroup group;
    group.collection = std::make_shared<GpuSync::PointCollectionNode>();
 
    group.physics_state.resize(vertices.size());
 
    for (size_t i = 0; i < vertices.size(); ++i) {
        group.physics_state[i] = PhysicsState {
            .velocity = glm::vec3(0.0F),
            .force = glm::vec3(0.0F),
            .mass = mass_multiplier
        };
    }
 
    group.collection->set_points(vertices);
    group.collection->compute_frame();
 
    uint32_t expected = 0;
    while (!m_access_token.compare_exchange_weak(expected, 1,
        std::memory_order_acquire, std::memory_order_relaxed)) {
        if (m_shutdown.load(std::memory_order_relaxed))
            return;
        expected = 0;
    }
 
    m_collections.push_back(std::move(group));
 
    m_access_token.store(0, std::memory_order_release);
 
    MF_DEBUG(Journal::Component::Nodes, Journal::Context::NodeProcessing,
        "Added collection #{} with {} points (mass_mult={:.2f})",
        m_collections.size(), vertices.size(), mass_multiplier);
}
 
void PhysicsOperator::seed_from_upstream(const GraphicsOperator* upstream)
{
    if (get_vertex_count() > 0 || !upstream)
        return;
 
    const auto data = upstream->get_vertex_data();
    if (data.empty())
        return;
 
    const size_t count = data.size() / sizeof(PointVertex);
    if (count == 0)
        return;
 
    std::vector<PointVertex> vertices(count);
    std::memcpy(vertices.data(), data.data(), data.size());
    initialize(vertices);
 
    MF_DEBUG(Journal::Component::Nodes, Journal::Context::NodeProcessing,
        "PhysicsOperator seeded {} vertices from upstream '{}'",
        count, upstream->get_type_name());
}
 
//-----------------------------------------------------------------------------
// Processing
//-----------------------------------------------------------------------------
 
void PhysicsOperator::process(float dt)
{
    uint32_t expected = 0;
    while (!m_access_token.compare_exchange_weak(expected, 1,
        std::memory_order_acquire, std::memory_order_relaxed)) {
        if (m_shutdown.load(std::memory_order_relaxed))
            return;
        expected = 0;
    }
 
    if (m_collections.empty()) {
        m_access_token.store(0, std::memory_order_release);
        return;
    }
 
    const float effective_dt = m_force_internal_dt ? m_internal_dt : dt;
    apply_forces();
    integrate(effective_dt);
    handle_boundary_conditions();
    sync_to_point_collection();
 
    for (auto& group : m_collections) {
        group.collection->compute_frame();
    }
 
    m_access_token.store(0, std::memory_order_release);
}
 
//-----------------------------------------------------------------------------
// Parameter System
//-----------------------------------------------------------------------------
 
std::optional<PhysicsParameter> PhysicsOperator::string_to_parameter(std::string_view param)
{
    return Reflect::string_to_enum_case_insensitive<PhysicsParameter>(param);
}
 
void PhysicsOperator::set_parameter(std::string_view param, double value)
{
    auto param_enum = string_to_parameter(param);
 
    if (!param_enum) {
        try {
            Reflect::string_to_enum_or_throw_case_insensitive<PhysicsParameter>(
                param, "PhysicsOperator parameter");
        } catch (const std::invalid_argument& e) {
            MF_ERROR(Journal::Component::Nodes, Journal::Context::NodeProcessing,
                "{}", e.what());
        }
        MF_WARN(Journal::Component::Nodes, Journal::Context::NodeProcessing,
            "Unknown physics parameter: '{}'", param);
        return;
    }
 
    switch (*param_enum) {
    case PhysicsParameter::GRAVITY_X:
        m_gravity.x = static_cast<float>(value);
        break;
    case PhysicsParameter::GRAVITY_Y:
        m_gravity.y = static_cast<float>(value);
        break;
    case PhysicsParameter::GRAVITY_Z:
        m_gravity.z = static_cast<float>(value);
        break;
    case PhysicsParameter::DRAG:
        m_drag = glm::clamp(static_cast<float>(value), 0.0F, 1.0F);
        break;
    case PhysicsParameter::INTERACTION_RADIUS:
        m_interaction_radius = static_cast<float>(value);
        break;
    case PhysicsParameter::SPRING_STIFFNESS:
        m_spring_stiffness = static_cast<float>(value);
        break;
    case PhysicsParameter::REPULSION_STRENGTH:
        m_repulsion_strength = static_cast<float>(value);
        break;
    case PhysicsParameter::SPATIAL_INTERACTIONS:
        m_spatial_interactions_enabled = (value > 0.5);
        break;
    case PhysicsParameter::POINT_SIZE:
        m_point_size = static_cast<float>(value);
        for (auto& group : m_collections) {
            auto& points = group.collection->get_points();
            for (auto& pt : points) {
                pt.size = m_point_size;
            }
        }
        break;
    case PhysicsParameter::ATTRACTION_STRENGTH:
        m_attraction_strength = static_cast<float>(value);
        break;
    case PhysicsParameter::TURBULENCE:
        m_turbulence_strength = static_cast<float>(value);
        break;
    }
}
 
std::optional<double> PhysicsOperator::query_state(std::string_view query) const
{
    if (query == "point_count") {
        return static_cast<double>(get_point_count());
    }
    if (query == "collection_count") {
        return static_cast<double>(m_collections.size());
    }
    if (query == "avg_velocity") {
        glm::vec3 avg(0.0F);
        size_t total_points = 0;
 
        for (const auto& group : m_collections) {
            for (const auto& state : group.physics_state) {
                avg += state.velocity;
                ++total_points;
            }
        }
 
        if (total_points > 0) {
            avg /= static_cast<float>(total_points);
        }
        return static_cast<double>(glm::length(avg));
    }
 
    return std::nullopt;
}
 
//-----------------------------------------------------------------------------
// GraphicsOperator Interface (Data Extraction)
//-----------------------------------------------------------------------------
 
std::vector<PointVertex> PhysicsOperator::extract_vertices() const
{
    std::vector<PointVertex> positions;
 
    for (const auto& group : m_collections) {
        const auto& points = group.collection->get_points();
        for (const auto& pt : points) {
            positions.push_back(pt);
        }
    }
 
    return positions;
}
 
std::span<const uint8_t> PhysicsOperator::get_vertex_data_for_collection(uint32_t idx) const
{
    if (m_collections.empty() || idx >= m_collections.size()) {
        return {};
    }
    return m_collections[idx].collection->get_vertex_data();
}
 
std::span<const uint8_t> PhysicsOperator::get_vertex_data() const
{
    m_vertex_data_aggregate.clear();
    for (const auto& group : m_collections) {
        auto span = group.collection->get_vertex_data();
        m_vertex_data_aggregate.insert(
            m_vertex_data_aggregate.end(),
            span.begin(), span.end());
    }
    return { m_vertex_data_aggregate.data(), m_vertex_data_aggregate.size() };
}
 
std::optional<double> PhysicsOperator::get_particle_velocity(size_t global_index) const
{
    size_t current_offset = 0;
 
    for (const auto& group : m_collections) {
        size_t group_size = group.physics_state.size();
 
        if (global_index < current_offset + group_size) {
            size_t local_index = global_index - current_offset;
            return static_cast<double>(glm::length(group.physics_state[local_index].velocity));
        }
 
        current_offset += group_size;
    }
 
    return std::nullopt;
}
 
Kakshya::VertexLayout PhysicsOperator::get_vertex_layout() const
{
    if (m_collections.empty()) {
        return {};
    }
 
    auto layout_opt = m_collections[0].collection->get_vertex_layout();
    if (!layout_opt.has_value()) {
        return {};
    }
 
    auto layout = *layout_opt;
    layout.vertex_count = static_cast<uint32_t>(get_vertex_count());
    return layout;
}
 
size_t PhysicsOperator::get_vertex_count() const
{
    size_t total = 0;
    for (const auto& group : m_collections) {
        total += group.collection->get_vertex_count();
    }
    return total;
}
 
bool PhysicsOperator::is_vertex_data_dirty() const
{
    return std::ranges::any_of(
        m_collections,
        [](const auto& group) { return group.collection->needs_gpu_update(); });
}
 
void PhysicsOperator::mark_vertex_data_clean()
{
    for (auto& group : m_collections) {
        group.collection->mark_vertex_data_dirty(false);
    }
}
 
size_t PhysicsOperator::get_point_count() const
{
    size_t total = 0;
    for (const auto& group : m_collections) {
        total += group.collection->get_point_count();
    }
    return total;
}
 
void PhysicsOperator::set_bounds(const glm::vec3& min, const glm::vec3& max)
{
    m_bounds.min = min;
    m_bounds.max = max;
}
 
void PhysicsOperator::set_attraction_point(const glm::vec3& point)
{
    m_attraction_point = point;
    m_has_attraction_point = true;
}
 
//-----------------------------------------------------------------------------
// ONE_TO_ONE Parameter Mapping
//-----------------------------------------------------------------------------
 
void PhysicsOperator::apply_one_to_one(
    std::string_view param,
    const std::shared_ptr<NodeNetwork>& source)
{
    size_t point_count = get_point_count();
 
    if (source->get_node_count() != point_count) {
        MF_WARN(Journal::Component::Nodes, Journal::Context::NodeProcessing,
            "ONE_TO_ONE size mismatch: {} particles vs {} source nodes",
            point_count, source->get_node_count());
        return;
    }
 
    if (param == "force_x" || param == "force_y" || param == "force_z") {
        apply_per_particle_force(param, source);
    } else if (param == "mass") {
        apply_per_particle_mass(source);
    } else {
        GraphicsOperator::apply_one_to_one(param, source);
    }
}
 
void PhysicsOperator::apply_per_particle_force(
    std::string_view param,
    const std::shared_ptr<NodeNetwork>& source)
{
    size_t global_index = 0;
 
    for (auto& group : m_collections) {
        for (auto& i : group.physics_state) {
            auto val = source->get_node_output(global_index++);
            if (!val)
                continue;
 
            auto force = static_cast<float>(*val);
 
            if (param == "force_x") {
                i.force.x += force;
            } else if (param == "force_y") {
                i.force.y += force;
            } else if (param == "force_z") {
                i.force.z += force;
            }
        }
    }
}
 
void PhysicsOperator::apply_per_particle_mass(
    const std::shared_ptr<NodeNetwork>& source)
{
    size_t global_index = 0;
 
    for (auto& group : m_collections) {
        for (auto& i : group.physics_state) {
            auto val = source->get_node_output(global_index++);
            if (!val)
                continue;
 
            i.mass = std::max(0.1F, static_cast<float>(*val));
        }
    }
}
 
//-----------------------------------------------------------------------------
// Physics Simulation
//-----------------------------------------------------------------------------
 
void PhysicsOperator::apply_forces()
{
    for (auto& group : m_collections) {
        for (auto& state : group.physics_state) {
            state.force = m_gravity * state.mass;
        }
    }
 
    if (m_has_attraction_point) {
        apply_attraction_forces();
    }
 
    if (m_turbulence_strength > 0.001F) {
        apply_turbulence();
    }
 
    if (m_spatial_interactions_enabled && m_interaction_radius > 0.0F) {
        apply_spatial_interactions();
    }
 
    if (!m_bond_root.empty()) {
        apply_bond_forces();
    }
 
    if (!m_force_fields.empty()) {
        for (auto& group : m_collections) {
            auto& points = group.collection->get_points();
 
            for (size_t i = 0; i < points.size(); ++i) {
                for (const auto& field : m_force_fields) {
                    group.physics_state[i].force += field(points[i].position);
                }
            }
        }
    }
}
 
void PhysicsOperator::apply_turbulence()
{
    auto field = Kinesis::ForceFields::turbulence(m_turbulence_strength, m_random_generator);
 
    for (auto& group : m_collections) {
        for (auto& state : group.physics_state) {
            state.force += field(glm::vec3(0.0F));
        }
    }
}
 
void PhysicsOperator::apply_spatial_interactions()
{
    for (size_t g1 = 0; g1 < m_collections.size(); ++g1) {
        auto& group1 = m_collections[g1];
        auto& points1 = group1.collection->get_points();
 
        for (size_t i = 0; i < points1.size(); ++i) {
            const auto& pos_i = points1[i].position;
            auto& state_i = group1.physics_state[i];
 
            for (size_t g2 = 0; g2 < m_collections.size(); ++g2) {
                auto& group2 = m_collections[g2];
                auto& points2 = group2.collection->get_points();
 
                size_t start_j = (g1 == g2) ? i + 1 : 0;
 
                for (size_t j = start_j; j < points2.size(); ++j) {
                    const auto& pos_j = points2[j].position;
                    auto& state_j = group2.physics_state[j];
 
                    glm::vec3 delta = pos_j - pos_i;
                    float distance = glm::length(delta);
 
                    if (distance < m_interaction_radius && distance > 0.001F) {
                        glm::vec3 direction = delta / distance;
 
                        float spring_force = m_spring_stiffness * (distance - m_interaction_radius * 0.5F);
 
                        float repulsion_force = 0.0F;
                        if (distance < m_interaction_radius * 0.3F) {
                            repulsion_force = m_repulsion_strength / (distance * distance);
                        }
 
                        glm::vec3 force = direction * (spring_force - repulsion_force);
 
                        state_i.force += force;
                        state_j.force -= force;
                    }
                }
            }
        }
    }
}
 
std::optional<PhysicsOperator::GroupIndex> PhysicsOperator::resolve_global_index(size_t global_index) const
{
    size_t offset = 0;
    for (size_t g = 0; g < m_collections.size(); ++g) {
        size_t count = m_collections[g].collection->get_point_count();
        if (global_index < offset + count) {
            return GroupIndex { .group = g, .local = global_index - offset };
        }
        offset += count;
    }
    return std::nullopt;
}
 
/**
 * @brief Pull each bonded particle toward its root with a spring.
 *
 * Rest length scales with cbrt(root's accreted mass): a grown body's
 * satellites settle at a wider spread than a fresh one's, the same way any
 * fixed-density accumulation of matter occupies more volume as it gains
 * mass. Growth is a wider swarm of small particles, not a bigger point.
 */
void PhysicsOperator::apply_bond_forces()
{
    for (size_t i = 0; i < m_bond_root.size(); ++i) {
        uint32_t root = m_bond_root[i];
        if (root == i) {
            continue;
        }
 
        auto self_idx = resolve_global_index(i);
        auto root_idx = resolve_global_index(root);
        if (!self_idx || !root_idx) {
            continue;
        }
 
        auto& self_state = m_collections[self_idx->group].physics_state[self_idx->local];
        const auto& self_pos = m_collections[self_idx->group].collection->get_points()[self_idx->local].position;
        const auto& root_pos = m_collections[root_idx->group].collection->get_points()[root_idx->local].position;
 
        glm::vec3 delta = root_pos - self_pos;
        float distance = glm::length(delta);
        if (distance < 0.001F) {
            continue;
        }
 
        const float rest_length = m_bond_rest_length * std::cbrt(get_accreted_mass(root));
 
        glm::vec3 direction = delta / distance;
        self_state.force += direction * (m_bond_stiffness * (distance - rest_length));
    }
}
 
void PhysicsOperator::sync_bonds_from_claims(std::span<const uint32_t> claimed_by)
{
    if (!m_bonds_enabled) {
        return;
    }
 
    if (m_accreted_mass.size() != claimed_by.size()) {
        m_accreted_mass.assign(claimed_by.size(), 1.0F);
    }
 
    for (size_t i = 0; i < claimed_by.size(); ++i) {
        uint32_t root = claimed_by[i];
        if (root == i || root >= m_accreted_mass.size()) {
            continue;
        }
        m_accreted_mass[root] += m_accreted_mass[i];
        m_accreted_mass[i] = 0.0F;
    }
 
    m_bond_root.assign(claimed_by.begin(), claimed_by.end());
}
 
void PhysicsOperator::clear_bonds()
{
    m_bond_root.clear();
}
 
void PhysicsOperator::enable_bonds(bool enable)
{
    m_bonds_enabled = enable;
    if (!enable) {
        clear_bonds();
    }
}
 
float PhysicsOperator::get_accreted_mass(size_t global_index) const
{
    if (global_index >= m_accreted_mass.size()) {
        return 1.0F;
    }
    return m_accreted_mass[global_index];
}
 
std::span<const float> PhysicsOperator::get_accreted_mass_span()
{
    const size_t count = get_point_count();
    if (count > 0 && m_accreted_mass.size() != count) {
        m_accreted_mass.assign(count, 1.0F);
    }
    return m_accreted_mass;
}
 
float PhysicsOperator::get_total_mass() const
{
    if (m_accreted_mass.empty()) {
        float count = 0.0F;
        for (const auto& group : m_collections) {
            count += static_cast<float>(group.collection->get_point_count());
        }
        return count;
    }
 
    float total = 0.0F;
    for (float mass : m_accreted_mass) {
        total += mass;
    }
    return total;
}
 
size_t PhysicsOperator::bond_count() const
{
    size_t count = 0;
    for (size_t i = 0; i < m_bond_root.size(); ++i) {
        if (m_bond_root[i] != i) {
            ++count;
        }
    }
    return count;
}
 
void PhysicsOperator::apply_attraction_forces()
{
    auto field = Kinesis::ForceFields::point_attractor(m_attraction_point, m_attraction_strength);
 
    for (auto& group : m_collections) {
        auto& points = group.collection->get_points();
 
        for (size_t i = 0; i < points.size(); ++i) {
            group.physics_state[i].force += field(points[i].position) * group.physics_state[i].mass;
        }
    }
}
 
void PhysicsOperator::integrate(float dt)
{
    for (auto& group : m_collections) {
        auto& points = group.collection->get_points();
 
        for (size_t i = 0; i < points.size(); ++i) {
            auto& state = group.physics_state[i];
            auto& vertex = points[i];
 
            glm::vec3 acceleration = state.force / state.mass;
            state.velocity += acceleration * dt;
            state.velocity *= (1.0F - m_drag);
            vertex.position += state.velocity * dt;
 
            state.force = glm::vec3(0.0F);
        }
    }
}
 
void PhysicsOperator::handle_boundary_conditions()
{
    if (m_bounds_mode == BoundsMode::NONE) {
        return;
    }
 
    constexpr float damping = 0.8F;
 
    for (auto& group : m_collections) {
        auto& points = group.collection->get_points();
 
        for (size_t i = 0; i < points.size(); ++i) {
            auto& vertex = points[i];
            auto& state = group.physics_state[i];
 
            for (int axis = 0; axis < 3; ++axis) {
                if (vertex.position[axis] < m_bounds.min[axis]) {
                    switch (m_bounds_mode) {
                    case BoundsMode::BOUNCE:
                        vertex.position[axis] = m_bounds.min[axis];
                        state.velocity[axis] *= -damping;
                        break;
                    case BoundsMode::WRAP:
                        vertex.position[axis] = m_bounds.max[axis];
                        break;
                    case BoundsMode::CLAMP:
                        vertex.position[axis] = m_bounds.min[axis];
                        state.velocity[axis] = 0.0F;
                        break;
                    case BoundsMode::NONE:
                        break;
                    }
                } else if (vertex.position[axis] > m_bounds.max[axis]) {
                    switch (m_bounds_mode) {
                    case BoundsMode::BOUNCE:
                        vertex.position[axis] = m_bounds.max[axis];
                        state.velocity[axis] *= -damping;
                        break;
                    case BoundsMode::WRAP:
                        vertex.position[axis] = m_bounds.min[axis];
                        break;
                    case BoundsMode::CLAMP:
                        vertex.position[axis] = m_bounds.max[axis];
                        state.velocity[axis] = 0.0F;
                        break;
                    case BoundsMode::NONE:
                        break;
                    }
                }
            }
        }
    }
}
 
void PhysicsOperator::sync_to_point_collection()
{
    for (auto& group : m_collections) {
        group.collection->mark_vertex_data_dirty(true);
    }
}
 
void* PhysicsOperator::get_data_at(size_t global_index)
{
    size_t offset = 0;
    for (auto& group : m_collections) {
        if (global_index < offset + group.collection->get_point_count()) {
            size_t local_index = global_index - offset;
            return &group.collection->get_points()[local_index];
        }
        offset += group.collection->get_point_count();
    }
    return nullptr;
}
 
void PhysicsOperator::apply_global_impulse(const glm::vec3& impulse)
{
    for (auto& group : m_collections) {
        for (auto& state : group.physics_state) {
            state.velocity += impulse / state.mass;
        }
    }
}
 
void PhysicsOperator::apply_impulse(size_t index, const glm::vec3& impulse)
{
    size_t offset = 0;
    for (auto& group : m_collections) {
        if (index < offset + group.collection->get_point_count()) {
            size_t local_index = index - offset;
            group.physics_state[local_index].velocity += impulse / group.physics_state[local_index].mass;
            return;
        }
        offset += group.collection->get_point_count();
    }
}
 
std::vector<uint32_t> PhysicsOperator::build_cluster_ids() const
{
    std::vector<uint32_t> ids(get_vertex_count(), 0U);
 
    if (m_collections.size() <= 1) {
        return ids;
    }
 
    size_t offset = 0;
    uint32_t cluster = 0;
    for (const auto& group : m_collections) {
        const size_t count = group.collection->get_vertex_count();
        for (size_t i = 0; i < count && offset + i < ids.size(); ++i) {
            ids[offset + i] = cluster;
        }
        offset += count;
        ++cluster;
    }
 
    return ids;
}
 
void PhysicsOperator::add_force_field(Kinesis::VectorField field)
{
    m_force_fields.push_back(std::move(field));
 
    MF_DEBUG(Journal::Component::Nodes, Journal::Context::NodeProcessing,
        "Added force field #{}", m_force_fields.size());
}
 
void PhysicsOperator::clear_force_fields()
{
    m_force_fields.clear();
 
    MF_DEBUG(Journal::Component::Nodes, Journal::Context::NodeProcessing,
        "Cleared all force fields");
}
 
} // namespace MayaFlux::Nodes::Network