Sandals  v0.0.0
A C++ library for ODEs/DAEs integration
Loading...
Searching...
No Matches
Linear.hh
Go to the documentation of this file.
1/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *\
2 * Copyright (c) 2025, Davide Stocco and Enrico Bertolazzi. *
3 * *
4 * The Sandals project is distributed under the BSD 2-Clause License. *
5 * *
6 * Davide Stocco Enrico Bertolazzi *
7 * University of Trento University of Trento *
8 * e-mail: davide.stocco@unitn.it e-mail: enrico.bertolazzi@unitn.it *
9\* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
10#pragma once
11
12#ifndef SANDALS_LINEAR_SYSTEM_HH
13#define SANDALS_LINEAR_SYSTEM_HH
14
15#include <Sandals.hh>
17
18namespace Sandals {
19
20 /*\
21 | _ _
22 | | | (_)_ __ ___ __ _ _ __
23 | | | | | '_ \ / _ \/ _` | '__|
24 | | |___| | | | | __/ (_| | |
25 | |_____|_|_| |_|\___|\__,_|_|
26 |
27 \*/
28
39 template <typename Real, Integer N, Integer M = 0>
40 class Linear : public Explicit<Real, N, M>
41 {
42 public:
43 using Pointer = std::shared_ptr<SemiExplicit<Real, N, M>>;
50
51 private:
52 Eigen::FullPivLU<MatrixE> m_lu;
53
54 public:
58 Linear() : Explicit<Real, N, M>(Type::LINEAR, "(missing name)") {}
59
64 Linear(std::string t_name) : Explicit<Real, N, M>(Type::LINEAR, t_name) {}
65
79 VectorF F(VectorF const & x, VectorF const & x_dot, Real const t) const override
80 {
81 return this->E(t)*x_dot - this->A(t)*x - this->b(t);
82 }
83
99 MatrixJF JF_x(VectorF const & /*x*/, VectorF const & /*x_dot*/, Real const t) const override {return -this->A(t);}
100
116 MatrixJF JF_x_dot(VectorF const & /*x*/, VectorF const & /*x_dot*/, Real const t) const override {return this->E(t);}
117
129 VectorF f(VectorF const & x, Real const t) const override
130 {
131 this->m_lu.compute(this->E(t));
132 SANDALS_ASSERT(this->m_lu.rank() == N, "Sandals:Linear::f(...): singular mass matrix E(t) detected.");
133 return this->m_lu.solve(this->A(t)*x + this->b(t));
134 }
135
149 MatrixJF Jf_x(VectorF const & /*x*/, Real const t) const override
150 {
151 this->m_lu.compute(this->E(t));
152 SANDALS_ASSERT(this->m_lu.rank() == N, "Sandals:Linear::Jf_x(...): singular mass matrix E(t) detected.");
153 return this->m_lu.solve(this->A(t));
154 }
155
161 virtual MatrixE E(Real const t) const = 0;
162
168 virtual MatrixA A(Real const t) const = 0;
169
175 virtual VectorB b(Real const t) const = 0;
176
177 }; // class Linear
178
179 /*\
180 | _ _ __ __
181 | | | (_)_ __ ___ __ _ _ _\ \ / / __ __ _ _ __ _ __ ___ _ __
182 | | | | | '_ \ / _ \/ _` | '__\ \ /\ / / '__/ _` | '_ \| '_ \ / _ \ '__|
183 | | |___| | | | | __/ (_| | | \ V V /| | | (_| | |_) | |_) | __/ |
184 | |_____|_|_| |_|\___|\__,_|_| \_/\_/ |_| \__,_| .__/| .__/ \___|_|
185 | |_| |_|
186 \*/
187
198 template <typename Real, Integer N, Integer M = 0>
199 class LinearWrapper : public Linear<Real, N, M>
200 {
201 public:
202 using Pointer = std::shared_ptr<LinearWrapper<Real, N, M>>;
203 using typename Linear<Real, N, M>::VectorF;
208 using typename Explicit<Real, N, M>::VectorH;
209 using typename Explicit<Real, N, M>::MatrixJH;
210 using FunctionE = std::function<MatrixE(Real const)>;
211 using FunctionA = std::function<MatrixA(Real const)>;
212 using FunctionB = std::function<VectorB(Real const)>;
213 using FunctionH = std::function<VectorH(VectorF const &, Real const)>;
214 using FunctionJH = std::function<MatrixJH(VectorF const &, Real const)>;
215 using FunctionID = std::function<bool(VectorF const &, Real const)>;
216
217 inline const static FunctionH DefaultH = [](VectorF const &, Real const) {return VectorH::Zero();};
218 inline const static FunctionJH DefaultJH = [](VectorF const &, Real const) {return MatrixJH::Zero();};
219 inline const static FunctionID DefaultID = [](VectorF const &, Real const) {return true;};
220
221 private:
222 FunctionE m_E{nullptr};
223 FunctionA m_A{nullptr};
224 FunctionB m_b{nullptr};
225 FunctionH m_h{nullptr};
228
229 public:
240 t_Jh_x = DefaultJH, FunctionID t_in_domain = DefaultID) : Linear<Real, N, M>("(missing name)"),
241 m_E(t_E), m_A(t_A), m_b(t_b), m_h(t_h), m_Jh_x(t_Jh_x), m_in_domain(t_in_domain)
242 {}
243
254 LinearWrapper(std::string t_name, FunctionE t_E, FunctionA t_A, FunctionB t_b, FunctionH t_h = DefaultH,
255 FunctionJH t_Jh_x = DefaultJH, FunctionID t_in_domain = DefaultID) : Linear<Real, N, M>(t_name),
256 m_E(t_E), m_A(t_A), m_b(t_b), m_h(t_h), m_Jh_x(t_Jh_x), m_in_domain(t_in_domain)
257 {}
258
263
268 FunctionE & E() {return this->m_E;}
269
274 FunctionA & A() {return this->m_A;}
275
280 FunctionB & b() {return this->m_b;}
281
286 FunctionH & h() {return this->m_h;}
287
292 FunctionJH & Jh_x() {return this->m_Jh_x;}
293
298 FunctionID & in_domain() {return this->m_in_domain;}
299
305 MatrixE E(Real const t) const override
306 {
307 return this->m_E(t);
308 }
309
315 MatrixA A(Real const t) const override
316 {
317 return this->m_A(t);
318 }
319
325 VectorB b(Real const t) const override
326 {
327 return this->m_b(t);
328 }
329
336 VectorH h(VectorF const & x, Real const t) const override
337 {
338 return this->m_h(x, t);
339 }
340
354 MatrixJH Jh_x(VectorF const & x, Real const t) const override
355 {
356 return this->m_Jh_x(x, t);
357 }
358
366 bool in_domain(VectorF const & x, Real const t) const override
367 {
368 return this->m_in_domain(x, t);
369 }
370
371 }; // class LinearWrapper
372
373} // namespace Sandals
374
375#endif // SANDALS_LINEAR_SYSTEM_HH
#define SANDALS_ASSERT(COND, MSG)
Definition Sandals.hh:44
Explicit(Type t_type, std::string t_name)
Definition Explicit.hh:55
typename Implicit< Real, N, M >::Type Type
Definition Explicit.hh:47
typename Implicit< Real, N, M >::MatrixJF MatrixJF
Definition Explicit.hh:46
typename Implicit< Real, N, M >::VectorF VectorF
Definition Explicit.hh:45
Eigen::Vector< Real, M > VectorH
Definition Implicit.hh:48
Eigen::Matrix< Real, M, N > MatrixJH
Definition Implicit.hh:49
typename Explicit< Real, N, M >::Type Type
Definition Linear.hh:49
virtual VectorB b(Real const t) const =0
Eigen::FullPivLU< MatrixE > m_lu
Definition Linear.hh:52
typename Explicit< Real, N, M >::VectorF VectorF
Definition Linear.hh:44
virtual MatrixE E(Real const t) const =0
typename Explicit< Real, N, M >::MatrixJF MatrixA
Definition Linear.hh:47
MatrixJF JF_x_dot(VectorF const &, VectorF const &, Real const t) const override
Definition Linear.hh:116
typename Explicit< Real, N, M >::MatrixJF MatrixE
Definition Linear.hh:46
virtual MatrixA A(Real const t) const =0
Linear()
Definition Linear.hh:58
typename Explicit< Real, N, M >::MatrixJF MatrixJF
Definition Linear.hh:45
VectorF F(VectorF const &x, VectorF const &x_dot, Real const t) const override
Definition Linear.hh:79
typename Explicit< Real, N, M >::VectorF VectorB
Definition Linear.hh:48
VectorF f(VectorF const &x, Real const t) const override
Definition Linear.hh:129
std::shared_ptr< SemiExplicit< Real, N, M > > Pointer
Definition Linear.hh:43
MatrixJF Jf_x(VectorF const &, Real const t) const override
Definition Linear.hh:149
Linear(std::string t_name)
Definition Linear.hh:64
MatrixJF JF_x(VectorF const &, VectorF const &, Real const t) const override
Definition Linear.hh:99
LinearWrapper(FunctionE t_E, FunctionA t_A, FunctionB t_b, FunctionH t_h=DefaultH, FunctionJH t_Jh_x=DefaultJH, FunctionID t_in_domain=DefaultID)
Definition Linear.hh:239
MatrixA A(Real const t) const override
Definition Linear.hh:315
FunctionH & h()
Definition Linear.hh:286
MatrixJH Jh_x(VectorF const &x, Real const t) const override
Definition Linear.hh:354
~LinearWrapper()
Definition Linear.hh:262
std::function< VectorB(Real const)> FunctionB
Definition Linear.hh:212
static const FunctionH DefaultH
Definition Linear.hh:217
VectorH h(VectorF const &x, Real const t) const override
Definition Linear.hh:336
std::function< MatrixE(Real const)> FunctionE
Definition Linear.hh:210
VectorB b(Real const t) const override
Definition Linear.hh:325
FunctionH m_h
Definition Linear.hh:225
FunctionB m_b
Definition Linear.hh:224
FunctionID & in_domain()
Definition Linear.hh:298
static const FunctionJH DefaultJH
Definition Linear.hh:218
LinearWrapper(std::string t_name, FunctionE t_E, FunctionA t_A, FunctionB t_b, FunctionH t_h=DefaultH, FunctionJH t_Jh_x=DefaultJH, FunctionID t_in_domain=DefaultID)
Definition Linear.hh:254
FunctionB & b()
Definition Linear.hh:280
FunctionE & E()
Definition Linear.hh:268
FunctionE m_E
Definition Linear.hh:222
FunctionA m_A
Definition Linear.hh:223
FunctionID m_in_domain
Definition Linear.hh:227
static const FunctionID DefaultID
Definition Linear.hh:219
std::function< MatrixJH(VectorF const &, Real const)> FunctionJH
Definition Linear.hh:214
std::function< bool(VectorF const &, Real const)> FunctionID
Definition Linear.hh:215
MatrixE E(Real const t) const override
Definition Linear.hh:305
typename Linear< Real, N, M >::MatrixJF MatrixJF
Definition Linear.hh:204
bool in_domain(VectorF const &x, Real const t) const override
Definition Linear.hh:366
FunctionJH & Jh_x()
Definition Linear.hh:292
typename Linear< Real, N, M >::MatrixJF MatrixA
Definition Linear.hh:206
std::function< MatrixA(Real const)> FunctionA
Definition Linear.hh:211
std::shared_ptr< LinearWrapper< Real, N, M > > Pointer
Definition Linear.hh:202
FunctionA & A()
Definition Linear.hh:274
std::function< VectorH(VectorF const &, Real const)> FunctionH
Definition Linear.hh:213
typename Linear< Real, N, M >::MatrixJF MatrixE
Definition Linear.hh:205
typename Linear< Real, N, M >::VectorF VectorB
Definition Linear.hh:207
FunctionJH m_Jh_x
Definition Linear.hh:226
The namespace for the Sandals library.
Definition Sandals.hh:89