API¶
Solver Interface¶
The main interface to the celerite solver is via the
celerite::solver::CholeskySolver
class but many of the key
interfaces are abstracted into the celerite::solver::Solver
.
-
template<typename
T
>
classcelerite::solver
::
Solver
¶ Subclassed by celerite::solver::CholeskySolver< T, SIZE >, celerite::solver::DirectSolver< T >
Public Functions
-
void
compute
(const T &jitter, const vector_t &a_real, const vector_t &c_real, const vector_t &a_comp, const vector_t &b_comp, const vector_t &c_comp, const vector_t &d_comp, const Eigen::VectorXd &A, const Eigen::MatrixXd &U, const Eigen::MatrixXd &V, const Eigen::VectorXd &x, const Eigen::VectorXd &diag) = 0¶ Compute the matrix and factorize
This method is overloaded by subclasses to provide the specific implementation.
- Parameters
jitter
: The jitter of the kernel.a_real
: The coefficients of the real terms.c_real
: The exponents of the real terms.a_comp
: The real part of the coefficients of the complex terms.b_comp
: The imaginary part of the of the complex terms.c_comp
: The real part of the exponents of the complex terms.d_comp
: The imaginary part of the exponents of the complex terms.x
: The sorted array of input coordinates.diag
: An array that should be added to the diagonal of the matrix. This often corresponds to measurement uncertainties and in that case,diag
should be the measurement variance (i.e. sigma^2).
-
T
dot_solve
(const Eigen::VectorXd &b) const¶ Solve the system
b^T . A^-1 . b
A previous call to
solver::Solver::compute
defines a matrixA
and this method solvesb^T . A^-1 . b
for a vectorb
.- Parameters
b
: The right hand side of the linear system.
-
bool
computed
() const¶ Flag indicating if
compute
was successfully executed.
-
void
compute
(const T &jitter, const vector_t &a_real, const vector_t &c_real, const vector_t &a_comp, const vector_t &b_comp, const vector_t &c_comp, const vector_t &d_comp, const Eigen::VectorXd &x, const Eigen::VectorXd &diag)¶ Compute the matrix and factorize for purely real alphas
- Parameters
jitter
: The jitter of the kernel.a_real
: The coefficients of the real terms.c_real
: The exponents of the real terms.a_comp
: The real part of the coefficients of the complex terms.b_comp
: The imaginary part of the coefficients of the complex terms.c_comp
: The real part of the exponents of the complex terms.d_comp
: The imaginary part of the exponents of the complex terms.x
: The sorted array of input coordinates.diag
: An array that should be added to the diagonal of the matrix. This often corresponds to measurement uncertainties and in that case,diag
should be the measurement variance (i.e. sigma^2).
-
void
compute
(const T &jitter, const vector_t &a_real, const vector_t &c_real, const vector_t &a_comp, const vector_t &c_comp, const vector_t &d_comp, const Eigen::VectorXd &x, const Eigen::VectorXd &diag)¶ Compute the matrix and factorize for purely real alphas
- Parameters
jitter
: The jitter of the kernel.a_real
: The coefficients of the real terms.c_real
: The exponents of the real terms.a_comp
: The real part of the coefficients of the complex terms.c_comp
: The real part of the exponents of the complex terms.d_comp
: The imaginary part of the exponents of the complex terms.x
: The sorted array of input coordinates.diag
: An array that should be added to the diagonal of the matrix. This often corresponds to measurement uncertainties and in that case,diag
should be the measurement variance (i.e. sigma^2).
-
void
compute
(const T &jitter, const vector_t &a_real, const vector_t &c_real, const Eigen::VectorXd &x, const Eigen::VectorXd &diag)¶ Compute the matrix and factorize for a set of purely real terms
- Parameters
jitter
: The jitter of the kernel.a_real
: The coefficients of the real terms.c_real
: The exponents of the real terms.x
: The sorted array of input coordinates.diag
: An array that should be added to the diagonal of the matrix. This often corresponds to measurement uncertainties and in that case,diag
should be the measurement variance (i.e. sigma^2).
-
void
compute
(const T &jitter, const vector_t &a_comp, const vector_t &c_comp, const vector_t &d_comp, const Eigen::VectorXd &x, const Eigen::VectorXd &diag)¶ Compute the matrix and factorize for purely complex terms with real alphas
- Parameters
jitter
: The jitter of the kernel.a_comp
: The real part of the coefficients of the complex terms.c_comp
: The real part of the exponents of the complex terms.d_comp
: The imaginary part of the exponents of the complex terms.x
: The sorted array of input coordinates.diag
: An array that should be added to the diagonal of the matrix. This often corresponds to measurement uncertainties and in that case,diag
should be the measurement variance (i.e. sigma^2).
-
void
compute
(const T &jitter, const vector_t &a_comp, const vector_t &b_comp, const vector_t &c_comp, const vector_t &d_comp, const Eigen::VectorXd &x, const Eigen::VectorXd &diag)¶ Compute the matrix and factorize for purely complex terms
- Parameters
jitter
: The jitter of the kernel.a_comp
: The real part of the coefficients of the complex terms.b_comp
: The imaginary part of the coefficients of the complex terms.c_comp
: The real part of the exponents of the complex terms.d_comp
: The imaginary part of the exponents of the complex terms.x
: The sorted array of input coordinates.diag
: An array that should be added to the diagonal of the matrix. This often corresponds to measurement uncertainties and in that case,diag
should be the measurement variance (i.e. sigma^2).
-
void
Fast Solver¶
The celerite::solver::CholeskySolver
is an \(\mathcal{O}(N)\)
solver that exploits the semi-separable structure of the matrix to solve and
compute determinants.
-
template<typename
T
, intSIZE
= Eigen::Dynamic>
classcelerite::solver
::
CholeskySolver
: public celerite::solver::Solver<T>¶ Public Functions
-
void
compute
(const T &jitter, const vector_t &a_real, const vector_t &c_real, const vector_t &a_comp, const vector_t &b_comp, const vector_t &c_comp, const vector_t &d_comp, const Eigen::VectorXd &A, const Eigen::MatrixXd &U, const Eigen::MatrixXd &V, const Eigen::VectorXd &x, const Eigen::VectorXd &diag)¶ Compute the Cholesky factorization of the matrix
- Parameters
jitter
: The jitter of the kernel.a_real
: The coefficients of the real terms.c_real
: The exponents of the real terms.a_comp
: The real part of the coefficients of the complex terms.b_comp
: The imaginary part of the of the complex terms.c_comp
: The real part of the exponents of the complex terms.d_comp
: The imaginary part of the exponents of the complex terms.x
: The sorted array of input coordinates.diag
: An array that should be added to the diagonal of the matrix. This often corresponds to measurement uncertainties and in that case,diag
should be the measurement variance (i.e. sigma^2).
-
matrix_t
solve
(const Eigen::MatrixXd &b) const¶ Solve the system
b^T . A^-1 . b
A previous call to
solver::CholeskySolver::compute
defines a matrixA
and this method solvesb^T . A^-1 . b
for a vectorb
.- Parameters
b
: The right hand side of the linear system.
-
T
dot_solve
(const Eigen::VectorXd &b) const¶ Solve the system
b^T . A^-1 . b
A previous call to
solver::CholeskySolver::compute
defines a matrixA
and this method solvesb^T . A^-1 . b
for a vectorb
.- Parameters
b
: The right hand side of the linear system.
-
matrix_t
dot_L
(const Eigen::MatrixXd &z) const¶ Compute the dot product of the square root of a celerite matrix
This method computes
L.z
whereA = L.L^T
is the matrix defined incompute
.- Parameters
z
: The matrix z from above.
-
matrix_t
dot
(const T &jitter, const vector_t &a_real, const vector_t &c_real, const vector_t &a_comp, const vector_t &b_comp, const vector_t &c_comp, const vector_t &d_comp, const vector_t &A, const matrix_t &U, const matrix_t &V, const Eigen::VectorXd &x, const Eigen::MatrixXd &z)¶ Compute the dot product of a celerite matrix with another matrix
- Parameters
jitter
: The jitter of the kernel.a_real
: The coefficients of the real terms.c_real
: The exponents of the real terms.a_comp
: The real part of the coefficients of the complex terms.b_comp
: The imaginary part of the of the complex terms.c_comp
: The real part of the exponents of the complex terms.d_comp
: The imaginary part of the exponents of the complex terms.x
: The sorted array of input coordinates.z
: The matrix that will be dotted with the celerite matrix.
-
vector_t
predict
(const Eigen::VectorXd &y, const Eigen::VectorXd &x) const¶ Compute the conditional mean for a GP in O(N)
This method computes
K_*.(K+sigma)^{-1}.y
whereK
is the matrix defined incompute
.- Parameters
y
: The data vector to condition onx
: The input coordinates to compute the predictions at
-
void