Class AbstractModel
Defined in File abstract_model.h
Inheritance Relationships
Derived Type
public amici::Model
(Class Model)
Class Documentation
-
class amici::AbstractModel
Abstract base class of amici::Model defining functions that need to be implemented in an AMICI model.
Some functions have empty default implementations or throw. This class shall not have any data members.
Subclassed by amici::Model
Public Functions
-
virtual ~AbstractModel() = default
-
virtual std::unique_ptr<Solver> getSolver() = 0
Retrieves the solver object.
- Returns
The Solver instance
-
virtual void froot(const realtype t, const AmiVector &x, const AmiVector &dx, gsl::span<realtype> root) = 0
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
-
virtual void fxdot(const realtype t, const AmiVector &x, const AmiVector &dx, AmiVector &xdot) = 0
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
-
virtual void fsxdot(const realtype t, const AmiVector &x, const AmiVector &dx, int ip, const AmiVector &sx, const AmiVector &sdx, AmiVector &sxdot) = 0
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
-
virtual void fxBdot_ss(const realtype t, const AmiVector &xB, const AmiVector &dxB, AmiVector &xBdot) = 0
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
-
virtual void fJSparseB_ss(SUNMatrix JB) = 0
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) = 0
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
-
virtual void fJ(const realtype t, realtype cj, const AmiVector &x, const AmiVector &dx, const AmiVector &xdot, SUNMatrix J) = 0
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
-
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) = 0
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
-
virtual void fJSparse(const realtype t, realtype cj, const AmiVector &x, const AmiVector &dx, const AmiVector &xdot, SUNMatrix J) = 0
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
-
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) = 0
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
-
virtual void fJDiag(const realtype t, AmiVector &Jdiag, realtype cj, const AmiVector &x, const AmiVector &dx) = 0
Diagonal Jacobian function.
- Parameters
t – time
Jdiag – array to which the diagonal of the Jacobian will be written
cj – scaling factor (inverse of timestep, DAE only)
x – state
dx – time derivative of state (DAE only)
-
virtual void fdxdotdp(const realtype t, const AmiVector &x, const AmiVector &dx) = 0
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 fJv(const realtype t, const AmiVector &x, const AmiVector &dx, const AmiVector &xdot, const AmiVector &v, AmiVector &nJv, realtype cj) = 0
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)
-
virtual std::string getAmiciVersion() const
Returns the AMICI version that was used to generate the model.
- Returns
AMICI version string
-
virtual std::string getAmiciCommit() const
Returns the AMICI commit that was used to generate the model.
- Returns
AMICI commit string
-
virtual void fx0(realtype *x0, const realtype t, const realtype *p, const realtype *k)
Model-specific implementation of fx0.
- Parameters
x0 – initial state
t – initial time
p – parameter vector
k – constant vector
-
virtual bool isFixedParameterStateReinitializationAllowed() const
Function indicating whether reinitialization of states depending on fixed parameters is permissible.
- Returns
flag indicating whether reinitialization of states depending on fixed parameters is permissible
-
virtual void fx0_fixedParameters(realtype *x0, const realtype t, const realtype *p, const realtype *k, gsl::span<const int> reinitialization_state_idxs)
Model-specific implementation of fx0_fixedParameters.
- Parameters
x0 – initial state
t – initial time
p – parameter vector
k – constant vector
reinitialization_state_idxs – Indices of states to be reinitialized based on provided constants / fixed parameters.
-
virtual void fsx0_fixedParameters(realtype *sx0, const realtype t, const realtype *x0, const realtype *p, const realtype *k, int ip, gsl::span<const int> reinitialization_state_idxs)
Model-specific implementation of fsx0_fixedParameters.
- Parameters
sx0 – initial state sensitivities
t – initial time
x0 – initial state
p – parameter vector
k – constant vector
ip – sensitivity index
reinitialization_state_idxs – Indices of states to be reinitialized based on provided constants / fixed parameters.
-
virtual void fsx0(realtype *sx0, const realtype t, const realtype *x0, const realtype *p, const realtype *k, int ip)
Model-specific implementation of fsx0.
- Parameters
sx0 – initial state sensitivities
t – initial time
x0 – initial state
p – parameter vector
k – constant vector
ip – sensitivity index
-
virtual void fdx0(AmiVector &x0, AmiVector &dx0)
Initial value for time derivative of states (only necessary for DAEs)
- Parameters
x0 – Vector with the initial states
dx0 – Vector to which the initial derivative states will be written (only DAE)
-
virtual void fstau(realtype *stau, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl, const realtype *sx, int ip, int ie)
Model-specific implementation of fstau.
- Parameters
stau – total derivative of event timepoint
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
tcl – total abundances for conservation laws
sx – current state sensitivity
ip – sensitivity index
ie – event index
-
virtual void fy(realtype *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w)
Model-specific implementation of fy.
- Parameters
y – model output at current timepoint
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
w – repeating elements vector
-
virtual void fdydp(realtype *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, int ip, const realtype *w, const realtype *dwdp)
Model-specific implementation of fdydp (MATLAB-only)
- Parameters
dydp – partial derivative of observables y w.r.t. model parameters p
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
ip – parameter index w.r.t. which the derivative is requested
w – repeating elements vector
dwdp – Recurring terms in xdot, parameter derivative
-
virtual void fdydp(realtype *dydp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, int ip, const realtype *w, const realtype *tcl, const realtype *dtcldp)
Model-specific implementation of fdydp (Python)
- Parameters
dydp – partial derivative of observables y w.r.t. model parameters p
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
ip – parameter index w.r.t. which the derivative is requested
w – repeating elements vector
tcl – total abundances for conservation laws
dtcldp – Sensitivities of total abundances for conservation laws
-
virtual void fdydx(realtype *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx)
Model-specific implementation of fdydx.
- Parameters
dydx – partial derivative of observables y w.r.t. model states x
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
w – repeating elements vector
dwdx – Recurring terms in xdot, state derivative
-
virtual void fz(realtype *z, int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h)
Model-specific implementation of fz.
- Parameters
z – value of event output
ie – event index
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
-
virtual void fsz(realtype *sz, int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *sx, int ip)
Model-specific implementation of fsz.
- Parameters
sz – Sensitivity of rz, total derivative
ie – event index
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
sx – current state sensitivity
ip – sensitivity index
-
virtual void frz(realtype *rz, int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h)
Model-specific implementation of frz.
- Parameters
rz – value of root function at current timepoint (non-output events not included)
ie – event index
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
-
virtual void fsrz(realtype *srz, int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *sx, int ip)
Model-specific implementation of fsrz.
- Parameters
srz – Sensitivity of rz, total derivative
ie – event index
t – current time
x – current state
p – parameter vector
k – constant vector
sx – current state sensitivity
h – Heaviside vector
ip – sensitivity index
-
virtual void fdzdp(realtype *dzdp, int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, int ip)
Model-specific implementation of fdzdp.
- Parameters
dzdp – partial derivative of event-resolved output z w.r.t. model parameters p
ie – event index
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
ip – parameter index w.r.t. which the derivative is requested
-
virtual void fdzdx(realtype *dzdx, int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h)
Model-specific implementation of fdzdx.
- Parameters
dzdx – partial derivative of event-resolved output z w.r.t. model states x
ie – event index
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
-
virtual void fdrzdp(realtype *drzdp, int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, int ip)
Model-specific implementation of fdrzdp.
- Parameters
drzdp – partial derivative of root output rz w.r.t. model parameters p
ie – event index
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
ip – parameter index w.r.t. which the derivative is requested
-
virtual void fdrzdx(realtype *drzdx, int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h)
Model-specific implementation of fdrzdx.
- Parameters
drzdx – partial derivative of root output rz w.r.t. model states x
ie – event index
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
-
virtual void fdeltax(realtype *deltax, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, int ie, const realtype *xdot, const realtype *xdot_old)
Model-specific implementation of fdeltax.
- Parameters
deltax – state update
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
ie – event index
xdot – new model right hand side
xdot_old – previous model right hand side
-
virtual void fdeltasx(realtype *deltasx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, int ip, int ie, const realtype *xdot, const realtype *xdot_old, const realtype *sx, const realtype *stau, const realtype *tcl)
Model-specific implementation of fdeltasx.
- Parameters
deltasx – sensitivity update
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
w – repeating elements vector
ip – sensitivity index
ie – event index
xdot – new model right hand side
xdot_old – previous model right hand side
sx – state sensitivity
stau – event-time sensitivity
tcl – total abundances for conservation laws
-
virtual void fdeltaxB(realtype *deltaxB, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, int ie, const realtype *xdot, const realtype *xdot_old, const realtype *xB)
Model-specific implementation of fdeltaxB.
- Parameters
deltaxB – adjoint state update
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
ie – event index
xdot – new model right hand side
xdot_old – previous model right hand side
xB – current adjoint state
-
virtual void fdeltaqB(realtype *deltaqB, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, int ip, int ie, const realtype *xdot, const realtype *xdot_old, const realtype *xB)
Model-specific implementation of fdeltaqB.
- Parameters
deltaqB – sensitivity update
t – current time
x – current state
p – parameter vector
k – constant vector
h – Heaviside vector
ip – sensitivity index
ie – event index
xdot – new model right hand side
xdot_old – previous model right hand side
xB – adjoint state
-
virtual void fsigmay(realtype *sigmay, const realtype t, const realtype *p, const realtype *k, const realtype *y)
Model-specific implementation of fsigmay.
- Parameters
sigmay – standard deviation of measurements
t – current time
p – parameter vector
k – constant vector
y – model output at timepoint t
-
virtual void fdsigmaydp(realtype *dsigmaydp, const realtype t, const realtype *p, const realtype *k, const realtype *y, int ip)
Model-specific implementation of fdsigmaydp.
- Parameters
dsigmaydp – partial derivative of standard deviation of measurements
t – current time
p – parameter vector
k – constant vector
y – model output at timepoint t
ip – sensitivity index
-
virtual void fdsigmaydy(realtype *dsigmaydy, const realtype t, const realtype *p, const realtype *k, const realtype *y)
Model-specific implementation of fsigmay.
- Parameters
dsigmaydy – partial derivative of standard deviation of measurements w.r.t. model outputs
t – current time
p – parameter vector
k – constant vector
y – model output at timepoint t
-
virtual void fsigmaz(realtype *sigmaz, const realtype t, const realtype *p, const realtype *k)
Model-specific implementation of fsigmaz.
- Parameters
sigmaz – standard deviation of event measurements
t – current time
p – parameter vector
k – constant vector
-
virtual void fdsigmazdp(realtype *dsigmazdp, const realtype t, const realtype *p, const realtype *k, int ip)
Model-specific implementation of fsigmaz.
- Parameters
dsigmazdp – partial derivative of standard deviation of event measurements
t – current time
p – parameter vector
k – constant vector
ip – sensitivity index
-
virtual void fJy(realtype *nllh, int iy, const realtype *p, const realtype *k, const realtype *y, const realtype *sigmay, const realtype *my)
Model-specific implementation of fJy.
- Parameters
nllh – negative log-likelihood for measurements y
iy – output index
p – parameter vector
k – constant vector
y – model output at timepoint
sigmay – measurement standard deviation at timepoint
my – measurements at timepoint
-
virtual void fJz(realtype *nllh, int iz, const realtype *p, const realtype *k, const realtype *z, const realtype *sigmaz, const realtype *mz)
Model-specific implementation of fJz.
- Parameters
nllh – negative log-likelihood for event measurements z
iz – event output index
p – parameter vector
k – constant vector
z – model event output at timepoint
sigmaz – event measurement standard deviation at timepoint
mz – event measurements at timepoint
-
virtual void fJrz(realtype *nllh, int iz, const realtype *p, const realtype *k, const realtype *z, const realtype *sigmaz)
Model-specific implementation of fJrz.
- Parameters
nllh – regularization for event measurements z
iz – event output index
p – parameter vector
k – constant vector
z – model event output at timepoint
sigmaz – event measurement standard deviation at timepoint
-
virtual void fdJydy(realtype *dJydy, int iy, const realtype *p, const realtype *k, const realtype *y, const realtype *sigmay, const realtype *my)
Model-specific implementation of fdJydy.
- Parameters
dJydy – partial derivative of time-resolved measurement negative log-likelihood Jy
iy – output index
p – parameter vector
k – constant vector
y – model output at timepoint
sigmay – measurement standard deviation at timepoint
my – measurement at timepoint
-
virtual void fdJydy_colptrs(SUNMatrixWrapper &dJydy, int index)
Model-specific implementation of fdJydy colptrs.
- Parameters
dJydy – sparse matrix to which colptrs will be written
index – ytrue index
-
virtual void fdJydy_rowvals(SUNMatrixWrapper &dJydy, int index)
Model-specific implementation of fdJydy rowvals.
- Parameters
dJydy – sparse matrix to which rowvals will be written
index –
ytrue
index
-
virtual void fdJydsigma(realtype *dJydsigma, int iy, const realtype *p, const realtype *k, const realtype *y, const realtype *sigmay, const realtype *my)
Model-specific implementation of fdJydsigma.
- Parameters
dJydsigma – Sensitivity of time-resolved measurement negative log-likelihood Jy w.r.t. standard deviation sigmay
iy – output index
p – parameter vector
k – constant vector
y – model output at timepoint
sigmay – measurement standard deviation at timepoint
my – measurement at timepoint
-
virtual void fdJzdz(realtype *dJzdz, int iz, const realtype *p, const realtype *k, const realtype *z, const realtype *sigmaz, const realtype *mz)
Model-specific implementation of fdJzdz.
- Parameters
dJzdz – partial derivative of event measurement negative log-likelihood Jz
iz – event output index
p – parameter vector
k – constant vector
z – model event output at timepoint
sigmaz – event measurement standard deviation at timepoint
mz – event measurement at timepoint
-
virtual void fdJzdsigma(realtype *dJzdsigma, int iz, const realtype *p, const realtype *k, const realtype *z, const realtype *sigmaz, const realtype *mz)
Model-specific implementation of fdJzdsigma.
- Parameters
dJzdsigma – Sensitivity of event measurement negative log-likelihood Jz w.r.t. standard deviation sigmaz
iz – event output index
p – parameter vector
k – constant vector
z – model event output at timepoint
sigmaz – event measurement standard deviation at timepoint
mz – event measurement at timepoint
-
virtual void fdJrzdz(realtype *dJrzdz, int iz, const realtype *p, const realtype *k, const realtype *rz, const realtype *sigmaz)
Model-specific implementation of fdJrzdz.
- Parameters
dJrzdz – partial derivative of event penalization Jrz
iz – event output index
p – parameter vector
k – constant vector
rz – model root output at timepoint
sigmaz – event measurement standard deviation at timepoint
-
virtual void fdJrzdsigma(realtype *dJrzdsigma, int iz, const realtype *p, const realtype *k, const realtype *rz, const realtype *sigmaz)
Model-specific implementation of fdJrzdsigma.
- Parameters
dJrzdsigma – Sensitivity of event penalization Jrz w.r.t. standard deviation sigmaz
iz – event output index
p – parameter vector
k – constant vector
rz – model root output at timepoint
sigmaz – event measurement standard deviation at timepoint
-
virtual void fw(realtype *w, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl)
Model-specific implementation of fw.
- Parameters
w – Recurring terms in xdot
t – timepoint
x – vector with the states
p – parameter vector
k – constants vector
h – Heaviside vector
tcl – total abundances for conservation laws
-
virtual void fdwdp(realtype *dwdp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *tcl, const realtype *stcl)
Model-specific sparse implementation of dwdp.
- Parameters
dwdp – Recurring terms in xdot, parameter derivative
t – timepoint
x – vector with the states
p – parameter vector
k – constants vector
h – Heaviside vector
w – vector with helper variables
tcl – total abundances for conservation laws
stcl – sensitivities of total abundances for conservation laws
-
virtual void fdwdp_colptrs(SUNMatrixWrapper &dwdp)
Model-specific implementation for dwdp, column pointers.
- Parameters
dwdp – sparse matrix to which colptrs will be written
-
virtual void fdwdp_rowvals(SUNMatrixWrapper &dwdp)
Model-specific implementation for dwdp, row values.
- Parameters
dwdp – sparse matrix to which rowvals will be written
-
virtual void fdwdp(realtype *dwdp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *tcl, const realtype *stcl, int ip)
Model-specific sensitivity implementation of dwdp.
- Parameters
dwdp – Recurring terms in xdot, parameter derivative
t – timepoint
x – vector with the states
p – parameter vector
k – constants vector
h – Heaviside vector
w – vector with helper variables
tcl – total abundances for conservation laws
stcl – sensitivities of total abundances for conservation laws
ip – sensitivity parameter index
-
virtual void fdwdx(realtype *dwdx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *tcl)
Model-specific implementation of dwdx, data part.
- Parameters
dwdx – Recurring terms in xdot, state derivative
t – timepoint
x – vector with the states
p – parameter vector
k – constants vector
h – Heaviside vector
w – vector with helper variables
tcl – total abundances for conservation laws
-
virtual void fdwdx_colptrs(SUNMatrixWrapper &dwdx)
Model-specific implementation for dwdx, column pointers.
- Parameters
dwdx – sparse matrix to which colptrs will be written
-
virtual void fdwdx_rowvals(SUNMatrixWrapper &dwdx)
Model-specific implementation for dwdx, row values.
- Parameters
dwdx – sparse matrix to which rowvals will be written
-
virtual void fdwdw(realtype *dwdw, realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *tcl)
Model-specific implementation of fdwdw, no w chainrule (Py)
- Parameters
dwdw – partial derivative w wrt w
t – timepoint
x – Vector with the states
p – parameter vector
k – constants vector
h – Heaviside vector
w – vector with helper variables
tcl – Total abundances for conservation laws
-
virtual void fdwdw_colptrs(SUNMatrixWrapper &dwdw)
Model-specific implementation of fdwdw, colptrs part.
- Parameters
dwdw – sparse matrix to which colptrs will be written
-
virtual void fdwdw_rowvals(SUNMatrixWrapper &dwdw)
Model-specific implementation of fdwdw, rowvals part.
- Parameters
dwdw – sparse matrix to which rowvals will be written
-
virtual void fdx_rdatadx_solver(realtype *dx_rdatadx_solver, const realtype *x, const realtype *tcl, const realtype *p, const realtype *k)
Compute dx_rdata / dx_solver.
- Parameters
dx_rdatadx_solver – dx_rdata / dx_solver
p – parameter vector
k – constant vector
x – State variables with conservation laws applied
tcl – Total abundances for conservation laws
-
virtual void fdx_rdatadx_solver_colptrs(SUNMatrixWrapper &dxrdatadxsolver)
Model-specific implementation of fdx_rdatadx_solver, colptrs part.
- Parameters
dxrdatadxsolver – sparse matrix to which colptrs will be written
-
virtual void fdx_rdatadx_solver_rowvals(SUNMatrixWrapper &dxrdatadxsolver)
Model-specific implementation of fdx_rdatadx_solver, rowvals part.
- Parameters
dxrdatadxsolver – sparse matrix to which rowvals will be written
-
virtual void fdx_rdatadp(realtype *dx_rdatadp, const realtype *x, const realtype *tcl, const realtype *p, const realtype *k, const int ip)
Compute dx_rdata / dp.
- Parameters
dx_rdatadp – dx_rdata / dp
p – parameter vector
k – constant vector
x – State variables with conservation laws applied
tcl – Total abundances for conservation laws
ip – Sensitivity index
-
virtual void fdx_rdatadtcl(realtype *dx_rdatadtcl, const realtype *x, const realtype *tcl, const realtype *p, const realtype *k)
Compute dx_rdata / dtcl.
- Parameters
dx_rdatadtcl – dx_rdata / dtcl
p – parameter vector
k – constant vector
x – State variables with conservation laws applied
tcl – Total abundances for conservation laws
-
virtual void fdx_rdatadtcl_colptrs(SUNMatrixWrapper &dx_rdatadtcl)
Model-specific implementation of fdx_rdatadtcl, colptrs part.
- Parameters
dx_rdatadtcl – sparse matrix to which colptrs will be written
-
virtual void fdx_rdatadtcl_rowvals(SUNMatrixWrapper &dx_rdatadtcl)
Model-specific implementation of fdx_rdatadtcl, rowvals part.
- Parameters
dx_rdatadtcl – sparse matrix to which rowvals will be written
-
virtual void fdtotal_cldp(realtype *dtotal_cldp, const realtype *x_rdata, const realtype *p, const realtype *k, const int ip)
Compute dtotal_cl / dp.
- Parameters
dtotal_cldp – dtotal_cl / dp
x_rdata – State variables with conservation laws applied
p – parameter vector
k – constant vector
ip – Sensitivity index
-
virtual void fdtotal_cldx_rdata(realtype *dtotal_cldx_rdata, const realtype *x_rdata, const realtype *p, const realtype *k, const realtype *tcl)
Compute dtotal_cl / dx_rdata.
- Parameters
dtotal_cldx_rdata – dtotal_cl / dx_rdata
x_rdata – State variables with conservation laws applied
p – parameter vector
k – constant vector
tcl – Total abundances for conservation laws
-
virtual void fdtotal_cldx_rdata_colptrs(SUNMatrixWrapper &dtotal_cldx_rdata)
Model-specific implementation of fdtotal_cldx_rdata, colptrs part.
- Parameters
dtotal_cldx_rdata – sparse matrix to which colptrs will be written
-
virtual void fdtotal_cldx_rdata_rowvals(SUNMatrixWrapper &dtotal_cldx_rdata)
Model-specific implementation of fdtotal_cldx_rdata, rowvals part.
- Parameters
dtotal_cldx_rdata – sparse matrix to which rowvals will be written
-
virtual ~AbstractModel() = default