Provides an abstract interface for defining vector-valued model components. More...
#include <ModPiece.h>
Provides an abstract interface for defining vector-valued model components.
Let \(x_i \in \mathbb{R}^{N_{xi}}\) be a collection of \(M_x\) vectors index by \(i\). Each vector could have a different size. Now let \( f_i(x_1,\ldots, x_{M_x}) \in \mathbb{R}^{N_{fi}}\) denote a collection of \(M_f\) model outputs, each of which can depend on all of the input vectors \(x_i\). The ModPiece class defines an interface for defining the mappings of the form \([x_1, \ldots, x_{M_x}] \rightarrow [f_1(x_1, \ldots, x_N), \ldots, f_{M_f}(x_1, \ldots, x_{M_{x}})] \). While this form is quite general, functions with a single input and a single output can easily be defined as a ModPiece with \(M_x = M_f = 1\).
The ModPiece class (and therefore all its children) contains two vectors, inputSizes
and outputSizes
, that are useful for investigating the number and size of the inputs and outputs.
Defining a new model component The ModPiece
class defines a general modeling interface that allows algorithms (such as MCMC) to work with arbitrary models, without needing to know the model details. Users can add create new models in two ways:
As an example of creating a new ModPiece, consider two vectors \(x_1\in\mathbb{R}^2\) and \(x_2\in\mathbb{R}^2\) as well as three functions \(f_1(x_1,x_2) = sin(x_1*x_2)\), \(f_2(x_1,x_2) = cos(x_1*x_2)\), and \(f_3(x_1,x_2) = tan(x_1*x_2)\). The following code implements this relationship with a ModPiece. Notice that the input sizes are \([2,2]\) and the output sizes (assuming componentwise products) are \([2,2,2]\).
If MUQ was compiled with Python bindings, it is also possible to create new children of ModPiece in Python. In Python, we use numpy arrays instead of Eigen::VectorXd. With this, the above class in Python might take the form
Once created, we can evaluate a ModPiece using the Evaluate
function. Notice that the user-level Evaluate
function and the EvaluateImpl
function discussed above are different. Direct calls to EvaluateImpl
should be avoided (hence the reason for making it protected) because the Evaluate
function performs additional checks to make sure the inputs and outputs are the correct size. The Evaluate
function also instruments the runtime and number of EvaluateImpl
calls.
For example, to evaluate the model defined in the MyTrigonometricPiece class, we might use something like this:
Similarly, in python
Note: in Python, a list of input vectors must be used. We cannot pass each vector individually like we can in c++.
Definition at line 147 of file ModPiece.h.
Public Member Functions | |
ModPiece (std::vector< int > const &inputSizes, std::vector< int > const &outputSizes) | |
ModPiece (Eigen::VectorXi const &inputSizes, Eigen::VectorXi const &outputSizes) | |
virtual | ~ModPiece ()=default |
virtual double | GetRunTime (const std::string &method="Evaluate") const override |
Get the average run time for one of the implemented methods. More... | |
virtual void | ResetCallTime () override |
Resets the number of call and times. More... | |
virtual unsigned long int | GetNumCalls (const std::string &method="Evaluate") const override |
get the number of times one of the implemented methods has been called. More... | |
virtual std::vector< Eigen::VectorXd > const & | Evaluate (std::vector< Eigen::VectorXd > const &input) |
Evaluate the ModPiece. More... | |
virtual std::vector< Eigen::VectorXd > const & | Evaluate (ref_vector< Eigen::VectorXd > const &input) |
VARIADIC_TO_REFVECTOR (Evaluate, Eigen::VectorXd, std::vector< Eigen::VectorXd > const &) | |
virtual Eigen::VectorXd const & | Gradient (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity) |
Compute the Gradient \(J^Tv\). More... | |
virtual Eigen::VectorXd const & | Gradient (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity) |
Eigen::VectorXd const & | Gradient (unsigned int outWrt, unsigned int inWrt, Eigen::VectorXd const &last, Eigen::VectorXd const &sens) |
template<typename... Args> | |
Eigen::VectorXd const & | Gradient (unsigned int wrtOut, unsigned int wrtIn, Args const &... args) |
Eigen::VectorXd const & | ApplyHessian (unsigned int outWrt, unsigned int inWrt1, unsigned int inWrt2, Eigen::VectorXd const &last, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec) |
template<typename... Args> | |
Eigen::VectorXd const & | ApplyHessian (unsigned int wrtOut, unsigned int wrtIn1, unsigned int wrtIn2, Args const &... args) |
virtual Eigen::MatrixXd const & | Jacobian (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input) |
Compute the Jacobian of this ModPiece. More... | |
virtual Eigen::MatrixXd const & | Jacobian (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input) |
template<typename... Args> | |
Eigen::MatrixXd const & | Jacobian (unsigned int outWrt, unsigned int inWrt, Args const &... args) |
template<typename... Args> | |
Eigen::MatrixXd | JacobianByFD (unsigned int outWrt, unsigned int inWrt, Args const &... args) |
template<typename... Args> | |
Eigen::MatrixXd | ApplyJacobianByFD (unsigned int outWrt, unsigned int inWrt, Args const &... args) |
virtual Eigen::VectorXd const & | ApplyJacobian (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec) |
Apply the Jacobian of this ModPiece to a vector. More... | |
virtual Eigen::VectorXd const & | ApplyJacobian (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec) |
Eigen::VectorXd const & | ApplyJacobian (unsigned int outWrt, unsigned int inWrt, Eigen::VectorXd const &last, Eigen::VectorXd const &sens) |
template<typename... Args> | |
Eigen::VectorXd const & | ApplyJacobian (unsigned int wrtOut, unsigned int wrtIn, Args const &... args) |
virtual Eigen::VectorXd | GradientByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity) |
virtual Eigen::VectorXd | GradientByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sensitivity) |
virtual Eigen::MatrixXd | JacobianByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input) |
virtual Eigen::MatrixXd | JacobianByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input) |
virtual Eigen::VectorXd | ApplyJacobianByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, std::vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec) |
virtual Eigen::VectorXd | ApplyJacobianByFD (unsigned int const outputDimWrt, unsigned int const inputDimWrt, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &vec) |
virtual Eigen::VectorXd const & | ApplyHessian (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) |
virtual Eigen::VectorXd const & | ApplyHessian (unsigned int const outWrt, unsigned int const inWrt1, unsigned int const inWrt2, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec) |
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) |
virtual Eigen::VectorXd | ApplyHessianByFD (unsigned int const outWrt, unsigned int const inWrt1, unsigned int const inWrt2, ref_vector< Eigen::VectorXd > const &input, Eigen::VectorXd const &sens, Eigen::VectorXd const &vec) |
void | EnableCache () |
void | DisableCache () |
bool | CacheStatus () const |
virtual void | SetWarnLevel (unsigned int newLevel) |
Public Member Functions inherited from muq::Modeling::WorkPiece | |
WorkPiece () | |
Create a muq::Modeling::WorkPiece with no fixed number of inputs and outputs and variable input/output types. More... | |
WorkPiece (int const num, WorkPiece::Fix const fix=WorkPiece::Fix::Inputs) | |
Create a muq::Modeling::WorkPiece with either a fixed number of inputs or outputs and variable input/output types. More... | |
WorkPiece (int const numIns, int const numOuts) | |
Create a muq::Modeling::WorkPiece with a fixed number of inputs and outputs but variable input/output types. More... | |
WorkPiece (std::vector< std::string > const &types, WorkPiece::Fix const fix=WorkPiece::Fix::Inputs) | |
Create a muq::Modeling::WorkPiece with either a fixed number of inputs with specified types or a fixed number of outputs with specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &types, WorkPiece::Fix const fix=WorkPiece::Fix::Inputs) | |
Create a muq::Modeling::WorkPiece where either some of the inputs have specified types or some of the outputs have specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &types, int const num, WorkPiece::Fix const fixTypes=WorkPiece::Fix::Inputs, WorkPiece::Fix const fixNum=WorkPiece::Fix::Inputs) | |
Create a muq::Modeling::WorkPiece where either some of the inputs have specified types or some of the outputs have specified types and either the number of inputs or the number of outputs is fixed. More... | |
WorkPiece (std::vector< std::string > const &types, int const num) | |
Create a muq::Modeling::WorkPiece with a fixed number of inputs with specified types and a fixed number of outputs (of uknown type) More... | |
WorkPiece (int const num, std::vector< std::string > const &types) | |
Create a muq::Modeling::WorkPiece with a fixed number of outputs with specified types and a fixed number of inputs (of uknown type) More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, int const numIns, int const numOuts) | |
Create a muq::Modeling::WorkPiece where some of the inputs are known and we know the input and output numbers. More... | |
WorkPiece (int const numIns, std::map< unsigned int, std::string > const &outTypes, int const numOuts) | |
Create a muq::Modeling::WorkPiece where some of the outputs are known and we know the input and output numbers. More... | |
WorkPiece (std::vector< std::string > const &inTypes, std::vector< std::string > const &outTypes) | |
Create a muq::Modeling::WorkPiece with a fixed number of inputs and outputs with specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, std::vector< std::string > const &outTypes) | |
Create a muq::Modeling::WorkPiece where some of the inputs are known and all of the outputs have specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, int const num, std::vector< std::string > const &outTypes) | |
Create a muq::Modeling::WorkPiece where some of the inputs are known with a known number of inputs and all of the outputs have specified types. More... | |
WorkPiece (std::vector< std::string > const &inTypes, std::map< unsigned int, std::string > const &outTypes) | |
Create a muq::Modeling::WorkPiece where some of the outputs and all of the inputs have specified types. More... | |
WorkPiece (std::vector< std::string > const &inTypes, std::map< unsigned int, std::string > const &outTypes, int const num) | |
Create a muq::Modeling::WorkPiece where some of the outputs with a known number of outputs and all of the inputs have specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, std::map< unsigned int, std::string > const &outTypes) | |
Create a muq::Mdoeling::WorkPiece where some of the inputs and some of the outputs have specified types. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, int const numIn, std::map< unsigned int, std::string > const &outTypes) | |
Create a muq::Mdoeling::WorkPiece where some of the inputs and some of the outputs have specified types with a fixed number of inputs. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, std::map< unsigned int, std::string > const &outTypes, int const numOut) | |
Create a muq::Mdoeling::WorkPiece where some of the inputs and some of the outputs have specified types with a fixed number of outputs. More... | |
WorkPiece (std::map< unsigned int, std::string > const &inTypes, int const numIn, std::map< unsigned int, std::string > const &outTypes, int const numOut) | |
Create a muq::Mdoeling::WorkPiece where some of the inputs and some of the outputs have specified types with a fixed number of inputs and outputs. More... | |
virtual | ~WorkPiece () |
Default destructor. More... | |
std::vector< boost::any > const & | Evaluate (std::vector< boost::any > const &ins) |
Evaluate this muq::Modeling::WorkPiece. More... | |
std::vector< boost::any > const & | Evaluate (ref_vector< boost::any > const &ins) |
Evaluate this muq::Modeling::WorkPiece using references to the inputs. More... | |
std::vector< boost::any > const & | Evaluate () |
Evaluate this muq::Modeling::WorkPiece in the case that there are no inputs. More... | |
template<typename... Args> | |
std::vector< boost::any > const & | Evaluate (Args... args) |
Evalaute this muq::Modeling::WorkPiece using multiple arguments. More... | |
std::string const & | Name () |
Get the (unique) name of this work piece. More... | |
void | SetName (std::string const &newName) |
Set the name of this work piece. More... | |
std::string | InputType (unsigned int inputNum, bool const demangle=true) const |
Get the input type (if we know it) for a specific input. More... | |
int | InputSize (unsigned int inputNum) const |
Get the length of a vector valued input with fixed size. More... | |
void | SetInputSize (unsigned int inputNum, int newSize) |
std::string | OutputType (unsigned int outputNum, bool const demangle=true) const |
Get the output type (if we know it) for a specific output. More... | |
std::map< unsigned int, std::string > | OutputTypes () const |
Get the output types. More... | |
std::map< unsigned int, std::string > | InputTypes () const |
Get the input types. More... | |
unsigned int | ID () const |
Get the unique ID number. More... | |
Public Attributes | |
const Eigen::VectorXi | inputSizes |
const Eigen::VectorXi | outputSizes |
Public Attributes inherited from muq::Modeling::WorkPiece | |
int | numInputs |
The number of inputs. More... | |
int | numOutputs |
The number of outputs. More... | |
Additional Inherited Members | |
Static Public Member Functions inherited from muq::Modeling::WorkPiece | |
static ref_vector< const boost::any > | ToRefVector (std::vector< boost::any > const &anyVec) |
Create vector of references from a vector of boost::any's. More... | |
static ref_vector< const Eigen::VectorXd > | ToRefVector (std::vector< Eigen::VectorXd > const &anyVec) |
ModPiece::ModPiece | ( | std::vector< int > const & | inputSizes, |
std::vector< int > const & | outputSizes | ||
) |
Definition at line 9 of file ModPiece.cpp.
ModPiece::ModPiece | ( | Eigen::VectorXi const & | inputSizes, |
Eigen::VectorXi const & | outputSizes | ||
) |
Definition at line 14 of file ModPiece.cpp.
|
virtualdefault |
|
virtual |
Definition at line 368 of file ModPiece.cpp.
References ApplyHessianImpl(), muq::Utilities::demangle(), hessAction, hessActTime, inputSizes, numHessActCalls, and outputSizes.
|
virtual |
Assume we have a simple function given by the composition of two functions \(h(x) = f(g(x))\), where \(g:\mathbb{R}^N\rightarrow\mathbb{M}\) and \(f:\mathbb{R}^M \rightarrow \mathbb{R}$. The gradient of \)@_fakenlh \( with respect to \)x \( is given by the chain rule as \f[ \nabla_x h = \nabla_g f J_x(g), \f] where \)@_fakenlJ_x(g) \( denotes the Jacobian matrix of \)g \( with respect to \)x
Definition at line 358 of file ModPiece.cpp.
References ApplyHessian(), and muq::Modeling::WorkPiece::ToRefVector().
|
inline |
Definition at line 245 of file ModPiece.h.
References nlohmann::detail::last.
Referenced by ApplyHessian(), and ApplyHessianRecurse().
|
inline |
Definition at line 257 of file ModPiece.h.
References ApplyHessianRecurse().
|
virtual |
Definition at line 458 of file ModPiece.cpp.
References nlohmann::detail::dtoa_impl::e, Gradient(), inputSizes, and numHessActFDCalls.
|
virtual |
Applies the Hessian to a vector using finite differences. Note that only two evaluation are the gradient function are needed to compute the action of the Hessian on a vector.
Definition at line 448 of file ModPiece.cpp.
References muq::Modeling::WorkPiece::ToRefVector().
Referenced by muq::Optimization::CostFunction::ApplyHessian(), muq::Modeling::PyModPiece::ApplyHessianImpl(), ApplyHessianImpl(), and muq::Modeling::UMBridgeModPiece::ApplyHessianImpl().
|
protectedvirtual |
Reimplemented in muq::Modeling::CwiseUnaryOperator< T1, T2, T3 >, muq::Modeling::UMBridgeModPiece, muq::Optimization::CostFunction, muq::Modeling::SplitVector, muq::Modeling::SumPiece, muq::Modeling::ModGraphPiece, muq::Modeling::DensityBase, and muq::Modeling::PyModPiece.
Definition at line 436 of file ModPiece.cpp.
References ApplyHessianByFD(), and hessAction.
Referenced by ApplyHessian().
|
inlineprivate |
Definition at line 590 of file ModPiece.h.
Referenced by ApplyHessian().
|
inlineprivate |
Definition at line 607 of file ModPiece.h.
References ApplyHessian(), nlohmann::detail::last, and nlohmann::detail::dtoa_impl::value.
|
virtual |
Definition at line 161 of file ModPiece.cpp.
References ApplyJacobianImpl(), CheckInputs(), jacActTime, jacobianAction, and numJacActCalls.
|
virtual |
Apply the Jacobian of this ModPiece to a vector.
This function computes the matrix-vector product \(J_{ij} v \), where \(J_{ij}\) is the Jacobian matrix of output \(f_i\) with respect to input \(x_j\).
Like the Evaluate
function, this function performs some checks on the inputs and then calls the ApplyJacobianImpl
function. It also keeps track of how many times ApplyJacobian
is called and the runtime of ApplyJacobianImpl
.
If a child class does not override the ApplyJacobianImpl
function, finite differences are used by call the ApplyJacobianByFD
function. Computing \(J_{ij} v \) with finite differences requires two evaluations regardless of the dimension of \(v\).
[in] | outputDimWrt | The index \(i\) describing the output of interest. (Think of the derivative \(\partial f_i / \partial x_j\)) |
[in] | inputDimWrt | The index \(j\) describing which input we want to take derivatives with respect to. |
[in] | input | A vector of inputs corresponding to \(x_1,\ldots, x_{M_x}\) |
[in] | vec | The vector \(v\) that we want to multiply with the Jacobian matrix. |
Definition at line 153 of file ModPiece.cpp.
References muq::Modeling::WorkPiece::ToRefVector().
Referenced by ApplyJacobian(), and ApplyJacobianRecurse().
|
inline |
Definition at line 339 of file ModPiece.h.
References ApplyJacobian(), and nlohmann::detail::last.
|
inline |
Definition at line 345 of file ModPiece.h.
References ApplyJacobianRecurse().
|
virtual |
Definition at line 334 of file ModPiece.cpp.
References nlohmann::detail::dtoa_impl::e, muq::Modeling::WorkPiece::Evaluate(), and numJacActFDCalls.
|
virtual |
Definition at line 327 of file ModPiece.cpp.
References ApplyJacobianByFD(), and muq::Modeling::WorkPiece::ToRefVector().
|
inline |
Definition at line 303 of file ModPiece.h.
Referenced by ApplyJacobianByFD(), ApplyJacobianImpl(), muq::Modeling::PyModPiece::ApplyJacobianImpl(), and muq::Modeling::UMBridgeModPiece::ApplyJacobianImpl().
|
inlineprivate |
Definition at line 655 of file ModPiece.h.
References ApplyJacobianByFD().
|
inlineprivate |
Definition at line 660 of file ModPiece.h.
References ApplyJacobianByFD(), and nlohmann::detail::last.
|
protectedvirtual |
Reimplemented in muq::Modeling::ODE, muq::Modeling::CwiseUnaryOperator< T1, T2, T3 >, muq::Modeling::UMBridgeModPiece, muq::Modeling::SplitVector, muq::Modeling::ScaleVector, muq::Modeling::CombineVectors, muq::Modeling::SumPiece, muq::Modeling::PyModPiece, muq::Modeling::MultiLogisticLikelihood, muq::Modeling::ModGraphPiece, muq::Modeling::GradientPiece, muq::Modeling::DensityBase, and muq::Modeling::AffineOperator.
Definition at line 262 of file ModPiece.cpp.
References ApplyJacobianByFD(), muq::Utilities::demangle(), fdWarnLevel, and jacobianAction.
Referenced by ApplyJacobian().
|
inlineprivate |
Definition at line 631 of file ModPiece.h.
References nlohmann::detail::dtoa_impl::value.
Referenced by ApplyJacobian().
|
inlineprivate |
Definition at line 638 of file ModPiece.h.
References ApplyJacobian(), nlohmann::detail::last, and nlohmann::detail::dtoa_impl::value.
|
inline |
Definition at line 461 of file ModPiece.h.
References cacheEnabled.
|
protected |
Definition at line 211 of file ModPiece.cpp.
References muq::Utilities::demangle(), inputSizes, and nlohmann::to_string().
Referenced by ApplyJacobian(), Evaluate(), Gradient(), and Jacobian().
|
inline |
Definition at line 460 of file ModPiece.h.
References cacheEnabled.
|
inline |
The ModPiece class has support for caching one evaluation. If more extensive caching is needed, the muq::Modeling::FlannCache class can be used. The one step caching in the ModPiece class can be useful when subsequent calls to the ModPiece might have the exact same inputs. If enabled, both the input and output of EvaluateImpl calls will be stored. If the next call to Evaluate has the same input, then EvaluateImpl will not be called, the cached value will just be returned.
Note that caching introduces some additional overhead and is not always a good idea. It is therefore disabled by default and should be turned on manually for expensive ModPieces.
Enabling the one-step cache is particularly useful in gradient based optimization and MCMC methods, where subsequent calls to Evaluate and Gradient will occur with the same inputs. Without using caching, a call to Gradient in a ModGraphPiece could result in duplicate calls to Evaluate.
Definition at line 459 of file ModPiece.h.
References cacheEnabled.
|
virtual |
Definition at line 54 of file ModPiece.cpp.
References cacheEnabled, cacheInput, CheckInputs(), muq::Modeling::WorkPiece::evalTime, EvaluateImpl(), ExistsInCache(), muq::Modeling::WorkPiece::numEvalCalls, and outputs.
|
virtual |
Evaluate the ModPiece.
This function takes in a vector of inputs and returns a vector of outputs. It calls the EvaluateImpl
function after checking the size of the inputs. It also measures the run time of EvaluateImpl
and keeps track of the total number of Evaluate
calls.
[in] | input | A vector of inputs corresponding to \(x_1,\ldots, x_{M_x}\) |
Definition at line 22 of file ModPiece.cpp.
References muq::Modeling::WorkPiece::Evaluate(), and muq::Modeling::WorkPiece::ToRefVector().
|
overrideprotectedvirtual |
User-implemented function that determines the behavior of this muq::Modeling::WorkPiece.
This function determines how the WorkPiece::inputs determine WorkPiece::outputs. Must be implemented by a child.
WorkPiece::Evaluate() calls this function after checking the inputs and storing them in WorkPiece::inputs. This function populates WorkPiece::outputs, the outputs of this muq::Modeling::WorkPiece. WorkPiece::Evaluate() checks the outputs after calling this function.
Implements muq::Modeling::WorkPiece.
Definition at line 179 of file ModPiece.cpp.
References outputs, and muq::Modeling::WorkPiece::outputs.
Referenced by Evaluate().
|
protectedpure virtual |
Implemented in muq::Modeling::SplitVector, muq::Modeling::ScaleVector, muq::Modeling::ODE, muq::Modeling::ModGraphPiece, muq::Modeling::FlannCache, muq::Modeling::RandomVariable, muq::Modeling::DensityBase, muq::Modeling::CombineVectors, muq::Modeling::SumPiece, muq::Modeling::PyModPiece, muq::Modeling::OneStepCachePiece, muq::Modeling::JacobianPiece, muq::Modeling::GradientPiece, muq::Modeling::CwiseUnaryOperator< T1, T2, T3 >, muq::Modeling::UMBridgeModPiece, muq::Modeling::ReplicateOperator, muq::Modeling::ProductPiece, muq::Modeling::MultiLogisticLikelihood, muq::Modeling::ConstantVector, muq::Approximation::LocalRegression, muq::Approximation::BasisExpansion, muq::Optimization::CostFunction, muq::Modeling::AffineOperator, and muq::Approximation::Regression::PoisednessConstraint.
|
protected |
Checks to see if the input matches the value stored in the one-step cache.
Definition at line 27 of file ModPiece.cpp.
References cacheInput.
Referenced by Evaluate().
|
overridevirtual |
get the number of times one of the implemented methods has been called.
This function returns the number of times the EvaluateImpl, GradientImpl, JacobianImpl, ApplyJacobianImpl, or ApplyHessianImpl functions have been called.
[in] | method | The implemented function of interest. Possible options are "Evaluate", "Gradient", "Jacobian", "JacobianAction", "HessianAction", "GradientFD", "JacobianActionFD", "HessianActionFD" |
Reimplemented from muq::Modeling::WorkPiece.
Definition at line 535 of file ModPiece.cpp.
References muq::Modeling::WorkPiece::numEvalCalls, numGradCalls, numGradFDCalls, numHessActCalls, numHessActFDCalls, numJacActCalls, numJacActFDCalls, numJacCalls, and numJacFDCalls.
|
overridevirtual |
Get the average run time for one of the implemented methods.
This function returns the average wall clock time (in milliseconds) for the EvaluateImpl, GradientImpl, JacobianImpl, JacobianActionImpl, or HessianImp functions. If the function was never called, -1 is returned.
[in] | method | The implemented function of interest. Possible options are "Evaluate", "Gradient", "Jacobian", "JacobianAction", or "HessianAction" |
Reimplemented from muq::Modeling::WorkPiece.
Definition at line 497 of file ModPiece.cpp.
References nlohmann::detail::dtoa_impl::e, muq::Modeling::WorkPiece::evalTime, gradTime, hessActTime, jacActTime, jacTime, muq::Modeling::WorkPiece::numEvalCalls, numGradCalls, numHessActCalls, numJacActCalls, and numJacCalls.
|
virtual |
Definition at line 92 of file ModPiece.cpp.
References CheckInputs(), gradient, GradientImpl(), gradTime, and numGradCalls.
|
virtual |
Compute the Gradient \(J^Tv\).
Consider some scalar-valued function \(h(f_i) : \mathbb{R}^{N_{fi}} \rightarrow \mathbb{R}\) and let \(v = \nabla_{f_i} h\) denote the gradient of \(h\) with respect to \(f_i\). Let \([x_1, \ldots, x_{M_x}]\) be the inputs to \(f_i\), i.e.,
\[ h(f_i) = h( f_i(x_1, \ldots, x_{M_x}) ) \]
Given \(v = \nabla_{f_i} h\), this function computes \( \nabla_{x_j} h\) for some user-specified index \(j\). Thus, this function computes one step of the chain rule. If \(J_{ij}\) is the Jacobian matrix of \(f_i\) with respect to \(x_j\). Then,
\[ \nabla_{x_j} h = J_{ij}^T \left[ \nabla_{f_i} h \right]. \]
Like the Evaluate
function, this function performs some checks on the inputs and then calls the GradientImpl
function. It also keeps track of how many times Gradient
is called and the runtime of GradientImpl
.
If the GradientImpl function is not overriden by a child class, it will default to using finite differences by calling the GradientByFD function.
[in] | outputDimWrt | The index \(i\) describing the output of interest. (Think of the derivative \(\partial f_i / \partial x_j\)) |
[in] | inputDimWrt | The index \(j\) describing which input we want to take derivatives with respect to. |
[in] | input | A vector of inputs corresponding to \(x_1,\ldots, x_{M_x}\) |
[in] | sensitivity | The vector \(\nabla_{f_i} h\) containing the gradient of some scalar function \(h\) with respect to output \(i\) of this ModPiece. |
Definition at line 83 of file ModPiece.cpp.
References muq::Modeling::WorkPiece::ToRefVector().
Referenced by ApplyHessianByFD(), Gradient(), GradientRecurse(), and muq::Modeling::MultiLogisticLikelihood::JacobianImpl().
|
inline |
Definition at line 227 of file ModPiece.h.
References Gradient(), and nlohmann::detail::last.
|
inline |
Definition at line 237 of file ModPiece.h.
References GradientRecurse().
|
virtual |
Definition at line 119 of file ModPiece.cpp.
References JacobianByFD(), and numGradFDCalls.
|
virtual |
Definition at line 111 of file ModPiece.cpp.
References muq::Modeling::WorkPiece::ToRefVector().
Referenced by muq::Optimization::CostFunction::Gradient(), muq::Modeling::PyModPiece::GradientImpl(), and muq::Modeling::UMBridgeModPiece::GradientImpl().
|
protectedvirtual |
Reimplemented in muq::Modeling::ODE, muq::Modeling::CwiseUnaryOperator< T1, T2, T3 >, muq::Modeling::UMBridgeModPiece, muq::Approximation::Regression::PoisednessCost, muq::Optimization::CostFunction, muq::Modeling::SplitVector, muq::Modeling::ScaleVector, muq::Modeling::CombineVectors, muq::Modeling::SumPiece, muq::Modeling::PyModPiece, muq::Modeling::MultiLogisticLikelihood, muq::Modeling::ModGraphPiece, muq::Modeling::DensityBase, and muq::Modeling::AffineOperator.
Definition at line 237 of file ModPiece.cpp.
References gradient, and Jacobian().
Referenced by Gradient().
|
inlineprivate |
Definition at line 563 of file ModPiece.h.
Referenced by Gradient().
|
inlineprivate |
Definition at line 579 of file ModPiece.h.
References Gradient(), nlohmann::detail::last, and nlohmann::detail::dtoa_impl::value.
|
virtual |
Definition at line 135 of file ModPiece.cpp.
References CheckInputs(), jacobian, JacobianImpl(), jacTime, and numJacCalls.
|
virtual |
Compute the Jacobian of this ModPiece.
This function computes the Jacobian matrix of some output \(f_i(x_1, \ldots, x_{M_x})\) with respect to one of the inputs \(x_j\).
Like the Evaluate
function, this function performs some checks on the inputs and then calls the JacobianImpl
function. It also keeps track of how many times Jacobian
is called and the runtime of JacobianImpl
.
If a child class does not override the JacobianImpl
function, finite differences are used by call the JacobianByFD
function.
[in] | outputDimWrt | The index \(i\) describing the output of interest. (Think of the derivative \(\partial f_i / \partial x_j\)) |
[in] | inputDimWrt | The index \(j\) describing which input we want to take derivatives with respect to. |
[in] | input | A vector of inputs corresponding to \(x_1,\ldots, x_{M_x}\) |
Definition at line 128 of file ModPiece.cpp.
References muq::Modeling::WorkPiece::ToRefVector().
Referenced by muq::Modeling::MultiLogisticLikelihood::ApplyJacobianImpl(), GradientImpl(), and Jacobian().
|
inline |
Definition at line 291 of file ModPiece.h.
References Jacobian().
|
inlineprivate |
Definition at line 621 of file ModPiece.h.
References Jacobian().
|
inlineprivate |
Definition at line 625 of file ModPiece.h.
References Jacobian(), and nlohmann::detail::last.
|
virtual |
Definition at line 297 of file ModPiece.cpp.
References nlohmann::detail::dtoa_impl::e, muq::Modeling::WorkPiece::Evaluate(), inputSizes, numJacFDCalls, and outputSizes.
|
virtual |
Definition at line 290 of file ModPiece.cpp.
References JacobianByFD(), and muq::Modeling::WorkPiece::ToRefVector().
|
inline |
Definition at line 297 of file ModPiece.h.
Referenced by GradientByFD(), JacobianByFD(), JacobianImpl(), and muq::Modeling::PyModPiece::JacobianImpl().
|
inlineprivate |
Definition at line 645 of file ModPiece.h.
References JacobianByFD().
|
inlineprivate |
Definition at line 649 of file ModPiece.h.
References JacobianByFD(), and nlohmann::detail::last.
|
protectedvirtual |
Reimplemented in muq::Modeling::ODE, muq::Modeling::CwiseUnaryOperator< T1, T2, T3 >, muq::Optimization::CostFunction, muq::Approximation::BasisExpansion, muq::Modeling::SplitVector, muq::Modeling::ScaleVector, muq::Modeling::CombineVectors, muq::Modeling::SumPiece, muq::Modeling::PyModPiece, muq::Modeling::MultiLogisticLikelihood, muq::Modeling::ModGraphPiece, muq::Modeling::DensityBase, muq::Modeling::AffineOperator, and muq::Approximation::Regression::PoisednessConstraint.
Definition at line 245 of file ModPiece.cpp.
References muq::Utilities::demangle(), fdWarnLevel, jacobian, and JacobianByFD().
Referenced by Jacobian().
|
overridevirtual |
Resets the number of call and times.
Sets the number of calls to each implemented function to zero and the recorded wall clock times to zero.
Reimplemented from muq::Modeling::WorkPiece.
Definition at line 520 of file ModPiece.cpp.
References muq::Modeling::WorkPiece::evalTime, gradTime, hessActTime, jacActTime, jacTime, muq::Modeling::WorkPiece::numEvalCalls, numGradCalls, numHessActCalls, numJacActCalls, and numJacCalls.
|
inlinevirtual |
Set the warning level for finite difference defaults. If newLevel==0, no warnings or errors are thrown. If newLevel==1, a warning is printed to std::cout. If newLevel==2, an exception is thrown.
Definition at line 469 of file ModPiece.h.
References fdWarnLevel.
|
protected |
muq::Modeling::ModPiece::VARIADIC_TO_REFVECTOR | ( | Evaluate | , |
Eigen::VectorXd | , | ||
std::vector< Eigen::VectorXd > const & | |||
) |
|
protected |
Definition at line 479 of file ModPiece.h.
Referenced by muq::Modeling::ODE::ApplyJacobianImpl(), CacheStatus(), DisableCache(), EnableCache(), Evaluate(), muq::Modeling::ODE::GradientImpl(), and muq::Modeling::ODE::JacobianImpl().
|
protected |
Definition at line 480 of file ModPiece.h.
Referenced by muq::Modeling::ODE::ApplyJacobianImpl(), Evaluate(), ExistsInCache(), muq::Modeling::ODE::GradientImpl(), and muq::Modeling::ODE::JacobianImpl().
|
protected |
Specifies is MUQ should warn users when derivative functions (e.g, Jacobian, ApplyJacobian, Gradient, ApplyHessian) default to finite differences. If fdWarnLevel==0, no warnings or errors are thrown. If fdWarnLevel==1, a warning is printed to std::cout. If fdWarnLevel==2, an exception is thrown.
Definition at line 516 of file ModPiece.h.
Referenced by ApplyJacobianImpl(), JacobianImpl(), and SetWarnLevel().
|
protected |
Definition at line 504 of file ModPiece.h.
Referenced by Gradient(), muq::Modeling::AffineOperator::GradientImpl(), GradientImpl(), muq::Modeling::DensityBase::GradientImpl(), muq::Modeling::ModGraphPiece::GradientImpl(), muq::Modeling::MultiLogisticLikelihood::GradientImpl(), muq::Modeling::SumPiece::GradientImpl(), muq::Modeling::PyModPiece::GradientImpl(), muq::Modeling::CombineVectors::GradientImpl(), muq::Modeling::ScaleVector::GradientImpl(), muq::Modeling::SplitVector::GradientImpl(), muq::Modeling::UMBridgeModPiece::GradientImpl(), and muq::Modeling::ODE::GradientImpl().
|
protected |
Definition at line 498 of file ModPiece.h.
Referenced by GetRunTime(), Gradient(), and ResetCallTime().
|
protected |
Definition at line 507 of file ModPiece.h.
Referenced by ApplyHessian(), muq::Modeling::PyModPiece::ApplyHessianImpl(), ApplyHessianImpl(), muq::Modeling::DensityBase::ApplyHessianImpl(), muq::Modeling::ModGraphPiece::ApplyHessianImpl(), muq::Modeling::SumPiece::ApplyHessianImpl(), muq::Modeling::SplitVector::ApplyHessianImpl(), and muq::Modeling::UMBridgeModPiece::ApplyHessianImpl().
|
protected |
Definition at line 501 of file ModPiece.h.
Referenced by ApplyHessian(), GetRunTime(), and ResetCallTime().
const Eigen::VectorXi muq::Modeling::ModPiece::inputSizes |
Definition at line 471 of file ModPiece.h.
Referenced by ApplyHessian(), ApplyHessianByFD(), muq::Modeling::DensityBase::ApplyHessianImpl(), muq::Modeling::SplitVector::ApplyHessianImpl(), muq::Modeling::CombineVectors::ApplyJacobianImpl(), muq::Modeling::SplitVector::ApplyJacobianImpl(), CheckInputs(), muq::Approximation::PolynomialChaosExpansion::Covariance(), muq::Modeling::FlannCache::FlannCache(), muq::Approximation::BasisExpansion::FromHDF5(), muq::Modeling::GradientPiece::GetInputSizes(), muq::Modeling::JacobianPiece::GetInputSizes(), muq::Modeling::ModGraphPiece::GradientGraph(), muq::Modeling::MultiLogisticLikelihood::GradientImpl(), muq::Modeling::CombineVectors::GradientImpl(), muq::Modeling::SplitVector::GradientImpl(), muq::Modeling::ODE::GradientImpl(), JacobianByFD(), muq::Modeling::ModGraphPiece::JacobianGraph(), muq::Modeling::ModGraphPiece::JacobianImpl(), muq::Modeling::SumPiece::JacobianImpl(), muq::Modeling::CombineVectors::JacobianImpl(), muq::Modeling::ScaleVector::JacobianImpl(), muq::Modeling::SplitVector::JacobianImpl(), muq::Modeling::ODE::JacobianImpl(), muq::Approximation::PolynomialChaosExpansion::MainSensitivity(), muq::Modeling::ModGraphPiece::ModGraphPiece(), muq::Approximation::LocalRegression::SetUp(), muq::Approximation::PolynomialChaosExpansion::SobolSensitivity(), muq::Approximation::BasisExpansion::ToHDF5(), muq::Approximation::PolynomialChaosExpansion::TotalSensitivity(), and muq::Approximation::PolynomialChaosExpansion::Variance().
|
protected |
Definition at line 500 of file ModPiece.h.
Referenced by ApplyJacobian(), GetRunTime(), and ResetCallTime().
|
protected |
Definition at line 506 of file ModPiece.h.
Referenced by Jacobian(), muq::Modeling::AffineOperator::JacobianImpl(), JacobianImpl(), muq::Modeling::DensityBase::JacobianImpl(), muq::Modeling::ModGraphPiece::JacobianImpl(), muq::Modeling::MultiLogisticLikelihood::JacobianImpl(), muq::Modeling::SumPiece::JacobianImpl(), muq::Modeling::PyModPiece::JacobianImpl(), muq::Modeling::CombineVectors::JacobianImpl(), muq::Modeling::ScaleVector::JacobianImpl(), muq::Modeling::SplitVector::JacobianImpl(), muq::Approximation::BasisExpansion::JacobianImpl(), and muq::Modeling::ODE::JacobianImpl().
|
protected |
Definition at line 505 of file ModPiece.h.
Referenced by ApplyJacobian(), muq::Modeling::AffineOperator::ApplyJacobianImpl(), ApplyJacobianImpl(), muq::Modeling::DensityBase::ApplyJacobianImpl(), muq::Modeling::GradientPiece::ApplyJacobianImpl(), muq::Modeling::ModGraphPiece::ApplyJacobianImpl(), muq::Modeling::MultiLogisticLikelihood::ApplyJacobianImpl(), muq::Modeling::SumPiece::ApplyJacobianImpl(), muq::Modeling::PyModPiece::ApplyJacobianImpl(), muq::Modeling::CombineVectors::ApplyJacobianImpl(), muq::Modeling::ScaleVector::ApplyJacobianImpl(), muq::Modeling::SplitVector::ApplyJacobianImpl(), muq::Modeling::UMBridgeModPiece::ApplyJacobianImpl(), and muq::Modeling::ODE::ApplyJacobianImpl().
|
protected |
Definition at line 499 of file ModPiece.h.
Referenced by GetRunTime(), Jacobian(), and ResetCallTime().
|
protected |
Definition at line 487 of file ModPiece.h.
Referenced by GetNumCalls(), GetRunTime(), Gradient(), and ResetCallTime().
|
protected |
Definition at line 491 of file ModPiece.h.
Referenced by GetNumCalls(), and GradientByFD().
|
protected |
Definition at line 490 of file ModPiece.h.
Referenced by ApplyHessian(), GetNumCalls(), GetRunTime(), and ResetCallTime().
|
protected |
Definition at line 494 of file ModPiece.h.
Referenced by ApplyHessianByFD(), and GetNumCalls().
|
protected |
Definition at line 489 of file ModPiece.h.
Referenced by ApplyJacobian(), GetNumCalls(), GetRunTime(), and ResetCallTime().
|
protected |
Definition at line 493 of file ModPiece.h.
Referenced by ApplyJacobianByFD(), and GetNumCalls().
|
protected |
Definition at line 488 of file ModPiece.h.
Referenced by GetNumCalls(), GetRunTime(), Jacobian(), and ResetCallTime().
|
protected |
Definition at line 492 of file ModPiece.h.
Referenced by GetNumCalls(), and JacobianByFD().
|
protected |
Definition at line 503 of file ModPiece.h.
Referenced by muq::Modeling::ODE::ApplyJacobianImpl(), muq::Modeling::ConstantVector::ConstantVector(), Evaluate(), muq::Modeling::AffineOperator::EvaluateImpl(), muq::Approximation::BasisExpansion::EvaluateImpl(), muq::Approximation::LocalRegression::EvaluateImpl(), muq::Modeling::MultiLogisticLikelihood::EvaluateImpl(), muq::Modeling::ProductPiece::EvaluateImpl(), muq::Modeling::ReplicateOperator::EvaluateImpl(), muq::Modeling::UMBridgeModPiece::EvaluateImpl(), EvaluateImpl(), muq::Modeling::GradientPiece::EvaluateImpl(), muq::Modeling::JacobianPiece::EvaluateImpl(), muq::Modeling::OneStepCachePiece::EvaluateImpl(), muq::Modeling::SumPiece::EvaluateImpl(), muq::Modeling::CombineVectors::EvaluateImpl(), muq::Modeling::DensityBase::EvaluateImpl(), muq::Modeling::RandomVariable::EvaluateImpl(), muq::Modeling::FlannCache::EvaluateImpl(), muq::Modeling::ModGraphPiece::EvaluateImpl(), muq::Modeling::ODE::EvaluateImpl(), muq::Modeling::ScaleVector::EvaluateImpl(), muq::Modeling::SplitVector::EvaluateImpl(), muq::Modeling::ODE::GradientImpl(), muq::Modeling::ODE::JacobianImpl(), muq::Modeling::ModGraphPiece::MatchInputs(), muq::Modeling::ConstantVector::SetValue(), and muq::Modeling::UMBridgeModPiece::UMBridgeModPiece().
const Eigen::VectorXi muq::Modeling::ModPiece::outputSizes |
Definition at line 472 of file ModPiece.h.
Referenced by ApplyHessian(), muq::Modeling::CombineVectors::ApplyJacobianImpl(), muq::Approximation::PolynomialChaosExpansion::Covariance(), muq::Modeling::ReplicateOperator::EvaluateImpl(), muq::Modeling::CombineVectors::EvaluateImpl(), muq::Modeling::ODE::EvaluateImpl(), muq::Modeling::GradientPiece::GetOutputSizes(), muq::Modeling::JacobianPiece::GetOutputSizes(), muq::Modeling::ModGraphPiece::GradientGraph(), muq::Modeling::ODE::GradientImpl(), JacobianByFD(), muq::Modeling::ModGraphPiece::JacobianGraph(), muq::Modeling::CombineVectors::JacobianImpl(), muq::Modeling::ScaleVector::JacobianImpl(), muq::Modeling::ODE::JacobianImpl(), muq::Approximation::LocalRegression::SetUp(), and muq::Modeling::ConstantVector::SetValue().