Helios++
Helios software for LiDAR simulations
DiffDesignMatrixInterpolator.h
1 #pragma once
2 
3 #include <fluxionum/DiffDesignMatrix.h>
4 #include <fluxionum/LinearPiecesFunction.h>
5 #include <fluxionum/ParametricLinearPiecesFunction.h>
6 #include <fluxionum/FixedIterativeEulerMethod.h>
7 #include <fluxionum/FixedParametricIterativeEulerMethod.h>
8 #include <fluxionum/ClosestLesserSampleFunction.h>
9 #include <fluxionum/ParametricClosestLesserSampleFunction.h>
10 #include <fluxionum/FluxionumException.h>
11 
12 #include <armadillo>
13 
14 namespace fluxionum {
15 
21 namespace DiffDesignMatrixInterpolator{
22 
23 // *** MAKE METHODS *** //
24 // ********************** //
36 template <typename A, typename B>
38  DiffDesignMatrix<A, B> const &ddm,
39  arma::Col<B> const &slope,
40  arma::Col<B> const &intercept
41 ){
43  ddm.getTimeVector(),
44  slope,
45  intercept
46  );
47 }
48 
54 template <typename A, typename B>
56  DiffDesignMatrix<A, B> const &ddm,
57  DesignMatrix<B> const &dm,
58  size_t const colIdx,
59  arma::Col<B> *intercept,
60  arma::Col<B> *slope
61 ){
62  switch(ddm.getDiffType()){
63  case DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES:{
64  *intercept = dm.getColumnCopy(colIdx).subvec(0, dm.getNumRows()-2);
65  *slope = ddm.getA().col(colIdx);
67  ddm,
68  *slope,
69  *intercept
70  );
71  }
72  case DiffDesignMatrixType::CENTRAL_FINITE_DIFFERENCES:{
73  throw FluxionumException(
74  "DiffDesignMatrixInterpolator::makeLinearPiecesFunction:\n"
75  "\tCentral finite differences not supported"
76  );
77  }
78  default:{
79  throw FluxionumException(
80  "DiffDesignMatrixInterpolator::makeLinearPiecesFunction:\n"
81  "\tUnexpected differential type"
82  );
83  }
84  }
85 }
86 
97 template <typename A, typename B>
99  DiffDesignMatrix<A, B> const &ddm,
100  arma::Mat<B> const &intercepts
101 ){
103  ddm.getTimeVector(),
104  ddm.getA(),
105  intercepts
106  );
107 }
108 
114 template <typename A, typename B>
116  DiffDesignMatrix<A, B> const &ddm,
117  DesignMatrix<B> const &dm
118 ){
119  switch(ddm.getDiffType()){
120  case DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES:{
122  ddm,
123  dm.getX()
124  );
125  }
126  case DiffDesignMatrixType::CENTRAL_FINITE_DIFFERENCES:{
127  throw FluxionumException(
128  "DiffDesignMatrixInterpolator::"
129  "makeParametricLinearPiecesFunction:\n"
130  "\tCentral finite differences not supported"
131  );
132  }
133  default:{
134  throw FluxionumException(
135  "DiffDesignMatrixInterpolator::"
136  "makeParametricLinearPiecesFunction:\n"
137  "\tUnexpected differential type"
138  );
139  }
140  }
141 }
142 
153 template <typename A, typename B>
155  DiffDesignMatrix<A, B> const &ddm,
156  arma::Col<B> const &y,
157  Function<A, B> &dydt
158 ){
160  dydt,
161  ddm.getTimeVector()(0),
162  y(0),
163  ddm.getTimeVector(),
164  y,
165  0
166  );
167 }
168 
175 template <typename A, typename B>
177  DiffDesignMatrix<A, B> const &ddm,
178  DesignMatrix<B> const &dm,
179  size_t const colIdx,
180  arma::Col<B> *y,
182  arma::Col<B> *dydtSamples
183 ){
184  *y = arma::Col<B>(dm.getColumn(colIdx).subvec(0, dm.getNumRows()-2));
185  *dydtSamples = ddm.getA().col(colIdx);
187  ddm.getTimeVector(),
188  *dydtSamples,
189  0
190  );
191  switch(ddm.getDiffType()){
192  case DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES:{
194  ddm,
195  *y,
196  *dydt
197  );
198  }
199  case DiffDesignMatrixType::CENTRAL_FINITE_DIFFERENCES:{
200  throw FluxionumException(
201  "DiffDesignMatrixInterpolator::"
202  "makeFixedIterativeEulerMethod:\n"
203  "\tCentral finite differences not supported"
204  );
205  }
206  default:{
207  throw FluxionumException(
208  "DiffDesignMatrixInterpolator::"
209  "makeFixedIterativeEulerMethod:\n"
210  "\tUnexpected differential type"
211  );
212  }
213  }
214 }
215 
228 template <typename A, typename B>
231  DiffDesignMatrix<A, B> const &ddm,
232  arma::Mat<B> const &y,
233  Function<A, arma::Col<B>> &dydt
234 ){
236  dydt,
237  ddm.getTimeVector()(0),
238  y.row(0).as_col(),
239  ddm.getTimeVector(),
240  y,
241  0
242  );
243 }
244 
253 template <typename A, typename B>
256  DiffDesignMatrix<A, B> const &ddm,
257  DesignMatrix<B> const &dm,
259 ){
261  ddm.getTimeVector(),
262  ddm.getA(),
263  0
264  );
265  switch(ddm.getDiffType()){
266  case DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES:{
268  ddm,
269  dm.getX(),
270  *dydt
271  );
272  }
273  case DiffDesignMatrixType::CENTRAL_FINITE_DIFFERENCES:{
274  throw FluxionumException(
275  "DiffDesignMatrixInterpolator::"
276  "makeFixedParametricIterativeEulerMethod:\n"
277  "\tCentral finite differences not supported"
278  );
279  }
280  default:{
281  throw FluxionumException(
282  "DiffDesignMatrixInterpolator::"
283  "makeFixedParametricIterativeEulerMethod:\n"
284  "\tUnexpected differential type"
285  );
286  }
287  }
288 }
289 
290 
291 }}
Closest lesser sample function.
Definition: ClosestLesserSampleFunction.h:28
This class represents a DesignMatrix .
Definition: DesignMatrix.h:41
arma::Col< T > getColumnCopy(size_t const j) const
Like DesignMatrix::getColumn(size_t const) but returning a copy by value instead of a view-like refer...
Definition: DesignMatrix.h:279
size_t getNumRows() const override
Obtain the number of rows of the DesignMatrix .
Definition: DesignMatrix.h:237
arma::Mat< T > const & getX() const
Obtain a constant/read reference to the matrix.
Definition: DesignMatrix.h:219
arma::subview_col< T > const getColumn(size_t const j) const
Obtain the -th column of the DesignMatrix .
Definition: DesignMatrix.h:272
The heart of a differential design matrix is the idea that the columns are defining the values of var...
Definition: DiffDesignMatrix.h:169
arma::Mat< VarType > const & getA() const
Obtain a constant/read reference to the matrix.
Definition: DiffDesignMatrix.h:481
arma::Col< TimeType > const & getTimeVector() const
Obtain the sorted series of time values DiffDesignMatrix::t.
Definition: DiffDesignMatrix.h:475
DiffDesignMatrixType getDiffType() const
Obtain the type of differential of the DiffDesignMatrix.
Definition: DiffDesignMatrix.h:463
Fixed iterative Euler method.
Definition: FixedIterativeEulerMethod.h:25
Fixed parametric iterative Euler method.
Definition: FixedParametricIterativeEulerMethod.h:27
Base class for fluxionum exceptions.
Definition: FluxionumException.h:16
Abstract class representing a function.
Definition: Function.h:27
Linear pieces function.
Definition: LinearPiecesFunction.h:28
Definition: ParametricClosestLesserSampleFunction.h:10
Parametric linear pieces function.
Definition: ParametricLinearPiecesFunction.h:35
FixedParametricIterativeEulerMethod< A, B > makeFixedParametricIterativeEulerMethod(DiffDesignMatrix< A, B > const &ddm, arma::Mat< B > const &y, Function< A, arma::Col< B >> &dydt)
Obtain a function representing a fixed parametric iterative Euler method from given DiffDesignMatrix,...
Definition: DiffDesignMatrixInterpolator.h:230
ParametricLinearPiecesFunction< A, B > makeParametricLinearPiecesFunction(DiffDesignMatrix< A, B > const &ddm, arma::Mat< B > const &intercepts)
Obtain a parametric linear pieces function from given DiffDesignMatrix and known values.
Definition: DiffDesignMatrixInterpolator.h:98
FixedIterativeEulerMethod< A, B > makeFixedIterativeEulerMethod(DiffDesignMatrix< A, B > const &ddm, arma::Col< B > const &y, Function< A, B > &dydt)
Obtain a function representing a fixed iterative Euler method from given DiffDesignMatrix,...
Definition: DiffDesignMatrixInterpolator.h:154
LinearPiecesFunction< A, B > makeLinearPiecesFunction(DiffDesignMatrix< A, B > const &ddm, arma::Col< B > const &slope, arma::Col< B > const &intercept)
Obtain a linear pieces function from given DiffDesignMatrix and known values.
Definition: DiffDesignMatrixInterpolator.h:37