Loading [MathJax]/extensions/tex2jax.js
PolyFEM
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
NLProblem.hpp
Go to the documentation of this file.
1#pragma once
2
6
8{
9 class Solver;
10}
11
12namespace polyfem::solver
13{
14 class NLProblem : public FullNLProblem
15 {
16 public:
17 using typename FullNLProblem::Scalar;
18 using typename FullNLProblem::THessian;
19 using typename FullNLProblem::TVector;
20
21 protected:
23 const int full_size,
24 const std::vector<std::shared_ptr<Form>> &forms,
25 const std::vector<std::shared_ptr<AugmentedLagrangianForm>> &penalty_forms,
26 const std::shared_ptr<polysolve::linear::Solver> &solver);
27
28 public:
29 NLProblem(const int full_size,
30 const std::shared_ptr<utils::PeriodicBoundary> &periodic_bc,
31 const double t,
32 const std::vector<std::shared_ptr<Form>> &forms,
33 const std::vector<std::shared_ptr<AugmentedLagrangianForm>> &penalty_forms,
34 const std::shared_ptr<polysolve::linear::Solver> &solver);
35 virtual ~NLProblem() = default;
36
37 virtual double value(const TVector &x) override;
38 virtual void gradient(const TVector &x, TVector &gradv) override;
39 virtual void hessian(const TVector &x, THessian &hessian) override;
40
41 virtual bool is_step_valid(const TVector &x0, const TVector &x1) override;
42 virtual bool is_step_collision_free(const TVector &x0, const TVector &x1) override;
43 virtual double max_step_size(const TVector &x0, const TVector &x1) override;
44 void line_search_begin(const TVector &x0, const TVector &x1) override;
45 virtual void post_step(const polysolve::nonlinear::PostStepData &data) override;
46
47 void solution_changed(const TVector &new_x) override;
48
49 void init_lagging(const TVector &x) override;
50 void update_lagging(const TVector &x, const int iter_num) override;
51
52 // --------------------------------------------------------------------
53
54 virtual void update_quantities(const double t, const TVector &x);
55
56 int full_size() const { return full_size_; }
57 int reduced_size() const { return reduced_size_; }
58
61
62 TVector full_to_reduced(const TVector &full) const;
63 virtual TVector full_to_reduced_grad(const TVector &full) const;
64 TVector reduced_to_full(const TVector &reduced) const;
65
67
68 double normalize_forms() override;
69
70 protected:
71 const int full_size_;
73
74 enum class CurrentSize
75 {
78 };
80 int current_size() const
81 {
83 }
84
85 double t_;
86
87 protected:
88 std::vector<std::shared_ptr<AugmentedLagrangianForm>> penalty_forms_;
89 // The decomposion comes from sec 1.3 of https://www.cs.cornell.edu/courses/cs6241/2021sp/meetings/nb-2021-03-11.pdf
94 Eigen::PermutationMatrix<Eigen::Dynamic, Eigen::Dynamic> P_;
95 TVector Q1R1iTb_;
96 std::shared_ptr<polysolve::linear::Solver> solver_;
97
98 std::shared_ptr<FullNLProblem> penalty_problem_;
100
101 void setup_constraints();
103 };
104} // namespace polyfem::solver
int x
std::vector< std::shared_ptr< Form > > & forms()
const int full_size_
Size of the full problem.
Definition NLProblem.hpp:71
StiffnessMatrix R1_
R1 block of the QR decomposition of the constraints matrix.
Definition NLProblem.hpp:93
void line_search_begin(const TVector &x0, const TVector &x1) override
Eigen::PermutationMatrix< Eigen::Dynamic, Eigen::Dynamic > P_
Permutation matrix of the QR decomposition of the constraints matrix.
Definition NLProblem.hpp:94
double normalize_forms() override
StiffnessMatrix Q2_
Q2 block of the QR decomposition of the constraints matrix.
Definition NLProblem.hpp:91
virtual bool is_step_valid(const TVector &x0, const TVector &x1) override
int reduced_size_
Size of the reduced problem.
Definition NLProblem.hpp:72
virtual void post_step(const polysolve::nonlinear::PostStepData &data) override
virtual TVector full_to_reduced_grad(const TVector &full) const
StiffnessMatrix Q2t_
Q2 transpose.
Definition NLProblem.hpp:92
TVector full_to_reduced(const TVector &full) const
virtual void update_quantities(const double t, const TVector &x)
virtual void gradient(const TVector &x, TVector &gradv) override
void init_lagging(const TVector &x) override
virtual bool is_step_collision_free(const TVector &x0, const TVector &x1) override
TVector reduced_to_full(const TVector &reduced) const
void update_lagging(const TVector &x, const int iter_num) override
std::vector< std::shared_ptr< AugmentedLagrangianForm > > penalty_forms_
Definition NLProblem.hpp:88
virtual double value(const TVector &x) override
StiffnessMatrix Q1_
Q1 block of the QR decomposition of the constraints matrix.
Definition NLProblem.hpp:90
virtual double max_step_size(const TVector &x0, const TVector &x1) override
virtual ~NLProblem()=default
virtual void hessian(const TVector &x, THessian &hessian) override
std::shared_ptr< polysolve::linear::Solver > solver_
Definition NLProblem.hpp:96
void solution_changed(const TVector &new_x) override
std::shared_ptr< FullNLProblem > penalty_problem_
Definition NLProblem.hpp:98
CurrentSize current_size_
Current size of the problem (either full or reduced size)
Definition NLProblem.hpp:79
TVector Q1R1iTb_
Q1_ * (R1_.transpose().triangularView<Eigen::Upper>().solve(constraint_values_))
Definition NLProblem.hpp:95
void full_hessian_to_reduced_hessian(StiffnessMatrix &hessian) const
Eigen::SparseMatrix< double, Eigen::ColMajor > StiffnessMatrix
Definition Types.hpp:22