3 #include <fluxionum/UnivariateNewtonRaphsonMinimizer.h>
4 #include <fluxionum/DesignMatrix.h>
5 #include <fluxionum/TemporalDesignMatrix.h>
6 #include <fluxionum/IndexedDesignMatrix.h>
7 #include <fluxionum/DiffDesignMatrix.h>
8 #include <fluxionum/DiffDesignMatrixInterpolator.h>
9 #include <fluxionum/LinearPiecesFunction.h>
10 #include <fluxionum/ParametricLinearPiecesFunction.h>
11 #include <fluxionum/FixedIterativeEulerMethod.h>
12 #include <fluxionum/FixedParametricIterativeEulerMethod.h>
13 #include <fluxionum/ClosestLesserSampleFunction.h>
14 #include <fluxionum/ParametricClosestLesserSampleFunction.h>
21 namespace HeliosTests{
23 using namespace fluxionum;
38 double const eps = 0.00001;
67 bool testUnivariateNewtonRaphsonMinimization();
72 bool testDesignMatrixBuilding();
77 bool testDesignMatrixMethods();
82 bool testDiffDesignMatrix();
87 bool testDesignFunctions();
94 if(!testUnivariateNewtonRaphsonMinimization())
return false;
95 if(!testDesignMatrixBuilding())
return false;
96 if(!testDesignMatrixMethods())
return false;
97 if(!testDiffDesignMatrix())
return false;
98 if(!testDesignFunctions())
return false;
105 double expected = -3.0;
107 [] (
double x) ->
double {
return std::pow(x, 2)/2.0 + 3.0*x;},
108 [] (
double x) ->
double {
return x+3.0;},
109 [] (
double x) ->
double {
return 1;}
111 double x = unrm.
argmin(0.0);
112 return std::fabs(x-expected) <= eps;
119 arma::Mat<double>
const &
120 )> validateDesignMatrix = [&] (
122 vector<string>
const &colNames,
123 arma::Mat<double>
const &X
128 if(X.n_rows != nRows || X.n_cols != nCols)
return false;
131 for(
size_t j = 0 ; j < nCols ; ++j){
135 else if(!colNames.empty())
return false;
136 for(
size_t i = 0 ; i < nRows ; ++i){
137 for(
size_t j = 0 ; j < nCols ; ++j){
138 if(std::fabs(dm(i, j)-X(i, j)) > eps)
return false;
145 vector<string>
const &,
146 arma::Col<double>
const &, arma::Mat<double>
const &
147 )> validateTemporalDesignMatrix = [&](
149 string const &timeName,
150 vector<string>
const &colNames,
151 arma::Col<double>
const &t,
152 arma::Mat<double>
const &X
155 if(!validateDesignMatrix(tdm, colNames, X))
return false;
158 for(
size_t i = 0 ; i < nRows ; ++i){
159 if(std::fabs(tdm[i]-t[i]) > eps)
return false;
165 vector<string>
const&,
166 vector<int>
const&, arma::Mat<double>
const&
167 )> validateIndexedDesignMatrix = [&](
169 string const &indexName,
170 vector<string>
const &colNames,
171 vector<int>
const &ids,
172 arma::Mat<double>
const &X
175 if(!validateDesignMatrix(idm, colNames, X))
return false;
177 if(ids.size() != nRows)
return false;
179 for(
size_t i = 0 ; i < nRows ; ++i){
180 if(idm[i] != ids[i])
return false;
186 vector<string> colNamesR2({
"x",
"y"});
187 arma::Mat<double> X1 = arma::randn(5, 2);
189 if(!validateDesignMatrix(dm1, colNamesR2, X1))
return false;
191 if(!validateDesignMatrix(dm2, colNamesR2, arma::Mat<double>()))
194 if(!validateDesignMatrix(dm3, vector<string>(0), arma::Mat<double>()))
197 if(!validateDesignMatrix(dm4, colNamesR2, X1))
return false;
199 if(!validateDesignMatrix(dm4, colNamesR2, X1))
return false;
202 vector<string> colNamesR2t({
"x",
"y",
"t"});
203 arma::Col<double> t1 = arma::randn(5, 1);
204 arma::Mat<double> X2 = arma::randn(5, 3);
206 if(!validateTemporalDesignMatrix(tdm1,
"time", colNamesR2, t1, X1))
209 if(!validateTemporalDesignMatrix(tdm2,
"time", colNamesR2, t1, X1))
212 if(!validateTemporalDesignMatrix(tdm2,
"time", colNamesR2, t1, X1))
215 if(!validateTemporalDesignMatrix(tdm3,
"TIME", colNamesR2, t1, X1))
218 if(!validateTemporalDesignMatrix(tdm4,
"time", colNamesR2, t1, X1))
221 if(!validateTemporalDesignMatrix(
222 tdm5,
"t", colNamesR2, X2.col(2), X2.cols(0, 1)
226 vector<int> ids1({0, 1, 2, 3, 4});
228 if(!validateIndexedDesignMatrix(idm1,
"index", colNamesR2, ids1, X1))
231 if(!validateIndexedDesignMatrix(idm2,
"index", colNamesR2, ids1, X1))
234 if(!validateIndexedDesignMatrix(idm2,
"index", colNamesR2, ids1, X1))
237 if(!validateIndexedDesignMatrix(idm3,
"INDEX", colNamesR2, ids1, X1))
240 if(!validateIndexedDesignMatrix(idm4,
"index", colNamesR2, ids1, X1))
243 if(!validateIndexedDesignMatrix(
244 idm5,
"idx", colNamesR2,
256 vector<string> hf1({
"x",
"y",
"z"});
257 arma::Mat<double> Xf1(
268 string const dmf1Path = testDir +
"design_matrix_header.txt";
270 if(!validateDesignMatrix(dmf1, hf1, Xf1))
return false;
271 string const dmf2Path = testDir +
"design_matrix_nonheader.txt";
273 if(!validateDesignMatrix(dmf2, vector<string>(0), Xf1))
return false;
276 vector<string> htf1({
"x",
"y"});
277 arma::Mat<double> Xtf1(
288 arma::Col<double> ttf1(
"0 0.1 0.2 0.3 0.4 0.5 0.6 0.8 1.0");
289 string const tdmf1Path = testDir +
"temporal_design_matrix_header.txt";
291 if(!validateTemporalDesignMatrix(tdmf1,
"t", htf1, ttf1, Xtf1))
293 string const tdmf2Path = testDir +
"temporal_design_matrix_nonheader.txt";
295 if(!validateTemporalDesignMatrix(
296 tdmf2,
"time", vector<string>(0), ttf1, Xtf1
300 vector<string> hif1({
"x",
"y"});
301 arma::Mat<double> Xif1(
312 vector<int> iif1({0, 1, 2, 3, 4, 6, 8, 7, 5});
313 string const idmf1Path = testDir +
"indexed_design_matrix_header.txt";
315 if(!validateIndexedDesignMatrix(idmf1,
"idx", hif1, iif1, Xif1))
317 string const idmf2Path = testDir +
"indexed_design_matrix_nonheader.txt";
319 if(!validateIndexedDesignMatrix(
320 idmf2,
"index", vector<string>(0), iif1, Xif1
349 arma::Col<double>(
"0.0 0.1 0.2")
357 arma::Col<double>(
"0.3 0.4 0.5")
365 arma::Col<double>(
"0.4 0.5 0.3")
381 arma::Col<double>(
"-5 -4 -3 -2 -1 0 1 2 3 4 5")
389 vector<int>({1, 3, 5})
397 vector<int>({6, 8, 10})
403 arma::Mat<double> em(
411 if(arma::any(arma::vectorise(arma::abs(dm.
getX()-em)) > eps))
return false;
414 em = arma::Mat<double>(
422 arma::Col<double> et({
423 0.0, 0.1, 0.2, 0.3, 0.4, 0.5
425 if(arma::any(arma::vectorise(arma::abs(tdm.
getX()-em)) > eps))
430 em = arma::Mat<double>(
441 if(arma::any(arma::vectorise(arma::abs(idm.
getX()-em)) > eps))
443 for(
size_t i = 0 ; i < ei.size() ; ++i){
444 if(std::fabs(ei[i]-idm.
getIndices()[i]) > eps)
return false;
450 em = arma::Mat<double>(
455 if(arma::any(arma::vectorise(arma::abs(dm.
getX()-em)) > eps))
return false;
458 em = arma::Mat<double>(
463 if(arma::any(arma::vectorise(arma::abs(tdm.
getX()-em)) > eps))
467 em = arma::Mat<double>(
472 if(arma::any(arma::vectorise(arma::abs(idm.
getX()-em)) > eps))
486 et = arma::Col<double>(
"-5 -2 2 5");
487 em = arma::Mat<double>(
494 if(arma::any(arma::abs(tdm4.
getTimeVector()-et) > eps))
return false;
495 if(arma::any(arma::abs(arma::vectorise(tdm4.
getX()-em)) > eps))
506 vector<string>
const&,
507 arma::Col<double>
const&, arma::Mat<double>
const&,
509 )> validateDiffDesignMatrix = [&] (
511 string const &timeName,
512 vector<string>
const &colNames,
513 arma::Col<double>
const &t,
514 arma::Mat<double>
const &A,
515 DiffDesignMatrixType diffType
521 if(A.n_rows != nRows || A.n_cols != nCols)
return false;
524 for(
size_t j = 0 ; j < nCols ; ++j){
528 else if(!colNames.empty())
return false;
529 for(
size_t i = 0 ; i < nRows ; ++i){
530 for(
size_t j = 0 ; j < nCols ; ++j){
531 if(std::fabs(ddm(i, j)-A(i, j)) > eps)
return false;
536 for(
size_t i = 0 ; i < nRows ; ++i){
537 if(std::fabs(ddm[i]-t[i]) > eps)
return false;
544 arma::Col<double> Et1(
"-3 -2 -1 0 1 2");
545 arma::Col<double> Et3(
"-2.5 -1.5 -0.5 0.5 1.5");
546 arma::Mat<double> EA1(
554 arma::Mat<double> EA3(
561 vector<string> EcolNames1({
"f1",
"f2"});
562 string EtimeName1 =
"t";
563 string EtimeName2 =
"time";
566 vector<string> colNames1({
"t",
"f1",
"f2"});
567 arma::Mat<double> X1 = arma::Mat<double>(
576 arma::Mat<double> X2 = arma::Mat<double>(
594 tdm1, DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES
597 tdm2, DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES
600 tdm2, DiffDesignMatrixType::CENTRAL_FINITE_DIFFERENCES
604 if(!validateDiffDesignMatrix(
605 ddm1, EtimeName1, EcolNames1, Et1, EA1,
606 DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES
609 if(!validateDiffDesignMatrix(
610 ddm2, EtimeName2, vector<string>(0), Et1, EA1,
611 DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES
613 if(!validateDiffDesignMatrix(
614 ddm3, EtimeName2, vector<string>(0), Et3, EA3,
615 DiffDesignMatrixType::CENTRAL_FINITE_DIFFERENCES
619 string const ddmf1Path = testDir +
"diff_design_matrix_header.txt";
621 if(!validateDiffDesignMatrix(
622 ddmf1, EtimeName1, EcolNames1, Et1, EA1,
623 DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES
625 string const ddmf2Path = testDir +
"diff_design_matrix_nonheader.txt";
627 if(!validateDiffDesignMatrix(
628 ddmf2, EtimeName2, vector<string>(0), Et3, EA3,
629 DiffDesignMatrixType::CENTRAL_FINITE_DIFFERENCES
659 DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES
662 DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES
666 arma::Col<double> intercept(
669 arma::Col<double> slope = ffd1.
getA().col(0);
671 DiffDesignMatrixInterpolator::makeLinearPiecesFunction(
677 DiffDesignMatrixInterpolator::makeLinearPiecesFunction(
678 ffd1, tdm1, 0, &intercept, &slope
680 arma::Col<double> lpf1t(
"-1 0 1 2 3 4 5 6 7.5 9");
681 arma::Col<double> lpf1E(
"-1 0 1 2 2.5 3 3 3 3.5 4");
682 for(
size_t i = 0 ; i < lpf1t.n_elem ; ++i){
683 if(std::fabs(lpf1(lpf1t.at(i)) - lpf1E.at(i)) > eps)
return false;
684 if(std::fabs(lpf2(lpf1t.at(i)) - lpf1E.at(i)) > eps)
return false;
689 DiffDesignMatrixInterpolator::makeParametricLinearPiecesFunction(
694 DiffDesignMatrixInterpolator::makeParametricLinearPiecesFunction(
697 arma::Mat<double> plpf1E(
709 for(
size_t i = 0 ; i < lpf1t.n_elem ; ++i){
710 arma::Col<double> plpf1y = plpf1(lpf1t.at(i));
711 arma::Col<double> plpf2y = plpf2(lpf1t.at(i));
712 for(
size_t j = 0 ; j < 2 ; ++j){
713 if(std::fabs(plpf1y.at(j) - plpf1E.at(i, j)) > eps)
return false;
714 if(std::fabs(plpf2y.at(j) - plpf1E.at(i, j)) > eps)
return false;
719 arma::Col<double> fiem1t(
"0 1 2 3 4 5 6 7.5 9");
720 arma::Col<double> fiem1E(
"0 1 2 2.5 3 3 3 3.5 4");
721 arma::Col<double> ySamples1(
724 arma::Col<double> ySamples2;
725 arma::Col<double> dydtSamples1(
728 arma::Col<double> dydtSamples2;
736 DiffDesignMatrixInterpolator::makeFixedIterativeEulerMethod(
742 DiffDesignMatrixInterpolator::makeFixedIterativeEulerMethod(
750 for(
size_t i = 0 ; i < fiem1t.n_elem ; ++i){
751 double const h = (i==0) ? 0 : fiem1t.at(i) - fiem1t.at(i-1);
752 if(std::fabs(fiem1(h) - fiem1E.at(i)) > eps)
return false;
753 if(std::fabs(fiem2(h) - fiem1E.at(i)) > eps)
return false;
757 arma::Mat<double> pfiem1E(
768 arma::Mat<double>
const &pySamples1 = tdm2.
getX();
769 arma::Mat<double>
const &pdydtSamples1 = ffd2.
getA();
776 DiffDesignMatrixInterpolator::makeFixedParametricIterativeEulerMethod(
783 DiffDesignMatrixInterpolator::makeFixedParametricIterativeEulerMethod(
789 for(
size_t i = 0 ; i < fiem1t.n_elem ; ++i){
790 double const h = (i==0) ? 0 : fiem1t.at(i) - fiem1t.at(i-1);
791 if(arma::any(arma::abs(fpiem1(h) - pfiem1E.row(i).as_col()) > eps))
793 if(arma::any(arma::abs(fpiem2(h) - pfiem1E.row(i).as_col()) > eps))
BaseTest class.
Definition: BaseTest.h:20
Fluxionum test.
Definition: FluxionumTest.h:31
bool run() override
Definition: FluxionumTest.h:92
std::string testDir
The directory where test files are located.
Definition: FluxionumTest.h:42
bool testDesignMatrixBuilding()
Test the building of different design matrices.
Definition: FluxionumTest.h:115
bool testDesignFunctions()
Test functions from design matrices.
Definition: FluxionumTest.h:636
bool testUnivariateNewtonRaphsonMinimization()
Test univariate Newton-Raphson minimization.
Definition: FluxionumTest.h:104
FluxionumTest(std::string testDir)
Fluxionum test constructor.
Definition: FluxionumTest.h:49
bool testDesignMatrixMethods()
Test the methods of DesignMatrix.
Definition: FluxionumTest.h:327
bool testDiffDesignMatrix()
Test the generation of differential design matrices.
Definition: FluxionumTest.h:502
string const & getColumnName(size_t const j) const
Obtain the name of the -th column.
Definition: AbstractDesignMatrix.h:111
bool hasColumnNames() const
Check whether there are available column names for the AbstractDesignMatrix (true) or not (false)
Definition: AbstractDesignMatrix.h:104
vector< string > const & getColumnNames() const
Obtain a constant/read reference to the column names.
Definition: AbstractDesignMatrix.h:126
Closest lesser sample function.
Definition: ClosestLesserSampleFunction.h:28
This class represents a DesignMatrix .
Definition: DesignMatrix.h:41
virtual void mergeInPlace(DesignMatrix const &dm)
Merge given DesignMatrix into this DesignMatrix.
Definition: DesignMatrix.h:119
size_t getNumColumns() const override
Obtain the number of columns of the DesignMatrix .
Definition: DesignMatrix.h:244
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
void swapColumns(vector< long unsigned int > const &indices)
Definition: DesignMatrix.h:125
The heart of a differential design matrix is the idea that the columns are defining the values of var...
Definition: DiffDesignMatrix.h:169
size_t getNumRows() const override
Obtain the number of rows of the DiffDesignMatrix .
Definition: DiffDesignMatrix.h:443
string const & getTimeName() const
Obtain the name of the time attribute.
Definition: DiffDesignMatrix.h:469
size_t getNumColumns() const override
Obtain the number of columns of the DiffDesignMatrix .
Definition: DiffDesignMatrix.h:450
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
Definition: IndexedDesignMatrix.h:43
void mergeInPlace(DesignMatrix< VarType > const &dm) override
Extend DesignMatrix::mergeInPlace method so the indices are also merged.
Definition: IndexedDesignMatrix.h:269
vector< IndexType > const & getIndices() const
Obtain a constant/read reference to the vector of indices.
Definition: IndexedDesignMatrix.h:287
string const getIndexName() const
Obtain the name of the index attribute.
Definition: IndexedDesignMatrix.h:293
Linear pieces function.
Definition: LinearPiecesFunction.h:28
Parametric linear pieces function.
Definition: ParametricLinearPiecesFunction.h:35
This class represents a DesignMatrix where each row satisfy that its elements can be placed into a co...
Definition: TemporalDesignMatrix.h:74
virtual size_t slopeFilter(VarType const tau)
Apply a filter based on the first derivative understood as the slope.
Definition: TemporalDesignMatrix.h:391
string const & getTimeName() const
Obtain the name of the time attribute.
Definition: TemporalDesignMatrix.h:424
DiffDesignMatrix< TimeType, VarType > toDiffDesignMatrix(DiffDesignMatrixType diffType=DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES, bool const sort=true) const
Build a DiffDesignMatrix from the TemporalDesignMatrix.
arma::Col< TimeType > const & getTimeVector() const
Obtain a constant/read reference to the time vector .
Definition: TemporalDesignMatrix.h:418
virtual void sortByTime()
Sort all rows of matrix and all components of vector so , considering the row is associated with ...
Definition: TemporalDesignMatrix.h:351
void mergeInPlace(DesignMatrix< VarType > const &dm) override
Extend DesignMatrix::mergeInPlace method so the time values are also merged.
Definition: TemporalDesignMatrix.h:315
virtual void shiftTime(TimeType const x)
Add given time to all time values in vector .
Definition: TemporalDesignMatrix.h:343
Implementation of univariate Newton-Raphson minimizer.
Definition: UnivariateNewtonRaphsonMinimizer.h:23
IT argmin(IT x) override
Implementation of the univariate Newton-Raphson minimization.