Conditional Map Base#

template<typename MemorySpace>
class ConditionalMapBase : public mpart::ParameterizedFunctionBase<MemorySpace>#

Provides an abstract base class for conditional transport maps where the input dimension might be larger than output dimension.

This class provides an interface for functions \(T:\mathbb{R}^{N+M}\rightarrow \mathbb{R}^M\), where $N\geq 0$. Let \(x_1\in \mathbb{R}^N\) denote the first block of inputs and \(x_2\in\mathbb{R}^M\) denote the second block of inputs. Let \(r\in\mathbb{R}^M\) denote the map output, \(r=T(x_2; x_1)\). The conditional maps implemented by children of this class guarantee that for fixed \(x_1\), the conditional mapping from \(x_1 \rightarrow r\) is invertible and the Jacobian matrix with respect to \(x_2\), \(\nabla_{x_2} T\), is positive definite.

Subclassed by mpart::AffineMap< MemorySpace >, mpart::ComposedMap< MemorySpace >, mpart::IdentityMap< MemorySpace >, mpart::InnerMarginalAffineMap< MemorySpace >, mpart::MonotoneComponent< ExpansionType, PosFuncType, QuadratureType, MemorySpace >, mpart::RectifiedMultivariateExpansion< MemorySpace, OffdiagEval, DiagEval, Rectifier >, mpart::SummarizedMap< MemorySpace >, mpart::TriangularMap< MemorySpace >, mpart::UnivariateExpansion< MemorySpace, BasisType >

Public Functions

inline ConditionalMapBase(unsigned int inDim, unsigned int outDim, unsigned int nCoeffs)#

Construct a new Conditional Map Base object.

Typically the ConditionalMapBase class is not constructed directly. The preferred way of creating a ConditionalMapBase object is through factory methods in the mpart::MapFactory namespace.

Parameters:
  • inDim – The dimension \(N\) of the input to this map.

  • outDim – The dimension \(M\) of the output from this map.

  • nCoeffs – The number of coefficients in the map parameterization.

inline ConditionalMapBase(unsigned int inDim, unsigned int outDim, unsigned int nCoeffs, Kokkos::View<const double*, MemorySpace> coeffsIn)#
virtual ~ConditionalMapBase() = default#
inline virtual std::shared_ptr<ParameterizedFunctionBase<MemorySpace>> GetBaseFunction()#

For Monotone parameterizations that are based on a non-monotone base function, this function will return the base function. If the monotone parameterization is not constructed from a non-monotone base, then this function will return a nullptr.

template<typename AnyMemorySpace>
Kokkos::View<double*, AnyMemorySpace> LogDeterminant(StridedMatrix<const double, AnyMemorySpace> const &pts)#

Computes the log determinant of the map Jacobian. For a map \(T:\mathbb{R}^N\rightarrow \mathbb{R}^M\) with \(M\leq N\) and components \(T_i(x_{1:N-M+i})\), this function computes the determinant of the Jacobian of \(T\) with respect to \(x_{N-M:N}\). While the map is rectangular, the Jacobian with respect to these inputs will be square. The fact that the map is lower triangular will then imply that the determinant is given by the product of diagonal derviatives.

\[ \det{\nabla_{x_{N-M:N}} T} = \prod_{i=1}^M \frac{\partial T_i}{\partial x_{N-M+i}}. \]

Parameters:

pts – The points where we want to evaluate the log determinant.

template<typename ...AllTraits>
inline Kokkos::View<double*, typename Kokkos::View<double**, AllTraits...>::memory_space> LogDeterminant(Kokkos::View<double**, AllTraits...> pts)#

LogDeterminant function with conversion between default view layout and const strided matrix.

template<typename AnyMemorySpace>
inline Kokkos::View<double*, AnyMemorySpace> LogDeterminant(StridedMatrix<double, AnyMemorySpace> const &pts)#

LogDeterminant function with conversion from regular strided matrix to const strided matrix.

Eigen::VectorXd LogDeterminant(Eigen::Ref<const Eigen::RowMatrixXd> const &pts)#

LogDeterminant function with conversion from Eigen to Kokkos (and possibly copy to/from device.)

virtual void LogDeterminantImpl(StridedMatrix<const double, MemorySpace> const &pts, StridedVector<double, MemorySpace> output) = 0#
template<typename AnyMemorySpace>
StridedMatrix<double, AnyMemorySpace> Inverse(StridedMatrix<const double, AnyMemorySpace> const &x1, StridedMatrix<const double, AnyMemorySpace> const &r)#

Computes the inverse of the map.

The ConditionalMapBase class defines invertible functions \(T:\mathbb{R}^{d_{in}}\rightarrow\mathbb{R}^{d_{out}}\) where the input dimension \(d_{in}\geq d_{out}\). An input \(x\in\mathbb{R}^{d_{in}}\) can then be split into two parts \(x=[x_1,x_2]\), where \(x_1\in\mathbb{R}^{d_{in}-d_{out}}\) represents the “extra” inputs and \(x_2\in\mathbb{R}^{d_{out}}\) has the same size as the map output. This function computes the second block \(x_2\) of the input given the first block \(x_1\) and a vector \(r=T(x_1,x_2)\).

Parameters:
  • x1 – A \(d_{in}-d_{out}\times N\) or \(d_{in}\times N\) matrix containing \(N\) values of the first input block: \(\{x_1^{(1)},\ldots x_1^{(N)}\}\). If \(d_{in}\) rows are given, then the last \(d_{out}\) rows may be used as an initial guess for \(x_2\) by certain solvers.

  • r – A \(d_{out}\times N\) matrix containing \(N\) values of the map output: \(\{r^{(1)},\ldots, r^{(N)}\}\).

Returns:

A \(d_{out} \times N\) matrix containing the computed values of \(\{x_2^{(1)},\ldots,x_2^{(N)}\}\).

template<typename ViewType1, typename ViewType2>
inline StridedMatrix<double, typename ViewType1::memory_space> Inverse(ViewType1 x, ViewType2 r)#

Inverse function with conversion between general view type and const strided matrix.

Eigen::RowMatrixXd Inverse(Eigen::Ref<const Eigen::RowMatrixXd> const &x1, Eigen::Ref<const Eigen::RowMatrixXd> const &r)#

Inverse function with conversion between eigen matrix and Kokkos view.

virtual void InverseImpl(StridedMatrix<const double, MemorySpace> const &x1, StridedMatrix<const double, MemorySpace> const &r, StridedMatrix<double, MemorySpace> output) = 0#

Pure abstract function overridden by child classes.

template<typename AnyMemorySpace>
StridedMatrix<double, AnyMemorySpace> LogDeterminantCoeffGrad(StridedMatrix<const double, AnyMemorySpace> const &pts)#

Computes the gradient of the log determinant with respect to the map coefficients.

For a map \(T(x; w) : \mathbb{R}^N \rightarrow \mathbb{R}^M\) parameterized by coefficients \(w\in\mathbb{R}^K\), this function computes

\[ \nabla_w \det{\nabla_x T(x_i; w)}, \]
at multiple points \(x_i\).

Parameters:

pts – A collection of points where we want to evaluate the gradient. Each column corresponds to a point.

Returns:

A matrix containing the coefficient gradient at each input point. The \(i^{th}\) column contains \(\nabla_w \det{\nabla_x T(x_i; w)}\).

template<typename ViewType>
inline StridedMatrix<double, typename ViewType::memory_space> LogDeterminantCoeffGrad(ViewType pts)#

Include conversion between general view type and Strided matrix.

Eigen::RowMatrixXd LogDeterminantCoeffGrad(Eigen::Ref<const Eigen::RowMatrixXd> const &pts)#

Evaluation with additional conversion of Eigen matrix to Kokkos unmanaged view.

virtual void LogDeterminantCoeffGradImpl(StridedMatrix<const double, MemorySpace> const &pts, StridedMatrix<double, MemorySpace> output) = 0#
template<typename AnyMemorySpace>
StridedMatrix<double, AnyMemorySpace> LogDeterminantInputGrad(StridedMatrix<const double, AnyMemorySpace> const &pts)#
template<typename ViewType>
inline StridedMatrix<double, typename ViewType::memory_space> LogDeterminantInputGrad(ViewType pts)#

Include conversion between general view type and Strided matrix.

Eigen::RowMatrixXd LogDeterminantInputGrad(Eigen::Ref<const Eigen::RowMatrixXd> const &pts)#

Evaluation with additional conversion of Eigen matrix to Kokkos unmanaged view.

virtual void LogDeterminantInputGradImpl(StridedMatrix<const double, MemorySpace> const &pts, StridedMatrix<double, MemorySpace> output) = 0#
virtual std::pair<Kokkos::View<double*, Kokkos::HostSpace>, Kokkos::View<double*, Kokkos::HostSpace>> CoeffBounds() const#

Returns pair with lower and upper bounds on the map parameters that are needed to ensure the map is monotone. Will return vectors with -inf and inf if no bounds exist.

Typical usage could be: Kokkos::View<double*,Kokkos::HostSpace> lb, ub; std::tie(lb,ub) = comp->GetParameterBounds();

virtual void FillCoeffBoundsImpl(Kokkos::View<double*, Kokkos::HostSpace> lb, Kokkos::View<double*, Kokkos::HostSpace> ub) const#