Class Solver¶
Defined in File solver.h
Inheritance Relationships¶
Derived Types¶
public amici::CVodeSolver
(Class CVodeSolver)public amici::IDASolver
(Class IDASolver)
Class Documentation¶
-
class
amici
::
Solver
¶ The Solver class provides a generic interface to CVODES and IDAS solvers, individual realizations are realized in the CVodeSolver and the IDASolver class. All transient private/protected members (CVODES/IDAS memory, interface variables and status flags) are specified as mutable and not included in serialization or equality checks. No solver setting parameter should be marked mutable.
NOTE: Any changes in data members here must be propagated to copy ctor, equality operator, serialization functions in serialization.h, and amici::hdf5::readSolverSettingsFromHDF5 in hdf5.cpp.
Subclassed by amici::CVodeSolver, amici::IDASolver
Public Functions
-
Solver
() = default¶
-
Solver
(AmiciApplication *app)¶ Constructor.
- Parameters
app
: AMICI application context
-
~Solver
() = default¶
-
int
run
(realtype tout) const¶ runs a forward simulation until the specified timepoint
- Return
status flag
- Parameters
tout
: next timepoint
-
int
step
(realtype tout) const¶ makes a single step in the simulation
- Return
status flag
- Parameters
tout
: next timepoint
-
void
runB
(realtype tout) const¶ runs a backward simulation until the specified timepoint
- Parameters
tout
: next timepoint
-
void
setup
(realtype t0, Model *model, const AmiVector &x0, const AmiVector &dx0, const AmiVectorArray &sx0, const AmiVectorArray &sdx0) const¶ Initializes the ami memory object and applies specified options.
- Parameters
t0
: initial timepointmodel
: pointer to the model instancex0
: initial statesdx0
: initial derivative statessx0
: initial state sensitivitiessdx0
: initial derivative state sensitivities
-
void
setupB
(int *which, realtype tf, Model *model, const AmiVector &xB0, const AmiVector &dxB0, const AmiVector &xQB0) const¶ Initializes the AMI memory object for the backwards problem.
- Parameters
which
: index of the backward problem, will be set by this routinetf
: final timepoint (initial timepoint for the bwd problem)model
: pointer to the model instancexB0
: initial adjoint statesdxB0
: initial adjoint derivative statesxQB0
: initial adjoint quadratures
-
void
setupSteadystate
(const realtype t0, Model *model, const AmiVector &x0, const AmiVector &dx0, const AmiVector &xB0, const AmiVector &dxB0, const AmiVector &xQ0) const¶ Initializes the ami memory for quadrature computation.
- Parameters
t0
: initial timepointmodel
: pointer to the model instancex0
: initial statesdx0
: initial derivative statesxB0
: initial adjoint statesdxB0
: initial derivative adjoint statesxQ0
: initial quadrature vector
-
void
updateAndReinitStatesAndSensitivities
(Model *model)¶ Reinitializes state and respective sensitivities (if necessary) according to changes in fixedParameters.
- Parameters
model
: pointer to the model instance
-
void
getRootInfo
(int *rootsfound) const = 0¶ getRootInfo extracts information which event occurred
- Parameters
rootsfound
: array with flags indicating whether the respective event occurred
-
void
calcIC
(realtype tout1) const = 0¶ Calculates consistent initial conditions, assumes initial states to be correct (DAE only)
- Parameters
tout1
: next timepoint to be computed (sets timescale)
-
void
calcICB
(int which, realtype tout1) const = 0¶ Calculates consistent initial conditions for the backwards problem, assumes initial states to be correct (DAE only)
- Parameters
which
: identifier of the backwards problemtout1
: next timepoint to be computed (sets timescale)
-
void
solveB
(realtype tBout, int itaskB) const = 0¶ Solves the backward problem until a predefined timepoint (adjoint only)
- Parameters
tBout
: timepoint until which simulation should be performeditaskB
: task identifier, can be CV_NORMAL or CV_ONE_STEP
-
void
turnOffRootFinding
() const = 0¶ Disable rootfinding.
-
SensitivityMethod
getSensitivityMethod
() const¶ Return current sensitivity method.
- Return
method enum
-
void
setSensitivityMethod
(SensitivityMethod sensi_meth)¶ Set sensitivity method.
- Parameters
sensi_meth
:
-
SensitivityMethod
getSensitivityMethodPreequilibration
() const¶ Return current sensitivity method during preequilibration.
- Return
method enum
-
void
setSensitivityMethodPreequilibration
(SensitivityMethod sensi_meth_preeq)¶ Set sensitivity method for preequilibration.
- Parameters
sensi_meth_preeq
:
-
void
switchForwardSensisOff
() const¶ Disable forward sensitivity integration (used in steady state sim)
-
int
getNewtonMaxSteps
() const¶ Get maximum number of allowed Newton steps for steady state computation.
- Return
-
void
setNewtonMaxSteps
(int newton_maxsteps)¶ Set maximum number of allowed Newton steps for steady state computation.
- Parameters
newton_maxsteps
:
-
bool
getPreequilibration
() const¶ Get if model preequilibration is enabled.
- Return
-
void
setPreequilibration
(bool require_preequilibration)¶ Enable/disable model preequilibration.
- Parameters
require_preequilibration
:
-
int
getNewtonMaxLinearSteps
() const¶ Get maximum number of allowed linear steps per Newton step for steady state computation.
- Return
-
void
setNewtonMaxLinearSteps
(int newton_maxlinsteps)¶ Set maximum number of allowed linear steps per Newton step for steady state computation.
- Parameters
newton_maxlinsteps
:
-
NewtonDampingFactorMode
getNewtonDampingFactorMode
() const¶ Get a state of the damping factor used in the Newton solver.
- Return
-
void
setNewtonDampingFactorMode
(NewtonDampingFactorMode dampingFactorMode)¶ Turn on/off a damping factor in the Newton method.
- Parameters
dampingFactorMode
:
-
double
getNewtonDampingFactorLowerBound
() const¶ Get a lower bound of the damping factor used in the Newton solver.
- Return
-
void
setNewtonDampingFactorLowerBound
(double dampingFactorLowerBound)¶ Set a lower bound of the damping factor in the Newton solver.
- Parameters
dampingFactorLowerBound
:
-
SensitivityOrder
getSensitivityOrder
() const¶ Get sensitivity order.
- Return
sensitivity order
-
void
setSensitivityOrder
(SensitivityOrder sensi)¶ Set the sensitivity order.
- Parameters
sensi
: sensitivity order
-
double
getRelativeTolerance
() const¶ Get the relative tolerances for the forward problem.
Same tolerance is used for the backward problem if not specified differently via setRelativeToleranceASA.
- Return
relative tolerances
-
void
setRelativeTolerance
(double rtol)¶ Sets the relative tolerances for the forward problem.
Same tolerance is used for the backward problem if not specified differently via setRelativeToleranceASA.
- Parameters
rtol
: relative tolerance (non-negative number)
-
double
getAbsoluteTolerance
() const¶ Get the absolute tolerances for the forward problem.
Same tolerance is used for the backward problem if not specified differently via setAbsoluteToleranceASA.
- Return
absolute tolerances
-
void
setAbsoluteTolerance
(double atol)¶ Sets the absolute tolerances for the forward problem.
Same tolerance is used for the backward problem if not specified differently via setAbsoluteToleranceASA.
- Parameters
atol
: absolute tolerance (non-negative number)
-
double
getRelativeToleranceFSA
() const¶ Returns the relative tolerances for the forward sensitivity problem.
- Return
relative tolerances
-
void
setRelativeToleranceFSA
(double rtol)¶ Sets the relative tolerances for the forward sensitivity problem.
- Parameters
rtol
: relative tolerance (non-negative number)
-
double
getAbsoluteToleranceFSA
() const¶ Returns the absolute tolerances for the forward sensitivity problem.
- Return
absolute tolerances
-
void
setAbsoluteToleranceFSA
(double atol)¶ Sets the absolute tolerances for the forward sensitivity problem.
- Parameters
atol
: absolute tolerance (non-negative number)
-
double
getRelativeToleranceB
() const¶ Returns the relative tolerances for the adjoint sensitivity problem.
- Return
relative tolerances
-
void
setRelativeToleranceB
(double rtol)¶ Sets the relative tolerances for the adjoint sensitivity problem.
- Parameters
rtol
: relative tolerance (non-negative number)
-
double
getAbsoluteToleranceB
() const¶ Returns the absolute tolerances for the backward problem for adjoint sensitivity analysis.
- Return
absolute tolerances
-
void
setAbsoluteToleranceB
(double atol)¶ Sets the absolute tolerances for the backward problem for adjoint sensitivity analysis.
- Parameters
atol
: absolute tolerance (non-negative number)
-
double
getRelativeToleranceQuadratures
() const¶ Returns the relative tolerance for the quadrature problem.
- Return
relative tolerance
-
void
setRelativeToleranceQuadratures
(double rtol)¶ sets the relative tolerance for the quadrature problem
- Parameters
rtol
: relative tolerance (non-negative number)
-
double
getAbsoluteToleranceQuadratures
() const¶ returns the absolute tolerance for the quadrature problem
- Return
absolute tolerance
-
void
setAbsoluteToleranceQuadratures
(double atol)¶ sets the absolute tolerance for the quadrature problem
- Parameters
atol
: absolute tolerance (non-negative number)
-
double
getRelativeToleranceSteadyState
() const¶ returns the relative tolerance for the steady state problem
- Return
relative tolerance
-
void
setRelativeToleranceSteadyState
(double rtol)¶ sets the relative tolerance for the steady state problem
- Parameters
rtol
: relative tolerance (non-negative number)
-
double
getAbsoluteToleranceSteadyState
() const¶ returns the absolute tolerance for the steady state problem
- Return
absolute tolerance
-
void
setAbsoluteToleranceSteadyState
(double atol)¶ sets the absolute tolerance for the steady state problem
- Parameters
atol
: absolute tolerance (non-negative number)
-
double
getRelativeToleranceSteadyStateSensi
() const¶ returns the relative tolerance for the sensitivities of the steady state problem
- Return
relative tolerance
-
void
setRelativeToleranceSteadyStateSensi
(double rtol)¶ sets the relative tolerance for the sensitivities of the steady state problem
- Parameters
rtol
: relative tolerance (non-negative number)
-
double
getAbsoluteToleranceSteadyStateSensi
() const¶ returns the absolute tolerance for the sensitivities of the steady state problem
- Return
absolute tolerance
-
void
setAbsoluteToleranceSteadyStateSensi
(double atol)¶ sets the absolute tolerance for the sensitivities of the steady state problem
- Parameters
atol
: absolute tolerance (non-negative number)
-
long int
getMaxSteps
() const¶ returns the maximum number of solver steps for the forward problem
- Return
maximum number of solver steps
-
void
setMaxSteps
(long int maxsteps)¶ sets the maximum number of solver steps for the forward problem
- Parameters
maxsteps
: maximum number of solver steps (positive number)
-
long int
getMaxStepsBackwardProblem
() const¶ returns the maximum number of solver steps for the backward problem
- Return
maximum number of solver steps
-
void
setMaxStepsBackwardProblem
(long int maxsteps)¶ sets the maximum number of solver steps for the backward problem
- Note
default behaviour (100 times the value for the forward problem) can be restored by passing maxsteps=0
- Parameters
maxsteps
: maximum number of solver steps (non-negative number)
-
LinearMultistepMethod
getLinearMultistepMethod
() const¶ returns the linear system multistep method
- Return
linear system multistep method
-
void
setLinearMultistepMethod
(LinearMultistepMethod lmm)¶ sets the linear system multistep method
- Parameters
lmm
: linear system multistep method
-
NonlinearSolverIteration
getNonlinearSolverIteration
() const¶ returns the nonlinear system solution method
- Return
-
void
setNonlinearSolverIteration
(NonlinearSolverIteration iter)¶ sets the nonlinear system solution method
- Parameters
iter
: nonlinear system solution method
-
InterpolationType
getInterpolationType
() const¶ getInterpolationType
- Return
-
void
setInterpolationType
(InterpolationType interpType)¶ sets the interpolation of the forward solution that is used for the backwards problem
- Parameters
interpType
: interpolation type
-
int
getStateOrdering
() const¶ Gets KLU / SuperLUMT state ordering mode.
- Return
State-ordering as integer according to SUNLinSolKLU::StateOrdering or SUNLinSolSuperLUMT::StateOrdering (which differ).
-
void
setStateOrdering
(int ordering)¶ Sets KLU / SuperLUMT state ordering mode.
This only applies when linsol is set to LinearSolver::KLU or LinearSolver::SuperLUMT. Mind the difference between SUNLinSolKLU::StateOrdering and SUNLinSolSuperLUMT::StateOrdering.
- Parameters
ordering
: state ordering
-
bool
getStabilityLimitFlag
() const¶ returns stability limit detection mode
- Return
stldet can be false (deactivated) or true (activated)
-
void
setStabilityLimitFlag
(bool stldet)¶ set stability limit detection mode
- Parameters
stldet
: can be false (deactivated) or true (activated)
-
LinearSolver
getLinearSolver
() const¶ getLinearSolver
- Return
-
void
setLinearSolver
(LinearSolver linsol)¶ setLinearSolver
- Parameters
linsol
:
-
InternalSensitivityMethod
getInternalSensitivityMethod
() const¶ returns the internal sensitivity method
- Return
internal sensitivity method
-
void
setInternalSensitivityMethod
(InternalSensitivityMethod ism)¶ sets the internal sensitivity method
- Parameters
ism
: internal sensitivity method
-
RDataReporting
getReturnDataReportingMode
() const¶ returns the ReturnData reporting mode
- Return
ReturnData reporting mode
-
void
setReturnDataReportingMode
(RDataReporting rdrm)¶ sets the ReturnData reporting mode
- Parameters
rdrm
: ReturnData reporting mode
-
void
writeSolution
(realtype *t, AmiVector &x, AmiVector &dx, AmiVectorArray &sx, AmiVector &xQ) const¶ write solution from forward simulation
- Parameters
t
: timex
: statedx
: derivative statesx
: state sensitivityxQ
: quadrature
-
void
writeSolutionB
(realtype *t, AmiVector &xB, AmiVector &dxB, AmiVector &xQB, int which) const¶ write solution from forward simulation
- Parameters
t
: timexB
: adjoint statedxB
: adjoint derivative statexQB
: adjoint quadraturewhich
: index of adjoint problem
-
const AmiVector &
getState
(realtype t) const¶ Access state solution at time t.
- Return
x or interpolated solution dky
- Parameters
t
: time
-
const AmiVector &
getDerivativeState
(realtype t) const¶ Access derivative state solution at time t.
- Return
dx or interpolated solution dky
- Parameters
t
: time
-
const AmiVectorArray &
getStateSensitivity
(realtype t) const¶ Access state sensitivity solution at time t.
- Return
(interpolated) solution sx
- Parameters
t
: time
-
const AmiVector &
getAdjointState
(int which, realtype t) const¶ Access adjoint solution at time t.
- Return
(interpolated) solution xB
- Parameters
which
: adjoint problem indext
: time
-
const AmiVector &
getAdjointDerivativeState
(int which, realtype t) const¶ Access adjoint derivative solution at time t.
- Return
(interpolated) solution dxB
- Parameters
which
: adjoint problem indext
: time
-
const AmiVector &
getAdjointQuadrature
(int which, realtype t) const¶ Access adjoint quadrature solution at time t.
- Return
(interpolated) solution xQB
- Parameters
which
: adjoint problem indext
: time
-
const AmiVector &
getQuadrature
(realtype t) const¶ Access quadrature solution at time t.
- Return
(interpolated) solution xQ
- Parameters
t
: time
-
void
reInit
(realtype t0, const AmiVector &yy0, const AmiVector &yp0) const = 0¶ Reinitializes the states in the solver after an event occurrence.
- Parameters
t0
: reinitialization timepointyy0
: initial state variablesyp0
: initial derivative state variables (DAE only)
-
void
sensReInit
(const AmiVectorArray &yyS0, const AmiVectorArray &ypS0) const = 0¶ Reinitializes the state sensitivities in the solver after an event occurrence.
- Parameters
yyS0
: new state sensitivityypS0
: new derivative state sensitivities (DAE only)
-
void
sensToggleOff
() const = 0¶ Switches off computation of state sensitivities without deallocating the memory for sensitivities.
-
void
reInitB
(int which, realtype tB0, const AmiVector &yyB0, const AmiVector &ypB0) const = 0¶ Reinitializes the adjoint states after an event occurrence.
- Parameters
which
: identifier of the backwards problemtB0
: reinitialization timepointyyB0
: new adjoint stateypB0
: new adjoint derivative state
-
void
quadReInitB
(int which, const AmiVector &yQB0) const = 0¶ Reinitialize the adjoint states after an event occurrence.
- Parameters
which
: identifier of the backwards problemyQB0
: new adjoint quadrature state
-
int
nx
() const¶ number of states with which the solver was initialized
- Return
x.getLength()
-
int
nplist
() const¶ number of parameters with which the solver was initialized
- Return
sx.getLength()
-
int
nquad
() const¶ number of quadratures with which the solver was initialized
- Return
xQB.getLength()
-
bool
computingFSA
() const¶ check if FSA is being computed
- Return
flag
-
bool
computingASA
() const¶ check if ASA is being computed
- Return
flag
-
void
resetDiagnosis
() const¶ Resets vectors containing diagnosis information.
-
void
storeDiagnosis
() const¶ Stores diagnosis information from solver memory block for forward problem.
-
void
storeDiagnosisB
(int which) const¶ Stores diagnosis information from solver memory block for backward problem.
- Parameters
which
: identifier of the backwards problem
-
std::vector<int> const &
getNumSteps
() const¶ Accessor ns.
- Return
ns
-
std::vector<int> const &
getNumStepsB
() const¶ Accessor nsB.
- Return
nsB
-
std::vector<int> const &
getNumRhsEvals
() const¶ Accessor nrhs.
- Return
nrhs
-
std::vector<int> const &
getNumRhsEvalsB
() const¶ Accessor nrhsB.
- Return
nrhsB
-
std::vector<int> const &
getNumErrTestFails
() const¶ Accessor netf.
- Return
netf
-
std::vector<int> const &
getNumErrTestFailsB
() const¶ Accessor netfB.
- Return
netfB
-
std::vector<int> const &
getNumNonlinSolvConvFails
() const¶ Accessor nnlscf.
- Return
nnlscf
-
std::vector<int> const &
getNumNonlinSolvConvFailsB
() const¶ Accessor nnlscfB.
- Return
nnlscfB
-
std::vector<int> const &
getLastOrder
() const¶ Accessor order.
- Return
order
Public Members
-
AmiciApplication *
app
= &defaultContext¶ AMICI context
Protected Functions
-
void
setStopTime
(realtype tstop) const = 0¶ Sets a timepoint at which the simulation will be stopped.
- Parameters
tstop
: timepoint until which simulation should be performed
-
int
solve
(realtype tout, int itask) const = 0¶ Solves the forward problem until a predefined timepoint.
- Return
status flag indicating success of execution
- Parameters
tout
: timepoint until which simulation should be performeditask
: task identifier, can be CV_NORMAL or CV_ONE_STEP
-
int
solveF
(realtype tout, int itask, int *ncheckPtr) const = 0¶ Solves the forward problem until a predefined timepoint (adjoint only)
- Return
status flag indicating success of execution
- Parameters
tout
: timepoint until which simulation should be performeditask
: task identifier, can be CV_NORMAL or CV_ONE_STEPncheckPtr
: pointer to a number that counts the internal checkpoints
-
void
reInitPostProcessF
(realtype tnext) const = 0¶ reInitPostProcessF postprocessing of the solver memory after a discontinuity in the forward problem
- Parameters
tnext
: next timepoint (defines integration direction)
-
void
reInitPostProcessB
(realtype tnext) const = 0¶ reInitPostProcessB postprocessing of the solver memory after a discontinuity in the backward problem
- Parameters
tnext
: next timepoint (defines integration direction)
-
void
getSens
() const = 0¶ extracts the state sensitivity at the current timepoint from solver memory and writes it to the sx member variable
-
void
getB
(int which) const = 0¶ extracts the adjoint state at the current timepoint from solver memory and writes it to the xB member variable
- Parameters
which
: index of the backwards problem
-
void
getQuadB
(int which) const = 0¶ extracts the adjoint quadrature state at the current timepoint from solver memory and writes it to the xQB member variable
- Parameters
which
: index of the backwards problem
-
void
getQuad
(realtype &t) const = 0¶ extracts the quadrature at the current timepoint from solver memory and writes it to the xQ member variable
- Parameters
t
: timepoint for quadrature extraction
-
void
init
(realtype t0, const AmiVector &x0, const AmiVector &dx0) const = 0¶ Initializes the states at the specified initial timepoint.
- Parameters
t0
: initial timepointx0
: initial statesdx0
: initial derivative states
-
void
initSteadystate
(realtype t0, const AmiVector &x0, const AmiVector &dx0) const = 0¶ Initializes the states at the specified initial timepoint.
- Parameters
t0
: initial timepointx0
: initial statesdx0
: initial derivative states
-
void
sensInit1
(const AmiVectorArray &sx0, const AmiVectorArray &sdx0) const = 0¶ Initializes the forward sensitivities.
- Parameters
sx0
: initial states sensitivitiessdx0
: initial derivative states sensitivities
-
void
binit
(int which, realtype tf, const AmiVector &xB0, const AmiVector &dxB0) const = 0¶ Initialize the adjoint states at the specified final timepoint.
- Parameters
which
: identifier of the backwards problemtf
: final timepointxB0
: initial adjoint statedxB0
: initial adjoint derivative state
-
void
qbinit
(int which, const AmiVector &xQB0) const = 0¶ Initialize the quadrature states at the specified final timepoint.
- Parameters
which
: identifier of the backwards problemxQB0
: initial adjoint quadrature state
-
void
rootInit
(int ne) const = 0¶ Initializes the rootfinding for events.
- Parameters
ne
: number of different events
-
void
initializeNonLinearSolverSens
(const Model *model) const¶ Initalize non-linear solver for sensitivities.
- Parameters
model
: Model instance
-
void
setDenseJacFn
() const = 0¶ Set the dense Jacobian function.
-
void
setSparseJacFn
() const = 0¶ sets the sparse Jacobian function
-
void
setBandJacFn
() const = 0¶ sets the banded Jacobian function
-
void
setJacTimesVecFn
() const = 0¶ sets the Jacobian vector multiplication function
-
void
setDenseJacFnB
(int which) const = 0¶ sets the dense Jacobian function
- Parameters
which
: identifier of the backwards problem
-
void
setSparseJacFnB
(int which) const = 0¶ sets the sparse Jacobian function
- Parameters
which
: identifier of the backwards problem
-
void
setBandJacFnB
(int which) const = 0¶ sets the banded Jacobian function
- Parameters
which
: identifier of the backwards problem
-
void
setJacTimesVecFnB
(int which) const = 0¶ sets the Jacobian vector multiplication function
- Parameters
which
: identifier of the backwards problem
-
void
setSparseJacFn_ss
() const = 0¶ sets the sparse Jacobian function for backward steady state case
-
void
allocateSolver
() const = 0¶ Create specifies solver method and initializes solver memory for the forward problem.
-
void
setSStolerances
(double rtol, double atol) const = 0¶ sets scalar relative and absolute tolerances for the forward problem
- Parameters
rtol
: relative tolerancesatol
: absolute tolerances
-
void
setSensSStolerances
(double rtol, const double *atol) const = 0¶ activates sets scalar relative and absolute tolerances for the sensitivity variables
- Parameters
rtol
: relative tolerancesatol
: array of absolute tolerances for every sensitivity variable
-
void
setSensErrCon
(bool error_corr) const = 0¶ SetSensErrCon specifies whether error control is also enforced for sensitivities for the forward problem
- Parameters
error_corr
: activation flag
-
void
setQuadErrConB
(int which, bool flag) const = 0¶ Specifies whether error control is also enforced for the backward quadrature problem.
- Parameters
which
: identifier of the backwards problemflag
: activation flag
-
void
setQuadErrCon
(bool flag) const = 0¶ Specifies whether error control is also enforced for the forward quadrature problem.
- Parameters
flag
: activation flag
-
void
setErrHandlerFn
() const = 0¶ Attaches the error handler function (errMsgIdAndTxt) to the solver.
-
void
setUserData
(Model *model) const = 0¶ Attaches the user data instance (here this is a Model) to the forward problem.
- Parameters
model
: Model instance
-
void
setUserDataB
(int which, Model *model) const = 0¶ attaches the user data instance (here this is a Model) to the backward problem
- Parameters
which
: identifier of the backwards problemmodel
: Model instance
-
void
setMaxNumSteps
(long int mxsteps) const = 0¶ specifies the maximum number of steps for the forward problem
- Note
in contrast to the SUNDIALS method, this sets the overall maximum, not the maximum between output times.
- Parameters
mxsteps
: number of steps
-
void
setMaxNumStepsB
(int which, long int mxstepsB) const = 0¶ specifies the maximum number of steps for the forward problem
- Note
in contrast to the SUNDIALS method, this sets the overall maximum, not the maximum between output times.
- Parameters
which
: identifier of the backwards problemmxstepsB
: number of steps
-
void
setStabLimDet
(int stldet) const = 0¶ activates stability limit detection for the forward problem
- Parameters
stldet
: flag for stability limit detection (TRUE or FALSE)
-
void
setStabLimDetB
(int which, int stldet) const = 0¶ activates stability limit detection for the backward problem
- Parameters
which
: identifier of the backwards problemstldet
: flag for stability limit detection (TRUE or FALSE)
-
void
setId
(const Model *model) const = 0¶ specify algebraic/differential components (DAE only)
- Parameters
model
: model specification
-
void
setSuppressAlg
(bool flag) const = 0¶ deactivates error control for algebraic components (DAE only)
- Parameters
flag
: deactivation flag
-
void
setSensParams
(const realtype *p, const realtype *pbar, const int *plist) const = 0¶ specifies the scaling and indexes for sensitivity computation
- Parameters
p
: parameterspbar
: parameter scaling constantsplist
: parameter index list
-
void
getDky
(realtype t, int k) const = 0¶ interpolates the (derivative of the) solution at the requested timepoint
- Parameters
t
: timepointk
: derivative order
-
void
getDkyB
(realtype t, int k, int which) const = 0¶ interpolates the (derivative of the) solution at the requested timepoint
- Parameters
t
: timepointk
: derivative orderwhich
: index of backward problem
-
void
getSensDky
(realtype t, int k) const = 0¶ interpolates the (derivative of the) solution at the requested timepoint
- Parameters
t
: timepointk
: derivative order
-
void
getQuadDkyB
(realtype t, int k, int which) const = 0¶ interpolates the (derivative of the) solution at the requested timepoint
- Parameters
t
: timepointk
: derivative orderwhich
: index of backward problem
-
void
getQuadDky
(realtype t, int k) const = 0¶ interpolates the (derivative of the) solution at the requested timepoint
- Parameters
t
: timepointk
: derivative order
-
void
adjInit
() const = 0¶ initializes the adjoint problem
-
void
quadInit
(const AmiVector &xQ0) const = 0¶ initializes the quadratures
- Parameters
xQ0
: vector with initial values for xQ
-
void
allocateSolverB
(int *which) const = 0¶ Specifies solver method and initializes solver memory for the backward problem.
- Parameters
which
: identifier of the backwards problem
-
void
setSStolerancesB
(int which, realtype relTolB, realtype absTolB) const = 0¶ sets relative and absolute tolerances for the backward problem
- Parameters
which
: identifier of the backwards problemrelTolB
: relative tolerancesabsTolB
: absolute tolerances
-
void
quadSStolerancesB
(int which, realtype reltolQB, realtype abstolQB) const = 0¶ sets relative and absolute tolerances for the quadrature backward problem
- Parameters
which
: identifier of the backwards problemreltolQB
: relative tolerancesabstolQB
: absolute tolerances
-
void
quadSStolerances
(realtype reltolQB, realtype abstolQB) const = 0¶ sets relative and absolute tolerances for the quadrature problem
- Parameters
reltolQB
: relative tolerancesabstolQB
: absolute tolerances
-
void
getNumSteps
(const void *ami_mem, long int *numsteps) const = 0¶ reports the number of solver steps
- Parameters
ami_mem
: pointer to the solver memory instance (can be from forward or backward problem)numsteps
: output array
-
void
getNumRhsEvals
(const void *ami_mem, long int *numrhsevals) const = 0¶ reports the number of right hand evaluations
- Parameters
ami_mem
: pointer to the solver memory instance (can be from forward or backward problem)numrhsevals
: output array
-
void
getNumErrTestFails
(const void *ami_mem, long int *numerrtestfails) const = 0¶ reports the number of local error test failures
- Parameters
ami_mem
: pointer to the solver memory instance (can be from forward or backward problem)numerrtestfails
: output array
-
void
getNumNonlinSolvConvFails
(const void *ami_mem, long int *numnonlinsolvconvfails) const = 0¶ reports the number of nonlinear convergence failures
- Parameters
ami_mem
: pointer to the solver memory instance (can be from forward or backward problem)numnonlinsolvconvfails
: output array
-
void
getLastOrder
(const void *ami_mem, int *order) const = 0¶ Reports the order of the integration method during the last internal step.
- Parameters
ami_mem
: pointer to the solver memory instance (can be from forward or backward problem)order
: output array
-
void
initializeLinearSolver
(const Model *model) const¶ Initializes and sets the linear solver for the forward problem.
- Parameters
model
: pointer to the model object
-
void
initializeNonLinearSolver
() const¶ Sets the non-linear solver.
-
void
setLinearSolver
() const = 0¶ Sets the linear solver for the forward problem.
-
void
setLinearSolverB
(int which) const = 0¶ Sets the linear solver for the backward problem.
- Parameters
which
: index of the backward problem
-
void
setNonLinearSolver
() const = 0¶ Set the non-linear solver for the forward problem.
-
void
setNonLinearSolverB
(int which) const = 0¶ Set the non-linear solver for the backward problem.
- Parameters
which
: index of the backward problem
-
void
setNonLinearSolverSens
() const = 0¶ Set the non-linear solver for sensitivities.
-
void
initializeLinearSolverB
(const Model *model, int which) const¶ Initializes the linear solver for the backward problem.
- Parameters
model
: pointer to the model objectwhich
: index of the backward problem
-
void
initializeNonLinearSolverB
(int which) const¶ Initializes the non-linear solver for the backward problem.
- Parameters
which
: index of the backward problem
-
const Model *
getModel
() const = 0¶ Accessor function to the model stored in the user data
- Return
user data model
-
bool
getInitDone
() const¶ checks whether memory for the forward problem has been allocated
- Return
proxy for solverMemory->(cv|ida)_MallocDone
-
bool
getSensInitDone
() const¶ checks whether memory for forward sensitivities has been allocated
- Return
proxy for solverMemory->(cv|ida)_SensMallocDone
-
bool
getAdjInitDone
() const¶ checks whether memory for forward interpolation has been allocated
- Return
proxy for solverMemory->(cv|ida)_adjMallocDone
-
bool
getInitDoneB
(int which) const¶ checks whether memory for the backward problem has been allocated
- Return
proxy for solverMemoryB->(cv|ida)_MallocDone
- Parameters
which
: adjoint problem index
-
bool
getQuadInitDoneB
(int which) const¶ checks whether memory for backward quadratures has been allocated
- Return
proxy for solverMemoryB->(cv|ida)_QuadMallocDone
- Parameters
which
: adjoint problem index
-
bool
getQuadInitDone
() const¶ checks whether memory for quadratures has been allocated
- Return
proxy for solverMemory->(cv|ida)_QuadMallocDone
-
void
diag
() const = 0¶ attaches a diagonal linear solver to the forward problem
-
void
diagB
(int which) const = 0¶ attaches a diagonal linear solver to the backward problem
- Parameters
which
: identifier of the backwards problem
-
void
resetMutableMemory
(int nx, int nplist, int nquad) const¶ resets solverMemory and solverMemoryB
- Parameters
nx
: new number of state variablesnplist
: new number of sensitivity parametersnquad
: new number of quadratures (only differs from nplist for higher order sensitivity computation)
-
void *
getAdjBmem
(void *ami_mem, int which) const = 0¶ Retrieves the solver memory instance for the backward problem.
- Return
A (void *) pointer to the CVODES memory allocated for the backward problem.
- Parameters
which
: identifier of the backwards problemami_mem
: pointer to the forward solver memory instance
-
void
applyTolerances
() const¶ updates solver tolerances according to the currently specified member variables
-
void
applyTolerancesFSA
() const¶ updates FSA solver tolerances according to the currently specified member variables
-
void
applyTolerancesASA
(int which) const¶ updates ASA solver tolerances according to the currently specified member variables
- Parameters
which
: identifier of the backwards problem
-
void
applyQuadTolerancesASA
(int which) const¶ updates ASA quadrature solver tolerances according to the currently specified member variables
- Parameters
which
: identifier of the backwards problem
-
void
applyQuadTolerances
() const¶ updates quadrature solver tolerances according to the currently specified member variables
-
void
applySensitivityTolerances
() const¶ updates all sensitivity solver tolerances according to the currently specified member variables
-
void
setInitDone
() const¶ sets that memory for the forward problem has been allocated
-
void
setSensInitDone
() const¶ sets that memory for forward sensitivities has been allocated
-
void
setSensInitOff
() const¶ sets that memory for forward sensitivities has not been allocated
-
void
setAdjInitDone
() const¶ sets that memory for forward interpolation has been allocated
-
void
setInitDoneB
(int which) const¶ sets that memory for the backward problem has been allocated
- Parameters
which
: adjoint problem index
-
void
setQuadInitDoneB
(int which) const¶ sets that memory for backward quadratures has been allocated
- Parameters
which
: adjoint problem index
-
void
setQuadInitDone
() const¶ sets that memory for quadratures has been allocated
-
void
checkSensitivityMethod
(const SensitivityMethod sensi_meth, bool preequilibration) const¶ Sets sensitivity method (for simulation or preequilibration)
- Parameters
sensi_meth
: new value for sensi_meth[_preeq]preequilibration
: flag indicating preequilibration or simulation
Protected Attributes
-
std::unique_ptr<void, std::function<void(void*)>>
solver_memory_
¶ pointer to solver memory block
-
std::vector<std::unique_ptr<void, std::function<void(void*)>>>
solver_memory_B_
¶ pointer to solver memory block
-
InternalSensitivityMethod
ism_
= {InternalSensitivityMethod::simultaneous}¶ internal sensitivity method flag used to select the sensitivity solution method. Only applies for Forward Sensitivities.
-
LinearMultistepMethod
lmm_
= {LinearMultistepMethod::BDF}¶ specifies the linear multistep method.
-
NonlinearSolverIteration
iter_
= {NonlinearSolverIteration::newton}¶ specifies the type of nonlinear solver iteration
-
InterpolationType
interp_type_
= {InterpolationType::hermite}¶ interpolation type for the forward problem solution which is then used for the backwards problem.
-
long int
maxsteps_
= {10000}¶ maximum number of allowed integration steps
-
std::unique_ptr<SUNLinSolWrapper>
linear_solver_
¶ linear solver for the forward problem
-
std::unique_ptr<SUNLinSolWrapper>
linear_solver_B_
¶ linear solver for the backward problem
-
std::unique_ptr<SUNNonLinSolWrapper>
non_linear_solver_
¶ non-linear solver for the forward problem
-
std::unique_ptr<SUNNonLinSolWrapper>
non_linear_solver_B_
¶ non-linear solver for the backward problem
-
std::unique_ptr<SUNNonLinSolWrapper>
non_linear_solver_sens_
¶ non-linear solver for the sensitivities
-
bool
solver_was_called_F_
= {false}¶ flag indicating whether the forward solver has been called
-
bool
solver_was_called_B_
= {false}¶ flag indicating whether the backward solver has been called
-
AmiVectorArray
sx_
= {0, 0}¶ state sensitivities interface variable (dimension: nx_solver x nplist)
-
AmiVectorArray
sdx_
= {0, 0}¶ state derivative sensitivities dummy (dimension: nx_solver x nplist)
-
bool
force_reinit_postprocess_F_
= {false}¶ flag to force reInitPostProcessF before next call to solve
-
bool
force_reinit_postprocess_B_
= {false}¶ flag to force reInitPostProcessB before next call to solveB
-