MUQ  0.4.3
PyModPiece.cpp
Go to the documentation of this file.
2 
3 using namespace muq::Modeling;
4 
5 PyModPiece::PyModPiece(Eigen::VectorXi const& inputSizes,
6  Eigen::VectorXi const& outputSizes)
7  : ModPiece(inputSizes, outputSizes) {}
8 
10  this->EvaluateImpl(ToStdVec(input));
11 }
12 
13 void PyModPiece::GradientImpl(unsigned int const outputDimWrt,
14  unsigned int const inputDimWrt,
15  ref_vector<Eigen::VectorXd> const& input,
16  Eigen::VectorXd const& sensitivity) {
17  this->GradientImpl(outputDimWrt, inputDimWrt, ToStdVec(input), sensitivity);
18 }
19 
20 void PyModPiece::GradientImpl(unsigned int const outputDimWrt,
21  unsigned int const inputDimWrt,
22  std::vector<Eigen::VectorXd> const& input,
23  Eigen::VectorXd const& sensitivity)
24 {
25  gradient = GradientByFD(outputDimWrt, inputDimWrt, input,sensitivity);
26 }
27 
28 void PyModPiece::JacobianImpl(unsigned int const outputDimWrt,
29  unsigned int const inputDimWrt,
30  ref_vector<Eigen::VectorXd> const& input) {
31  this->JacobianImpl(outputDimWrt, inputDimWrt, ToStdVec(input));
32 }
33 
34 void PyModPiece::JacobianImpl(unsigned int const outputDimWrt,
35  unsigned int const inputDimWrt,
36  std::vector<Eigen::VectorXd> const& input)
37 {
38  jacobian = JacobianByFD(outputDimWrt, inputDimWrt, input);
39 }
40 
41 
42 void PyModPiece::ApplyJacobianImpl(unsigned int const outputDimWrt,
43  unsigned int const inputDimWrt,
44  ref_vector<Eigen::VectorXd> const& input,
45  Eigen::VectorXd const& vec) {
46  this->ApplyJacobianImpl(outputDimWrt, inputDimWrt, ToStdVec(input), vec);
47 }
48 
49 void PyModPiece::ApplyJacobianImpl(unsigned int const outputDimWrt,
50  unsigned int const inputDimWrt,
51  std::vector<Eigen::VectorXd> const& input,
52  Eigen::VectorXd const& vec)
53 {
54  jacobianAction = ApplyJacobianByFD(outputDimWrt, inputDimWrt, input, vec);
55 }
56 
57 void PyModPiece::ApplyHessianImpl(unsigned int const outputDimWrt,
58  unsigned int const inputDimWrt1,
59  unsigned int const inputDimWrt2,
60  ref_vector<Eigen::VectorXd> const& input,
61  Eigen::VectorXd const& sens,
62  Eigen::VectorXd const& vec)
63 {
64  this->ApplyHessianImpl(outputDimWrt, inputDimWrt1, inputDimWrt2, ToStdVec(input), sens, vec);
65 }
66 
67 void PyModPiece::ApplyHessianImpl(unsigned int const outputDimWrt,
68  unsigned int const inputDimWrt1,
69  unsigned int const inputDimWrt2,
70  std::vector<Eigen::VectorXd> const& input,
71  Eigen::VectorXd const& sens,
72  Eigen::VectorXd const& vec)
73 {
74  hessAction = ApplyHessianByFD(outputDimWrt, inputDimWrt1, inputDimWrt2, input, sens, vec);
75 }
Provides an abstract interface for defining vector-valued model components.
Definition: ModPiece.h:148
Eigen::MatrixXd jacobian
Definition: ModPiece.h:506
virtual Eigen::VectorXd GradientByFD(unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity)
Definition: ModPiece.cpp:111
virtual Eigen::VectorXd ApplyHessianByFD(unsigned int const outWrt, unsigned int const inWrt1, unsigned int const inWrt2, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec)
Definition: ModPiece.cpp:448
Eigen::MatrixXd ApplyJacobianByFD(unsigned int outWrt, unsigned int inWrt, Args const &... args)
Definition: ModPiece.h:303
Eigen::VectorXd gradient
Definition: ModPiece.h:504
Eigen::MatrixXd JacobianByFD(unsigned int outWrt, unsigned int inWrt, Args const &... args)
Definition: ModPiece.h:297
std::vector< Eigen::VectorXd > ToStdVec(ref_vector< Eigen::VectorXd > const &input)
Definition: ModPiece.cpp:586
Eigen::VectorXd hessAction
Definition: ModPiece.h:507
Eigen::VectorXd jacobianAction
Definition: ModPiece.h:505
virtual void GradientImpl(unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity) override
Definition: PyModPiece.cpp:13
virtual void ApplyJacobianImpl(unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec) override
Definition: PyModPiece.cpp:42
virtual void EvaluateImpl(ref_vector< Eigen::VectorXd > const &input) override
Definition: PyModPiece.cpp:9
virtual void ApplyHessianImpl(unsigned int const outputDimWrt, unsigned int const inputDimWrt1, unsigned int const inputDimWrt2, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec)
Definition: PyModPiece.cpp:67
virtual void JacobianImpl(unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input) override
Definition: PyModPiece.cpp:28
PyModPiece(Eigen::VectorXi const &inputSizes, Eigen::VectorXi const &outputSizes)
Definition: PyModPiece.cpp:5
std::vector< std::reference_wrapper< const T > > ref_vector
A vector of references to something ...
Definition: WorkPiece.h:37