Loading [MathJax]/jax/output/HTML-CSS/config.js
PolyFEM
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
BCLagrangianForm.hpp
Go to the documentation of this file.
1#pragma once
2
5
8
9namespace polyfem::solver
10{
13 {
14 public:
26 BCLagrangianForm(const int ndof,
27 const std::vector<int> &boundary_nodes,
28 const std::vector<mesh::LocalBoundary> &local_boundary,
29 const std::vector<mesh::LocalBoundary> &local_neumann_boundary,
30 const int n_boundary_samples,
31 const StiffnessMatrix &mass,
32 const assembler::RhsAssembler &rhs_assembler,
33 const size_t obstacle_ndof,
34 const bool is_time_dependent,
35 const double t);
36
37 std::string name() const override
38 {
39 return "bc-alagrangian";
40 }
41
48 BCLagrangianForm(const int ndof,
49 const std::vector<int> &boundary_nodes,
50 const StiffnessMatrix &mass,
51 const size_t obstacle_ndof,
52 const Eigen::MatrixXd &target_x);
53
57 double value_unweighted(const Eigen::VectorXd &x) const override;
58
62 void first_derivative_unweighted(const Eigen::VectorXd &x, Eigen::VectorXd &gradv) const override;
63
67 void second_derivative_unweighted(const Eigen::VectorXd &x, StiffnessMatrix &hessian) const override;
68
69 public:
73 void update_quantities(const double t, const Eigen::VectorXd &x) override;
74
75 void update_lagrangian(const Eigen::VectorXd &x, const double k_al) override;
76
77 double compute_error(const Eigen::VectorXd &x) const override;
78
79 virtual bool can_project() const override;
80 virtual void project_gradient(Eigen::VectorXd &grad) const override;
81 virtual void project_hessian(StiffnessMatrix &hessian) const override;
82
83 private:
84 const std::vector<int> &boundary_nodes_;
85 const std::vector<mesh::LocalBoundary> *local_boundary_;
86 const std::vector<mesh::LocalBoundary> *local_neumann_boundary_;
88 const int n_dofs_;
89 Eigen::VectorXi constraints_;
90 Eigen::VectorXi not_constraints_;
91
94
97
102 const StiffnessMatrix &mass,
103 const size_t obstacle_ndof);
104
107 void update_target(const double t);
108 };
109} // namespace polyfem::solver
int x
Form of the augmented lagrangian for bc constraints.
Eigen::VectorXi not_constraints_
Not Constraints.
void first_derivative_unweighted(const Eigen::VectorXd &x, Eigen::VectorXd &gradv) const override
Compute the first derivative of the value wrt x.
std::string name() const override
void update_quantities(const double t, const Eigen::VectorXd &x) override
Update time dependent quantities.
void update_target(const double t)
Update target x to the Dirichlet boundary values at time t.
virtual bool can_project() const override
void update_lagrangian(const Eigen::VectorXd &x, const double k_al) override
double compute_error(const Eigen::VectorXd &x) const override
const assembler::RhsAssembler * rhs_assembler_
Reference to the RHS assembler.
void second_derivative_unweighted(const Eigen::VectorXd &x, StiffnessMatrix &hessian) const override
Compute the second derivative of the value wrt x.
virtual void project_gradient(Eigen::VectorXd &grad) const override
void init_masked_lumped_mass(const StiffnessMatrix &mass, const size_t obstacle_ndof)
Initialize the masked lumped mass matrix.
const std::vector< mesh::LocalBoundary > * local_neumann_boundary_
const std::vector< mesh::LocalBoundary > * local_boundary_
StiffnessMatrix masked_lumped_mass_sqrt_
sqrt mass matrix masked by the AL dofs
double value_unweighted(const Eigen::VectorXd &x) const override
Compute the value of the form.
virtual void project_hessian(StiffnessMatrix &hessian) const override
StiffnessMatrix masked_lumped_mass_
mass matrix masked by the AL dofs
const std::vector< int > & boundary_nodes_
Eigen::VectorXi constraints_
Constraints.
Eigen::SparseMatrix< double, Eigen::ColMajor > StiffnessMatrix
Definition Types.hpp:22