Class AbstractModel
Defined in File abstract_model.h
Inheritance Relationships
Derived Type
public amici::Model
(Class Model)
Class Documentation
-
class 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(realtype const t, AmiVector const &x, AmiVector const &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(realtype const t, AmiVector const &x, AmiVector const &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(realtype const t, AmiVector const &x, AmiVector const &dx, int ip, AmiVector const &sx, AmiVector const &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(realtype const t, AmiVector const &xB, AmiVector const &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(realtype const t, realtype cj, AmiVector const &x, AmiVector const &dx, AmiVector const &xB, AmiVector const &dxB, AmiVector const &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(realtype const t, realtype cj, AmiVector const &x, AmiVector const &dx, AmiVector const &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(realtype const t, realtype cj, AmiVector const &x, AmiVector const &dx, AmiVector const &xB, AmiVector const &dxB, AmiVector const &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(realtype const t, realtype cj, AmiVector const &x, AmiVector const &dx, AmiVector const &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(realtype const t, realtype cj, AmiVector const &x, AmiVector const &dx, AmiVector const &xB, AmiVector const &dxB, AmiVector const &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(realtype const t, AmiVector &Jdiag, realtype cj, AmiVector const &x, AmiVector const &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(realtype const t, AmiVector const &x, AmiVector const &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(realtype const t, AmiVector const &x, AmiVector const &dx, AmiVector const &xdot, AmiVector const &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, realtype const t, realtype const *p, realtype const *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, realtype const t, realtype const *p, realtype const *k, gsl::span<int const> 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, realtype const t, realtype const *x0, realtype const *p, realtype const *k, int ip, gsl::span<int const> 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, realtype const t, realtype const *x0, realtype const *p, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *tcl, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, int ip, realtype const *w, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, int ip, realtype const *w, realtype const *tcl, realtype const *dtcldp, realtype const *spl, realtype const *sspl)
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
spl – spline value vector
sspl – sensitivities of spline values vector w.r.t. parameters \( p \)
-
virtual void fdydx(realtype *dydx, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *w, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, int ie, realtype const *xdot, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *w, int ip, int ie, realtype const *xdot, realtype const *xdot_old, realtype const *sx, realtype const *stau, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, int ie, realtype const *xdot, realtype const *xdot_old, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, int ip, int ie, realtype const *xdot, realtype const *xdot_old, realtype const *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, realtype const t, realtype const *p, realtype const *k, realtype const *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, realtype const t, realtype const *p, realtype const *k, realtype const *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, realtype const t, realtype const *p, realtype const *k, realtype const *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, realtype const t, realtype const *p, realtype const *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, realtype const t, realtype const *p, realtype const *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, realtype const *p, realtype const *k, realtype const *y, realtype const *sigmay, realtype const *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, realtype const *p, realtype const *k, realtype const *z, realtype const *sigmaz, realtype const *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, realtype const *p, realtype const *k, realtype const *z, realtype const *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, realtype const *p, realtype const *k, realtype const *y, realtype const *sigmay, realtype const *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, realtype const *p, realtype const *k, realtype const *y, realtype const *sigmay, realtype const *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, realtype const *p, realtype const *k, realtype const *z, realtype const *sigmaz, realtype const *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, realtype const *p, realtype const *k, realtype const *z, realtype const *sigmaz, realtype const *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, realtype const *p, realtype const *k, realtype const *rz, realtype const *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, realtype const *p, realtype const *k, realtype const *rz, realtype const *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, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *tcl, realtype const *spl, bool include_static = true)
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
spl – spline value vector
include_static – Whether to (re-)evaluate only dynamic expressions (false) or also static expressions (true). Dynamic expressions are those that depend directly or indirectly on time, static expressions are those that don’t.
-
virtual void fdwdp(realtype *dwdp, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *w, realtype const *tcl, realtype const *stcl, realtype const *spl, realtype const *sspl, bool include_static = true)
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
spl – spline value vector
sspl – sensitivities of spline values vector w.r.t. parameters \( p \)
include_static – Whether to (re-)evaluate only dynamic expressions (false) or also static expressions (true). Dynamic expressions are those that depend directly or indirectly on time, static expressions are those that don’t.
-
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 fdwdx(realtype *dwdx, realtype const t, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *w, realtype const *tcl, realtype const *spl, bool include_static = true)
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
spl – spline value vector
include_static – Whether to (re-)evaluate only dynamic expressions (false) or also static expressions (true). Dynamic expressions are those that depend directly or indirectly on time, static expressions are those that don’t.
-
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, realtype const *x, realtype const *p, realtype const *k, realtype const *h, realtype const *w, realtype const *tcl, bool include_static = true)
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
include_static – Whether to (re-)evaluate only dynamic expressions (false) or also static expressions (true). Dynamic expressions are those that depend directly or indirectly on time, static expressions are those that don’t.
-
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, realtype const *x, realtype const *tcl, realtype const *p, realtype const *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, realtype const *x, realtype const *tcl, realtype const *p, realtype const *k, int const 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, realtype const *x, realtype const *tcl, realtype const *p, realtype const *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, realtype const *x_rdata, realtype const *p, realtype const *k, int const 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, realtype const *x_rdata, realtype const *p, realtype const *k, realtype const *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 std::vector<HermiteSpline> fcreate_splines(realtype const *p, realtype const *k)
Model-specific implementation of spline creation.
- Parameters:
p – parameter vector
k – constants vector
- Returns:
Vector of splines used in the model
-
virtual void fdspline_valuesdp(realtype *dspline_valuesdp, realtype const *p, realtype const *k, int const ip)
Model-specific implementation the parametric derivatives of spline node values.
- Parameters:
dspline_valuesdp – vector to which derivatives will be written
p – parameter vector
k – constants vector
ip – Sensitivity index
-
virtual void fdspline_slopesdp(realtype *dspline_slopesdp, realtype const *p, realtype const *k, int const ip)
Model-specific implementation the parametric derivatives of slopevalues at spline nodes.
- Parameters:
dspline_slopesdp – vector to which derivatives will be written
p – parameter vector
k – constants vector
ip – Sensitivity index
-
virtual ~AbstractModel() = default