Chaste
Release::3.4
|
#include <AbstractCardiacCellInterface.hpp>
Public Member Functions | |
AbstractCardiacCellInterface (boost::shared_ptr< AbstractIvpOdeSolver > pOdeSolver, unsigned voltageIndex, boost::shared_ptr< AbstractStimulusFunction > pIntracellularStimulus) | |
virtual | ~AbstractCardiacCellInterface () |
virtual void | SetTimestep (double dt)=0 |
virtual unsigned | GetNumberOfStateVariables () const =0 |
virtual unsigned | GetNumberOfParameters () const =0 |
virtual std::vector< double > | GetStdVecStateVariables ()=0 |
virtual const std::vector < std::string > & | rGetStateVariableNames () const =0 |
virtual void | SetStateVariables (const std::vector< double > &rVariables)=0 |
virtual void | SetStateVariable (unsigned index, double newValue)=0 |
virtual void | SetStateVariable (const std::string &rName, double newValue)=0 |
virtual double | GetAnyVariable (const std::string &rName, double time)=0 |
virtual double | GetParameter (const std::string &rParameterName)=0 |
virtual double | GetParameter (unsigned parameterIndex)=0 |
virtual void | SetParameter (const std::string &rParameterName, double value)=0 |
virtual void | SolveAndUpdateState (double tStart, double tEnd)=0 |
virtual OdeSolution | Compute (double tStart, double tEnd, double tSamp=0.0)=0 |
virtual void | ComputeExceptVoltage (double tStart, double tEnd)=0 |
virtual double | GetIIonic (const std::vector< double > *pStateVariables=NULL)=0 |
virtual void | SetVoltage (double voltage)=0 |
virtual double | GetVoltage ()=0 |
unsigned | GetVoltageIndex () |
void | SetStimulusFunction (boost::shared_ptr< AbstractStimulusFunction > pStimulus) |
double | GetStimulus (double time) |
void | SetIntracellularStimulusFunction (boost::shared_ptr< AbstractStimulusFunction > pStimulus) |
double | GetIntracellularStimulus (double time) |
double | GetIntracellularAreaStimulus (double time) |
void | SetUsedInTissueSimulation (bool tissue=true) |
virtual boost::shared_ptr < RegularStimulus > | UseCellMLDefaultStimulus () |
bool | HasCellMLDefaultStimulus () |
virtual AbstractLookupTableCollection * | GetLookupTableCollection () |
boost::shared_ptr < AbstractStimulusFunction > | GetStimulusFunction () |
const boost::shared_ptr < AbstractStimulusFunction > | GetStimulusFunction () const |
const boost::shared_ptr < AbstractIvpOdeSolver > | GetSolver () const |
void | SetSolver (boost::shared_ptr< AbstractIvpOdeSolver > pSolver) |
virtual void | SetVoltageDerivativeToZero (bool clamp=true) |
void | SetFixedVoltage (double voltage) |
virtual void | SetStretch (double stretch) |
virtual double | GetIntracellularCalciumConcentration () |
Protected Attributes | |
unsigned | mVoltageIndex |
boost::shared_ptr < AbstractIvpOdeSolver > | mpOdeSolver |
boost::shared_ptr < AbstractStimulusFunction > | mpIntracellularStimulus |
bool | mSetVoltageDerivativeToZero |
bool | mIsUsedInTissue |
bool | mHasDefaultStimulusFromCellML |
double | mFixedVoltage |
Private Member Functions | |
template<class Archive > | |
void | serialize (Archive &archive, const unsigned int version) |
Friends | |
class | boost::serialization::access |
This class defines a common interface to AbstractCardiacCell and AbstractCvodeCell, primarily for the benefit of single-cell cardiac simulations, so that calling code doesn't need to know which solver is being used.
Strictly speaking this isn't an interface, since some methods have implementations defined. But the name AbstractCardiacCell was already taken.
Definition at line 58 of file AbstractCardiacCellInterface.hpp.
AbstractCardiacCellInterface::AbstractCardiacCellInterface | ( | boost::shared_ptr< AbstractIvpOdeSolver > | pOdeSolver, |
unsigned | voltageIndex, | ||
boost::shared_ptr< AbstractStimulusFunction > | pIntracellularStimulus | ||
) |
Create a new cardiac cell. The state variables of the cell will be set to AbstractOdeSystemInformation::GetInitialConditions(). Note that calls to SetDefaultInitialConditions() on a particular instance of this class will not modify its state variables. You can modify them directly with rGetStateVariables().
pOdeSolver | the ODE solver to use when simulating this cell (ignored for some subclasses; can be an empty pointer in these cases) |
voltageIndex | the index of the transmembrane potential within the system; see mVoltageIndex |
pIntracellularStimulus | the intracellular stimulus current |
Definition at line 53 of file AbstractCardiacCellInterface.cpp.
References Citations::Register().
|
virtual |
Virtual destructor
Definition at line 71 of file AbstractCardiacCellInterface.cpp.
|
pure virtual |
Simulates this cell's behaviour between the time interval [tStart, tEnd], and return state variable values. The timestep used will depend on the subclass implementation.
tStart | beginning of the time interval to simulate |
tEnd | end of the time interval to simulate |
tSamp | sampling interval for returned results (defaults to dt) |
Implemented in AbstractBackwardEulerCardiacCell< 0u >, AbstractCardiacCell, AbstractCvodeCell, AbstractBackwardEulerCardiacCell< SIZE >, AbstractGeneralizedRushLarsenCardiacCell, and AbstractRushLarsenCardiacCell.
|
pure virtual |
Simulates this cell's behaviour between the time interval [tStart, tEnd], but does not update the voltage. The timestep used will depend on the subclass implementation.
tStart | beginning of the time interval to simulate |
tEnd | end of the time interval to simulate |
Implemented in AbstractBackwardEulerCardiacCell< 0u >, AbstractCardiacCell, AbstractCvodeCell, AbstractBackwardEulerCardiacCell< SIZE >, AbstractGeneralizedRushLarsenCardiacCell, AbstractRushLarsenCardiacCell, and FakeBathCell.
|
pure virtual |
All subclasses must implement a method that returns a state variable value
This needs to be declared here as AbstractCardiacProblem uses it.
rName | variable name |
time | the time at which to get the variable (only needed when evaluating derived quantities). |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
|
pure virtual |
Computes the total current flowing through the cell membrane, using the current values of the state variables.
Chaste's value for C_m can be obtained from HeartConfig::Instance()->GetCapacitance() and is measured in uF/cm^2.
For state-variable interpolation (SVI) we need to interpolate state variables at nodes onto quadrature points and then pass these into this method (see optional argument). Otherwise for ionic-current interpolation (ICI) we use the cell's internal state variables.
pStateVariables | optionally can be supplied to evaluate the ionic current at the given state; by default the cell's internal state will be used. |
Implemented in CorriasBuistICCModified, CorriasBuistSMCModified, FakeBathCell, CML_noble_varghese_kohl_noble_1998_basic_with_sac, and FitzHughNagumo1961OdeSystem.
Referenced by MonodomainCorrectionTermAssembler< ELEM_DIM, SPACE_DIM >::ComputeVectorTerm(), and BidomainCorrectionTermAssembler< ELEM_DIM, SPACE_DIM >::ComputeVectorTerm().
time | the time at which to evaluate the stimulus |
Definition at line 106 of file AbstractCardiacCellInterface.cpp.
References GetIntracellularStimulus(), HeartConfig::GetSurfaceAreaToVolumeRatio(), HeartConfig::Instance(), and mIsUsedInTissue.
Referenced by FitzHughNagumo1961OdeSystem::EvaluateYDerivatives().
|
virtual |
[Ca_i] is needed for mechanics, so we explicitly have a Get method (rather than use a get by name type method, to avoid inefficiency when using different types of cells). This method by default throws an exception, so should be implemented in the concrete class if intracellular (cytosolic) calcium concentration is one of the state variables.
Reimplemented in CML_noble_varghese_kohl_noble_1998_basic_with_sac, and FakeBathCell.
Definition at line 173 of file AbstractCardiacCellInterface.cpp.
References EXCEPTION.
time | the time at which to evaluate the stimulus |
Definition at line 100 of file AbstractCardiacCellInterface.cpp.
References mpIntracellularStimulus.
Referenced by GetIntracellularAreaStimulus(), and GetStimulus().
|
inlinevirtual |
Must be implemented by subclasses iff they support lookup tables.
Definition at line 346 of file AbstractCardiacCellInterface.hpp.
Referenced by HeartConfigRelatedCellFactory< SPACE_DIM >::CreateCellWithIntracellularStimulus().
|
pure virtual |
All subclasses must implement this method to get the number of parameters
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
|
pure virtual |
All subclasses must implement this method to get the number of state variables.
This needs to be declared here as AbstractCorrectionTermAssembler uses it.
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
Referenced by AbstractCardiacTissue< ELEMENT_DIM, SPACE_DIM >::SolveCellSystems().
|
pure virtual |
All subclasses must implement a method that returns a parameter value.
This needs to be declared here as HeartConfigCellFactory uses it.
rParameterName | the name of a parameter to get the value of, |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
Referenced by HeartConfigRelatedCellFactory< SPACE_DIM >::SetCellParameters().
All subclasses must implement a method that returns a parameter value.
This needs to be declared here as HeartConfigCellFactory uses it.
parameterIndex | the index of a parameter to get the value of, |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
const boost::shared_ptr< AbstractIvpOdeSolver > AbstractCardiacCellInterface::GetSolver | ( | ) | const |
For boost archiving use only (that's why the 'consts' are required)
Definition at line 149 of file AbstractCardiacCellInterface.cpp.
References mpOdeSolver.
|
pure virtual |
All subclasses must implement a method that returns state variables as a std::vector (even if they work with another format), for compatibility with the PDE solvers.
This needs to be declared here as AbstractCorrectionTermAssembler uses it.
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
Referenced by AbstractCardiacTissue< ELEMENT_DIM, SPACE_DIM >::SolveCellSystems().
time | the time at which to evaluate the stimulus |
Definition at line 88 of file AbstractCardiacCellInterface.cpp.
References GetIntracellularStimulus().
Referenced by CML_noble_varghese_kohl_noble_1998_basic_with_sac::EvaluateYDerivatives(), CorriasBuistSMCModified::EvaluateYDerivatives(), and CorriasBuistICCModified::EvaluateYDerivatives().
boost::shared_ptr< AbstractStimulusFunction > AbstractCardiacCellInterface::GetStimulusFunction | ( | ) |
Definition at line 138 of file AbstractCardiacCellInterface.cpp.
References mpIntracellularStimulus.
Referenced by AbstractPurkinjeCellFactory< ELEMENT_DIM, SPACE_DIM >::CreateJunction().
const boost::shared_ptr< AbstractStimulusFunction > AbstractCardiacCellInterface::GetStimulusFunction | ( | ) | const |
For boost archiving use only (that's why the 'consts' are required)
Definition at line 144 of file AbstractCardiacCellInterface.cpp.
References mpIntracellularStimulus.
|
pure virtual |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
Referenced by PurkinjeVentricularJunctionStimulus::GetStimulus(), and SetVoltageDerivativeToZero().
unsigned AbstractCardiacCellInterface::GetVoltageIndex | ( | ) |
Definition at line 76 of file AbstractCardiacCellInterface.cpp.
References mVoltageIndex.
Referenced by AbstractGeneralizedRushLarsenCardiacCell::Compute(), AbstractGeneralizedRushLarsenCardiacCell::SolveAndUpdateState(), and AbstractRushLarsenCardiacCell::UpdateTransmembranePotential().
bool AbstractCardiacCellInterface::HasCellMLDefaultStimulus | ( | ) |
Definition at line 133 of file AbstractCardiacCellInterface.cpp.
References mHasDefaultStimulusFromCellML.
Referenced by CellMLLoader::LoadCellMLFile().
|
pure virtual |
All subclasses must implement this method to get variable names.
Needs to be declared here as AdaptiveBidomainProblem uses it.
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
|
inlineprivate |
Main boost serialization method, some member variables are handled by Save/Load constructs.
archive | the archive file. |
version | the version of archiving, defined in the macro at the bottom of the class. |
Definition at line 458 of file AbstractCardiacCellInterface.hpp.
References mHasDefaultStimulusFromCellML, mIsUsedInTissue, and mSetVoltageDerivativeToZero.
void AbstractCardiacCellInterface::SetFixedVoltage | ( | double | voltage | ) |
When the voltage derivative has been set to zero by SetVoltageDerivativeToZero, this method sets the transmembrane potential to use when computing the other derivatives and total ionic current.
voltage | the value of the transmembrane potential |
Definition at line 168 of file AbstractCardiacCellInterface.cpp.
References mFixedVoltage.
Referenced by AbstractCvodeCell::SetVoltage(), and AbstractCardiacCell::SetVoltage().
void AbstractCardiacCellInterface::SetIntracellularStimulusFunction | ( | boost::shared_ptr< AbstractStimulusFunction > | pStimulus | ) |
Set the intracellular stimulus. This should have units of uA/cm^2 for single-cell problems, or uA/cm^3 in a tissue simulation.
pStimulus | new stimulus function |
Definition at line 94 of file AbstractCardiacCellInterface.cpp.
References mpIntracellularStimulus.
Referenced by HeartConfigRelatedCellFactory< SPACE_DIM >::SetCellIntracellularStimulus(), and SetStimulusFunction().
|
pure virtual |
All subclasses must implement a method that sets a parameter value.
This needs to be declared here as HeartConfigCellFactory uses it.
rParameterName | the parameter name to set the value of, |
value | value to set it to. |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
Referenced by HeartConfigRelatedCellFactory< SPACE_DIM >::SetCellParameters().
void AbstractCardiacCellInterface::SetSolver | ( | boost::shared_ptr< AbstractIvpOdeSolver > | pSolver | ) |
Change the ODE solver used by this cell. Note that this will only work for cell models which do not have a solver built-in (e.g. it will NOT work for AbstractCvodeCell derivatives).
pSolver | the new solver to use |
Definition at line 154 of file AbstractCardiacCellInterface.cpp.
References mpOdeSolver.
|
pure virtual |
Set the value of a single state variable in the ODE system.
This needs to be declared here as AdaptiveBidomainProblem uses it.
index | index of the state variable to be set |
newValue | new value of the state variable |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
|
pure virtual |
Set the value of a single state variable in the ODE system.
This needs to be declared here as AdaptiveBidomainProblem uses the above, to avoid compiler confusion...!
rName | name of the state variable to be set |
newValue | new value of the state variable |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
|
pure virtual |
All subclasses must implement this method to set the state variables.
This needs to be declared here as AbstractCardiacTissue uses it.
rVariables | the state variables to take a copy of. |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
Referenced by AbstractCardiacTissue< ELEMENT_DIM, SPACE_DIM >::SolveCellSystems().
void AbstractCardiacCellInterface::SetStimulusFunction | ( | boost::shared_ptr< AbstractStimulusFunction > | pStimulus | ) |
Set the intracellular stimulus. Shorthand for SetIntracellularStimulusFunction.
pStimulus | new stimulus function |
Definition at line 82 of file AbstractCardiacCellInterface.cpp.
References SetIntracellularStimulusFunction().
Referenced by AbstractPurkinjeCellFactory< ELEMENT_DIM, SPACE_DIM >::CreateJunction().
|
inlinevirtual |
In electromechanics problems, the stretch is passed back to cell-model in case mechano-electric feedback has been modelled. We define an empty method here. Stretch-dependent cell models should overload this method to use the input stretch accordingly.
stretch | the stretch of the cell in the axial direction |
Reimplemented in CML_noble_varghese_kohl_noble_1998_basic_with_sac.
Definition at line 403 of file AbstractCardiacCellInterface.hpp.
|
pure virtual |
Set the timestep (or maximum timestep when using CVODE) to use for simulating this cell.
dt | the timestep |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
void AbstractCardiacCellInterface::SetUsedInTissueSimulation | ( | bool | tissue = true | ) |
Set whether this cell object exists in the context of a tissue simulation, or can be used for single cell simulations. This affects the units of the intracellular stimulus (see GetIntracellularStimulus) and so is used by GetIntracellularAreaStimulus to perform a units conversion if necessary.
tissue | true if cell is in a tissue |
Definition at line 121 of file AbstractCardiacCellInterface.cpp.
References mIsUsedInTissue.
Referenced by AbstractCardiacTissue< ELEMENT_DIM, SPACE_DIM >::SetUpHaloCells().
|
pure virtual |
Set the cellular transmembrane potential.
voltage | new value |
Implemented in AbstractCardiacCell, and AbstractCvodeCell.
|
virtual |
Set whether to clamp the voltage by setting its derivative to zero.
clamp | whether to clamp |
Reimplemented in AbstractCvodeCell.
Definition at line 159 of file AbstractCardiacCellInterface.cpp.
References GetVoltage(), mFixedVoltage, and mSetVoltageDerivativeToZero.
Referenced by AbstractGeneralizedRushLarsenCardiacCell::ComputeExceptVoltage(), AbstractCvodeCell::ComputeExceptVoltage(), AbstractCardiacCell::ComputeExceptVoltage(), and AbstractCvodeCell::SetVoltageDerivativeToZero().
|
pure virtual |
Simulate this cell's behaviour between the time interval [tStart, tEnd], updating the internal state variable values. The timestep used will depend on the subclass implementation.
tStart | beginning of the time interval to simulate |
tEnd | end of the time interval to simulate |
Implemented in AbstractBackwardEulerCardiacCell< 0u >, AbstractCardiacCell, AbstractCvodeCell, AbstractBackwardEulerCardiacCell< SIZE >, AbstractGeneralizedRushLarsenCardiacCell, and AbstractRushLarsenCardiacCell.
|
virtual |
Use CellML metadata to set up the default stimulus for this cell. By default this method will always throw an exception. For suitably annotated models, PyCml will override this to provide a RegularStimulus as defined in the CellML.
Definition at line 126 of file AbstractCardiacCellInterface.cpp.
References EXCEPTION, and mHasDefaultStimulusFromCellML.
Referenced by CellMLLoader::LoadCellMLFile().
|
friend |
Needed for serialization.
Definition at line 450 of file AbstractCardiacCellInterface.hpp.
|
protected |
The value of the fixed voltage if mSetVoltageDerivativeToZero is set.
Definition at line 446 of file AbstractCardiacCellInterface.hpp.
Referenced by SetFixedVoltage(), and SetVoltageDerivativeToZero().
|
protected |
Whether this cell has a default stimulus specified by CellML metadata.
Definition at line 443 of file AbstractCardiacCellInterface.hpp.
Referenced by HasCellMLDefaultStimulus(), AbstractCardiacCell::serialize(), serialize(), and UseCellMLDefaultStimulus().
|
protected |
Whether this cell exists in a tissue, or is an isolated cell.
Definition at line 440 of file AbstractCardiacCellInterface.hpp.
Referenced by GetIntracellularAreaStimulus(), AbstractCardiacCell::serialize(), serialize(), and SetUsedInTissueSimulation().
|
protected |
The intracellular stimulus current.
Definition at line 430 of file AbstractCardiacCellInterface.hpp.
Referenced by GetIntracellularStimulus(), GetStimulusFunction(), and SetIntracellularStimulusFunction().
|
protected |
Pointer to the solver used to simulate this cell.
Definition at line 427 of file AbstractCardiacCellInterface.hpp.
Referenced by AbstractCardiacCell::Compute(), AbstractCardiacCell::ComputeExceptVoltage(), GetSolver(), SetSolver(), and AbstractCardiacCell::SolveAndUpdateState().
|
protected |
Flag set to true if ComputeExceptVoltage is called, to indicate to subclass EvaluateYDerivatives methods that V should be considered fixed, and hence dV/dt set to zero.
Definition at line 437 of file AbstractCardiacCellInterface.hpp.
Referenced by AbstractRushLarsenCardiacCell::ComputeExceptVoltage(), FitzHughNagumo1961OdeSystem::EvaluateYDerivatives(), CML_noble_varghese_kohl_noble_1998_basic_with_sac::EvaluateYDerivatives(), CorriasBuistSMCModified::EvaluateYDerivatives(), CorriasBuistICCModified::EvaluateYDerivatives(), AbstractCardiacCell::serialize(), serialize(), AbstractCvodeCell::SetVoltageDerivativeToZero(), and SetVoltageDerivativeToZero().
|
protected |
The index of the voltage within this system. Usually it will be an index into the state variable vector, but this is not guaranteed. It will however always be suitable for use with AbstractParameterisedSystem::GetAnyVariable and OdeSolution::GetVariableAtIndex.
Definition at line 424 of file AbstractCardiacCellInterface.hpp.
Referenced by FitzHughNagumo1961OdeSystem::GetIIonic(), AbstractCvodeCell::GetVoltage(), AbstractCardiacCell::GetVoltage(), GetVoltageIndex(), AbstractCvodeCell::SetVoltage(), and AbstractCardiacCell::SetVoltage().