Public Types | Public Member Functions | Protected Attributes | Private Member Functions | Friends | List of all members
gtsam::HessianFactor Class Reference

A Gaussian factor using the canonical parameters (information form) More...

#include <HessianFactor.h>

Inheritance diagram for gtsam::HessianFactor:
Inheritance graph
[legend]

Public Types

typedef GaussianFactor Base
 Typedef to base class. More...
 
typedef SymmetricBlockMatrix::Block Block
 A block from the Hessian matrix. More...
 
typedef SymmetricBlockMatrix::constBlock constBlock
 A block from the Hessian matrix (const version) More...
 
typedef boost::shared_ptr< Thisshared_ptr
 A shared_ptr to this class. More...
 
typedef HessianFactor This
 Typedef to this class. More...
 
- Public Types inherited from gtsam::GaussianFactor
typedef Factor Base
 Our base class. More...
 
typedef boost::shared_ptr< Thisshared_ptr
 shared_ptr to this class More...
 
typedef GaussianFactor This
 This class. More...
 
- Public Types inherited from gtsam::Factor
typedef KeyVector::const_iterator const_iterator
 Const iterator over keys. More...
 
typedef KeyVector::iterator iterator
 Iterator over keys. More...
 

Public Member Functions

Matrix augmentedInformation () const override
 
Matrix augmentedJacobian () const override
 
GaussianFactor::shared_ptr clone () const override
 
double constantTerm () const
 
double & constantTerm ()
 
boost::shared_ptr< GaussianConditionaleliminateCholesky (const Ordering &keys)
 
bool empty () const override
 
bool equals (const GaussianFactor &lf, double tol=1e-9) const override
 
double error (const VectorValues &c) const override
 
DenseIndex getDim (const_iterator variable) const override
 
Vector gradient (Key key, const VectorValues &x) const override
 
VectorValues gradientAtZero () const override
 eta for Hessian More...
 
void gradientAtZero (double *d) const override
 Raw memory access version of gradientAtZero. More...
 
std::map< Key, MatrixhessianBlockDiagonal () const override
 Return the block diagonal of the Hessian for this factor. More...
 
void hessianDiagonal (double *d) const override
 Raw memory access version of hessianDiagonal. More...
 
void hessianDiagonalAdd (VectorValues &d) const override
 Add the current diagonal to a VectorValues instance. More...
 
 HessianFactor ()
 
 HessianFactor (Key j, const Matrix &G, const Vector &g, double f)
 
 HessianFactor (Key j, const Vector &mu, const Matrix &Sigma)
 
 HessianFactor (Key j1, Key j2, const Matrix &G11, const Matrix &G12, const Vector &g1, const Matrix &G22, const Vector &g2, double f)
 
 HessianFactor (Key j1, Key j2, Key j3, const Matrix &G11, const Matrix &G12, const Matrix &G13, const Vector &g1, const Matrix &G22, const Matrix &G23, const Vector &g2, const Matrix &G33, const Vector &g3, double f)
 
 HessianFactor (const KeyVector &js, const std::vector< Matrix > &Gs, const std::vector< Vector > &gs, double f)
 
template<typename KEYS >
 HessianFactor (const KEYS &keys, const SymmetricBlockMatrix &augmentedInformation)
 
 HessianFactor (const JacobianFactor &cg)
 
 HessianFactor (const GaussianFactor &factor)
 
 HessianFactor (const GaussianFactorGraph &factors, const Scatter &scatter)
 
 HessianFactor (const GaussianFactorGraph &factors)
 
const SymmetricBlockMatrixinfo () const
 Return underlying information matrix. More...
 
SymmetricBlockMatrixinfo ()
 
Matrix information () const override
 
Eigen::SelfAdjointView< SymmetricBlockMatrix::constBlock, Eigen::UpperinformationView () const
 Return self-adjoint view onto the information matrix (NOT augmented). More...
 
std::pair< Matrix, Vectorjacobian () const override
 Return (dense) matrix associated with factor. More...
 
SymmetricBlockMatrix::constBlock linearTerm (const_iterator j) const
 
SymmetricBlockMatrix::constBlock linearTerm () const
 
SymmetricBlockMatrix::Block linearTerm ()
 
void multiplyHessianAdd (double alpha, const VectorValues &x, VectorValues &y) const override
 
GaussianFactor::shared_ptr negate () const override
 
void print (const std::string &s="", const KeyFormatter &formatter=DefaultKeyFormatter) const override
 
size_t rows () const
 
VectorValues solve ()
 Solve the system A'*A delta = A'*b in-place, return delta as VectorValues. More...
 
void updateHessian (const KeyVector &keys, SymmetricBlockMatrix *info) const override
 
void updateHessian (HessianFactor *other) const
 
 ~HessianFactor () override
 
- Public Member Functions inherited from gtsam::GaussianFactor
 GaussianFactor ()
 
template<typename CONTAINER >
 GaussianFactor (const CONTAINER &keys)
 
VectorValues hessianDiagonal () const
 Return the diagonal of the Hessian for this factor. More...
 
virtual ~GaussianFactor ()
 
- Public Member Functions inherited from gtsam::Factor
virtual ~Factor ()=default
 Default destructor. More...
 
Key front () const
 First key. More...
 
Key back () const
 Last key. More...
 
const_iterator find (Key key) const
 find More...
 
const KeyVectorkeys () const
 Access the factor's involved variable keys. More...
 
const_iterator begin () const
 
const_iterator end () const
 
size_t size () const
 
virtual void printKeys (const std::string &s="Factor", const KeyFormatter &formatter=DefaultKeyFormatter) const
 print only keys More...
 
KeyVectorkeys ()
 
iterator begin ()
 
iterator end ()
 

Protected Attributes

SymmetricBlockMatrix info_
 The full augmented information matrix, s.t. the quadratic error is 0.5*[x -1]'H[x -1]. More...
 
- Protected Attributes inherited from gtsam::Factor
KeyVector keys_
 The keys involved in this factor. More...
 

Private Member Functions

void Allocate (const Scatter &scatter)
 Allocate for given scatter pattern. More...
 
 HessianFactor (const Scatter &scatter)
 Constructor with given scatter pattern, allocating but not initializing storage. More...
 
template<class ARCHIVE >
void serialize (ARCHIVE &ar, const unsigned int)
 

Friends

class boost::serialization::access
 
class NonlinearClusterTree
 
class NonlinearFactorGraph
 

Additional Inherited Members

- Static Public Member Functions inherited from gtsam::GaussianFactor
template<typename CONTAINER >
static DenseIndex Slot (const CONTAINER &keys, Key key)
 
- Protected Member Functions inherited from gtsam::Factor
 Factor ()
 
template<typename CONTAINER >
 Factor (const CONTAINER &keys)
 
template<typename ITERATOR >
 Factor (ITERATOR first, ITERATOR last)
 
bool equals (const This &other, double tol=1e-9) const
 check equality More...
 
- Static Protected Member Functions inherited from gtsam::Factor
template<typename CONTAINER >
static Factor FromKeys (const CONTAINER &keys)
 
template<typename ITERATOR >
static Factor FromIterators (ITERATOR first, ITERATOR last)
 

Detailed Description

A Gaussian factor using the canonical parameters (information form)

HessianFactor implements a general quadratic factor of the form

\[ E(x) = 0.5 x^T G x - x^T g + 0.5 f \]

that stores the matrix $ G $, the vector $ g $, and the constant term $ f $.

When $ G $ is positive semidefinite, this factor represents a Gaussian, in which case $ G $ is the information matrix $ \Lambda $, $ g $ is the information vector $ \eta $, and $ f $ is the residual sum-square-error at the mean, when $ x = \mu $.

Indeed, the negative log-likelihood of a Gaussian is (up to a constant) $ E(x) = 0.5(x-\mu)^T P^{-1} (x-\mu) $ with $ \mu $ the mean and $ P $ the covariance matrix. Expanding the product we get

\[ E(x) = 0.5 x^T P^{-1} x - x^T P^{-1} \mu + 0.5 \mu^T P^{-1} \mu \]

We define the Information matrix (or Hessian) $ \Lambda = P^{-1} $ and the information vector $ \eta = P^{-1} \mu = \Lambda \mu $ to arrive at the canonical form of the Gaussian:

\[ E(x) = 0.5 x^T \Lambda x - x^T \eta + 0.5 \mu^T \Lambda \mu \]

This factor is one of the factors that can be in a GaussianFactorGraph. It may be returned from NonlinearFactor::linearize(), but is also used internally to store the Hessian during Cholesky elimination.

This can represent a quadratic factor with characteristics that cannot be represented using a JacobianFactor (which has the form $ E(x) = \Vert Ax - b \Vert^2 $ and stores the Jacobian $ A $ and error vector $ b $, i.e. is a sum-of-squares factor). For example, a HessianFactor need not be positive semidefinite, it can be indefinite or even negative semidefinite.

If a HessianFactor is indefinite or negative semi-definite, then in order for solving the linear system to be possible, the Hessian of the full system must be positive definite (i.e. when all small Hessians are combined, the result must be positive definite). If this is not the case, an error will occur during elimination.

This class stores G, g, and f as an augmented matrix HessianFactor::matrix_. The upper-left n x n blocks of HessianFactor::matrix_ store the upper-right triangle of G, the upper-right-most column of length n of HessianFactor::matrix_ stores g, and the lower-right entry of HessianFactor::matrix_ stores f, i.e.

HessianFactor::matrix_ = [ G11 G12 G13 ... g1
0 G22 G23 ... g2
0 0 G33 ... g3
: : : :
0 0 0 ... f ]

Blocks can be accessed as follows:

G11 = info(begin(), begin());
G12 = info(begin(), begin()+1);
G23 = info(begin()+1, begin()+2);
.......

Definition at line 101 of file HessianFactor.h.

Member Typedef Documentation

Typedef to base class.

Definition at line 108 of file HessianFactor.h.

A block from the Hessian matrix.

Definition at line 111 of file HessianFactor.h.

A block from the Hessian matrix (const version)

Definition at line 112 of file HessianFactor.h.

typedef boost::shared_ptr<This> gtsam::HessianFactor::shared_ptr

A shared_ptr to this class.

Definition at line 110 of file HessianFactor.h.

Typedef to this class.

Definition at line 109 of file HessianFactor.h.

Constructor & Destructor Documentation

gtsam::HessianFactor::HessianFactor ( )

default constructor for I/O

Definition at line 84 of file HessianFactor.cpp.

gtsam::HessianFactor::HessianFactor ( Key  j,
const Matrix G,
const Vector g,
double  f 
)

Construct a unary factor. G is the quadratic term (Hessian matrix), g the linear term (a vector), and f the constant term. The quadratic error is: 0.5*(f - 2*x'*g + x'*G*x)

Definition at line 91 of file HessianFactor.cpp.

gtsam::HessianFactor::HessianFactor ( Key  j,
const Vector mu,
const Matrix Sigma 
)

Construct a unary factor, given a mean and covariance matrix. error is 0.5*(x-mu)'inv(Sigma)(x-mu)

Definition at line 104 of file HessianFactor.cpp.

gtsam::HessianFactor::HessianFactor ( Key  j1,
Key  j2,
const Matrix G11,
const Matrix G12,
const Vector g1,
const Matrix G22,
const Vector g2,
double  f 
)

Construct a binary factor. Gxx are the upper-triangle blocks of the quadratic term (the Hessian matrix), gx the pieces of the linear vector term, and f the constant term. JacobianFactor error is

\[ 0.5* (Ax-b)' M (Ax-b) = 0.5*x'A'MAx - x'A'Mb + 0.5*b'Mb \]

HessianFactor error is

\[ 0.5*(x'Gx - 2x'g + f) = 0.5*x'Gx - x'*g + 0.5*f \]

So, with $ A = [A1 A2] $ and $ G=A*'M*A = [A1';A2']*M*[A1 A2] $ we have

n1*n1 G11 = A1'*M*A1
n1*n2 G12 = A1'*M*A2
n2*n2 G22 = A2'*M*A2
n1*1 g1 = A1'*M*b
n2*1 g2 = A2'*M*b
1*1 f = b'*M*b

Definition at line 115 of file HessianFactor.cpp.

gtsam::HessianFactor::HessianFactor ( Key  j1,
Key  j2,
Key  j3,
const Matrix G11,
const Matrix G12,
const Matrix G13,
const Vector g1,
const Matrix G22,
const Matrix G23,
const Vector g2,
const Matrix G33,
const Vector g3,
double  f 
)

Construct a ternary factor. Gxx are the upper-triangle blocks of the quadratic term (the Hessian matrix), gx the pieces of the linear vector term, and f the constant term.

Definition at line 128 of file HessianFactor.cpp.

gtsam::HessianFactor::HessianFactor ( const KeyVector js,
const std::vector< Matrix > &  Gs,
const std::vector< Vector > &  gs,
double  f 
)

Construct an n-way factor. Gs contains the upper-triangle blocks of the quadratic term (the Hessian matrix) provided in row-order, gs the pieces of the linear vector term, and f the constant term.

Definition at line 158 of file HessianFactor.cpp.

template<typename KEYS >
gtsam::HessianFactor::HessianFactor ( const KEYS &  keys,
const SymmetricBlockMatrix augmentedInformation 
)

Constructor with an arbitrary number of keys and with the augmented information matrix specified as a block matrix.

Definition at line 25 of file HessianFactor-inl.h.

gtsam::HessianFactor::HessianFactor ( const JacobianFactor cg)
explicit

Construct from a JacobianFactor (or from a GaussianConditional since it derives from it)

Definition at line 230 of file HessianFactor.cpp.

gtsam::HessianFactor::HessianFactor ( const GaussianFactor factor)
explicit

Attempt to construct from any GaussianFactor - currently supports JacobianFactor, HessianFactor, GaussianConditional, or any derived classes.

Definition at line 237 of file HessianFactor.cpp.

gtsam::HessianFactor::HessianFactor ( const GaussianFactorGraph factors,
const Scatter scatter 
)
explicit

Combine a set of factors into a single dense HessianFactor

Definition at line 252 of file HessianFactor.cpp.

gtsam::HessianFactor::HessianFactor ( const GaussianFactorGraph factors)
inlineexplicit

Combine a set of factors into a single dense HessianFactor

Definition at line 182 of file HessianFactor.h.

gtsam::HessianFactor::~HessianFactor ( )
inlineoverride

Destructor

Definition at line 186 of file HessianFactor.h.

gtsam::HessianFactor::HessianFactor ( const Scatter scatter)
private

Constructor with given scatter pattern, allocating but not initializing storage.

Definition at line 79 of file HessianFactor.cpp.

Member Function Documentation

void gtsam::HessianFactor::Allocate ( const Scatter scatter)
private

Allocate for given scatter pattern.

Definition at line 61 of file HessianFactor.cpp.

Matrix gtsam::HessianFactor::augmentedInformation ( ) const
overridevirtual

Return the augmented information matrix represented by this GaussianFactor. The augmented information matrix contains the information matrix with an additional column holding the information vector, and an additional row holding the transpose of the information vector. The lower-right entry contains the constant error term (when $ \delta x = 0 $). The augmented information matrix is described in more detail in HessianFactor, which in fact stores an augmented information matrix.

For HessianFactor, this is the same as info() except that this function returns a complete symmetric matrix whereas info() returns a matrix where only the upper triangle is valid, but should be interpreted as symmetric. This is because info() returns only a reference to the internal representation of the augmented information matrix, which stores only the upper triangle.

Implements gtsam::GaussianFactor.

Definition at line 289 of file HessianFactor.cpp.

Matrix gtsam::HessianFactor::augmentedJacobian ( ) const
overridevirtual

Return (dense) matrix associated with factor The returned system is an augmented matrix: [A b]

Parameters
setweight to use whitening to bake in weights

Implements gtsam::GaussianFactor.

Definition at line 334 of file HessianFactor.cpp.

GaussianFactor::shared_ptr gtsam::HessianFactor::clone ( ) const
inlineoverridevirtual

Clone this HessianFactor

Implements gtsam::GaussianFactor.

Definition at line 189 of file HessianFactor.h.

double gtsam::HessianFactor::constantTerm ( ) const
inline

Return the constant term $ f $ as described above

Returns
The constant term $ f $

Definition at line 230 of file HessianFactor.h.

double& gtsam::HessianFactor::constantTerm ( )
inline

Return the constant term $ f $ as described above

Returns
The constant term $ f $

Definition at line 238 of file HessianFactor.h.

boost::shared_ptr< GaussianConditional > gtsam::HessianFactor::eliminateCholesky ( const Ordering keys)

In-place elimination that returns a conditional on (ordered) keys specified, and leaves this factor to be on the remaining keys (separator) only. Does dense partial Cholesky.

Definition at line 474 of file HessianFactor.cpp.

bool gtsam::HessianFactor::empty ( ) const
inlineoverridevirtual

Check if the factor is empty. TODO: How should this be defined?

Implements gtsam::GaussianFactor.

Definition at line 225 of file HessianFactor.h.

bool gtsam::HessianFactor::equals ( const GaussianFactor lf,
double  tol = 1e-9 
) const
overridevirtual

Compare to another factor for testing (implementing Testable)

Implements gtsam::GaussianFactor.

Definition at line 280 of file HessianFactor.cpp.

double gtsam::HessianFactor::error ( const VectorValues c) const
overridevirtual

Evaluate the factor error f(x). returns 0.5*[x -1]'H[x -1] (also see constructor documentation)

Implements gtsam::GaussianFactor.

Definition at line 344 of file HessianFactor.cpp.

DenseIndex gtsam::HessianFactor::getDim ( const_iterator  variable) const
inlineoverridevirtual

Return the dimension of the variable pointed to by the given key iterator todo: Remove this in favor of keeping track of dimensions with variables?

Parameters
variableAn iterator pointing to the slot in this factor. You can use, for example, begin() + 2 to get the 3rd variable in this factor.

Implements gtsam::GaussianFactor.

Definition at line 210 of file HessianFactor.h.

Vector gtsam::HessianFactor::gradient ( Key  key,
const VectorValues x 
) const
overridevirtual

Compute the gradient at a key: f(x_i) = G_ij*x_j - g_i

Implements gtsam::GaussianFactor.

Definition at line 453 of file HessianFactor.cpp.

VectorValues gtsam::HessianFactor::gradientAtZero ( ) const
overridevirtual

eta for Hessian

Implements gtsam::GaussianFactor.

Definition at line 437 of file HessianFactor.cpp.

void gtsam::HessianFactor::gradientAtZero ( double *  d) const
overridevirtual

Raw memory access version of gradientAtZero.

Implements gtsam::GaussianFactor.

Reimplemented in gtsam::RegularHessianFactor< D >.

Definition at line 447 of file HessianFactor.cpp.

map< Key, Matrix > gtsam::HessianFactor::hessianBlockDiagonal ( ) const
overridevirtual

Return the block diagonal of the Hessian for this factor.

Implements gtsam::GaussianFactor.

Definition at line 323 of file HessianFactor.cpp.

void gtsam::HessianFactor::hessianDiagonal ( double *  d) const
overridevirtual

Raw memory access version of hessianDiagonal.

Implements gtsam::GaussianFactor.

Reimplemented in gtsam::RegularHessianFactor< D >.

Definition at line 317 of file HessianFactor.cpp.

void gtsam::HessianFactor::hessianDiagonalAdd ( VectorValues d) const
overridevirtual

Add the current diagonal to a VectorValues instance.

Implements gtsam::GaussianFactor.

Definition at line 305 of file HessianFactor.cpp.

const SymmetricBlockMatrix& gtsam::HessianFactor::info ( ) const
inline

Return underlying information matrix.

Definition at line 265 of file HessianFactor.h.

SymmetricBlockMatrix& gtsam::HessianFactor::info ( )
inline

Return non-const information matrix. TODO(gareth): Review the sanity of having non-const access to this.

Definition at line 269 of file HessianFactor.h.

Matrix gtsam::HessianFactor::information ( ) const
overridevirtual

Return the non-augmented information matrix represented by this GaussianFactor.

Implements gtsam::GaussianFactor.

Definition at line 300 of file HessianFactor.cpp.

Eigen::SelfAdjointView< SymmetricBlockMatrix::constBlock, Eigen::Upper > gtsam::HessianFactor::informationView ( ) const

Return self-adjoint view onto the information matrix (NOT augmented).

Definition at line 295 of file HessianFactor.cpp.

std::pair< Matrix, Vector > gtsam::HessianFactor::jacobian ( ) const
overridevirtual

Return (dense) matrix associated with factor.

Implements gtsam::GaussianFactor.

Definition at line 339 of file HessianFactor.cpp.

SymmetricBlockMatrix::constBlock gtsam::HessianFactor::linearTerm ( const_iterator  j) const
inline

Return the part of linear term $ g $ as described above corresponding to the requested variable.

Parameters
jWhich block row to get, as an iterator pointing to the slot in this factor. You can use, for example, begin() + 2 to get the 3rd variable in this factor.
Returns
The linear term $ g $

Definition at line 244 of file HessianFactor.h.

SymmetricBlockMatrix::constBlock gtsam::HessianFactor::linearTerm ( ) const
inline

Return the complete linear term $ g $ as described above.

Returns
The linear term $ g $

Definition at line 251 of file HessianFactor.h.

SymmetricBlockMatrix::Block gtsam::HessianFactor::linearTerm ( )
inline

Return the complete linear term $ g $ as described above.

Returns
The linear term $ g $

Definition at line 259 of file HessianFactor.h.

void gtsam::HessianFactor::multiplyHessianAdd ( double  alpha,
const VectorValues x,
VectorValues y 
) const
overridevirtual

y += alpha * A'*A*x

Implements gtsam::GaussianFactor.

Reimplemented in gtsam::RegularHessianFactor< D >.

Definition at line 398 of file HessianFactor.cpp.

GaussianFactor::shared_ptr gtsam::HessianFactor::negate ( ) const
overridevirtual

Construct the corresponding anti-factor to negate information stored stored in this factor.

Returns
a HessianFactor with negated Hessian matrices

Implements gtsam::GaussianFactor.

Definition at line 390 of file HessianFactor.cpp.

void gtsam::HessianFactor::print ( const std::string &  s = "",
const KeyFormatter formatter = DefaultKeyFormatter 
) const
overridevirtual

Print the factor for debugging and testing (implementing Testable)

Implements gtsam::GaussianFactor.

Definition at line 268 of file HessianFactor.cpp.

size_t gtsam::HessianFactor::rows ( void  ) const
inline

Return the number of columns and rows of the Hessian matrix, including the information vector.

Definition at line 215 of file HessianFactor.h.

template<class ARCHIVE >
void gtsam::HessianFactor::serialize ( ARCHIVE &  ar,
const unsigned  int 
)
inlineprivate

Definition at line 370 of file HessianFactor.h.

VectorValues gtsam::HessianFactor::solve ( )

Solve the system A'*A delta = A'*b in-place, return delta as VectorValues.

Definition at line 505 of file HessianFactor.cpp.

void gtsam::HessianFactor::updateHessian ( const KeyVector keys,
SymmetricBlockMatrix info 
) const
overridevirtual

Update an information matrix by adding the information corresponding to this factor (used internally during elimination).

Parameters
keysTHe ordered vector of keys for the information matrix to be updated
infoThe information matrix to be updated

Implements gtsam::GaussianFactor.

Definition at line 361 of file HessianFactor.cpp.

void gtsam::HessianFactor::updateHessian ( HessianFactor other) const
inline

Update another Hessian factor

Parameters
otherthe HessianFactor to be updated

Definition at line 328 of file HessianFactor.h.

Friends And Related Function Documentation

friend class boost::serialization::access
friend

Serialization function

Definition at line 368 of file HessianFactor.h.

friend class NonlinearClusterTree
friend

Definition at line 365 of file HessianFactor.h.

friend class NonlinearFactorGraph
friend

Definition at line 364 of file HessianFactor.h.

Member Data Documentation

SymmetricBlockMatrix gtsam::HessianFactor::info_
protected

The full augmented information matrix, s.t. the quadratic error is 0.5*[x -1]'H[x -1].

Definition at line 104 of file HessianFactor.h.


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


gtsam
Author(s):
autogenerated on Sat May 8 2021 02:58:12