Class Model

Inheritance Relationships

Base Type

Derived Types

Class Documentation

class amici::Model : public amici::AbstractModel

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(int nx_rdata, int nxtrue_rdata, int nx_solver, int nxtrue_solver, int nx_solver_reinit, int ny, int nytrue, int nz, int nztrue, int ne, int nJ, int nw, int ndwdx, int ndwdp, int ndwdw, int ndxdotdw, std::vector<int> ndJydy, int nnz, int ubw, int lbw, amici::SecondOrderMode o2mode, const std::vector<amici::realtype> &p, std::vector<amici::realtype> k, const std::vector<int> &plist, 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
  • nx_rdata: Number of state variables

  • nxtrue_rdata: Number of state variables of the non-augmented model

  • nx_solver: Number of state variables with conservation laws applied

  • nxtrue_solver: Number of state variables of the non-augmented model with conservation laws applied

  • nx_solver_reinit: Number of state variables with conservation laws subject to reinitialization

  • ny: Number of observables

  • nytrue: Number of observables of the non-augmented model

  • nz: Number of event observables

  • nztrue: Number of event observables of the non-augmented model

  • ne: Number of events

  • nJ: Number of objective functions

  • nw: Number of repeating elements

  • ndwdx: Number of nonzero elements in the x derivative of the repeating elements

  • ndwdp: Number of nonzero elements in the p derivative of the repeating elements

  • ndwdw: Number of nonzero elements in the w derivative of the repeating elements

  • ndxdotdw: Number of nonzero elements in the \(w\) derivative of \(xdot\)

  • ndJydy: Number of nonzero elements in the \(y\) derivative of \(dJy\) (dimension nytrue)

  • nnz: Number of nonzero elements in Jacobian

  • ubw: Upper matrix bandwidth in the Jacobian

  • lbw: Lower matrix bandwidth in the Jacobian

  • o2mode: Second order sensitivity mode

  • p: Parameters

  • k: Constants

  • plist: Indexes wrt to which sensitivities are to be computed

  • 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.

Return

Parameters
  • other: Object to copy from

Model *clone() const = 0

Clone this instance.

Return

The clone

void initialize(AmiVector &x, AmiVector &dx, AmiVectorArray &sx, AmiVectorArray &sdx, bool computeSensitivities)

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

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 initHeaviside(const AmiVector &x, const AmiVector &dx)

Initialize the Heaviside variables h at the initial time t0.

Heaviside variables activate/deactivate on event occurrences.

Parameters
  • x: Reference to state variables

  • dx: 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 and nx_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 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.

Return

Number of parameter IDs that matched the regex

Parameters
  • par_id_regex: Parameter ID regex

  • value: 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 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.

Return

Number of fixed parameter names that matched the regex

Parameters
  • par_name_regex: Parameter name regex

  • value: 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 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.

Return

Number of fixed parameter IDs that matched the regex

Parameters
  • par_id_regex: Fixed parameter name regex

  • value: 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 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.

Return

Number of fixed parameter names that matched the regex

Parameters
  • par_name_regex: Fixed parameter name regex

  • value: Fixed parameter value

std::string getName() const

Get the model name.

Return

Model name

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 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 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

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

std::vector<realtype> getInitialStates()

Get the initial states.

Return

Initial state vector

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, otherwise false

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

bool hasCustomInitialStateSensitivities() const

Return whether custom initial state sensitivities have been set.

Return

true if has custom initial state sensitivities, otherwise false.

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.

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 (dimension: 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 (dimension: ny)

  • t: Current timepoint

  • x: 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 (dimension: ny x nplist, 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 (dimension: ny)

  • it: Timepoint index

  • edata: Pointer to experimental data instance (optional, pass nullptr 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 (dimension: ny x nplist, row-major)

  • 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 (dimension: 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 (dimension: nplist)

  • s2llh: Second-order buffer (dimension: nJ - 1 x nplist, 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 (dimension: nplist)

  • s2llh: Second order output buffer (dimension: nJ - 1 x nplist, 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 (dimension: nJ x nx_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 (dimension: 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 (dimension: 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 (dimension: 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 (dimension: 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 (dimension: 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 (dimension: 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 (dimension: 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 (dimension: 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 (dimension: 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 (dimension: nplist)

  • s2llh: Second order buffer (dimension: (nJ-1) x nplist, 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 (dimension: nplist)

  • s2llh: Second order buffer (dimension: (nJ-1) x nplist, 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 (dimension: nJ x nx_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 (dimension: 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, 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 check

  • fun: 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(AmiVector &x)

Compute/get initial states.

Parameters
  • x: Output buffer.

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

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)

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)

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

int nx_rdata = {0}

Number of states

int nxtrue_rdata = {0}

Number of states in the unaugmented system

int nx_solver = {0}

Number of states with conservation laws applied

int nxtrue_solver = {0}

Number of states in the unaugmented system with conservation laws applied

int nx_solver_reinit = {0}

Number of solver states subject to reinitialization

int ny = {0}

Number of observables

int nytrue = {0}

Number of observables in the unaugmented system

int nz = {0}

Number of event outputs

int nztrue = {0}

Number of event outputs in the unaugmented system

int ne = {0}

Number of events

int nw = {0}

Number of common expressions

int nnz = {0}

Number of nonzero entries in Jacobian

int nJ = {0}

Dimension of the augmented objective function for 2nd order ASA

int ubw = {0}

Upper bandwidth of the Jacobian

int lbw = {0}

Lower bandwidth of the Jacobian

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

std::vector<realtype> idlist

Flag array for DAE equations

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 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 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 parameters p

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 states x.

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

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 expanded

  • x_solver: State variables with conservation laws applied

  • tcl: 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 expanded

  • sx_solver: State sensitivity variables with conservation laws applied

  • stcl: Sensitivities of total abundances for conservation laws

  • ip: 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 applied

  • x_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 expanded

  • sx_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 laws

  • x_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 laws

  • sx_rdata: State sensitivity variables with conservation laws expanded

  • ip: Sensitivity index

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 are false, this function directly returns x, otherwise all entries of x are copied in to amici::Model::x_pos_tmp_ and negative values are replaced by 0 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

SUNMatrixWrapper J_

Sparse Jacobian (dimension: amici::Model::nnz)

SUNMatrixWrapper JB_

Sparse Backwards Jacobian (dimension: amici::Model::nnz)

SUNMatrixWrapper dxdotdw_

Sparse dxdotdw temporary storage (dimension: ndxdotdw)

SUNMatrixWrapper dwdx_

Sparse dwdx temporary storage (dimension: ndwdx)

SUNMatrixWrapper dwdp_

Sparse dwdp temporary storage (dimension: ndwdp)

SUNMatrixWrapper M_

Dense Mass matrix (dimension: nx_solver x nx_solver)

SUNMatrixWrapper dxdotdp_full

Temporary storage of dxdotdp_full data across functions (Python only) (dimension: nplist x nx_solver, nnz: dynamic, type CSC_MAT)

SUNMatrixWrapper dxdotdp_explicit

Temporary storage of dxdotdp_explicit data across functions (Python only) (dimension: nplist x nx_solver, nnz: ndxdotdp_explicit, type CSC_MAT)

SUNMatrixWrapper dxdotdp_implicit

Temporary storage of dxdotdp_implicit data across functions, Python-only (dimension: nplist x nx_solver, nnz: dynamic, type CSC_MAT)

SUNMatrixWrapper dxdotdx_explicit

Temporary storage of dxdotdx_explicit data across functions (Python only) (dimension: nplist x nx_solver, nnz: ‘nxdotdotdx_explicit’, type CSC_MAT)

SUNMatrixWrapper dxdotdx_implicit

Temporary storage of dxdotdx_implicit data across functions, Python-only (dimension: nplist x nx_solver, nnz: dynamic, type CSC_MAT)

AmiVectorArray dxdotdp = {0, 0}

Temporary storage of dxdotdp data across functions, Matlab only (dimension: nplist x nx_solver, row-major)

std::vector<realtype> my_

Current observable (dimension: nytrue)

std::vector<realtype> mz_

Current event measurement (dimension: nztrue)

std::vector<SUNMatrixWrapper> dJydy_

Sparse observable derivative of data likelihood, only used if pythonGenerated == true (dimension nytrue, nJ x ny, row-major)

std::vector<realtype> dJydy_matlab_

Observable derivative of data likelihood, only used if pythonGenerated == false (dimension nJ x ny x nytrue, row-major)

std::vector<realtype> dJydsigma_

Observable sigma derivative of data likelihood (dimension nJ x ny x nytrue, row-major)

std::vector<realtype> dJydx_

State derivative of data likelihood (dimension nJ x nx_solver, row-major)

std::vector<realtype> dJydp_

Parameter derivative of data likelihood for current timepoint (dimension: nJ x nplist, row-major)

std::vector<realtype> dJzdz_

event output derivative of event likelihood (dimension nJ x nz x nztrue, row-major)

std::vector<realtype> dJzdsigma_

event sigma derivative of event likelihood (dimension nJ x nz x nztrue, row-major)

std::vector<realtype> dJrzdz_

event output derivative of event likelihood at final timepoint (dimension nJ x nz x nztrue, row-major)

std::vector<realtype> dJrzdsigma_

event sigma derivative of event likelihood at final timepoint (dimension nJ x nz x nztrue, row-major)

std::vector<realtype> dJzdx_

state derivative of event likelihood (dimension nJ x nx_solver, row-major)

std::vector<realtype> dJzdp_

parameter derivative of event likelihood for current timepoint (dimension: nJ x nplist x, row-major)

std::vector<realtype> dzdx_

state derivative of event output (dimension: nz x nx_solver, row-major)

std::vector<realtype> dzdp_

parameter derivative of event output (dimension: nz x nplist, row-major)

std::vector<realtype> drzdx_

state derivative of event regularization variable (dimension: nz x nx_solver, row-major)

std::vector<realtype> drzdp_

parameter derivative of event regularization variable (dimension: nz x nplist, row-major)

std::vector<realtype> dydp_

parameter derivative of observable (dimension: ny x nplist, row-major)

std::vector<realtype> dydx_

state derivative of time-resolved observable (dimension: nx_solver x ny, row-major)

std::vector<realtype> w_

temporary storage of w data across functions (dimension: nw)

std::vector<realtype> sx_

temporary storage for flattened sx, (dimension: nx_solver x nplist, row-major)

std::vector<realtype> x_rdata_

temporary storage for x_rdata (dimension: nx_rdata)

std::vector<realtype> sx_rdata_

temporary storage for sx_rdata slice (dimension: nx_rdata)

std::vector<realtype> y_

temporary storage for time-resolved observable (dimension: ny)

std::vector<realtype> sigmay_

data standard deviation for current timepoint (dimension: ny)

std::vector<realtype> dsigmaydp_

temporary storage for parameter derivative of data standard deviation, (dimension: ny x nplist, row-major)

std::vector<realtype> z_

temporary storage for event-resolved observable (dimension: nz)

std::vector<realtype> rz_

temporary storage for event regularization (dimension: nz)

std::vector<realtype> sigmaz_

temporary storage for event standard deviation (dimension: nz)

std::vector<realtype> dsigmazdp_

temporary storage for parameter derivative of event standard deviation, (dimension: nz x nplist, row-major)

std::vector<realtype> deltax_

temporary storage for change in x after event (dimension: nx_solver)

std::vector<realtype> deltasx_

temporary storage for change in sx after event (dimension: nx_solver x nplist, row-major)

std::vector<realtype> deltaxB_

temporary storage for change in xB after event (dimension: nx_solver)

std::vector<realtype> deltaqB_

temporary storage for change in qB after event (dimension: nJ x nplist, row-major)

AmiVector x_pos_tmp_ = {0}

temporary storage of positified state variables according to stateIsNonNegative (dimension: nx_solver)

std::vector<realtype> original_parameters_

original user-provided, possibly scaled parameter array (dimension: np)

std::vector<int> z2event_

index indicating to which event an event output belongs

std::vector<realtype> x0data_

state initialization (size nx_solver)

std::vector<realtype> sx0data_

sensitivity initialization (size nx_rdata x nplist, row-major)

std::vector<realtype> ts_

timepoints (size nt)

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

std::vector<ParameterScaling> pscale_

parameter transformation of originalParameters (dimension np)

realtype tstart_ = {0.0}

starting time

SteadyStateSensitivityMode steadystate_sensitivity_mode_ = {SteadyStateSensitivityMode::newtonOnly}

flag indicating whether steadystate sensitivities are to be computed via FSA when steadyStateSimulation is used

bool reinitialize_fixed_parameter_initial_states_ = {false}

flag indicating whether reinitialization of states depending on fixed parameters is activated

bool always_check_finite_ = {false}

Indicates whether the result of every call to Model::f* should be checked for finiteness

Friends

template<class Archive>
friend void serialize(Archive &ar, Model &m, unsigned int version)

Serialize Model (see boost::serialization::serialize).

Parameters
  • ar: Archive to serialize to

  • m: Data to serialize

  • version: Version number

friend bool operator==(const Model &a, const Model &b)

Check equality of data members.

Return

Equality

Parameters
  • a: First model instance

  • b: Second model instance