ContactForm.hpp

Go to the documentation of this file.

#pragma once
 
#include "Form.hpp"
 
#include <polyfem/Common.hpp>
#include <polyfem/utils/Types.hpp>
 
#include <ipc/collisions/collisions.hpp>
#include <ipc/collision_mesh.hpp>
#include <ipc/broad_phase/broad_phase.hpp>
#include <ipc/potentials/barrier_potential.hpp>
 
// map BroadPhaseMethod values to JSON as strings
namespace ipc
{
    NLOHMANN_JSON_SERIALIZE_ENUM(
        ipc::BroadPhaseMethod,
        {{ipc::BroadPhaseMethod::HASH_GRID, "hash_grid"}, // also default
         {ipc::BroadPhaseMethod::HASH_GRID, "HG"},
         {ipc::BroadPhaseMethod::BRUTE_FORCE, "brute_force"},
         {ipc::BroadPhaseMethod::BRUTE_FORCE, "BF"},
         {ipc::BroadPhaseMethod::SPATIAL_HASH, "spatial_hash"},
         {ipc::BroadPhaseMethod::SPATIAL_HASH, "SH"},
         {ipc::BroadPhaseMethod::BVH, "bvh"},
         {ipc::BroadPhaseMethod::BVH, "BVH"},
         {ipc::BroadPhaseMethod::SWEEP_AND_PRUNE, "sweep_and_prune"},
         {ipc::BroadPhaseMethod::SWEEP_AND_PRUNE, "SAP"},
         {ipc::BroadPhaseMethod::SWEEP_AND_TINIEST_QUEUE, "sweep_and_tiniest_queue"},
         {ipc::BroadPhaseMethod::SWEEP_AND_TINIEST_QUEUE, "STQ"}})
} // namespace ipc
 
namespace polyfem::solver
{
    class ContactForm : public Form
    {
    public:
        ContactForm(const ipc::CollisionMesh &collision_mesh,
                    const double dhat,
                    const double avg_mass,
                    const bool use_convergent_formulation,
                    const bool use_adaptive_barrier_stiffness,
                    const bool is_time_dependent,
                    const bool enable_shape_derivatives,
                    const ipc::BroadPhaseMethod broad_phase_method,
                    const double ccd_tolerance,
                    const int ccd_max_iterations);
        virtual ~ContactForm() = default;
 
        std::string name() const override { return "contact"; }
 
        void init(const Eigen::VectorXd &x) override;
 
        virtual void force_shape_derivative(const ipc::Collisions &collision_set, const Eigen::MatrixXd &solution, const Eigen::VectorXd &adjoint_sol, Eigen::VectorXd &term);
 
    protected:
        virtual double value_unweighted(const Eigen::VectorXd &x) const override;
 
        Eigen::VectorXd value_per_element_unweighted(const Eigen::VectorXd &x) const override;
 
        virtual void first_derivative_unweighted(const Eigen::VectorXd &x, Eigen::VectorXd &gradv) const override;
 
        virtual void second_derivative_unweighted(const Eigen::VectorXd &x, StiffnessMatrix &hessian) const override;
 
    public:
        void update_quantities(const double t, const Eigen::VectorXd &x) override;
 
        double max_step_size(const Eigen::VectorXd &x0, const Eigen::VectorXd &x1) const override;
 
        void line_search_begin(const Eigen::VectorXd &x0, const Eigen::VectorXd &x1) override;
 
        void line_search_end() override;
 
        void solution_changed(const Eigen::VectorXd &new_x) override;
 
        void post_step(const polysolve::nonlinear::PostStepData &data) override;
 
        bool is_step_collision_free(const Eigen::VectorXd &x0, const Eigen::VectorXd &x1) const override;
 
        virtual void update_barrier_stiffness(const Eigen::VectorXd &x, const Eigen::MatrixXd &grad_energy);
 
        Eigen::MatrixXd compute_displaced_surface(const Eigen::VectorXd &x) const;
 
        double barrier_stiffness() const { return barrier_stiffness_; }
        void set_barrier_stiffness(const double barrier_stiffness) { barrier_stiffness_ = barrier_stiffness; }
        bool use_adaptive_barrier_stiffness() const { return use_adaptive_barrier_stiffness_; }
        bool use_convergent_formulation() const { return collision_set_.use_convergent_formulation(); }
 
        bool enable_shape_derivatives() const { return enable_shape_derivatives_; }
 
        double weight() const override { return weight_ * barrier_stiffness_; }
 
        bool save_ccd_debug_meshes = false;
 
        double dhat() const { return dhat_; }
        const ipc::Collisions &collision_set() const { return collision_set_; }
        const ipc::BarrierPotential &barrier_potential() const { return barrier_potential_; }
 
    protected:
        void update_collision_set(const Eigen::MatrixXd &displaced_surface);
 
        const ipc::CollisionMesh &collision_mesh_;
 
        const double dhat_;
 
        const double dmin_ = 0;
 
        const bool use_adaptive_barrier_stiffness_;
        double barrier_stiffness_;
        double max_barrier_stiffness_;
 
        const double avg_mass_;
 
        const bool is_time_dependent_;
 
        const bool enable_shape_derivatives_;
 
        const ipc::BroadPhaseMethod broad_phase_method_;
        const double ccd_tolerance_;
        const int ccd_max_iterations_;
 
        double prev_distance_;
 
        bool use_cached_candidates_ = false;
        ipc::Collisions collision_set_;
        ipc::Candidates candidates_;
 
        const ipc::BarrierPotential barrier_potential_;
    };
} // namespace polyfem::solver