Program Listing for File steadystateproblem.h
↰ Return to documentation for file (include/amici/steadystateproblem.h
)
#ifndef AMICI_STEADYSTATEPROBLEM_H
#define AMICI_STEADYSTATEPROBLEM_H
#include "amici/defines.h"
#include "amici/forwardproblem.h"
#include "amici/solver_cvodes.h"
#include "amici/vector.h"
#include <amici/newton_solver.h>
#include <nvector/nvector_serial.h>
#include <functional>
#include <memory>
namespace amici {
class ExpData;
class Solver;
class Model;
class SteadystateProblem {
public:
explicit SteadystateProblem(const Solver &solver, const Model &model);
void workSteadyStateProblem(const Solver &solver, Model &model, int it);
void workSteadyStateBackwardProblem(const Solver &solver, Model &model,
const BackwardProblem *bwd);
const SimulationState &getFinalSimulationState() const { return state_; };
const AmiVector &getEquilibrationQuadratures() const { return xQB_; }
const AmiVector &getState() const { return state_.x; };
const AmiVectorArray &getStateSensitivity() const { return state_.sx; };
std::vector<realtype> const &getDJydx() const { return dJydx_; }
double getCPUTime() const { return cpu_time_; }
double getCPUTimeB() const { return cpu_timeB_; }
std::vector<SteadyStateStatus> const &getSteadyStateStatus() const {
return steady_state_status_;
}
realtype getSteadyStateTime() const { return state_.t; }
realtype getResidualNorm() const { return wrms_; }
const std::vector<int> &getNumSteps() const { return numsteps_; }
int getNumStepsB() const { return numstepsB_; }
void getAdjointUpdates(Model &model, const ExpData &edata);
AmiVector const &getAdjointState() const { return xB_; }
AmiVector const &getAdjointQuadrature() const { return xQB_; }
bool hasQuadrature() const { return hasQuadrature_; }
bool checkSteadyStateSuccess() const;
private:
void findSteadyState(const Solver &solver, Model &model, int it);
void findSteadyStateByNewtonsMethod(Model &model, bool newton_retry);
void findSteadyStateBySimulation(const Solver &solver, Model &model,
int it);
void computeSteadyStateQuadrature(const Solver &solver, Model &model);
void getQuadratureByLinSolve(Model &model);
void getQuadratureBySimulation(const Solver &solver, Model &model);
[[noreturn]] void handleSteadyStateFailure();
void writeErrorString(std::string *errorString,
SteadyStateStatus status) const;
bool getSensitivityFlag(const Model &model, const Solver &solver, int it,
SteadyStateContext context);
realtype getWrmsNorm(AmiVector const &x, AmiVector const &xdot,
realtype atol, realtype rtol, AmiVector &ewt) const;
realtype getWrms(Model &model, SensitivityMethod sensi_method);
realtype getWrmsFSA(Model &model);
void applyNewtonsMethod(Model &model, bool newton_retry);
void runSteadystateSimulation(const Solver &solver, Model &model,
bool backward);
std::unique_ptr<Solver> createSteadystateSimSolver(const Solver &solver,
Model &model,
bool forwardSensis,
bool backward) const;
void initializeForwardProblem(int it, const Solver &solver, Model &model);
bool initializeBackwardProblem(const Solver &solver, Model &model,
const BackwardProblem *bwd);
void computeQBfromQ(Model &model, const AmiVector &yQ,
AmiVector &yQB) const;
bool makePositiveAndCheckConvergence(Model &model);
bool updateDampingFactor(bool step_successful);
void flagUpdatedState();
void updateSensiSimulation(const Solver &solver);
void updateRightHandSide(Model &model);
void getNewtonStep(Model &model);
AmiVector delta_;
AmiVector delta_old_;
AmiVector ewt_;
AmiVector ewtQB_;
AmiVector x_old_;
AmiVector xdot_;
AmiVectorArray sdx_;
AmiVector xB_;
AmiVector xQ_;
AmiVector xQB_;
AmiVector xQBdot_;
int max_steps_{0};
realtype wrms_{NAN};
std::vector<realtype> dJydx_;
SimulationState state_;
std::vector<int> numsteps_{std::vector<int>(3, 0)};
int numstepsB_{0};
double cpu_time_{0.0};
double cpu_timeB_{0.0};
bool hasQuadrature_{false};
double gamma_{1.0};
std::vector<SteadyStateStatus> steady_state_status_;
realtype atol_{NAN};
realtype rtol_{NAN};
realtype atol_sensi_{NAN};
realtype rtol_sensi_{NAN};
realtype atol_quad_{NAN};
realtype rtol_quad_{NAN};
std::unique_ptr<NewtonSolver> newton_solver_{nullptr};
NewtonDampingFactorMode damping_factor_mode_{NewtonDampingFactorMode::on};
realtype damping_factor_lower_bound_{1e-8};
bool newton_step_conv_ {false};
bool check_sensi_conv_ {true};
bool xdot_updated_ {false};
bool delta_updated_ {false};
bool sensis_updated_ {false};
};
} // namespace amici
#endif // STEADYSTATEPROBLEM_H