API

Solver Interface

The main interface to the celerite solver is via the celerite::solver::BandSolver class but many of the key interfaces are abstracted into the celerite::solver::Solver.

template <typename T>
class celerite::solver::Solver

Subclassed by celerite::solver::BandSolver< T >, celerite::solver::DirectSolver< T >, celerite::solver::SparseSolver< T >

Public Functions

virtual int compute(const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_imag, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_imag, const Eigen::VectorXd &x, const Eigen::Matrix<T, Eigen::Dynamic, 1> &diag) = 0

Compute the matrix and factorize

This method is overloaded by subclasses to provide the specific implementation.

Return
0 on success. 1 for mismatched dimensions.
Parameters
  • alpha_real: The coefficients of the real terms.
  • beta_real: The exponents of the real terms.
  • alpha_complex_real: The real part of the coefficients of the complex terms.
  • alpha_complex_imag: The imaginary part of the of the complex terms.
  • beta_complex_real: The real part of the exponents of the complex terms.
  • beta_complex_imag: 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).

Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic> solve(const Eigen::MatrixXd &b) const

Solve a linear system for the matrix defined in compute

A previous call to solver::Solver::compute defines a matrix A and this method solves for x in the matrix equation A.x = b.

Parameters
  • b: The right hand side of the linear system.
  • x: A pointer that will be overwritten with the result.

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 matrix A and this method solves b^T . A^-1 . b for a vector b.

Parameters
  • b: The right hand side of the linear system.

T log_determinant() const

Get the log determinant of the matrix.

bool computed() const

Flag indicating if compute was successfully executed.

int compute(const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_imag, const Eigen::VectorXd &x, const Eigen::Matrix<T, Eigen::Dynamic, 1> &diag)

Compute the matrix and factorize for purely real alphas

Return
0 on success. 1 for mismatched dimensions.
Parameters
  • alpha_real: The coefficients of the real terms.
  • beta_real: The exponents of the real terms.
  • alpha_complex_real: The real part of the coefficients of the complex terms.
  • beta_complex_real: The real part of the exponents of the complex terms.
  • beta_complex_imag: 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).

int compute(const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_real, const Eigen::VectorXd &x, const Eigen::Matrix<T, Eigen::Dynamic, 1> &diag)

Compute the matrix and factorize for a set of purely real terms

Return
0 on success. 1 for mismatched dimensions.
Parameters
  • alpha_real: The coefficients of the real terms.
  • beta_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).

int compute(const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_imag, const Eigen::VectorXd &x, const Eigen::Matrix<T, Eigen::Dynamic, 1> &diag)

Compute the matrix and factorize for purely complex terms with real alphas

Return
0 on success. 1 for mismatched dimensions.
Parameters
  • alpha_complex_real: The real part of the coefficients of the complex terms.
  • beta_complex_real: The real part of the exponents of the complex terms.
  • beta_complex_imag: 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).

int compute(const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_imag, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_imag, const Eigen::VectorXd &x, const Eigen::Matrix<T, Eigen::Dynamic, 1> &diag)

Compute the matrix and factorize for purely complex terms

Return
0 on success. 1 for mismatched dimensions.
Parameters
  • alpha_complex_real: The real part of the coefficients of the complex terms.
  • alpha_complex_imag: The imaginary part of the coefficients of the complex terms.
  • beta_complex_real: The real part of the exponents of the complex terms.
  • beta_complex_imag: 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).

Fast Solver

The celerite::solver::BandSolver is an implementation of celerite::solver::Solver that exploits the “extended matrix” formalism to solve and compute determinants in \(\mathcal{O}(N)\) operations.

template <typename T>
class celerite::solver::BandSolver

The class provides all of the computation power for the fast GP solver

The celerite::solver::BandSolver::compute method must be called before most of the other methods.

Inherits from celerite::solver::Solver< T >

Public Functions

BandSolver(bool use_lapack = false)

You can decide to use LAPACK for solving if it is available

Parameters
  • use_lapack: If true, LAPACK will be used for solving the band

int compute(const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_imag, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_imag, const Eigen::VectorXd &x, const Eigen::Matrix<T, Eigen::Dynamic, 1> &diag)

Compute the extended matrix and perform the banded LU decomposition

Return
0 on success. 1 for mismatched dimensions.
Parameters
  • alpha_real: The coefficients of the real terms.
  • beta_real: The exponents of the real terms.
  • alpha_complex_real: The real part of the coefficients of the complex terms.
  • alpha_complex_imag: The imaginary part of the of the complex terms.
  • beta_complex_real: The real part of the exponents of the complex terms.
  • beta_complex_imag: 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 solve(const Eigen::MatrixXd &b, T *x) const

Solve a linear system for the matrix defined in compute

A previous call to solver.Solver.compute defines a matrix A and this method solves for x in the matrix equation A.x = b.

Parameters
  • b: The right hand side of the linear system.
  • x: A pointer that will be overwritten with the result.

Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic> dot(const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &alpha_complex_imag, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_real, const Eigen::Matrix<T, Eigen::Dynamic, 1> &beta_complex_imag, const Eigen::VectorXd &x, const Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic> &b_in)

Compute the dot product of a celerite matrix and another arbitrary matrix

This method computes A.b where A is defined by the parameters and b is an arbitrary matrix of the correct shape.

Return
The matrix A.b.
Parameters
  • alpha_real: The coefficients of the real terms.
  • beta_real: The exponents of the real terms.
  • alpha_complex_real: The real part of the coefficients of the complex terms.
  • alpha_complex_imag: The imaginary part of the of the complex terms.
  • beta_complex_real: The real part of the exponents of the complex terms.
  • beta_complex_imag: The imaginary part of the exponents of the complex terms.
  • x: The sorted array of input coordinates.
  • b_in: The matrix b described above.