autogenu-jupyter
An automatic code generator and the continuation/GMRES (C/GMRES) based numerical solvers for nonlinear MPC
Loading...
Searching...
No Matches
newton_gmres.hpp
Go to the documentation of this file.
1#ifndef CGMRES__NEWTON_GMRES_HPP_
2#define CGMRES__NEWTON_GMRES_HPP_
3
4#include <stdexcept>
5
6#include "cgmres/types.hpp"
8
9
10namespace cgmres {
11namespace detail {
12
13template <class NLP>
15public:
16 static constexpr int nx = NLP::nx;
17 static constexpr int dim = NLP::dim;
18
19 NewtonGMRES(const NLP& nlp, const Scalar finite_difference_epsilon)
20 : nlp_(nlp),
21 finite_difference_epsilon_(finite_difference_epsilon),
22 updated_solution_(Vector<dim>::Zero()),
23 fonc_(Vector<dim>::Zero()),
24 fonc_1_(Vector<dim>::Zero()) {
25 if (finite_difference_epsilon <= 0.0) {
26 throw std::invalid_argument("[NewtonGMRES]: 'finite_difference_epsilon' must be positive!");
27 }
28 }
29
30 NewtonGMRES() = default;
31
32 ~NewtonGMRES() = default;
33
34 Scalar optError() const {
35 return fonc_.template lpNorm<2>();
36 }
37
38 template <typename VectorType>
39 void eval_fonc(const Scalar t, const MatrixBase<VectorType>& x, const Vector<dim>& solution) {
40 nlp_.eval_fonc_hu(t, x, solution, fonc_);
41 }
42
43 template <typename VectorType1, typename VectorType2, typename VectorType3, typename VectorType4>
44 void eval_b(const Scalar t, const MatrixBase<VectorType1>& x,
45 const MatrixBase<VectorType2>& solution,
46 const MatrixBase<VectorType3>& solution_update,
47 const MatrixBase<VectorType4>& b_vec) {
48 assert(solution.size() == dim);
49 assert(solution_update.size() == dim);
50 assert(b_vec.size() == dim);
51 updated_solution_ = solution + finite_difference_epsilon_ * solution_update;
52 nlp_.eval_fonc_hu(t, x, solution, fonc_);
53 nlp_.eval_fonc_hu(t, x, updated_solution_, fonc_1_);
54 CGMRES_EIGEN_CONST_CAST(VectorType4, b_vec) = - fonc_ - (fonc_1_ - fonc_) / finite_difference_epsilon_;
55 }
56
57 template <typename VectorType1, typename VectorType2, typename VectorType3, typename VectorType4>
59 const MatrixBase<VectorType2>& solution,
60 const MatrixBase<VectorType3>& solution_update,
61 const MatrixBase<VectorType4>& ax_vec) {
62 assert(solution.size() == dim);
63 assert(solution_update.size() == dim);
64 assert(ax_vec.size() == dim);
65 updated_solution_ = solution + finite_difference_epsilon_ * solution_update;
66 nlp_.eval_fonc_hu(t, x, updated_solution_, fonc_1_);
67 CGMRES_EIGEN_CONST_CAST(VectorType4, ax_vec) = (fonc_1_ - fonc_) / finite_difference_epsilon_;
68 }
69
70 void retrive_dummy(Vector<dim>& solution, const Scalar min_dummy) {
71 fonc_1_.setZero();
72 nlp_.retrive_dummy(solution, fonc_1_, min_dummy);
73 }
74
75 void retrive_mu(Vector<dim>& solution) {
76 fonc_1_.setZero();
77 nlp_.retrive_mu(solution, fonc_1_);
78 }
79
80 decltype(auto) x() const { return nlp_.x(); }
81
82 decltype(auto) lmd() const { return nlp_.lmd(); }
83
84 const NLP& get_nlp() const { return nlp_; }
85
86 void synchronize_ocp() { nlp_.synchronize_ocp(); }
87
88 EIGEN_MAKE_ALIGNED_OPERATOR_NEW
89
90private:
91 NLP nlp_;
92 Scalar finite_difference_epsilon_;
93 Vector<dim> updated_solution_, fonc_, fonc_1_;
94};
95
96} // namespace detail
97} // namespace cgmres
98
99#endif // CGMRES__NEWTON_GMRES_HPP_
Definition: newton_gmres.hpp:14
decltype(auto) x() const
Definition: newton_gmres.hpp:80
void eval_Ax(const Scalar t, const MatrixBase< VectorType1 > &x, const MatrixBase< VectorType2 > &solution, const MatrixBase< VectorType3 > &solution_update, const MatrixBase< VectorType4 > &ax_vec)
Definition: newton_gmres.hpp:58
void retrive_dummy(Vector< dim > &solution, const Scalar min_dummy)
Definition: newton_gmres.hpp:70
void eval_b(const Scalar t, const MatrixBase< VectorType1 > &x, const MatrixBase< VectorType2 > &solution, const MatrixBase< VectorType3 > &solution_update, const MatrixBase< VectorType4 > &b_vec)
Definition: newton_gmres.hpp:44
const NLP & get_nlp() const
Definition: newton_gmres.hpp:84
NewtonGMRES(const NLP &nlp, const Scalar finite_difference_epsilon)
Definition: newton_gmres.hpp:19
static constexpr int nx
Definition: newton_gmres.hpp:16
void synchronize_ocp()
Definition: newton_gmres.hpp:86
Scalar optError() const
Definition: newton_gmres.hpp:34
void retrive_mu(Vector< dim > &solution)
Definition: newton_gmres.hpp:75
void eval_fonc(const Scalar t, const MatrixBase< VectorType > &x, const Vector< dim > &solution)
Definition: newton_gmres.hpp:39
decltype(auto) lmd() const
Definition: newton_gmres.hpp:82
static constexpr int dim
Definition: newton_gmres.hpp:17
#define CGMRES_EIGEN_CONST_CAST(TYPE, OBJ)
Definition: macros.hpp:7
Definition: continuation_gmres.hpp:11
Eigen::Matrix< Scalar, size, 1 > Vector
Alias of Eigen::Vector.
Definition: types.hpp:23
double Scalar
Alias of double.
Definition: types.hpp:11
Eigen::MatrixBase< MatrixType > MatrixBase
Alias of Eigen::MatrixBase.
Definition: types.hpp:29