Class Model_DAE

Inheritance Relationships

Base Type

Class Documentation

class amici::Model_DAE : public amici::Model

The Model class represents an AMICI DAE model.

The model does not contain any data, but represents the state of the model at a specific time t. The states must not always be in sync, but may be updated asynchronously.

Public Functions

Model_DAE() = default

default constructor

inline Model_DAE(const ModelDimensions &model_dimensions, SimulationParameters simulation_parameters, const SecondOrderMode o2mode, std::vector<realtype> const &idlist, std::vector<int> const &z2event, const bool pythonGenerated = false, const int ndxdotdp_explicit = 0)

Constructor with model dimensions.

Parameters
  • model_dimensionsModel dimensions

  • simulation_parameters – Simulation parameters

  • o2mode – second order sensitivity mode

  • idlist – indexes indicating algebraic components (DAE only)

  • z2event – mapping of event outputs to events

  • pythonGenerated – flag indicating matlab or python wrapping

  • ndxdotdp_explicit – number of nonzero elements dxdotdp_explicit

virtual void fJ(realtype t, realtype cj, const AmiVector &x, const AmiVector &dx, const AmiVector &xdot, SUNMatrix J) override

Dense Jacobian function.

Parameters
  • t – time

  • cj – scaling factor (inverse of timestep, DAE only)

  • x – state

  • dx – time derivative of state (DAE only)

  • xdot – values of residual function (unused)

  • J – dense matrix to which values of the jacobian will be written

void fJ(realtype t, realtype cj, const_N_Vector x, const_N_Vector dx, const_N_Vector xdot, SUNMatrix J)

Jacobian of xdot with respect to states x.

Parameters
  • t – timepoint

  • cj – scaling factor, inverse of the step size

  • x – Vector with the states

  • dx – Vector with the derivative states

  • xdot – Vector with the right hand side

  • J – Matrix to which the Jacobian will be written

virtual void fJB(const realtype t, realtype cj, const AmiVector &x, const AmiVector &dx, const AmiVector &xB, const AmiVector &dxB, const AmiVector &xBdot, SUNMatrix JB) override

Dense Jacobian function.

Parameters
  • t – time

  • cj – scaling factor (inverse of timestep, DAE only)

  • x – state

  • dx – time derivative of state (DAE only)

  • xB – Vector with the adjoint states

  • dxB – Vector with the adjoint derivative states

  • xBdot – Vector with the adjoint right hand side (unused)

  • JB – dense matrix to which values of the jacobian will be written

void fJB(realtype t, realtype cj, const_N_Vector x, const_N_Vector dx, const_N_Vector xB, const_N_Vector dxB, SUNMatrix JB)

Jacobian of xBdot with respect to adjoint state xB.

Parameters
  • t – timepoint

  • cj – scaling factor, inverse of the step size

  • x – Vector with the states

  • dx – Vector with the derivative states

  • xB – Vector with the adjoint states

  • dxB – Vector with the adjoint derivative states

  • JB – Matrix to which the Jacobian will be written

virtual void fJSparse(realtype t, realtype cj, const AmiVector &x, const AmiVector &dx, const AmiVector &xdot, SUNMatrix J) override

Sparse Jacobian function.

Parameters
  • t – time

  • cj – scaling factor (inverse of timestep, DAE only)

  • x – state

  • dx – time derivative of state (DAE only)

  • xdot – values of residual function (unused)

  • J – sparse matrix to which values of the Jacobian will be written

void fJSparse(realtype t, realtype cj, const_N_Vector x, const_N_Vector dx, SUNMatrix J)

J in sparse form (for sparse solvers from the SuiteSparse Package)

Parameters
  • t – timepoint

  • cj – scalar in Jacobian (inverse stepsize)

  • x – Vector with the states

  • dx – Vector with the derivative states

  • J – Matrix to which the Jacobian will be written

virtual void fJSparseB(const realtype t, realtype cj, const AmiVector &x, const AmiVector &dx, const AmiVector &xB, const AmiVector &dxB, const AmiVector &xBdot, SUNMatrix JB) override

Sparse Jacobian function.

Parameters
  • t – time

  • cj – scaling factor (inverse of timestep, DAE only)

  • x – state

  • dx – time derivative of state (DAE only)

  • xB – Vector with the adjoint states

  • dxB – Vector with the adjoint derivative states

  • xBdot – Vector with the adjoint right hand side (unused)

  • JB – dense matrix to which values of the jacobian will be written

void fJSparseB(realtype t, realtype cj, const_N_Vector x, const_N_Vector dx, const_N_Vector xB, const_N_Vector dxB, SUNMatrix JB)

JB in sparse form (for sparse solvers from the SuiteSparse Package)

Parameters
  • t – timepoint

  • cj – scalar in Jacobian

  • x – Vector with the states

  • dx – Vector with the derivative states

  • xB – Vector with the adjoint states

  • dxB – Vector with the adjoint derivative states

  • JB – Matrix to which the Jacobian will be written

virtual void fJDiag(realtype t, AmiVector &JDiag, realtype cj, const AmiVector &x, const AmiVector &dx) override

Diagonal of the Jacobian (for preconditioning)

Parameters
  • t – timepoint

  • JDiag – Vector to which the Jacobian diagonal will be written

  • cj – scaling factor, inverse of the step size

  • x – Vector with the states

  • dx – Vector with the derivative states

virtual void fJv(realtype t, const AmiVector &x, const AmiVector &dx, const AmiVector &xdot, const AmiVector &v, AmiVector &nJv, realtype cj) override

Jacobian multiply function.

Parameters
  • t – time

  • x – state

  • dx – time derivative of state (DAE only)

  • xdot – values of residual function (unused)

  • v – multiplication vector (unused)

  • nJv – array to which result of multiplication will be written

  • cj – scaling factor (inverse of timestep, DAE only)

void fJv(realtype t, const_N_Vector x, const_N_Vector dx, const_N_Vector v, N_Vector Jv, realtype cj)

Matrix vector product of J with a vector v (for iterative solvers)

Parameters
  • t – timepoint

  • cj – scaling factor, inverse of the step size

  • x – Vector with the states

  • dx – Vector with the derivative states

  • v – Vector with which the Jacobian is multiplied

  • Jv – Vector to which the Jacobian vector product will be written

void fJvB(realtype t, const_N_Vector x, const_N_Vector dx, const_N_Vector xB, const_N_Vector dxB, const_N_Vector vB, N_Vector JvB, realtype cj)

Matrix vector product of JB with a vector v (for iterative solvers)

Parameters
  • t – timepoint

  • x – Vector with the states

  • dx – Vector with the derivative states

  • xB – Vector with the adjoint states

  • dxB – Vector with the adjoint derivative states

  • vB – Vector with which the Jacobian is multiplied

  • JvB – Vector to which the Jacobian vector product will be written

  • cj – scalar in Jacobian (inverse stepsize)

virtual void froot(realtype t, const AmiVector &x, const AmiVector &dx, gsl::span<realtype> root) override

Root function.

Parameters
  • t – time

  • x – state

  • dx – time derivative of state (DAE only)

  • root – array to which values of the root function will be written

void froot(realtype t, const_N_Vector x, const_N_Vector dx, gsl::span<realtype> root)

Event trigger function for events.

Parameters
  • t – timepoint

  • x – Vector with the states

  • dx – Vector with the derivative states

  • root – array with root function values

virtual void fxdot(realtype t, const AmiVector &x, const AmiVector &dx, AmiVector &xdot) override

Residual function.

Parameters
  • t – time

  • x – state

  • dx – time derivative of state (DAE only)

  • xdot – array to which values of the residual function will be written

void fxdot(realtype t, const_N_Vector x, const_N_Vector dx, N_Vector xdot)

Residual function of the DAE.

Parameters
  • t – timepoint

  • x – Vector with the states

  • dx – Vector with the derivative states

  • xdot – Vector with the right hand side

void fxBdot(realtype t, const_N_Vector x, const_N_Vector dx, const_N_Vector xB, const_N_Vector dxB, N_Vector xBdot)

Right hand side of differential equation for adjoint state xB.

Parameters
  • t – timepoint

  • x – Vector with the states

  • dx – Vector with the derivative states

  • xB – Vector with the adjoint states

  • dxB – Vector with the adjoint derivative states

  • xBdot – Vector with the adjoint right hand side

void fqBdot(realtype t, const_N_Vector x, const_N_Vector dx, const_N_Vector xB, const_N_Vector dxB, N_Vector qBdot)

Right hand side of integral equation for quadrature states qB.

Parameters
  • t – timepoint

  • x – Vector with the states

  • dx – Vector with the derivative states

  • xB – Vector with the adjoint states

  • dxB – Vector with the adjoint derivative states

  • qBdot – Vector with the adjoint quadrature right hand side

virtual void fxBdot_ss(const realtype t, const AmiVector &xB, const AmiVector &dxB, AmiVector &xBdot) override

Residual function backward when running in steady state mode.

Parameters
  • t – time

  • xB – adjoint state

  • dxB – time derivative of state (DAE only)

  • xBdot – array to which values of the residual function will be written

void fxBdot_ss(realtype t, const_N_Vector xB, const_N_Vector dxB, N_Vector xBdot) const

Implementation of fxBdot for steady state case at the N_Vector level.

Parameters
  • t – timepoint

  • xB – Vector with the adjoint state

  • dxB – Vector with the adjoint derivative states

  • xBdot – Vector with the adjoint right hand side

void fqBdot_ss(realtype t, const_N_Vector xB, const_N_Vector dxB, N_Vector qBdot) const

Implementation of fqBdot for steady state at the N_Vector level.

Parameters
  • t – timepoint

  • xB – Vector with the adjoint states

  • dxB – Vector with the adjoint derivative states

  • qBdot – Vector with the adjoint quadrature right hand side

virtual void fJSparseB_ss(SUNMatrix JB) override

Sparse Jacobian function backward, steady state case.

Parameters

JB – sparse matrix to which values of the Jacobian will be written

virtual void writeSteadystateJB(const realtype t, realtype cj, const AmiVector &x, const AmiVector &dx, const AmiVector &xB, const AmiVector &dxB, const AmiVector &xBdot) override

Computes the sparse backward Jacobian for steadystate integration and writes it to the model member.

Parameters
  • t – timepoint

  • cj – scalar in Jacobian

  • x – Vector with the states

  • dx – Vector with the derivative states

  • xB – Vector with the adjoint states

  • dxB – Vector with the adjoint derivative states

  • xBdot – Vector with the adjoint state right hand side

void fdxdotdp(realtype t, const const_N_Vector x, const const_N_Vector dx)

Sensitivity of dx/dt wrt model parameters p.

Parameters
  • t – timepoint

  • x – Vector with the states

  • dx – Vector with the derivative states

inline virtual void fdxdotdp(const realtype t, const AmiVector &x, const AmiVector &dx) override

Model-specific sparse implementation of explicit parameter derivative of right hand side.

Parameters
  • t – time

  • x – state

  • dx – time derivative of state (DAE only)

virtual void fsxdot(realtype t, const AmiVector &x, const AmiVector &dx, int ip, const AmiVector &sx, const AmiVector &sdx, AmiVector &sxdot) override

Sensitivity Residual function.

Parameters
  • t – time

  • x – state

  • dx – time derivative of state (DAE only)

  • ip – parameter index

  • sx – sensitivity state

  • sdx – time derivative of sensitivity state (DAE only)

  • sxdot – array to which values of the sensitivity residual function will be written

void fsxdot(realtype t, const_N_Vector x, const_N_Vector dx, int ip, const_N_Vector sx, const_N_Vector sdx, N_Vector sxdot)

Right hand side of differential equation for state sensitivities sx.

Parameters
  • t – timepoint

  • x – Vector with the states

  • dx – Vector with the derivative states

  • ip – parameter index

  • sx – Vector with the state sensitivities

  • sdx – Vector with the derivative state sensitivities

  • sxdot – Vector with the sensitivity right hand side

void fM(realtype t, const_N_Vector x)

Mass matrix for DAE systems.

Parameters
  • t – timepoint

  • x – Vector with the states

virtual std::unique_ptr<Solver> getSolver() override

Retrieves the solver object.

Returns

The Solver instance

Protected Functions

virtual void fJSparse(SUNMatrixContent_Sparse JSparse, realtype t, const realtype *x, const double *p, const double *k, const realtype *h, realtype cj, const realtype *dx, const realtype *w, const realtype *dwdx) = 0

Model specific implementation for fJSparse.

Parameters
  • JSparse – Matrix to which the Jacobian will be written

  • t – timepoint

  • x – Vector with the states

  • p – parameter vector

  • k – constants vector

  • h – Heaviside vector

  • cj – scaling factor, inverse of the step size

  • dx – Vector with the derivative states

  • w – vector with helper variables

  • dwdx – derivative of w wrt x

virtual void froot(realtype *root, realtype t, const realtype *x, const double *p, const double *k, const realtype *h, const realtype *dx)

Model specific implementation for froot.

Parameters
  • root – values of the trigger function

  • t – timepoint

  • x – Vector with the states

  • p – parameter vector

  • k – constants vector

  • h – Heaviside vector

  • dx – Vector with the derivative states

virtual void fxdot(realtype *xdot, realtype t, const realtype *x, const double *p, const double *k, const realtype *h, const realtype *dx, const realtype *w) = 0

Model specific implementation for fxdot.

Parameters
  • xdot – residual function

  • t – timepoint

  • x – Vector with the states

  • p – parameter vector

  • k – constants vector

  • h – Heaviside vector

  • w – vector with helper variables

  • dx – Vector with the derivative states

virtual void fdxdotdp(realtype *dxdotdp, realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, int ip, const realtype *dx, const realtype *w, const realtype *dwdp)

Model specific implementation of fdxdotdp.

Parameters
  • dxdotdp – partial derivative xdot wrt p

  • t – timepoint

  • x – Vector with the states

  • p – parameter vector

  • k – constants vector

  • h – Heaviside vector

  • ip – parameter index

  • dx – Vector with the derivative states

  • w – vector with helper variables

  • dwdp – derivative of w wrt p

virtual void fM(realtype *M, const realtype t, const realtype *x, const realtype *p, const realtype *k)

Model specific implementation of fM.

Parameters
  • M – mass matrix

  • t – timepoint

  • x – Vector with the states

  • p – parameter vector

  • k – constants vector