Chaste  Release::2018.1
NhsModelWithBackwardSolver Class Reference

#include <NhsModelWithBackwardSolver.hpp>

+ Inheritance diagram for NhsModelWithBackwardSolver:
+ Collaboration diagram for NhsModelWithBackwardSolver:

Public Member Functions

 NhsModelWithBackwardSolver ()
 
void RunDoNotUpdate (double startTime, double endTime, double timestep)
 
double GetNextActiveTension ()
 
void RunAndUpdate (double startTime, double endTime, double timestep)
 
- Public Member Functions inherited from NhsContractionModel
 NhsContractionModel ()
 
void SetStretchAndStretchRate (double lambda, double dlambdaDt)
 
void SetInputParameters (ContractionModelInputParameters &rInputParameters)
 
void SetIntracellularCalciumConcentration (double calciumConcentration)
 
double GetCalciumTroponinValue ()
 
void EvaluateYDerivatives (double time, const std::vector< double > &rY, std::vector< double > &rDY)
 
double GetActiveTension ()
 
bool IsStretchDependent ()
 
bool IsStretchRateDependent ()
 
- Public Member Functions inherited from AbstractOdeBasedContractionModel
 AbstractOdeBasedContractionModel (unsigned numStateVariables)
 
void UpdateStateVariables ()
 
- Public Member Functions inherited from AbstractOdeSystem
 AbstractOdeSystem (unsigned numberOfStateVariables)
 
virtual ~AbstractOdeSystem ()
 
virtual bool CalculateStoppingEvent (double time, const std::vector< double > &rY)
 
virtual double CalculateRootFunction (double time, const std::vector< double > &rY)
 
bool GetUseAnalyticJacobian ()
 
const std::vector< double > & rGetConstStateVariables () const
 
- Public Member Functions inherited from AbstractParameterisedSystem< std::vector< double > >
 AbstractParameterisedSystem (unsigned numberOfStateVariables)
 
std::vector< double > & rGetStateVariables ()
 
std::vector< doubleGetStateVariables ()
 
void SetStateVariables (const std::vector< double > &rStateVariables)
 
double GetStateVariable (unsigned index) const
 
double GetStateVariable (const std::string &rName) const
 
void SetStateVariable (unsigned index, double newValue)
 
void SetStateVariable (const std::string &rName, double newValue)
 
virtual void VerifyStateVariables ()
 
void SetDefaultInitialConditions (const std::vector< double > &rInitialConditions)
 
void SetDefaultInitialCondition (unsigned index, double initialCondition)
 
std::vector< doubleGetInitialConditions () const
 
void ResetToInitialConditions ()
 
double GetParameter (unsigned index) const
 
double GetParameter (const std::string &rName) const
 
void SetParameter (const std::string &rName, double value)
 
void SetParameter (unsigned index, double value)
 
double GetAnyVariable (unsigned index, double time=0.0, std::vector< double > *pDerivedQuantities=NULL)
 
double GetAnyVariable (const std::string &rName, double time=0.0, std::vector< double > *pDerivedQuantities=NULL)
 
void SetAnyVariable (unsigned index, double value)
 
void SetAnyVariable (const std::string &rName, double value)
 
virtual std::vector< doubleComputeDerivedQuantities (double time, const std::vector< double > &rState)
 
std::vector< doubleComputeDerivedQuantitiesFromCurrentState (double time)
 
- Public Member Functions inherited from AbstractUntemplatedParameterisedSystem
 AbstractUntemplatedParameterisedSystem (unsigned numberOfStateVariables)
 
virtual ~AbstractUntemplatedParameterisedSystem ()
 
boost::shared_ptr< const
AbstractOdeSystemInformation
GetSystemInformation () const
 
std::string GetSystemName () const
 
unsigned GetNumberOfAttributes () const
 
bool HasAttribute (const std::string &rName) const
 
double GetAttribute (const std::string &rName) const
 
unsigned GetNumberOfStateVariables () const
 
const std::vector< std::string > & rGetStateVariableNames () const
 
const std::vector< std::string > & rGetStateVariableUnits () const
 
unsigned GetStateVariableIndex (const std::string &rName) const
 
bool HasStateVariable (const std::string &rName) const
 
std::string GetStateVariableUnits (unsigned index) const
 
unsigned GetNumberOfParameters () const
 
const std::vector< std::string > & rGetParameterNames () const
 
const std::vector< std::string > & rGetParameterUnits () const
 
unsigned GetParameterIndex (const std::string &rName) const
 
bool HasParameter (const std::string &rName) const
 
std::string GetParameterUnits (unsigned index) const
 
unsigned GetNumberOfDerivedQuantities () const
 
const std::vector< std::string > & rGetDerivedQuantityNames () const
 
const std::vector< std::string > & rGetDerivedQuantityUnits () const
 
unsigned GetDerivedQuantityIndex (const std::string &rName) const
 
bool HasDerivedQuantity (const std::string &rName) const
 
std::string GetDerivedQuantityUnits (unsigned index) const
 
unsigned GetAnyVariableIndex (const std::string &rName) const
 
bool HasAnyVariable (const std::string &rName) const
 
std::string GetAnyVariableUnits (unsigned index) const
 
std::string GetAnyVariableUnits (const std::string &rName) const
 
- Public Member Functions inherited from AbstractContractionModel
 AbstractContractionModel ()
 
virtual ~AbstractContractionModel ()
 
void SetStretch (double stretch)
 

Private Member Functions

double ImplicitSolveForQ ()
 
void CalculateCaTropAndZDerivatives (double calciumTroponin, double z, double Q, double &dCaTrop, double &dz)
 
void CalculateBackwardEulerResidual (double calciumTroponin, double z, double Q, double &residualComponent1, double &residualComponent2)
 

Private Attributes

double mDt
 

Static Private Attributes

static const double mTolerance = 1e-10
 

Additional Inherited Members

- Protected Member Functions inherited from NhsContractionModel
void CalculateCalciumTrop50 ()
 
double CalculateT0 (double z)
 
- Protected Member Functions inherited from AbstractParameterisedSystem< std::vector< double > >
std::string DumpState (const std::string &rMessage)
 
std::string DumpState (const std::string &rMessage, std::vector< double >Y)
 
std::string DumpState (const std::string &rMessage, std::vector< double >Y, double time)
 
void CheckParametersOnLoad (const std::vector< double > &rParameters, const std::vector< std::string > &rParameterNames)
 
- Protected Attributes inherited from NhsContractionModel
double mLambda
 
double mDLambdaDt
 
double mCalciumI
 
double mCalciumTrop50
 
double mK1
 
double mK2
 
- Protected Attributes inherited from AbstractOdeBasedContractionModel
std::vector< doublemTemporaryStateVariables
 
double mTime
 
- Protected Attributes inherited from AbstractOdeSystem
bool mUseAnalyticJacobian
 
- Protected Attributes inherited from AbstractParameterisedSystem< std::vector< double > >
std::vector< doublemStateVariables
 
std::vector< doublemParameters
 
- Protected Attributes inherited from AbstractUntemplatedParameterisedSystem
unsigned mNumberOfStateVariables
 
boost::shared_ptr
< AbstractOdeSystemInformation
mpSystemInfo
 
- Static Protected Attributes inherited from NhsContractionModel
static const double mKon = 100
 
static const double mKrefoff = 0.2
 
static const double mGamma = 2
 
static const double mCalciumTroponinMax = 0.07
 
static const double mAlphaR1 = 0.002
 
static const double mAlphaR2 = 0.0017
 
static const double mKZ = 0.15
 
static const unsigned mNr = 3u
 
static const double mBeta1 = -4
 
static const double mAlpha0 = 0.008
 
static const unsigned mN = 3u
 
static const double mZp = 0.85
 
static const double mCalcium50ref = 0.00105
 
static const double mTref = 56.2
 
static const double mBeta0 = 4.9
 
static const double mA = 0.35
 
static const double mA1 = -29
 
static const double mA2 = 138
 
static const double mA3 = 129
 
static const double mAlpha1 = 0.03
 
static const double mAlpha2 = 0.130
 
static const double mAlpha3 = 0.625
 

Detailed Description

The full backward Euler method on the NHS system compute the solution w^{n+1} in R^5 by: w^{n+1} - dt g(w^{n+1}) = w^n where the ODE system is dw/dt = g(w).

This would involve solving a 5d nonlinear system at each timestep. However, only g1 and g2 (ie dCatrop/dt and dz/dt) are nonlinear, g3, g4 and g5 (corresponding to dQi/dt) are linear and uncoupled. Therefore the backward euler solutions of Q1,Q2,Q3 can be computed immediately, leaving a 2D nonlinear system to be solved using newton's method

Definition at line 56 of file NhsModelWithBackwardSolver.hpp.

Constructor & Destructor Documentation

NhsModelWithBackwardSolver::NhsModelWithBackwardSolver ( )

Member Function Documentation

void NhsModelWithBackwardSolver::CalculateBackwardEulerResidual ( double  calciumTroponin,
double  z,
double  Q,
double residualComponent1,
double residualComponent2 
)
private

Compute the residual function for the 2D nonlinear system when Backward Euler is used The backward Euler discretisation is: w^{n+1} - dt g(w^{n+1}) = w^n where: w = (ca_trop, z) and the ODE system is: dw/dt = g(w)

so the residual is: f = w^{n+1} - dt g(w^{n+1}) - w^n

Parameters
calciumTroponinCurrent guess for Ca_trop value
zcurrent guess for z
Q= Q1+Q2+Q3, where Qi already computed at next timestep
residualComponent1Returned value - first component of residual
residualComponent2Returned value - second component of residual

Definition at line 95 of file NhsModelWithBackwardSolver.cpp.

References CalculateCaTropAndZDerivatives(), mDt, and AbstractOdeBasedContractionModel::mTemporaryStateVariables.

Referenced by RunDoNotUpdate().

void NhsModelWithBackwardSolver::CalculateCaTropAndZDerivatives ( double  calciumTroponin,
double  z,
double  Q,
double dCaTrop,
double dz 
)
private

The same as EvaluateYDerivatives in NhsContractionModel, but doesn't use std::vectors (for efficiency, and because this class is hardcoded and hand-optimised for backward euler), and just returns the derivatives of the first two components

Parameters
calciumTroponinCalcium troponin
zz
QQ=Q1+Q2+Q3
dCaTropthe returned value of dCaTrop/dt
dzthe returned value of dz/dt

Definition at line 54 of file NhsModelWithBackwardSolver.cpp.

References NhsContractionModel::CalculateT0(), EXCEPTION, NhsContractionModel::mA, NhsContractionModel::mAlpha0, NhsContractionModel::mAlphaR1, NhsContractionModel::mAlphaR2, NhsContractionModel::mCalciumI, NhsContractionModel::mCalciumTrop50, NhsContractionModel::mCalciumTroponinMax, NhsContractionModel::mGamma, NhsContractionModel::mKon, NhsContractionModel::mKrefoff, NhsContractionModel::mKZ, NhsContractionModel::mN, NhsContractionModel::mNr, and NhsContractionModel::mTref.

Referenced by CalculateBackwardEulerResidual().

double NhsModelWithBackwardSolver::GetNextActiveTension ( )
virtual
Returns
the active tension corresponding to the stored state variables computed from the last RunDoNotUpdate(), ie the active tension at the next time. Note that calling GetActiveTension() on the base class will use the internal state variables and return the active tension at the last time, if RunDoNotUpdate() has been called but UpdateStateVariables() has not

Reimplemented from NhsContractionModel.

Definition at line 180 of file NhsModelWithBackwardSolver.cpp.

References NhsContractionModel::CalculateT0(), NhsContractionModel::mA, and AbstractOdeBasedContractionModel::mTemporaryStateVariables.

double NhsModelWithBackwardSolver::ImplicitSolveForQ ( )
private

Solve for Q1,Q2,Q3 (and therefore Q) implicitly using backward euler. These can be done directly as the rhs is linear in Qi

Returns
Q=Q1+Q2+Q3

Definition at line 45 of file NhsModelWithBackwardSolver.cpp.

References NhsContractionModel::mA1, NhsContractionModel::mA2, NhsContractionModel::mA3, NhsContractionModel::mAlpha1, NhsContractionModel::mAlpha2, NhsContractionModel::mAlpha3, NhsContractionModel::mDLambdaDt, mDt, and AbstractOdeBasedContractionModel::mTemporaryStateVariables.

Referenced by RunDoNotUpdate().

void NhsModelWithBackwardSolver::RunAndUpdate ( double  startTime,
double  endTime,
double  timestep 
)
virtual

Overload the RunAndUpdate() method too, as that would use the base class's default (euler) solver

Parameters
startTime
endTime
timestep

Reimplemented from AbstractOdeBasedContractionModel.

Definition at line 195 of file NhsModelWithBackwardSolver.cpp.

References RunDoNotUpdate(), and AbstractOdeBasedContractionModel::UpdateStateVariables().

void NhsModelWithBackwardSolver::RunDoNotUpdate ( double  startTime,
double  endTime,
double  timestep 
)
virtual

Solves for the new state variables at the given end time using the implicit method. Note that the internal state variables are not altered, the solution is saved instead. Call UpdateStateVariables() to update, and GetNextActiveTension() to get the solved active tension

The state variables are not updated because this solve will be called as part of the newton iteration (ie guess stretch, see what the new active tension is) in a fully implicit method.

Note: overloaded from the method in AbstractOdeBasedContractionModel, which just does a simple Euler solve

Parameters
startTime
endTime
timestep

Reimplemented from AbstractOdeBasedContractionModel.

Definition at line 111 of file NhsModelWithBackwardSolver.cpp.

References TimeStepper::AdvanceOneTimeStep(), CalculateBackwardEulerResidual(), ImplicitSolveForQ(), TimeStepper::IsTimeAtEnd(), mDt, AbstractParameterisedSystem< std::vector< double > >::mStateVariables, AbstractOdeBasedContractionModel::mTemporaryStateVariables, and mTolerance.

Referenced by RunAndUpdate().

Member Data Documentation

double NhsModelWithBackwardSolver::mDt
private

Timestep for the ODEs solving

Definition at line 63 of file NhsModelWithBackwardSolver.hpp.

Referenced by CalculateBackwardEulerResidual(), ImplicitSolveForQ(), and RunDoNotUpdate().

const double NhsModelWithBackwardSolver::mTolerance = 1e-10
staticprivate

Tolerance for solving nonlinear system which require newton's method

Definition at line 60 of file NhsModelWithBackwardSolver.hpp.

Referenced by RunDoNotUpdate().


The documentation for this class was generated from the following files: