VariableToSimulation.hpp

Go to the documentation of this file.

#pragma once
 
#include <polyfem/Common.hpp>
#include <polyfem/optimization/parametrization/Parametrization.hpp>
#include <polyfem/optimization/AdjointTools.hpp>
#include <polyfem/optimization/DiffCache.hpp>
#include <polyfem/optimization/parametrization/PeriodicMeshToMesh.hpp>
 
#include <Eigen/Core>
 
#include <cassert>
#include <functional>
#include <memory>
#include <string>
#include <utility>
#include <vector>
 
namespace polyfem::solver
{
    class VariableToSimulation
    {
    public:
        VariableToSimulation(std::vector<std::shared_ptr<State>> states,
                             std::vector<std::shared_ptr<DiffCache>> diff_caches,
                             CompositeParametrization parametrization) : states(std::move(states)),
                                                                         diff_caches(std::move(diff_caches)),
                                                                         parametrization(std::move(parametrization))
        {
            assert(this->states.size() > 0);
        }
 
        virtual ~VariableToSimulation() = default;
 
        inline virtual void update(const Eigen::VectorXd &x)
        {
            update_state(parametrization.eval(x), get_output_indexing(x));
        }
 
        virtual std::string name() const = 0;
 
        virtual ParameterType get_parameter_type() const = 0;
        virtual Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const = 0;
        virtual Eigen::VectorXd inverse_eval();
 
        virtual void set_output_indexing(const json &args);
        Eigen::VectorXi get_output_indexing(const Eigen::VectorXd &x) const;
 
        virtual Eigen::VectorXd apply_parametrization_jacobian(const Eigen::VectorXd &term, const Eigen::VectorXd &x) const;
 
        const std::vector<std::shared_ptr<State>> states;
        const std::vector<std::shared_ptr<DiffCache>> diff_caches;
        CompositeParametrization parametrization;
 
    protected:
        virtual void update_state(const Eigen::VectorXd &state_variable, const Eigen::VectorXi &indices);
 
        Eigen::VectorXi output_indexing_; // if a derived class overrides apply_parametrization_jacobian(term, x), this is not necessarily used.
    };
 
    class VariableToSimulationGroup
    {
    public:
        virtual ~VariableToSimulationGroup() = default;
 
        inline void update(const Eigen::VectorXd &x)
        {
            for (auto &v2s : data)
                v2s->update(x);
        }
 
        void compute_state_variable(const ParameterType type, const State *state_ptr, const Eigen::VectorXd &x, Eigen::VectorXd &state_variable) const;
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const;
 
        virtual Eigen::VectorXd apply_parametrization_jacobian(const ParameterType type, const State *state_ptr, const Eigen::VectorXd &x, const std::function<Eigen::VectorXd()> &grad) const;
 
        std::vector<std::shared_ptr<VariableToSimulation>> data;
    };
 
    // state variable dof = dim * n_vertices
    class ShapeVariableToSimulation : public VariableToSimulation
    {
    public:
        using VariableToSimulation::VariableToSimulation;
        virtual ~ShapeVariableToSimulation() {}
 
        std::string name() const override { return "shape"; }
 
        ParameterType get_parameter_type() const override { return ParameterType::Shape; }
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const override;
        Eigen::VectorXd inverse_eval() override;
 
        void set_output_indexing(const json &args) override;
 
    protected:
        virtual void update_state(const Eigen::VectorXd &state_variable, const Eigen::VectorXi &indices) override;
    };
 
    // To optimize per element elastic parameters
    // state variable dof = 2 * n_elements
    class ElasticVariableToSimulation : public VariableToSimulation
    {
    public:
        using VariableToSimulation::VariableToSimulation;
        virtual ~ElasticVariableToSimulation() {}
 
        std::string name() const override { return "elastic"; }
 
        ParameterType get_parameter_type() const override { return ParameterType::LameParameter; }
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const override;
        Eigen::VectorXd inverse_eval() override;
 
    protected:
        void update_state(const Eigen::VectorXd &state_variable, const Eigen::VectorXi &indices) override;
    };
 
    // To optimize the friction constant
    // state variable dof = 1
    class FrictionCoeffientVariableToSimulation : public VariableToSimulation
    {
    public:
        using VariableToSimulation::VariableToSimulation;
        virtual ~FrictionCoeffientVariableToSimulation() {}
 
        std::string name() const override { return "friction"; }
 
        ParameterType get_parameter_type() const override { return ParameterType::FrictionCoefficient; }
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const override;
        Eigen::VectorXd inverse_eval() override;
 
    protected:
        void update_state(const Eigen::VectorXd &state_variable, const Eigen::VectorXi &indices) override;
    };
 
    // To optimize the damping constant psi and phi
    // state variable dof = 2
    class DampingCoeffientVariableToSimulation : public VariableToSimulation
    {
    public:
        using VariableToSimulation::VariableToSimulation;
        virtual ~DampingCoeffientVariableToSimulation() {}
 
        std::string name() const override { return "damping"; }
 
        ParameterType get_parameter_type() const override { return ParameterType::DampingCoefficient; }
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const override;
        Eigen::VectorXd inverse_eval() override;
 
    protected:
        void update_state(const Eigen::VectorXd &state_variable, const Eigen::VectorXi &indices) override;
    };
 
    // To optimize the per node initial condition
    // state variable dof = dim * n_bases
    class InitialConditionVariableToSimulation : public VariableToSimulation
    {
    public:
        using VariableToSimulation::VariableToSimulation;
        virtual ~InitialConditionVariableToSimulation() {}
 
        std::string name() const override { return "initial"; }
 
        ParameterType get_parameter_type() const override { return ParameterType::InitialCondition; }
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const override;
        Eigen::VectorXd inverse_eval() override;
 
    protected:
        void update_state(const Eigen::VectorXd &state_variable, const Eigen::VectorXi &indices) override;
    };
 
    // To optimize the per node Dirichlet values, only work in transient simulations
    // state variable dof = dim * n_time_steps * n_dirichlet_nodes
    class DirichletVariableToSimulation : public VariableToSimulation
    {
    public:
        using VariableToSimulation::VariableToSimulation;
        virtual ~DirichletVariableToSimulation() {}
 
        std::string name() const override { return "dirichlet"; }
 
        void set_dirichlet_boundaries(const std::vector<int> &dirichlet_boundaries)
        {
            dirichlet_boundaries_ = dirichlet_boundaries;
        }
 
        ParameterType get_parameter_type() const override { return ParameterType::DirichletBC; }
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const override;
        Eigen::VectorXd inverse_eval() override;
 
        void set_output_indexing(const json &args) override;
 
    protected:
        void update_state(const Eigen::VectorXd &state_variable, const Eigen::VectorXi &indices) override;
 
    private:
        std::string variable_to_string(const Eigen::VectorXd &variable);
 
        std::vector<int> dirichlet_boundaries_;
    };
 
    class DirichletNodesVariableToSimulation : public VariableToSimulation
    {
    public:
        using VariableToSimulation::VariableToSimulation;
        virtual ~DirichletNodesVariableToSimulation() {}
 
        std::string name() const override { return "dirichlet-nodes"; }
 
        void set_dirichlet_nodes(const Eigen::VectorXi &dirichlet_nodes)
        {
            dirichlet_nodes_ = dirichlet_nodes;
        }
 
        ParameterType get_parameter_type() const override { return ParameterType::DirichletBC; }
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const override;
        Eigen::VectorXd inverse_eval() override;
 
        void set_output_indexing(const json &args) override;
 
    protected:
        void update_state(const Eigen::VectorXd &state_variable, const Eigen::VectorXi &indices) override;
 
    private:
        std::string variable_to_string(const Eigen::VectorXd &variable);
 
        Eigen::VectorXi dirichlet_nodes_;
    };
 
    // To optimize the per node pressure boundaries
    // Each pressure boundary will have the same value
    // state variable dof = dim * n_time_steps
    class PressureVariableToSimulation : public VariableToSimulation
    {
    public:
        using VariableToSimulation::VariableToSimulation;
        virtual ~PressureVariableToSimulation() {}
 
        std::string name() const override { return "pressure"; }
 
        void set_pressure_boundaries(const std::vector<int> &pressure_boundaries)
        {
            pressure_boundaries_ = pressure_boundaries;
        }
 
        ParameterType get_parameter_type() const override { return ParameterType::PressureBC; }
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const override;
        Eigen::VectorXd inverse_eval() override;
 
        void set_output_indexing(const json &args) override;
 
    protected:
        void update_state(const Eigen::VectorXd &state_variable, const Eigen::VectorXi &indices) override;
 
    private:
        std::string variable_to_string(const Eigen::VectorXd &variable);
 
        std::vector<int> pressure_boundaries_;
    };
 
    // state variable dof = dim * n_vertices - periodic dof
    class PeriodicShapeVariableToSimulation : public VariableToSimulation
    {
    public:
        using VariableToSimulation::VariableToSimulation;
        virtual ~PeriodicShapeVariableToSimulation() {}
 
        std::string name() const override { return "periodic-shape"; }
 
        ParameterType get_parameter_type() const override { return ParameterType::PeriodicShape; }
 
        Eigen::VectorXd compute_adjoint_term(const Eigen::VectorXd &x) const override;
        Eigen::VectorXd inverse_eval() override;
 
        void update(const Eigen::VectorXd &x) override;
 
        Eigen::VectorXd apply_parametrization_jacobian(const Eigen::VectorXd &term, const Eigen::VectorXd &x) const override;
 
    protected:
        std::unique_ptr<PeriodicMeshToMesh> periodic_mesh_map;
        Eigen::VectorXd periodic_mesh_representation;
    };
} // namespace polyfem::solver