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
amici
::
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 dimensionssimulation_parameters
: Simulation parameterso2mode
: Second order sensitivity modeidlist
: Indexes indicating algebraic components (DAE only)z2event
: Mapping of event outputs to eventspythonGenerated
: Flag indicating matlab or python wrappingndxdotdp_explicit
: Number of nonzero elements indxdotdp_explicit
ndxdotdx_explicit
: Number of nonzero elements indxdotdx_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.
- Return
- Parameters
other
: Object to copy from
-
void
initialize
(AmiVector &x, AmiVector &dx, AmiVectorArray &sx, AmiVectorArray &sdx, bool computeSensitivities)¶ Initialize model properties.
- Parameters
x
: Reference to state variablesdx
: Reference to time derivative of states (DAE only)sx
: Reference to state variable sensitivitiessdx
: Reference to time derivative of state sensitivities (DAE only)computeSensitivities
: Flag indicating whether sensitivities are to be computed
-
void
initializeB
(AmiVector &xB, AmiVector &dxB, AmiVector &xQB, bool posteq) const¶ Initialize model properties.
- Parameters
xB
: Adjoint state variablesdxB
: Time derivative of adjoint states (DAE only)xQB
: Adjoint quadraturesposteq
: 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 sensitivitiesx
: Reference to state variables
-
void
initHeaviside
(const AmiVector &x, const AmiVector &dx)¶ Initialize the Heaviside variables
h
at the initial timet0
.Heaviside variables activate/deactivate on event occurrences.
- Parameters
x
: Reference to state variablesdx
: Reference to time derivative of states (DAE only)
-
int
nplist
() const¶ Get number of parameters wrt to which sensitivities are computed.
- Return
Length of sensitivity index vector
-
int
np
() const¶ Get total number of model parameters.
- Return
Length of parameter vector
-
int
nk
() const¶ Get number of constants.
- Return
Length of constant vector
-
int
ncl
() const¶ Get number of conservation laws.
- Return
Number of conservation laws (i.e., difference between
nx_rdata
andnx_solver
).
-
int
nx_reinit
() const¶ Get number of solver states subject to reinitialization.
- Return
Model member
nx_solver_reinit
-
const double *
k
() const¶ Get fixed parameters.
- Return
Pointer to constants array
-
int
nMaxEvent
() const¶ Get maximum number of events that may occur for each type.
- Return
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.
- Return
Number of timepoints
-
std::vector<ParameterScaling> const &
getParameterScale
() const¶ Get parameter scale for each parameter.
- Return
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.
- Return
Unscaled parameters
-
std::vector<realtype> const &
getParameters
() const¶ Get parameter vector.
- Return
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.
- Return
Parameter value
- Parameters
par_id
: Parameter ID
-
realtype
getParameterByName
(std::string const &par_name) const¶ Get value of first model parameter with the specified name.
- Return
Parameter value
- Parameters
par_name
: Parameter name
-
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 valuesignoreErrors
: 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 IDvalue
: 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.
- Return
Number of parameter IDs that matched the regex
- Parameters
par_id_regex
: Parameter ID regexvalue
: Parameter value
-
void
setParameterByName
(std::string const &par_name, realtype value)¶ Set value of first model parameter with the specified name.
- Parameters
par_name
: Parameter namevalue
: 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 valuesignoreErrors
: 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.
- Return
Number of fixed parameter names that matched the regex
- Parameters
par_name_regex
: Parameter name regexvalue
: Parameter value
-
std::vector<realtype> const &
getFixedParameters
() const¶ Get values of fixed parameters.
- Return
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.
- Return
Parameter value
- Parameters
par_id
: Parameter ID
-
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.
- Return
Parameter value
- Parameters
par_name
: Parameter name
-
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 idvalue
: 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.
- Return
Number of fixed parameter IDs that matched the regex
- Parameters
par_id_regex
: Fixed parameter name regexvalue
: Fixed parameter value
-
void
setFixedParameterByName
(std::string const &par_name, realtype value)¶ Set value of first fixed parameter with the specified name.
- Parameters
par_name
: Fixed parameter IDvalue
: 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.
- Return
Number of fixed parameter names that matched the regex
- Parameters
par_name_regex
: Fixed parameter name regexvalue
: Fixed parameter value
-
bool
hasParameterNames
() const¶ Report whether the model has parameter names set.
- Return
Boolean indicating whether parameter names were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getParameterNames
() const¶ Get names of the model parameters.
- Return
The parameter names
-
bool
hasStateNames
() const¶ Report whether the model has state names set.
- Return
Boolean indicating whether state names were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getStateNames
() const¶ Get names of the model states.
- Return
State names
-
bool
hasFixedParameterNames
() const¶ Report whether the model has fixed parameter names set.
- Return
Boolean indicating whether fixed parameter names were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getFixedParameterNames
() const¶ Get names of the fixed model parameters.
- Return
Fixed parameter names
-
bool
hasObservableNames
() const¶ Report whether the model has observable names set.
- Return
Boolean indicating whether observable names were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getObservableNames
() const¶ Get names of the observables.
- Return
Observable names
-
bool
hasExpressionNames
() const¶ Report whether the model has expression names set.
- Return
Boolean indicating whether expression names were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getExpressionNames
() const¶ Get names of the expressions.
- Return
Expression names
-
bool
hasParameterIds
() const¶ Report whether the model has parameter IDs set.
- Return
Boolean indicating whether parameter IDs were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getParameterIds
() const¶ Get IDs of the model parameters.
- Return
Parameter IDs
-
bool
hasStateIds
() const¶ Report whether the model has state IDs set.
- Return
Boolean indicating whether state IDs were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getStateIds
() const¶ Get IDs of the model states.
- Return
Sate IDs
-
bool
hasFixedParameterIds
() const¶ Report whether the model has fixed parameter IDs set.
- Return
Boolean indicating whether fixed parameter IDs were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getFixedParameterIds
() const¶ Get IDs of the fixed model parameters.
- Return
Fixed parameter IDs
-
bool
hasObservableIds
() const¶ Report whether the model has observable IDs set.
- Return
Boolean indicating whether observable ids were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getObservableIds
() const¶ Get IDs of the observables.
- Return
Observable IDs
-
bool
hasExpressionIds
() const¶ Report whether the model has expression IDs set.
- Return
Boolean indicating whether expression ids were set. Also returns
true
if the number of corresponding variables is just zero.
-
std::vector<std::string>
getExpressionIds
() const¶ Get IDs of the expression.
- Return
Expression IDs
-
bool
hasQuadraticLLH
() const¶ Checks whether the defined noise model is gaussian, i.e., the nllh is quadratic.
- Return
boolean flag
-
std::vector<realtype> const &
getTimepoints
() const¶ Get the timepoint vector.
- Return
Timepoint vector
-
realtype
getTimepoint
(int it) const¶ Get simulation timepoint for time index
it
.- Return
Timepoint
- Parameters
it
: Time index
-
void
setTimepoints
(std::vector<realtype> const &ts)¶ Set the timepoint vector.
- Parameters
ts
: New timepoint vector
-
double
t0
() const¶ Get simulation start time.
- Return
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.
- Return
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.
-
ModelState const &
getModelState
() const¶ Get the current model state.
- Return
Current model state
-
void
setModelState
(ModelState const &state)¶ Set the current model state.
- Parameters
state
: Model state
-
std::vector<int> const &
getParameterList
() const¶ Get the list of parameters for which sensitivities are computed.
- Return
List of parameter indices
-
int
plist
(int pos) const¶ Get entry in parameter list by index.
- Return
Index in parameter list
- Parameters
pos
: Index in sensitivity 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.
- Return
true
if has custom initial states, otherwisefalse
-
std::vector<realtype>
getInitialStateSensitivities
()¶ Get the initial states sensitivities.
- Return
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.
- Return
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 amodel.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.
- Return
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.
- Return
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 (shapenw
)t
: Current timepointx
: Current state
-
void
getObservable
(gsl::span<realtype> y, const realtype t, const AmiVector &x)¶ Get time-resolved observables.
- Parameters
y
: Buffer (shapeny
)t
: Current timepointx
: Current state
-
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 (shapeny
xnplist
, row-major)t
: Timpointx
: State variablessx
: State sensitivities
-
void
getObservableSigma
(gsl::span<realtype> sigmay, const int it, const ExpData *edata)¶ Get time-resolved observable standard deviations.
- Parameters
sigmay
: Buffer (shapeny
)it
: Timepoint indexedata
: Pointer to experimental data instance (optional, passnullptr
to ignore)
-
void
getObservableSigmaSensitivity
(gsl::span<realtype> ssigmay, 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 (shapeny
xnplist
, row-major)it
: Timepoint indexedata
: Pointer to experimental data instance (optional, passnullptr
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 indexx
: State variablesedata
: 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 (shapenplist
)s2llh
: Second-order buffer (shapenJ - 1
xnplist
, row-major)it
: Timepoint indexx
: State variablessx
: State sensitivitiesedata
: 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 (shapenplist
)s2llh
: Second order output buffer (shapenJ - 1
xnplist
, row-major)it
: Timepoint indexx
: State variablesedata
: 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 (shapenJ
xnx_solver
, row-major)it
: Timepoint indexx
: State variablesedata
: 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 (shapenz
)ie
: Event indext
: Timepointx
: 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 (shapenz x nplist
, row-major)ie
: Event indext
: Timepointx
: State variablessx
: 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 (shapenz 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 (shapenz
)ie
: Event indext
: Timepointx
: 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 (shapenz x nplist
, row-major)ie
: Event indext
: Timepointx
: State variablessx
: 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 (shapenz
)ie
: Event indexnroots
: Event occurrencet
: Timepointedata
: Pointer to experimental data (optional, passnullptr
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 (shapenz x nplist
, row-major)ie
: Event indexnroots
: Event occurrencet
: Timepointedata
: Pointer to experimental data (optional, passnullptr
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 indexnroots
: Event occurrencet
: Timepointx
: State variablesedata
: 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 indexnroots
: Event occurrencet
: Timepointx
: State variablesedata
: 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 (shapenplist
)s2llh
: Second order buffer (shapenJ-1
xnplist
, row-major)ie
: Event indexnroots
: Event occurrencet
: Timepointx
: State variablessx
: State sensitivitiesedata
: 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 (shapenplist
)s2llh
: Second order buffer (shape(nJ-1)
xnplist
, row-major)ie
: Event indexnroots
: Event occurrencet
: Timepointx
: State variablesedata
: 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 (shapenJ
xnx_solver
, row-major)ie
: Event indexnroots
: Event occurrencet
: Timepointx
: State variablesedata
: 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 (shapenplist
)t
: Timepointie
: Event indexx
: State variablessx
: 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 indext
: Current timepointxdot
: Current residual function valuesxdot_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 indext
: Current timepointx_old
: Current statexdot
: Current residual function valuesxdot_old
: Value of residual function before eventstau
: Timepoint sensitivity, to be computed withModel::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 indext
: Current timepointx
: Current statexdot
: Current residual function valuesxdot_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 indext
: Current timepointx
: Current statexB
: Current adjoint statexdot
: Current residual function valuesxdot_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, const char *fun) const¶ Check if the given array has only finite elements.
If not, try to give hints by which other fields this could be caused.
- Return
amici::AMICI_RECOVERABLE_ERROR
if a NaN/Inf value was found,amici::AMICI_SUCCESS
otherwise- Parameters
array
: Array to checkfun
: Name of the function that generated the values (for more informative messages).
-
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.- Return
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 sensitivitiesx
: 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 sensitivitiesx
: State variables
-
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 inamici::ReturnData
).x_solver
: State variables with conservation laws applied (solver returns this)
-
void
fsx_rdata
(AmiVectorArray &sx_rdata, const AmiVectorArray &sx_solver)¶ Expand conservation law for state sensitivities.
- Parameters
sx_rdata
: Output buffer for state variables sensitivities with conservation laws expanded (stored inamici::ReturnData
).sx_solver
: State variables sensitivities 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.
- Return
Those indices.
-
const AmiVectorArray &
get_dxdotdp
() const¶ getter for dxdotdp (matlab generated)
- Return
dxdotdp
-
const SUNMatrixWrapper &
get_dxdotdp_full
() const¶ getter for dxdotdp (python generated)
- Return
dxdotdp
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
-
AmiciApplication *
app
= &defaultContext¶ AMICI application context
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 slicebuffer
: Output data sliceie
: 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 slicebuffer
: output data sliceie
: 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 informationsllh
: First order buffers2llh
: 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 buffers2llh
: 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 timepointx
: Current state
-
void
fdydp
(realtype t, const AmiVector &x)¶ Compute partial derivative of observables \(y\) w.r.t. model parameters
p
.- Parameters
t
: Current timepointx
: Current state
-
void
fdydx
(realtype t, const AmiVector &x)¶ Compute partial derivative of observables \(y\) w.r.t. state variables
x
.- Parameters
t
: Current timepointx
: Current state
-
void
fsigmay
(int it, const ExpData *edata)¶ Compute standard deviation of measurements.
- Parameters
it
: Timepoint indexedata
: 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 indexedata
: pointer toamici::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 addedit
: Timepoint indexy
: Simulated observableedata
: 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 indexx
: state variablesedata
: 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 indexx
: state variablesedata
: 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 indexx
: state variablesedata
: 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 indexx
: State variablesedata
: Pointer to experimental data instance
-
void
fz
(int ie, realtype t, const AmiVector &x)¶ Compute event-resolved output.
- Parameters
ie
: Event indext
: Current timepointx
: 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 indext
: current timepointx
: 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 indext
: Current timepointx
: 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 indext
: Current timepointx
: 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 indext
: Current timepointx
: 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 indext
: Current timepointx
: Current state
-
void
fsigmaz
(const int ie, const int nroots, const realtype t, const ExpData *edata)¶ Compute standard deviation of events.
- Parameters
ie
: Event indexnroots
: Event indext
: Current timepointedata
: 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 indexnroots
: Event occurrencet
: Current timepointedata
: 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 addednroots
: Event indexz
: Simulated eventedata
: 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 indexnroots
: Event indext
: Current timepointx
: State variablesedata
: 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 indexnroots
: Event indext
: Current timepointx
: State variablesedata
: 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 indexnroots
: Event indext
: Current timepointx
: State variablesedata
: 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 indexnroots
: Event indext
: Current timepointx
: State variablesedata
: 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 addednroots
: Event indexrz
: Regularization variableedata
: 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 indexnroots
: Event indext
: Current timepointx
: State variablesedata
: 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 indexnroots
: event indext
: current timepointx
: state variablesedata
: pointer to experimental data instance
-
void
fw
(realtype t, const realtype *x)¶ Compute recurring terms in xdot.
- Parameters
t
: Timepointx
: Array with the states
-
void
fdwdp
(realtype t, const realtype *x)¶ Compute parameter derivative for recurring terms in xdot.
- Parameters
t
: Timepointx
: Array with the states
-
void
fdwdx
(realtype t, const realtype *x)¶ Compute state derivative for recurring terms in xdot.
- Parameters
t
: Timepointx
: Array with the states
-
void
fdwdw
(realtype t, const realtype *x)¶ Compute self derivative for recurring terms in xdot.
- Parameters
t
: Timepointx
: Array with the states
-
void
fx_rdata
(realtype *x_rdata, const realtype *x_solver, const realtype *tcl)¶ Compute fx_rdata.
To be implemented by derived class if applicable.
- Parameters
x_rdata
: State variables with conservation laws expandedx_solver
: State variables with conservation laws appliedtcl
: Total abundances for conservation laws
-
void
fsx_rdata
(realtype *sx_rdata, const realtype *sx_solver, const realtype *stcl, int ip)¶ Compute fsx_solver.
To be implemented by derived class if applicable.
- Parameters
sx_rdata
: State sensitivity variables with conservation laws expandedsx_solver
: State sensitivity variables with conservation laws appliedstcl
: Sensitivities of total abundances for conservation lawsip
: Sensitivity index
-
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 appliedx_rdata
: State variables with conservation laws expanded
-
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 expandedsx_solver
: State sensitivity variables with conservation laws applied
-
void
ftotal_cl
(realtype *total_cl, const realtype *x_rdata)¶ Compute ftotal_cl.
To be implemented by derived class if applicable.
- Parameters
total_cl
: Total abundances of conservation lawsx_rdata
: State variables with conservation laws expanded
-
void
fstotal_cl
(realtype *stotal_cl, const realtype *sx_rdata, int ip)¶ Compute fstotal_cl.
To be implemented by derived class if applicable.
- Parameters
stotal_cl
: Sensitivities for the total abundances of conservation lawssx_rdata
: State sensitivity variables with conservation laws expandedip
: Sensitivity index
-
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.- Return
State vector with negative values replaced by
0
according to stateIsNonNegative- Parameters
x
: State vector possibly containing negative values
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
-
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
-