MUQ  0.4.3
muq::Modeling::ModPiece Class Referenceabstract

Provides an abstract interface for defining vector-valued model components. More...

#include <ModPiece.h>

Inheritance diagram for muq::Modeling::ModPiece:

Detailed Description

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.

std::shared_ptr<ModPiece> m;
// ... create a child of ModPiece and store in m
int m_x = m->inputSizes.size(); // number of inputs
int m_f = m->outputSizes.size(); // number of outputs
int n_x1 = m->intputSizes(0); // The number of components in the first input vector
int n_f1 = m->outputSizes(0); // The number of components in the first output vector

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:

  • Combining existing ModPiece classes in a WorkGraph to create more complicated relationships.
  • Creating a new ModPiece by defining a class that inherits from ModPiece.

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]\).

class MyTrigonometricPiece : public muq::Modeling::ModPiece
{
public:
MyTrigonometricPiece() : muq::Modeling::ModPiece(2*Eigen::VectorXi::Ones(2), // inputSizes = [2,2]
2*Eigen::VectorXi::Ones(3)){}; // outputSizes = [2,2,2]
protected:
virtual void EvaluateImpl(muq::Modeling::ref_vector<Eigen::VectorXd> const& inputs) override
{
Eigen::VectorXd const& x1 = inputs.at(0);
Eigen::VectorXd const& x2 = inputs.at(1);
Eigen::VectorXd xprod = x1.array()*x2.array();
Eigen::VectorXd sinOut, cosOut, tanOut;
sinOut = xprod.array().sin();
cosOut = xprod.array().cos();
tanOut = xprod.array().tan();
// Resize the outputs vector (which lives in the ModPiece base class)
outputs.resize(3);
outputs.at(0) = sinOut;
outputs.at(1) = cosOut;
outputs.at(2) = tanOut;
}
};
Provides an abstract interface for defining vector-valued model components.
Definition: ModPiece.h:148
std::vector< Eigen::VectorXd > outputs
Definition: ModPiece.h:503
ModPiece(std::vector< int > const &inputSizes, std::vector< int > const &outputSizes)
Definition: ModPiece.cpp:9
virtual void EvaluateImpl(ref_vector< boost::any > const &inputs) override
User-implemented function that determines the behavior of this muq::Modeling::WorkPiece.
Definition: ModPiece.cpp:179
std::vector< std::reference_wrapper< const T > > ref_vector
A vector of references to something ...
Definition: WorkPiece.h:37

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

import pymuqModeling as mm
import numpy as np
class MyTrigonometricPiece(mm.PyModPiece):
def __init__(self):
mm.PyModPiece.__init__(self, [2,2], [2,2,2])
def EvaluateImpl(self, inputs):
x1 = inputs[0]
x2 = inputs[1]
xprod = x1*x2
sinOut = np.sin(xprod)
cosOut = np.cos(xprod)
tanOut = np.tan(xprod)
self.outputs = [sinOut, cosOut, tanOut]

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:

int main()
{
auto piece = std::make_shared<MyTrigonometricPiece>();
std::vector<Eigen::VectorXd> inputs(2);
inputs.at(0) = Eigen::VectorXd::Random(2);
inputs.at(1) = Eigen::VectorXd::Random(2);
std::vector<Eigen::VectorXd> outputs = piece->Evaluate(inputs);
// analogously, we can pass in individual Eigen::VectorXd's instead of the std::vector
// outputs = piece->Evaluate(inputs.at(0), inputs.at(1));
return 0;
};
int main()

Similarly, in python

piece = MyTrigonometricPiece()
inputs = [ np.random.rand(2), np.random.rand(2) ]
outputs = piece.Evaluate(inputs)

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)
 

Constructor & Destructor Documentation

◆ ModPiece() [1/2]

ModPiece::ModPiece ( std::vector< int > const &  inputSizes,
std::vector< int > const &  outputSizes 
)

Definition at line 9 of file ModPiece.cpp.

◆ ModPiece() [2/2]

ModPiece::ModPiece ( Eigen::VectorXi const &  inputSizes,
Eigen::VectorXi const &  outputSizes 
)

Definition at line 14 of file ModPiece.cpp.

◆ ~ModPiece()

virtual muq::Modeling::ModPiece::~ModPiece ( )
virtualdefault

Member Function Documentation

◆ ApplyHessian() [1/4]

Eigen::VectorXd const & ModPiece::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

◆ ApplyHessian() [2/4]

Eigen::VectorXd const & ModPiece::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

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().

◆ ApplyHessian() [3/4]

Eigen::VectorXd const& muq::Modeling::ModPiece::ApplyHessian ( unsigned int  outWrt,
unsigned int  inWrt1,
unsigned int  inWrt2,
Eigen::VectorXd const &  last,
Eigen::VectorXd const &  sens,
Eigen::VectorXd const &  vec 
)
inline

Definition at line 245 of file ModPiece.h.

References nlohmann::detail::last.

Referenced by ApplyHessian(), and ApplyHessianRecurse().

◆ ApplyHessian() [4/4]

template<typename... Args>
Eigen::VectorXd const& muq::Modeling::ModPiece::ApplyHessian ( unsigned int  wrtOut,
unsigned int  wrtIn1,
unsigned int  wrtIn2,
Args const &...  args 
)
inline

Definition at line 257 of file ModPiece.h.

References ApplyHessianRecurse().

◆ ApplyHessianByFD() [1/2]

Eigen::VectorXd ModPiece::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 
)
virtual

◆ ApplyHessianByFD() [2/2]

Eigen::VectorXd ModPiece::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

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().

◆ ApplyHessianImpl()

void ModPiece::ApplyHessianImpl ( 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 
)
protectedvirtual

◆ ApplyHessianRecurse() [1/2]

template<typename NextType , typename... Args>
Eigen::VectorXd const& muq::Modeling::ModPiece::ApplyHessianRecurse ( unsigned int  outWrt,
unsigned int  inWrt1,
unsigned int  inWrt2,
ref_vector< Eigen::VectorXd > &  invec,
NextType const &  ith,
Args const &...  args 
)
inlineprivate

Definition at line 590 of file ModPiece.h.

Referenced by ApplyHessian().

◆ ApplyHessianRecurse() [2/2]

template<typename NextType >
Eigen::VectorXd const& muq::Modeling::ModPiece::ApplyHessianRecurse ( unsigned int  outWrt,
unsigned int  inWrt1,
unsigned int  inWrt2,
ref_vector< Eigen::VectorXd > &  invec,
NextType const &  last,
Eigen::VectorXd const &  sens,
Eigen::VectorXd const &  vec 
)
inlineprivate

◆ ApplyJacobian() [1/4]

Eigen::VectorXd const & ModPiece::ApplyJacobian ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
ref_vector< Eigen::VectorXd > const &  input,
Eigen::VectorXd const &  vec 
)
virtual

◆ ApplyJacobian() [2/4]

Eigen::VectorXd const & ModPiece::ApplyJacobian ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
std::vector< Eigen::VectorXd > const &  input,
Eigen::VectorXd const &  vec 
)
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\).

Parameters
[in]outputDimWrtThe index \(i\) describing the output of interest. (Think of the derivative \(\partial f_i / \partial x_j\))
[in]inputDimWrtThe index \(j\) describing which input we want to take derivatives with respect to.
[in]inputA vector of inputs corresponding to \(x_1,\ldots, x_{M_x}\)
[in]vecThe vector \(v\) that we want to multiply with the Jacobian matrix.
Returns
An Eigen::VectorXd containing the matrix-vector product \( J_{ij} v\).

Definition at line 153 of file ModPiece.cpp.

References muq::Modeling::WorkPiece::ToRefVector().

Referenced by ApplyJacobian(), and ApplyJacobianRecurse().

◆ ApplyJacobian() [3/4]

Eigen::VectorXd const& muq::Modeling::ModPiece::ApplyJacobian ( unsigned int  outWrt,
unsigned int  inWrt,
Eigen::VectorXd const &  last,
Eigen::VectorXd const &  sens 
)
inline

Definition at line 339 of file ModPiece.h.

References ApplyJacobian(), and nlohmann::detail::last.

◆ ApplyJacobian() [4/4]

template<typename... Args>
Eigen::VectorXd const& muq::Modeling::ModPiece::ApplyJacobian ( unsigned int  wrtOut,
unsigned int  wrtIn,
Args const &...  args 
)
inline

Definition at line 345 of file ModPiece.h.

References ApplyJacobianRecurse().

◆ ApplyJacobianByFD() [1/5]

Eigen::VectorXd ModPiece::ApplyJacobianByFD ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
ref_vector< Eigen::VectorXd > const &  input,
Eigen::VectorXd const &  vec 
)
virtual

◆ ApplyJacobianByFD() [2/5]

Eigen::VectorXd ModPiece::ApplyJacobianByFD ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
std::vector< Eigen::VectorXd > const &  input,
Eigen::VectorXd const &  vec 
)
virtual

Definition at line 327 of file ModPiece.cpp.

References ApplyJacobianByFD(), and muq::Modeling::WorkPiece::ToRefVector().

◆ ApplyJacobianByFD() [3/5]

template<typename... Args>
Eigen::MatrixXd muq::Modeling::ModPiece::ApplyJacobianByFD ( unsigned int  outWrt,
unsigned int  inWrt,
Args const &...  args 
)
inline

◆ ApplyJacobianByFD() [4/5]

template<typename NextType , typename... Args>
Eigen::MatrixXd muq::Modeling::ModPiece::ApplyJacobianByFD ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
NextType const &  ith,
Args const &...  args 
)
inlineprivate

Definition at line 655 of file ModPiece.h.

References ApplyJacobianByFD().

◆ ApplyJacobianByFD() [5/5]

template<typename NextType >
Eigen::MatrixXd muq::Modeling::ModPiece::ApplyJacobianByFD ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
NextType const &  last,
Eigen::VectorXd const &  sens 
)
inlineprivate

Definition at line 660 of file ModPiece.h.

References ApplyJacobianByFD(), and nlohmann::detail::last.

◆ ApplyJacobianImpl()

◆ ApplyJacobianRecurse() [1/2]

template<typename NextType , typename... Args>
Eigen::VectorXd const& muq::Modeling::ModPiece::ApplyJacobianRecurse ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
NextType const &  ith,
Args const &...  args 
)
inlineprivate

Definition at line 631 of file ModPiece.h.

References nlohmann::detail::dtoa_impl::value.

Referenced by ApplyJacobian().

◆ ApplyJacobianRecurse() [2/2]

template<typename NextType >
Eigen::VectorXd const& muq::Modeling::ModPiece::ApplyJacobianRecurse ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
NextType const &  last,
Eigen::VectorXd const &  sens 
)
inlineprivate

◆ CacheStatus()

bool muq::Modeling::ModPiece::CacheStatus ( ) const
inline

Definition at line 461 of file ModPiece.h.

References cacheEnabled.

◆ CheckInputs()

void ModPiece::CheckInputs ( ref_vector< Eigen::VectorXd > const &  input,
std::string const &  funcName 
)
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().

◆ DisableCache()

void muq::Modeling::ModPiece::DisableCache ( )
inline

Definition at line 460 of file ModPiece.h.

References cacheEnabled.

◆ EnableCache()

void muq::Modeling::ModPiece::EnableCache ( )
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.

◆ Evaluate() [1/2]

std::vector< Eigen::VectorXd > const & ModPiece::Evaluate ( ref_vector< Eigen::VectorXd > const &  input)
virtual

◆ Evaluate() [2/2]

std::vector< Eigen::VectorXd > const & ModPiece::Evaluate ( std::vector< Eigen::VectorXd > const &  input)
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.

Parameters
[in]inputA vector of inputs corresponding to \(x_1,\ldots, x_{M_x}\)
Returns
A vector of outputs corresponding to \(f_1(x_1, \ldots, x_N), \ldots, f_{M_f}(x_1, \ldots, x_{M_{x}})\)

Definition at line 22 of file ModPiece.cpp.

References muq::Modeling::WorkPiece::Evaluate(), and muq::Modeling::WorkPiece::ToRefVector().

◆ EvaluateImpl() [1/2]

void ModPiece::EvaluateImpl ( ref_vector< boost::any > const &  inputs)
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().

◆ EvaluateImpl() [2/2]

◆ ExistsInCache()

bool ModPiece::ExistsInCache ( ref_vector< Eigen::VectorXd > const &  input) const
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().

◆ GetNumCalls()

unsigned long int ModPiece::GetNumCalls ( const std::string &  method = "Evaluate") const
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.

Parameters
[in]methodThe implemented function of interest. Possible options are "Evaluate", "Gradient", "Jacobian", "JacobianAction", "HessianAction", "GradientFD", "JacobianActionFD", "HessianActionFD"
Returns
An integer with the number of calls.

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.

◆ GetRunTime()

double ModPiece::GetRunTime ( const std::string &  method = "Evaluate") const
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.

Parameters
[in]methodThe implemented function of interest. Possible options are "Evaluate", "Gradient", "Jacobian", "JacobianAction", or "HessianAction"
Returns
The average wall clock time in milli-seconds.

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.

◆ Gradient() [1/4]

Eigen::VectorXd const & ModPiece::Gradient ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
ref_vector< Eigen::VectorXd > const &  input,
Eigen::VectorXd const &  sensitivity 
)
virtual

Definition at line 92 of file ModPiece.cpp.

References CheckInputs(), gradient, GradientImpl(), gradTime, and numGradCalls.

◆ Gradient() [2/4]

Eigen::VectorXd const & ModPiece::Gradient ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
std::vector< Eigen::VectorXd > const &  input,
Eigen::VectorXd const &  sensitivity 
)
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.

Parameters
[in]outputDimWrtThe index \(i\) describing the output of interest. (Think of the derivative \(\partial f_i / \partial x_j\))
[in]inputDimWrtThe index \(j\) describing which input we want to take derivatives with respect to.
[in]inputA vector of inputs corresponding to \(x_1,\ldots, x_{M_x}\)
[in]sensitivityThe vector \(\nabla_{f_i} h\) containing the gradient of some scalar function \(h\) with respect to output \(i\) of this ModPiece.
Returns
An Eigen::VectorXd containing the new gradient \(\nabla_{x_j} h\).

Definition at line 83 of file ModPiece.cpp.

References muq::Modeling::WorkPiece::ToRefVector().

Referenced by ApplyHessianByFD(), Gradient(), GradientRecurse(), and muq::Modeling::MultiLogisticLikelihood::JacobianImpl().

◆ Gradient() [3/4]

Eigen::VectorXd const& muq::Modeling::ModPiece::Gradient ( unsigned int  outWrt,
unsigned int  inWrt,
Eigen::VectorXd const &  last,
Eigen::VectorXd const &  sens 
)
inline

Definition at line 227 of file ModPiece.h.

References Gradient(), and nlohmann::detail::last.

◆ Gradient() [4/4]

template<typename... Args>
Eigen::VectorXd const& muq::Modeling::ModPiece::Gradient ( unsigned int  wrtOut,
unsigned int  wrtIn,
Args const &...  args 
)
inline

Definition at line 237 of file ModPiece.h.

References GradientRecurse().

◆ GradientByFD() [1/2]

Eigen::VectorXd ModPiece::GradientByFD ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
ref_vector< Eigen::VectorXd > const &  input,
Eigen::VectorXd const &  sensitivity 
)
virtual

Definition at line 119 of file ModPiece.cpp.

References JacobianByFD(), and numGradFDCalls.

◆ GradientByFD() [2/2]

Eigen::VectorXd ModPiece::GradientByFD ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
std::vector< Eigen::VectorXd > const &  input,
Eigen::VectorXd const &  sensitivity 
)
virtual

◆ GradientImpl()

◆ GradientRecurse() [1/2]

template<typename NextType , typename... Args>
Eigen::VectorXd const& muq::Modeling::ModPiece::GradientRecurse ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
NextType const &  ith,
Args const &...  args 
)
inlineprivate

Definition at line 563 of file ModPiece.h.

Referenced by Gradient().

◆ GradientRecurse() [2/2]

template<typename NextType >
Eigen::VectorXd const& muq::Modeling::ModPiece::GradientRecurse ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
NextType const &  last,
Eigen::VectorXd const &  sens 
)
inlineprivate

◆ Jacobian() [1/5]

Eigen::MatrixXd const & ModPiece::Jacobian ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
ref_vector< Eigen::VectorXd > const &  input 
)
virtual

Definition at line 135 of file ModPiece.cpp.

References CheckInputs(), jacobian, JacobianImpl(), jacTime, and numJacCalls.

◆ Jacobian() [2/5]

Eigen::MatrixXd const & ModPiece::Jacobian ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
std::vector< Eigen::VectorXd > const &  input 
)
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.

Parameters
[in]outputDimWrtThe index \(i\) describing the output of interest. (Think of the derivative \(\partial f_i / \partial x_j\))
[in]inputDimWrtThe index \(j\) describing which input we want to take derivatives with respect to.
[in]inputA vector of inputs corresponding to \(x_1,\ldots, x_{M_x}\)
Returns
An Eigen::MatrixXd containing the new Jacobian matrix \( J_{ij}\).

Definition at line 128 of file ModPiece.cpp.

References muq::Modeling::WorkPiece::ToRefVector().

Referenced by muq::Modeling::MultiLogisticLikelihood::ApplyJacobianImpl(), GradientImpl(), and Jacobian().

◆ Jacobian() [3/5]

template<typename... Args>
Eigen::MatrixXd const& muq::Modeling::ModPiece::Jacobian ( unsigned int  outWrt,
unsigned int  inWrt,
Args const &...  args 
)
inline

Definition at line 291 of file ModPiece.h.

References Jacobian().

◆ Jacobian() [4/5]

template<typename... Args>
Eigen::MatrixXd const& muq::Modeling::ModPiece::Jacobian ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
Eigen::VectorXd const &  ith,
Args const &...  args 
)
inlineprivate

Definition at line 621 of file ModPiece.h.

References Jacobian().

◆ Jacobian() [5/5]

Eigen::MatrixXd const& muq::Modeling::ModPiece::Jacobian ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
Eigen::VectorXd const &  last 
)
inlineprivate

Definition at line 625 of file ModPiece.h.

References Jacobian(), and nlohmann::detail::last.

◆ JacobianByFD() [1/5]

Eigen::MatrixXd ModPiece::JacobianByFD ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
ref_vector< Eigen::VectorXd > const &  input 
)
virtual

◆ JacobianByFD() [2/5]

Eigen::MatrixXd ModPiece::JacobianByFD ( unsigned int const  outputDimWrt,
unsigned int const  inputDimWrt,
std::vector< Eigen::VectorXd > const &  input 
)
virtual

Definition at line 290 of file ModPiece.cpp.

References JacobianByFD(), and muq::Modeling::WorkPiece::ToRefVector().

◆ JacobianByFD() [3/5]

template<typename... Args>
Eigen::MatrixXd muq::Modeling::ModPiece::JacobianByFD ( unsigned int  outWrt,
unsigned int  inWrt,
Args const &...  args 
)
inline

◆ JacobianByFD() [4/5]

template<typename... Args>
Eigen::MatrixXd muq::Modeling::ModPiece::JacobianByFD ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
Eigen::VectorXd const &  ith,
Args const &...  args 
)
inlineprivate

Definition at line 645 of file ModPiece.h.

References JacobianByFD().

◆ JacobianByFD() [5/5]

Eigen::MatrixXd muq::Modeling::ModPiece::JacobianByFD ( unsigned int  outWrt,
unsigned int  inWrt,
ref_vector< Eigen::VectorXd > &  vec,
Eigen::VectorXd const &  last 
)
inlineprivate

Definition at line 649 of file ModPiece.h.

References JacobianByFD(), and nlohmann::detail::last.

◆ JacobianImpl()

◆ ResetCallTime()

void ModPiece::ResetCallTime ( )
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.

◆ SetWarnLevel()

virtual void muq::Modeling::ModPiece::SetWarnLevel ( unsigned int  newLevel)
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.

◆ ToStdVec()

◆ VARIADIC_TO_REFVECTOR()

muq::Modeling::ModPiece::VARIADIC_TO_REFVECTOR ( Evaluate  ,
Eigen::VectorXd  ,
std::vector< Eigen::VectorXd > const &   
)

Member Data Documentation

◆ cacheEnabled

bool muq::Modeling::ModPiece::cacheEnabled = false
protected

◆ cacheInput

std::vector<Eigen::VectorXd> muq::Modeling::ModPiece::cacheInput
protected

◆ fdWarnLevel

unsigned int muq::Modeling::ModPiece::fdWarnLevel = 0
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().

◆ gradient

◆ gradTime

double muq::Modeling::ModPiece::gradTime = 0
protected

Definition at line 498 of file ModPiece.h.

Referenced by GetRunTime(), Gradient(), and ResetCallTime().

◆ hessAction

◆ hessActTime

double muq::Modeling::ModPiece::hessActTime = 0
protected

Definition at line 501 of file ModPiece.h.

Referenced by ApplyHessian(), GetRunTime(), and ResetCallTime().

◆ inputSizes

const Eigen::VectorXi muq::Modeling::ModPiece::inputSizes

◆ jacActTime

double muq::Modeling::ModPiece::jacActTime = 0
protected

Definition at line 500 of file ModPiece.h.

Referenced by ApplyJacobian(), GetRunTime(), and ResetCallTime().

◆ jacobian

◆ jacobianAction

◆ jacTime

double muq::Modeling::ModPiece::jacTime = 0
protected

Definition at line 499 of file ModPiece.h.

Referenced by GetRunTime(), Jacobian(), and ResetCallTime().

◆ numGradCalls

unsigned long int muq::Modeling::ModPiece::numGradCalls = 0
protected

Definition at line 487 of file ModPiece.h.

Referenced by GetNumCalls(), GetRunTime(), Gradient(), and ResetCallTime().

◆ numGradFDCalls

unsigned long int muq::Modeling::ModPiece::numGradFDCalls = 0
protected

Definition at line 491 of file ModPiece.h.

Referenced by GetNumCalls(), and GradientByFD().

◆ numHessActCalls

unsigned long int muq::Modeling::ModPiece::numHessActCalls = 0
protected

Definition at line 490 of file ModPiece.h.

Referenced by ApplyHessian(), GetNumCalls(), GetRunTime(), and ResetCallTime().

◆ numHessActFDCalls

unsigned long int muq::Modeling::ModPiece::numHessActFDCalls = 0
protected

Definition at line 494 of file ModPiece.h.

Referenced by ApplyHessianByFD(), and GetNumCalls().

◆ numJacActCalls

unsigned long int muq::Modeling::ModPiece::numJacActCalls = 0
protected

Definition at line 489 of file ModPiece.h.

Referenced by ApplyJacobian(), GetNumCalls(), GetRunTime(), and ResetCallTime().

◆ numJacActFDCalls

unsigned long int muq::Modeling::ModPiece::numJacActFDCalls = 0
protected

Definition at line 493 of file ModPiece.h.

Referenced by ApplyJacobianByFD(), and GetNumCalls().

◆ numJacCalls

unsigned long int muq::Modeling::ModPiece::numJacCalls = 0
protected

Definition at line 488 of file ModPiece.h.

Referenced by GetNumCalls(), GetRunTime(), Jacobian(), and ResetCallTime().

◆ numJacFDCalls

unsigned long int muq::Modeling::ModPiece::numJacFDCalls = 0
protected

Definition at line 492 of file ModPiece.h.

Referenced by GetNumCalls(), and JacobianByFD().

◆ outputs

◆ outputSizes


The documentation for this class was generated from the following files: