![]() |
Stan Math Library
2.15.0
reverse mode automatic differentiation
|
#include <cholesky_decompose.hpp>
Public Member Functions | |
cholesky_scalar (const Eigen::Matrix< var, -1, -1 > &A, const Eigen::Matrix< double, -1, -1 > &L_A) | |
Constructor for cholesky function. More... | |
virtual void | chain () |
Reverse mode differentiation algorithm refernce: More... | |
![]() | |
vari (double x) | |
Construct a variable implementation from a value. More... | |
vari (double x, bool stacked) | |
virtual | ~vari () |
Throw an illegal argument exception. More... | |
void | init_dependent () |
Initialize the adjoint for this (dependent) variable to 1. More... | |
void | set_zero_adjoint () |
Set the adjoint value of this variable to 0. More... | |
Public Attributes | |
int | M_ |
vari ** | variRefA_ |
vari ** | variRefL_ |
![]() | |
const double | val_ |
The value of this variable. More... | |
double | adj_ |
The adjoint of this variable, which is the partial derivative of this variable with respect to the root variable. More... | |
Additional Inherited Members | |
![]() | |
static void * | operator new (size_t nbytes) |
Allocate memory from the underlying memory pool. More... | |
static void | operator delete (void *) |
Delete a pointer from the underlying memory pool. More... | |
Definition at line 137 of file cholesky_decompose.hpp.
|
inline |
Constructor for cholesky function.
Stores varis for A Instantiates and stores varis for L Instantiates and stores dummy vari for upper triangular part of var result returned in cholesky_decompose function call
variRefL aren't on the chainable autodiff stack, only used for storage and computation. Note that varis for L are constructed externally in cholesky_decompose.
A | matrix |
L_A | matrix, cholesky factor of A |
Definition at line 157 of file cholesky_decompose.hpp.
|
inlinevirtual |
Reverse mode differentiation algorithm refernce:
Mike Giles. An extended collection of matrix derivative results for forward and reverse mode AD. Jan. 2008.
Note algorithm as laid out in Giles is row-major, so Eigen::Matrices are explicitly storage order RowMajor, whereas Eigen defaults to ColumnMajor. Also note algorithm starts by calculating the adjoint for A(M_ - 1, M_ - 1), hence pos on line 94 is decremented to start at pos = M_ * (M_ + 1) / 2.
Reimplemented from stan::math::vari.
Definition at line 191 of file cholesky_decompose.hpp.
int stan::math::cholesky_scalar::M_ |
Definition at line 139 of file cholesky_decompose.hpp.
vari** stan::math::cholesky_scalar::variRefA_ |
Definition at line 140 of file cholesky_decompose.hpp.
vari** stan::math::cholesky_scalar::variRefL_ |
Definition at line 141 of file cholesky_decompose.hpp.