Class Model
Defined in File model.h
Inheritance Relationships
Base Types
public amici::AbstractModel
(Class AbstractModel)public amici::ModelDimensions
(Struct ModelDimensions)
Derived Types
public amici::Model_DAE
(Class Model_DAE)public amici::Model_ODE
(Class Model_ODE)
Class Documentation
-
class Model : public amici::AbstractModel, public amici::ModelDimensions
The Model class represents an AMICI ODE/DAE model.
The model can compute various model related quantities based on symbolically generated code.
Subclassed by amici::Model_DAE, amici::Model_ODE
Public Functions
-
Model() = default
Default constructor
-
Model(ModelDimensions const &model_dimensions, SimulationParameters simulation_parameters, amici::SecondOrderMode o2mode, std::vector<amici::realtype> idlist, std::vector<int> z2event, bool pythonGenerated = false, int ndxdotdp_explicit = 0, int ndxdotdx_explicit = 0, int w_recursion_depth = 0)
Constructor with model dimensions.
- Parameters
model_dimensions – Model 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 in
dxdotdp_explicit
ndxdotdx_explicit – Number of nonzero elements in
dxdotdx_explicit
w_recursion_depth – Recursion depth of fw
-
~Model() override = default
Destructor.
-
Model &operator=(Model const &other) = delete
Copy assignment is disabled until const members are removed.
- Parameters
other – Object to copy from
- Returns
-
void initialize(AmiVector &x, AmiVector &dx, AmiVectorArray &sx, AmiVectorArray &sdx, bool computeSensitivities, std::vector<int> &roots_found)
Initialize model properties.
- Parameters
x – Reference to state variables
dx – Reference to time derivative of states (DAE only)
sx – Reference to state variable sensitivities
sdx – Reference to time derivative of state sensitivities (DAE only)
computeSensitivities – Flag indicating whether sensitivities are to be computed
roots_found – boolean indicators indicating whether roots were found at t0 by this fun
-
void initializeB(AmiVector &xB, AmiVector &dxB, AmiVector &xQB, bool posteq) const
Initialize model properties.
- Parameters
xB – Adjoint state variables
dxB – Time derivative of adjoint states (DAE only)
xQB – Adjoint quadratures
posteq – Flag indicating whether postequilibration was performed
-
void initializeStates(AmiVector &x)
Initialize initial states.
- Parameters
x – State vector to be initialized
-
void initializeStateSensitivities(AmiVectorArray &sx, const AmiVector &x)
Initialize initial state sensitivities.
- Parameters
sx – Reference to state variable sensitivities
x – Reference to state variables
-
void initEvents(const AmiVector &x, const AmiVector &dx, std::vector<int> &roots_found)
Initialize the Heaviside variables
h
at the initial timet0
.Heaviside variables activate/deactivate on event occurrences.
- Parameters
x – Reference to state variables
dx – Reference to time derivative of states (DAE only)
roots_found – boolean indicators indicating whether roots were found at t0 by this fun
-
int nplist() const
Get number of parameters wrt to which sensitivities are computed.
- Returns
Length of sensitivity index vector
-
int np() const
Get total number of model parameters.
- Returns
Length of parameter vector
-
int nk() const
Get number of constants.
- Returns
Length of constant vector
-
int ncl() const
Get number of conservation laws.
- Returns
Number of conservation laws (i.e., difference between
nx_rdata
andnx_solver
).
-
int nx_reinit() const
Get number of solver states subject to reinitialization.
- Returns
Model member
nx_solver_reinit
-
const double *k() const
Get fixed parameters.
- Returns
Pointer to constants array
-
int nMaxEvent() const
Get maximum number of events that may occur for each type.
- Returns
Maximum number of events that may occur for each type
-
void setNMaxEvent(int nmaxevent)
Set maximum number of events that may occur for each type.
- Parameters
nmaxevent – Maximum number of events that may occur for each type
-
int nt() const
Get number of timepoints.
- Returns
Number of timepoints
-
std::vector<ParameterScaling> const &getParameterScale() const
Get parameter scale for each parameter.
- Returns
Vector of parameter scales
-
void setParameterScale(ParameterScaling pscale)
Set parameter scale for each parameter.
NOTE: Resets initial state sensitivities.
- Parameters
pscale – Scalar parameter scale to be set for all parameters
-
void setParameterScale(const std::vector<ParameterScaling> &pscaleVec)
Set parameter scale for each parameter.
NOTE: Resets initial state sensitivities.
- Parameters
pscaleVec – Vector of parameter scales
-
std::vector<realtype> const &getUnscaledParameters() const
Get parameters with transformation according to parameter scale applied.
- Returns
Unscaled parameters
-
std::vector<realtype> const &getParameters() const
Get parameter vector.
- Returns
The user-set parameters (see also
Model::getUnscaledParameters
)
-
realtype getParameterById(std::string const &par_id) const
Get value of first model parameter with the specified ID.
- Parameters
par_id – Parameter ID
- Returns
Parameter value
-
realtype getParameterByName(std::string const &par_name) const
Get value of first model parameter with the specified name.
- Parameters
par_name – Parameter name
- Returns
Parameter value
-
void setParameters(std::vector<realtype> const &p)
Set the parameter vector.
- Parameters
p – Vector of parameters
-
void setParameterById(std::map<std::string, realtype> const &p, bool ignoreErrors = false)
Set model parameters according to the parameter IDs and mapped values.
- Parameters
p – Map of parameters IDs and values
ignoreErrors – Ignore errors such as parameter IDs in p which are not model parameters
-
void setParameterById(std::string const &par_id, realtype value)
Set value of first model parameter with the specified ID.
- Parameters
par_id – Parameter ID
value – Parameter value
-
int setParametersByIdRegex(std::string const &par_id_regex, realtype value)
Set all values of model parameters with IDs matching the specified regular expression.
- Parameters
par_id_regex – Parameter ID regex
value – Parameter value
- Returns
Number of parameter IDs that matched the regex
-
void setParameterByName(std::string const &par_name, realtype value)
Set value of first model parameter with the specified name.
- Parameters
par_name – Parameter name
value – Parameter value
-
void setParameterByName(std::map<std::string, realtype> const &p, bool ignoreErrors = false)
Set model parameters according to the parameter name and mapped values.
- Parameters
p – Map of parameters names and values
ignoreErrors – Ignore errors such as parameter names in p which are not model parameters
-
int setParametersByNameRegex(std::string const &par_name_regex, realtype value)
Set all values of all model parameters with names matching the specified regex.
- Parameters
par_name_regex – Parameter name regex
value – Parameter value
- Returns
Number of fixed parameter names that matched the regex
-
std::vector<realtype> const &getFixedParameters() const
Get values of fixed parameters.
- Returns
Vector of fixed parameters with same ordering as in Model::getFixedParameterIds
-
realtype getFixedParameterById(std::string const &par_id) const
Get value of fixed parameter with the specified ID.
- Parameters
par_id – Parameter ID
- Returns
Parameter value
-
realtype getFixedParameterByName(std::string const &par_name) const
Get value of fixed parameter with the specified name.
If multiple parameters have the same name, the first parameter with matching name is returned.
- Parameters
par_name – Parameter name
- Returns
Parameter value
-
void setFixedParameters(std::vector<realtype> const &k)
Set values for constants.
- Parameters
k – Vector of fixed parameters
-
void setFixedParameterById(std::string const &par_id, realtype value)
Set value of first fixed parameter with the specified ID.
- Parameters
par_id – Fixed parameter id
value – Fixed parameter value
-
int setFixedParametersByIdRegex(std::string const &par_id_regex, realtype value)
Set values of all fixed parameters with the ID matching the specified regex.
- Parameters
par_id_regex – Fixed parameter name regex
value – Fixed parameter value
- Returns
Number of fixed parameter IDs that matched the regex
-
void setFixedParameterByName(std::string const &par_name, realtype value)
Set value of first fixed parameter with the specified name.
- Parameters
par_name – Fixed parameter ID
value – Fixed parameter value
-
int setFixedParametersByNameRegex(std::string const &par_name_regex, realtype value)
Set value of all fixed parameters with name matching the specified regex.
- Parameters
par_name_regex – Fixed parameter name regex
value – Fixed parameter value
- Returns
Number of fixed parameter names that matched the regex
-
virtual bool hasParameterNames() const
Report whether the model has parameter names set.
- Returns
Boolean indicating whether parameter names were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getParameterNames() const
Get names of the model parameters.
- Returns
The parameter names
-
virtual bool hasStateNames() const
Report whether the model has state names set.
- Returns
Boolean indicating whether state names were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getStateNames() const
Get names of the model states.
- Returns
State names
-
virtual std::vector<std::string> getStateNamesSolver() const
Get names of the solver states.
- Returns
State names
-
virtual bool hasFixedParameterNames() const
Report whether the model has fixed parameter names set.
- Returns
Boolean indicating whether fixed parameter names were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getFixedParameterNames() const
Get names of the fixed model parameters.
- Returns
Fixed parameter names
-
virtual bool hasObservableNames() const
Report whether the model has observable names set.
- Returns
Boolean indicating whether observable names were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getObservableNames() const
Get names of the observables.
- Returns
Observable names
-
virtual bool hasExpressionNames() const
Report whether the model has expression names set.
- Returns
Boolean indicating whether expression names were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getExpressionNames() const
Get names of the expressions.
- Returns
Expression names
-
virtual bool hasParameterIds() const
Report whether the model has parameter IDs set.
- Returns
Boolean indicating whether parameter IDs were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getParameterIds() const
Get IDs of the model parameters.
- Returns
Parameter IDs
-
virtual bool hasStateIds() const
Report whether the model has state IDs set.
- Returns
Boolean indicating whether state IDs were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getStateIds() const
Get IDs of the model states.
- Returns
State IDs
-
virtual std::vector<std::string> getStateIdsSolver() const
Get IDs of the solver states.
- Returns
State IDs
-
virtual bool hasFixedParameterIds() const
Report whether the model has fixed parameter IDs set.
- Returns
Boolean indicating whether fixed parameter IDs were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getFixedParameterIds() const
Get IDs of the fixed model parameters.
- Returns
Fixed parameter IDs
-
virtual bool hasObservableIds() const
Report whether the model has observable IDs set.
- Returns
Boolean indicating whether observable ids were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getObservableIds() const
Get IDs of the observables.
- Returns
Observable IDs
-
virtual bool hasExpressionIds() const
Report whether the model has expression IDs set.
- Returns
Boolean indicating whether expression ids were set. Also returns
true
if the number of corresponding variables is just zero.
-
virtual std::vector<std::string> getExpressionIds() const
Get IDs of the expression.
- Returns
Expression IDs
-
virtual bool hasQuadraticLLH() const
Checks whether the defined noise model is gaussian, i.e., the nllh is quadratic.
- Returns
boolean flag
-
std::vector<realtype> const &getTimepoints() const
Get the timepoint vector.
- Returns
Timepoint vector
-
realtype getTimepoint(int it) const
Get simulation timepoint for time index
it
.- Parameters
it – Time index
- Returns
Timepoint
-
void setTimepoints(std::vector<realtype> const &ts)
Set the timepoint vector.
- Parameters
ts – New timepoint vector
-
double t0() const
Get simulation start time.
- Returns
Simulation start time
-
void setT0(double t0)
Set simulation start time.
- Parameters
t0 – Simulation start time
-
std::vector<bool> const &getStateIsNonNegative() const
Get flags indicating whether states should be treated as non-negative.
- Returns
Vector of flags
-
void setStateIsNonNegative(std::vector<bool> const &stateIsNonNegative)
Set flags indicating whether states should be treated as non-negative.
- Parameters
stateIsNonNegative – Vector of flags
-
void setAllStatesNonNegative()
Set flags indicating that all states should be treated as non-negative.
-
inline ModelState const &getModelState() const
Get the current model state.
- Returns
Current model state
-
inline void setModelState(ModelState const &state)
Set the current model state.
- Parameters
state – Model state
-
inline void setMinimumSigmaResiduals(double min_sigma)
Sets the estimated lower boundary for sigma_y. When :meth:
setAddSigmaResiduals
is activated, this lower boundary must ensure that log(sigma) + min_sigma > 0.- Parameters
min_sigma – lower boundary
-
inline realtype getMinimumSigmaResiduals() const
Gets the specified estimated lower boundary for sigma_y.
- Returns
lower boundary
-
inline void setAddSigmaResiduals(bool sigma_res)
Specifies whether residuals should be added to account for parameter dependent sigma.
If set to true, additional residuals of the form \( \sqrt{\log(\sigma) + C} \) will be added. This enables least-squares optimization for variables with Gaussian noise assumption and parameter dependent standard deviation sigma. The constant \( C \) can be set via :meth:
setMinimumSigmaResiduals
.- Parameters
sigma_res – if true, additional residuals are added
-
inline bool getAddSigmaResiduals() const
Checks whether residuals should be added to account for parameter dependent sigma.
- Returns
sigma_res
-
std::vector<int> const &getParameterList() const
Get the list of parameters for which sensitivities are computed.
- Returns
List of parameter indices
-
int plist(int pos) const
Get entry in parameter list by index.
- Parameters
pos – Index in sensitivity parameter list
- Returns
Index in parameter list
-
void setParameterList(std::vector<int> const &plist)
Set the list of parameters for which sensitivities are to be computed.
NOTE: Resets initial state sensitivities.
- Parameters
plist – List of parameter indices
-
void setInitialStates(std::vector<realtype> const &x0)
Set the initial states.
- Parameters
x0 – Initial state vector
-
bool hasCustomInitialStates() const
Return whether custom initial states have been set.
- Returns
true
if has custom initial states, otherwisefalse
-
std::vector<realtype> getInitialStateSensitivities()
Get the initial states sensitivities.
- Returns
vector of initial state sensitivities
-
void setInitialStateSensitivities(std::vector<realtype> const &sx0)
Set the initial state sensitivities.
- Parameters
sx0 – vector of initial state sensitivities with chainrule applied. This could be a slice of ReturnData::sx or ReturnData::sx0
-
bool hasCustomInitialStateSensitivities() const
Return whether custom initial state sensitivities have been set.
- Returns
true
if has custom initial state sensitivities, otherwisefalse
.
-
void setUnscaledInitialStateSensitivities(std::vector<realtype> const &sx0)
Set the initial state sensitivities.
- Parameters
sx0 – Vector of initial state sensitivities without chainrule applied. This could be the readin from a
model.sx0data
saved to HDF5.
-
void setSteadyStateSensitivityMode(SteadyStateSensitivityMode mode)
Set the mode how sensitivities are computed in the steadystate simulation.
- Parameters
mode – Steadystate sensitivity mode
-
SteadyStateSensitivityMode getSteadyStateSensitivityMode() const
Gets the mode how sensitivities are computed in the steadystate simulation.
- Returns
Mode
-
void setReinitializeFixedParameterInitialStates(bool flag)
Set whether initial states depending on fixed parameters are to be reinitialized after preequilibration and presimulation.
- Parameters
flag – Fixed parameters reinitialized?
-
bool getReinitializeFixedParameterInitialStates() const
Get whether initial states depending on fixedParameters are to be reinitialized after preequilibration and presimulation.
- Returns
flag
true
/false
-
void requireSensitivitiesForAllParameters()
Require computation of sensitivities for all parameters p [0..np[ in natural order.
NOTE: Resets initial state sensitivities.
-
void getExpression(gsl::span<realtype> w, const realtype t, const AmiVector &x)
Get time-resolved
w
.- Parameters
w – Buffer (shape
nw
)t – Current timepoint
x – Current state
-
void getObservable(gsl::span<realtype> y, const realtype t, const AmiVector &x)
Get time-resolved observables.
- Parameters
y – Buffer (shape
ny
)t – Current timepoint
x – Current state
-
virtual ObservableScaling getObservableScaling(int iy) const
Get scaling type for observable.
- Parameters
iy – observable index
- Returns
scaling type
-
void getObservableSensitivity(gsl::span<realtype> sy, const realtype t, const AmiVector &x, const AmiVectorArray &sx)
Get sensitivity of time-resolved observables.
Total derivative \( sy = dydx * sx + dydp\) (only for forward sensitivities).
- Parameters
sy – buffer (shape
ny
xnplist
, row-major)t – Timpoint
x – State variables
sx – State sensitivities
-
void getObservableSigma(gsl::span<realtype> sigmay, const int it, const ExpData *edata)
Get time-resolved observable standard deviations.
- Parameters
sigmay – Buffer (shape
ny
)it – Timepoint index
edata – Pointer to experimental data instance (optional, pass
nullptr
to ignore)
-
void getObservableSigmaSensitivity(gsl::span<realtype> ssigmay, gsl::span<const realtype> sy, const int it, const ExpData *edata)
Sensitivity of time-resolved observable standard deviation.
Total derivative (can be used with both adjoint and forward sensitivity).
- Parameters
ssigmay – Buffer (shape
ny
xnplist
, row-major)sy – Sensitivity of time-resolved observables for current timepoint
it – Timepoint index
edata – Pointer to experimental data instance (optional, pass
nullptr
to ignore)
-
void addObservableObjective(realtype &Jy, const int it, const AmiVector &x, const ExpData &edata)
Add time-resolved measurement negative log-likelihood \( Jy \).
- Parameters
Jy – Buffer (shape 1)
it – Timepoint index
x – State variables
edata – Experimental data
-
void addObservableObjectiveSensitivity(std::vector<realtype> &sllh, std::vector<realtype> &s2llh, const int it, const AmiVector &x, const AmiVectorArray &sx, const ExpData &edata)
Add sensitivity of time-resolved measurement negative log-likelihood \( Jy \).
- Parameters
sllh – First-order buffer (shape
nplist
)s2llh – Second-order buffer (shape
nJ - 1
xnplist
, row-major)it – Timepoint index
x – State variables
sx – State sensitivities
edata – Experimental data
-
void addPartialObservableObjectiveSensitivity(std::vector<realtype> &sllh, std::vector<realtype> &s2llh, const int it, const AmiVector &x, const ExpData &edata)
Add sensitivity of time-resolved measurement negative log-likelihood \( Jy \).
Partial derivative (to be used with adjoint sensitivities).
- Parameters
sllh – First order output buffer (shape
nplist
)s2llh – Second order output buffer (shape
nJ - 1
xnplist
, row-major)it – Timepoint index
x – State variables
edata – Experimental data
-
void getAdjointStateObservableUpdate(gsl::span<realtype> dJydx, const int it, const AmiVector &x, const ExpData &edata)
Get state sensitivity of the negative loglikelihood \( Jy \), partial derivative (to be used with adjoint sensitivities).
- Parameters
dJydx – Output buffer (shape
nJ
xnx_solver
, row-major)it – Timepoint index
x – State variables
edata – Experimental data instance
-
void getEvent(gsl::span<realtype> z, const int ie, const realtype t, const AmiVector &x)
Get event-resolved observables.
- Parameters
z – Output buffer (shape
nz
)ie – Event index
t – Timepoint
x – State variables
-
void getEventSensitivity(gsl::span<realtype> sz, const int ie, const realtype t, const AmiVector &x, const AmiVectorArray &sx)
Get sensitivities of event-resolved observables.
Total derivative (only forward sensitivities).
- Parameters
sz – Output buffer (shape
nz x nplist
, row-major)ie – Event index
t – Timepoint
x – State variables
sx – State sensitivities
-
void getUnobservedEventSensitivity(gsl::span<realtype> sz, const int ie)
Get sensitivity of
z
at final timepoint.Ignores sensitivity of timepoint. Total derivative.
- Parameters
sz – Output buffer (shape
nz x nplist
, row-major)ie – Event index
-
void getEventRegularization(gsl::span<realtype> rz, const int ie, const realtype t, const AmiVector &x)
Get regularization for event-resolved observables.
- Parameters
rz – Output buffer (shape
nz
)ie – Event index
t – Timepoint
x – State variables
-
void getEventRegularizationSensitivity(gsl::span<realtype> srz, const int ie, const realtype t, const AmiVector &x, const AmiVectorArray &sx)
Get sensitivities of regularization for event-resolved observables.
Total derivative. Only forward sensitivities.
- Parameters
srz – Output buffer (shape
nz x nplist
, row-major)ie – Event index
t – Timepoint
x – State variables
sx – State sensitivities
-
void getEventSigma(gsl::span<realtype> sigmaz, const int ie, const int nroots, const realtype t, const ExpData *edata)
Get event-resolved observable standard deviations.
- Parameters
sigmaz – Output buffer (shape
nz
)ie – Event index
nroots – Event occurrence
t – Timepoint
edata – Pointer to experimental data (optional, pass
nullptr
to ignore)
-
void getEventSigmaSensitivity(gsl::span<realtype> ssigmaz, const int ie, const int nroots, const realtype t, const ExpData *edata)
Get sensitivities of event-resolved observable standard deviations.
Total derivative (only forward sensitivities).
- Parameters
ssigmaz – Output buffer (shape
nz x nplist
, row-major)ie – Event index
nroots – Event occurrence
t – Timepoint
edata – Pointer to experimental data (optional, pass
nullptr
to ignore)
-
void addEventObjective(realtype &Jz, const int ie, const int nroots, const realtype t, const AmiVector &x, const ExpData &edata)
Add event-resolved observable negative log-likelihood.
- Parameters
Jz – Output buffer (shape 1)
ie – Event index
nroots – Event occurrence
t – Timepoint
x – State variables
edata – Experimental data
-
void addEventObjectiveRegularization(realtype &Jrz, const int ie, const int nroots, const realtype t, const AmiVector &x, const ExpData &edata)
Add event-resolved observable negative log-likelihood.
- Parameters
Jrz – Output buffer (shape 1)
ie – Event index
nroots – Event occurrence
t – Timepoint
x – State variables
edata – Experimental data
-
void addEventObjectiveSensitivity(std::vector<realtype> &sllh, std::vector<realtype> &s2llh, const int ie, const int nroots, const realtype t, const AmiVector &x, const AmiVectorArray &sx, const ExpData &edata)
Add sensitivity of time-resolved measurement negative log-likelihood \( Jy \).
Total derivative (to be used with forward sensitivities).
- Parameters
sllh – First order buffer (shape
nplist
)s2llh – Second order buffer (shape
nJ-1
xnplist
, row-major)ie – Event index
nroots – Event occurrence
t – Timepoint
x – State variables
sx – State sensitivities
edata – Experimental data
-
void addPartialEventObjectiveSensitivity(std::vector<realtype> &sllh, std::vector<realtype> &s2llh, const int ie, const int nroots, const realtype t, const AmiVector &x, const ExpData &edata)
Add sensitivity of time-resolved measurement negative log-likelihood \( Jy \).
Partial derivative (to be used with adjoint sensitivities).
- Parameters
sllh – First order buffer (shape
nplist
)s2llh – Second order buffer (shape
(nJ-1)
xnplist
, row-major)ie – Event index
nroots – Event occurrence
t – Timepoint
x – State variables
edata – Experimental data
-
void getAdjointStateEventUpdate(gsl::span<realtype> dJzdx, const int ie, const int nroots, const realtype t, const AmiVector &x, const ExpData &edata)
State sensitivity of the negative loglikelihood \( Jz \).
Partial derivative (to be used with adjoint sensitivities).
- Parameters
dJzdx – Output buffer (shape
nJ
xnx_solver
, row-major)ie – Event index
nroots – Event occurrence
t – Timepoint
x – State variables
edata – Experimental data
-
void getEventTimeSensitivity(std::vector<realtype> &stau, const realtype t, const int ie, const AmiVector &x, const AmiVectorArray &sx)
Sensitivity of event timepoint, total derivative.
Only forward sensitivities.
- Parameters
stau – Timepoint sensitivity (shape
nplist
)t – Timepoint
ie – Event index
x – State variables
sx – State sensitivities
-
void addStateEventUpdate(AmiVector &x, const int ie, const realtype t, const AmiVector &xdot, const AmiVector &xdot_old)
Update state variables after event.
- Parameters
x – Current state (will be overwritten)
ie – Event index
t – Current timepoint
xdot – Current residual function values
xdot_old – Value of residual function before event
-
void addStateSensitivityEventUpdate(AmiVectorArray &sx, const int ie, const realtype t, const AmiVector &x_old, const AmiVector &xdot, const AmiVector &xdot_old, const std::vector<realtype> &stau)
Update state sensitivity after event.
- Parameters
sx – Current state sensitivity (will be overwritten)
ie – Event index
t – Current timepoint
x_old – Current state
xdot – Current residual function values
xdot_old – Value of residual function before event
stau – Timepoint sensitivity, to be computed with
Model::getEventTimeSensitivity
-
void addAdjointStateEventUpdate(AmiVector &xB, const int ie, const realtype t, const AmiVector &x, const AmiVector &xdot, const AmiVector &xdot_old)
Update adjoint state after event.
- Parameters
xB – Current adjoint state (will be overwritten)
ie – Event index
t – Current timepoint
x – Current state
xdot – Current residual function values
xdot_old – Value of residual function before event
-
void addAdjointQuadratureEventUpdate(AmiVector xQB, const int ie, const realtype t, const AmiVector &x, const AmiVector &xB, const AmiVector &xdot, const AmiVector &xdot_old)
Update adjoint quadratures after event.
- Parameters
xQB – Current quadrature state (will be overwritten)
ie – Event index
t – Current timepoint
x – Current state
xB – Current adjoint state
xdot – Current residual function values
xdot_old – Value of residual function before event
-
void updateHeaviside(const std::vector<int> &rootsfound)
Update the Heaviside variables
h
on event occurrences.- Parameters
rootsfound – Provides the direction of the zero-crossing, so adding it will give the right update to the Heaviside variables (zero if no root was found)
-
void updateHeavisideB(const int *rootsfound)
Updates the Heaviside variables
h
on event occurrences in the backward problem.- Parameters
rootsfound – Provides the direction of the zero-crossing, so adding it will give the right update to the Heaviside variables (zero if no root was found)
-
int checkFinite(gsl::span<const realtype> array, ModelQuantity model_quantity) const
Check if the given array has only finite elements.
For (1D) spans.
- Parameters
array –
model_quantity – The model quantity
array
corresponds to
- Returns
-
int checkFinite(gsl::span<const realtype> array, ModelQuantity model_quantity, size_t num_cols) const
Check if the given array has only finite elements.
For flattened 2D arrays.
- Parameters
array – Flattened matrix
model_quantity – The model quantity
array
corresponds tonum_cols – Number of columns of the non-flattened matrix
- Returns
-
int checkFinite(SUNMatrix m, ModelQuantity model_quantity, realtype t) const
Check if the given array has only finite elements.
For SUNMatrix.
- Parameters
m – Matrix to check
model_quantity – The model quantity
m
corresponds tot – current timepoint
- Returns
-
void setAlwaysCheckFinite(bool alwaysCheck)
Set whether the result of every call to
Model::f*
should be checked for finiteness.- Parameters
alwaysCheck –
-
bool getAlwaysCheckFinite() const
Get setting of whether the result of every call to
Model::f*
should be checked for finiteness.- Returns
that
-
void fx0_fixedParameters(AmiVector &x)
Set only those initial states that are specified via fixed parameters.
- Parameters
x – Output buffer.
-
void fsx0(AmiVectorArray &sx, const AmiVector &x)
Compute/get initial value for initial state sensitivities.
- Parameters
sx – Output buffer for state sensitivities
x – State variables
-
void fsx0_fixedParameters(AmiVectorArray &sx, const AmiVector &x)
Get only those initial states sensitivities that are affected from
amici::Model::fx0_fixedParameters
.- Parameters
sx – Output buffer for state sensitivities
x – State variables
-
virtual void fsdx0()
Compute sensitivity of derivative initial states sensitivities
sdx0
.Only necessary for DAEs.
-
void fx_rdata(AmiVector &x_rdata, const AmiVector &x_solver)
Expand conservation law for states.
- Parameters
x_rdata – Output buffer for state variables with conservation laws expanded (stored in
amici::ReturnData
).x_solver – State variables with conservation laws applied (solver returns this)
-
void fsx_rdata(AmiVectorArray &sx_rdata, const AmiVectorArray &sx_solver, const AmiVector &x_solver)
Expand conservation law for state sensitivities.
- Parameters
sx_rdata – Output buffer for state variables sensitivities with conservation laws expanded (stored in
amici::ReturnData
).sx_solver – State variables sensitivities with conservation laws applied (solver returns this)
x_solver – State variables with conservation laws applied (solver returns this)
-
void setReinitializationStateIdxs(const std::vector<int> &idxs)
Set indices of states to be reinitialized based on provided constants / fixed parameters.
- Parameters
idxs – Array of state indices
-
std::vector<int> const &getReinitializationStateIdxs() const
Return indices of states to be reinitialized based on provided constants / fixed parameters.
- Returns
Those indices.
-
const AmiVectorArray &get_dxdotdp() const
getter for dxdotdp (matlab generated)
- Returns
dxdotdp
-
const SUNMatrixWrapper &get_dxdotdp_full() const
getter for dxdotdp (python generated)
- Returns
dxdotdp
-
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 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 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 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 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 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 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 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 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 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 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 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 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 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(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 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, 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 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 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 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 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 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 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 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 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 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 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 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 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 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 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 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 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 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 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_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 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
Public Members
-
bool pythonGenerated
Flag indicating Matlab- or Python-based model generation
-
SecondOrderMode o2mode = {SecondOrderMode::none}
Flag indicating whether for
amici::Solver::sensi_
==amici::SensitivityOrder::second
directional or full second order derivative will be computed
Protected Functions
-
void writeSliceEvent(gsl::span<const realtype> slice, gsl::span<realtype> buffer, const int ie)
Write part of a slice to a buffer according to indices specified in z2event.
- Parameters
slice – Input data slice
buffer – Output data slice
ie – Event index
-
void writeSensitivitySliceEvent(gsl::span<const realtype> slice, gsl::span<realtype> buffer, const int ie)
Write part of a sensitivity slice to a buffer according to indices specified in z2event.
- Parameters
slice – source data slice
buffer – output data slice
ie – event index
-
void writeLLHSensitivitySlice(const std::vector<realtype> &dLLhdp, std::vector<realtype> &sllh, std::vector<realtype> &s2llh)
Separate first and second order objective sensitivity information and write them into the respective buffers.
- Parameters
dLLhdp – Data with mangled first- and second-order information
sllh – First order buffer
s2llh – Second order buffer
-
void checkLLHBufferSize(const std::vector<realtype> &sllh, const std::vector<realtype> &s2llh) const
Verify that the provided buffers have the expected size.
- Parameters
sllh – first order buffer
s2llh – second order buffer
-
void initializeVectors()
Set the nplist-dependent vectors to their proper sizes.
-
void fy(realtype t, const AmiVector &x)
Compute observables / measurements.
- Parameters
t – Current timepoint
x – Current state
-
void fdydp(realtype t, const AmiVector &x)
Compute partial derivative of observables \( y \) w.r.t. model parameters
p
.- Parameters
t – Current timepoint
x – Current state
-
void fdydx(realtype t, const AmiVector &x)
Compute partial derivative of observables \( y \) w.r.t. state variables
x
.- Parameters
t – Current timepoint
x – Current state
-
void fsigmay(int it, const ExpData *edata)
Compute standard deviation of measurements.
- Parameters
it – Timepoint index
edata – Experimental data
-
void fdsigmaydp(int it, const ExpData *edata)
Compute partial derivative of standard deviation of measurements w.r.t. model parameters.
- Parameters
it – Timepoint index
edata – pointer to
amici::ExpData
data instance holding sigma values
-
void fdsigmaydy(int it, const ExpData *edata)
Compute partial derivative of standard deviation of measurements w.r.t. model outputs.
- Parameters
it – Timepoint index
edata – pointer to
amici::ExpData
data instance holding sigma values
-
void fJy(realtype &Jy, int it, const AmiVector &y, const ExpData &edata)
Compute negative log-likelihood of measurements \( y \).
- Parameters
Jy – Variable to which llh will be added
it – Timepoint index
y – Simulated observable
edata – Pointer to experimental data instance
-
void fdJydy(int it, const AmiVector &x, const ExpData &edata)
Compute partial derivative of time-resolved measurement negative log-likelihood \( Jy \).
- Parameters
it – timepoint index
x – state variables
edata – Pointer to experimental data
-
void fdJydsigma(int it, const AmiVector &x, const ExpData &edata)
Sensitivity of time-resolved measurement negative log-likelihood Jy w.r.t. standard deviation sigma.
- Parameters
it – timepoint index
x – state variables
edata – pointer to experimental data instance
-
void fdJydp(const int it, const AmiVector &x, const ExpData &edata)
Compute sensitivity of time-resolved measurement negative log-likelihood \( Jy \) w.r.t. parameters for the given timepoint.
- Parameters
it – timepoint index
x – state variables
edata – pointer to experimental data instance
-
void fdJydx(const int it, const AmiVector &x, const ExpData &edata)
Sensitivity of time-resolved measurement negative log-likelihood \( Jy \) w.r.t. state variables.
- Parameters
it – Timepoint index
x – State variables
edata – Pointer to experimental data instance
-
void fz(int ie, realtype t, const AmiVector &x)
Compute event-resolved output.
- Parameters
ie – Event index
t – Current timepoint
x – Current state
-
void fdzdp(int ie, realtype t, const AmiVector &x)
Compute partial derivative of event-resolved output
z
w.r.t. model parametersp
- Parameters
ie – event index
t – current timepoint
x – current state
-
void fdzdx(int ie, realtype t, const AmiVector &x)
Compute partial derivative of event-resolved output
z
w.r.t. model statesx
.- Parameters
ie – Event index
t – Current timepoint
x – Current state
-
void frz(int ie, realtype t, const AmiVector &x)
Compute event root function of events.
Equal to
Model::froot
but does not include non-output events.- Parameters
ie – Event index
t – Current timepoint
x – Current state
-
void fdrzdp(int ie, realtype t, const AmiVector &x)
Compute sensitivity of event-resolved root output w.r.t. model parameters
p
.- Parameters
ie – Event index
t – Current timepoint
x – Current state
-
void fdrzdx(int ie, realtype t, const AmiVector &x)
Compute sensitivity of event-resolved measurements \( rz \) w.r.t. model states
x
.- Parameters
ie – Event index
t – Current timepoint
x – Current state
-
void fsigmaz(const int ie, const int nroots, const realtype t, const ExpData *edata)
Compute standard deviation of events.
- Parameters
ie – Event index
nroots – Event index
t – Current timepoint
edata – Experimental data
-
void fdsigmazdp(int ie, int nroots, realtype t, const ExpData *edata)
Compute sensitivity of standard deviation of events measurements w.r.t. model parameters
p
.- Parameters
ie – Event index
nroots – Event occurrence
t – Current timepoint
edata – Pointer to experimental data instance
-
void fJz(realtype &Jz, int nroots, const AmiVector &z, const ExpData &edata)
Compute negative log-likelihood of event-resolved measurements
z
.- Parameters
Jz – Variable to which llh will be added
nroots – Event index
z – Simulated event
edata – Experimental data
-
void fdJzdz(const int ie, const int nroots, const realtype t, const AmiVector &x, const ExpData &edata)
Compute partial derivative of event measurement negative log-likelihood \( Jz \).
- Parameters
ie – Event index
nroots – Event index
t – Current timepoint
x – State variables
edata – Experimental data
-
void fdJzdsigma(const int ie, const int nroots, const realtype t, const AmiVector &x, const ExpData &edata)
Compute sensitivity of event measurement negative log-likelihood \( Jz \) w.r.t. standard deviation sigmaz.
- Parameters
ie – Event index
nroots – Event index
t – Current timepoint
x – State variables
edata – Pointer to experimental data instance
-
void fdJzdp(const int ie, const int nroots, realtype t, const AmiVector &x, const ExpData &edata)
Compute sensitivity of event-resolved measurement negative log-likelihood Jz w.r.t. parameters.
- Parameters
ie – Event index
nroots – Event index
t – Current timepoint
x – State variables
edata – Pointer to experimental data instance
-
void fdJzdx(const int ie, const int nroots, realtype t, const AmiVector &x, const ExpData &edata)
Compute sensitivity of event-resolved measurement negative log-likelihood Jz w.r.t. state variables.
- Parameters
ie – Event index
nroots – Event index
t – Current timepoint
x – State variables
edata – Experimental data
-
void fJrz(realtype &Jrz, int nroots, const AmiVector &rz, const ExpData &edata)
Compute regularization of negative log-likelihood with roots of event-resolved measurements rz.
- Parameters
Jrz – Variable to which regularization will be added
nroots – Event index
rz – Regularization variable
edata – Experimental data
-
void fdJrzdz(const int ie, const int nroots, const realtype t, const AmiVector &x, const ExpData &edata)
Compute partial derivative of event measurement negative log-likelihood J.
- Parameters
ie – Event index
nroots – Event index
t – Current timepoint
x – State variables
edata – Experimental data
-
void fdJrzdsigma(const int ie, const int nroots, const realtype t, const AmiVector &x, const ExpData &edata)
Compute sensitivity of event measurement negative log-likelihood Jz w.r.t. standard deviation sigmaz.
- Parameters
ie – event index
nroots – event index
t – current timepoint
x – state variables
edata – pointer to experimental data instance
-
void fw(realtype t, const realtype *x)
Compute recurring terms in xdot.
- Parameters
t – Timepoint
x – Array with the states
-
void fdwdp(realtype t, const realtype *x)
Compute parameter derivative for recurring terms in xdot.
- Parameters
t – Timepoint
x – Array with the states
-
void fdwdx(realtype t, const realtype *x)
Compute state derivative for recurring terms in xdot.
- Parameters
t – Timepoint
x – Array with the states
-
void fdwdw(realtype t, const realtype *x)
Compute self derivative for recurring terms in xdot.
- Parameters
t – Timepoint
x – Array with the states
-
virtual void fx_rdata(realtype *x_rdata, const realtype *x_solver, const realtype *tcl, const realtype *p, const realtype *k)
Compute fx_rdata.
To be implemented by derived class if applicable.
- Parameters
x_rdata – State variables with conservation laws expanded
x_solver – State variables with conservation laws applied
tcl – Total abundances for conservation laws
p – parameter vector
k – constant vector
-
virtual void fsx_rdata(realtype *sx_rdata, const realtype *sx_solver, const realtype *stcl, const realtype *p, const realtype *k, const realtype *x_solver, const realtype *tcl, const int ip)
Compute fsx_solver.
To be implemented by derived class if applicable.
- Parameters
sx_rdata – State sensitivity variables with conservation laws expanded
sx_solver – State sensitivity variables with conservation laws applied
stcl – Sensitivities of total abundances for conservation laws
p – parameter vector
k – constant vector
x_solver – State variables with conservation laws applied
tcl – Total abundances for conservation laws
ip – Sensitivity index
-
virtual void fx_solver(realtype *x_solver, const realtype *x_rdata)
Compute fx_solver.
To be implemented by derived class if applicable.
- Parameters
x_solver – State variables with conservation laws applied
x_rdata – State variables with conservation laws expanded
-
virtual void fsx_solver(realtype *sx_solver, const realtype *sx_rdata)
Compute fsx_solver.
To be implemented by derived class if applicable.
- Parameters
sx_rdata – State sensitivity variables with conservation laws expanded
sx_solver – State sensitivity variables with conservation laws applied
-
virtual void ftotal_cl(realtype *total_cl, const realtype *x_rdata, const realtype *p, const realtype *k)
Compute ftotal_cl.
To be implemented by derived class if applicable.
- Parameters
total_cl – Total abundances of conservation laws
x_rdata – State variables with conservation laws expanded
p – parameter vector
k – constant vector
-
virtual void fstotal_cl(realtype *stotal_cl, const realtype *sx_rdata, const int ip, const realtype *x_rdata, const realtype *p, const realtype *k, const realtype *tcl)
Compute fstotal_cl.
To be implemented by derived class if applicable.
- Parameters
stotal_cl – Sensitivities for the total abundances of conservation laws
sx_rdata – State sensitivity variables with conservation laws expanded
ip – Sensitivity index
x_rdata – State variables with conservation laws expanded
p – parameter vector
k – constant vector
tcl – Total abundances for conservation laws
-
const_N_Vector computeX_pos(const_N_Vector x)
Compute non-negative state vector.
Compute non-negative state vector according to stateIsNonNegative. If anyStateNonNegative is set to
false
, i.e., all entries in stateIsNonNegative arefalse
, this function directly returnsx
, otherwise all entries of x are copied in toamici::Model::x_pos_tmp_
and negative values are replaced by0
where applicable.- Parameters
x – State vector possibly containing negative values
- Returns
State vector with negative values replaced by
0
according to stateIsNonNegative
-
const realtype *computeX_pos(AmiVector const &x)
Compute non-negative state vector.
Compute non-negative state vector according to stateIsNonNegative. If anyStateNonNegative is set to
false
, i.e., all entries in stateIsNonNegative arefalse
, this function directly returnsx
, otherwise all entries of x are copied in toamici::Model::x_pos_tmp_
and negative values are replaced by0
where applicable.- Parameters
x – State vector possibly containing negative values
- Returns
State vector with negative values replaced by
0
according to stateIsNonNegative
Protected Attributes
-
ModelState state_
All variables necessary for function evaluation
-
ModelStateDerived derived_state_
Storage for model quantities beyond ModelState for the current timepoint
-
std::vector<int> z2event_
index indicating to which event an event output belongs
-
std::vector<bool> state_is_non_negative_
vector of bools indicating whether state variables are to be assumed to be positive
-
std::vector<bool> root_initial_values_
Vector of booleans indicating the initial boolean value for every event trigger function. Events at t0 can only trigger if the initial value is set to
false
. Must be specified during model compilation by setting theinitialValue
attribute of an event trigger.
-
bool any_state_non_negative_ = {false}
boolean indicating whether any entry in stateIsNonNegative is
true
-
int nmaxevent_ = {10}
maximal number of events to track
-
SteadyStateSensitivityMode steadystate_sensitivity_mode_ = {SteadyStateSensitivityMode::newtonOnly}
flag indicating whether steadystate sensitivities are to be computed via FSA when steadyStateSimulation is used
-
bool always_check_finite_ = {false}
Indicates whether the result of every call to
Model::f*
should be checked for finiteness
-
bool sigma_res_ = {false}
indicates whether sigma residuals are to be added for every datapoint
-
Model() = default