AbstractCardiacCellInterface Class Reference

#include <AbstractCardiacCellInterface.hpp>

Inheritance diagram for AbstractCardiacCellInterface:

Inheritance graph
[legend]
Collaboration diagram for AbstractCardiacCellInterface:

Collaboration graph
[legend]

List of all members.

Public Member Functions

 AbstractCardiacCellInterface (boost::shared_ptr< AbstractIvpOdeSolver > pOdeSolver, unsigned voltageIndex, boost::shared_ptr< AbstractStimulusFunction > pIntracellularStimulus)
virtual ~AbstractCardiacCellInterface ()
virtual void ResetToInitialConditions ()=0
virtual void SetTimestep (double dt)=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 void UseCellMLDefaultStimulus ()
bool HasCellMLDefaultStimulus ()
virtual void VerifyStateVariables ()
virtual
AbstractLookupTableCollection
GetLookupTableCollection ()
boost::shared_ptr
< AbstractStimulusFunction
GetStimulusFunction ()
const boost::shared_ptr
< AbstractStimulusFunction
GetStimulusFunction () const
const boost::shared_ptr
< AbstractIvpOdeSolver
GetSolver () const

Protected Attributes

unsigned mVoltageIndex
boost::shared_ptr
< AbstractIvpOdeSolver
mpOdeSolver
boost::shared_ptr
< AbstractStimulusFunction
mpIntracellularStimulus
bool mSetVoltageDerivativeToZero
bool mIsUsedInTissue
bool mHasDefaultStimulusFromCellML

Private Member Functions

template<class Archive>
void serialize (Archive &archive, const unsigned int version)

Friends

class boost::serialization::access


Detailed Description

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.

Note that serialization is not defined for this class. AbstractCvodeCell does not have (or need) serialization at present, so adding serialization here would break archive backwards compatibility for little gain.

Definition at line 54 of file AbstractCardiacCellInterface.hpp.


Constructor & Destructor Documentation

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().

Parameters:
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 36 of file AbstractCardiacCellInterface.cpp.

AbstractCardiacCellInterface::~AbstractCardiacCellInterface (  )  [virtual]

Virtual destructor

Definition at line 50 of file AbstractCardiacCellInterface.cpp.


Member Function Documentation

virtual void AbstractCardiacCellInterface::ResetToInitialConditions (  )  [pure virtual]

Reset the model's state variables to the default initial conditions.

Implemented in AbstractCardiacCell, and AbstractCvodeCell.

virtual void AbstractCardiacCellInterface::SetTimestep ( double  dt  )  [pure virtual]

Set the timestep (or maximum timestep) to use for simulating this cell.

Parameters:
dt the timestep

Implemented in AbstractCardiacCell, and AbstractCvodeCell.

virtual void AbstractCardiacCellInterface::SolveAndUpdateState ( double  tStart,
double  tEnd 
) [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.

Parameters:
tStart beginning of the time interval to simulate
tEnd end of the time interval to simulate

Implemented in AbstractBackwardEulerCardiacCell< SIZE >, AbstractCardiacCell, and AbstractCvodeCell.

virtual OdeSolution AbstractCardiacCellInterface::Compute ( double  tStart,
double  tEnd,
double  tSamp = 0.0 
) [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.

Parameters:
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< SIZE >, AbstractCardiacCell, and AbstractCvodeCell.

virtual void AbstractCardiacCellInterface::ComputeExceptVoltage ( double  tStart,
double  tEnd 
) [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.

Parameters:
tStart beginning of the time interval to simulate
tEnd end of the time interval to simulate

Implemented in AbstractBackwardEulerCardiacCell< SIZE >, AbstractCardiacCell, AbstractCvodeCell, and FakeBathCell.

virtual double AbstractCardiacCellInterface::GetIIonic ( const std::vector< double > *  pStateVariables = NULL  )  [pure virtual]

Computes the total current flowing through the cell membrane, using the current values of the state variables.

Should return a value in units of microamps/cm^2. Note that many cell models do not use these dimensions (let alone these units) and so a complex conversion is required. There are 2 main cases:

  • Cell model uses amps per unit capacitance. Often in this case the units used for the cell capacitance don't make sense (e.g. uF/cm^2 is used, and dV/dt=I/C_m). Hence we suggest examining the equation for dV/dt given in the model to determine what the cell model really considers the value for C_m to be, and scaling by Chaste's C_m / cell model C_m (the latter implicitly being dimensionless).
  • Cell model uses amps. In this case you need to divide by an estimate of the cell surface area. Assuming the model represents a single cell, and gives C_m in farads, then scaling by Chaste's C_m / model C_m seems reasonable. If the 'cell model' doesn't actually represent a single whole cell, then you'll have to be more careful. In both cases additional scaling may be required to obtain correct units once the dimensions have been sorted out.

Chaste's value for C_m can be obtained from HeartConfig::Instance()->GetCapacitance() and is measured in uF/cm^2.

Parameters:
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 FakeBathCell, FitzHughNagumo1961OdeSystem, and CML_noble_varghese_kohl_noble_1998_basic_with_sac.

Referenced by MonodomainCorrectionTermAssembler< ELEM_DIM, SPACE_DIM >::ComputeVectorTerm(), and BidomainCorrectionTermAssembler< ELEM_DIM, SPACE_DIM >::ComputeVectorTerm().

virtual void AbstractCardiacCellInterface::SetVoltage ( double  voltage  )  [pure virtual]

Set the transmembrane potential

Parameters:
voltage new value

Implemented in AbstractCardiacCell, and AbstractCvodeCell.

virtual double AbstractCardiacCellInterface::GetVoltage (  )  [pure virtual]

Get the current value of the transmembrane potential, as given in our state variable vector.

Implemented in AbstractCardiacCell, and AbstractCvodeCell.

unsigned AbstractCardiacCellInterface::GetVoltageIndex (  ) 

Get the index of the transmembrane potential 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 55 of file AbstractCardiacCellInterface.cpp.

References mVoltageIndex.

void AbstractCardiacCellInterface::SetStimulusFunction ( boost::shared_ptr< AbstractStimulusFunction pStimulus  ) 

Set the intracellular stimulus. Shorthand for SetIntracellularStimulusFunction.

Parameters:
pStimulus new stimulus function

Definition at line 61 of file AbstractCardiacCellInterface.cpp.

References SetIntracellularStimulusFunction().

double AbstractCardiacCellInterface::GetStimulus ( double  time  ) 

Get the value of the intracellular stimulus. Shorthand for GetIntracellularStimulus.

Parameters:
time the time at which to evaluate the stimulus

Definition at line 67 of file AbstractCardiacCellInterface.cpp.

References GetIntracellularStimulus().

Referenced by CML_noble_varghese_kohl_noble_1998_basic_with_sac::EvaluateYDerivatives(), and FitzHughNagumo1961OdeSystem::EvaluateYDerivatives().

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.

Parameters:
pStimulus new stimulus function

Definition at line 73 of file AbstractCardiacCellInterface.cpp.

References mpIntracellularStimulus.

Referenced by HeartConfigRelatedCellFactory< SPACE_DIM >::SetCellIntracellularStimulus(), and SetStimulusFunction().

double AbstractCardiacCellInterface::GetIntracellularStimulus ( double  time  ) 

Get the value of the intracellular stimulus. This will have units of uA/cm^2 for single-cell problems, or uA/cm^3 in a tissue simulation.

Parameters:
time the time at which to evaluate the stimulus

Definition at line 79 of file AbstractCardiacCellInterface.cpp.

References mpIntracellularStimulus.

Referenced by GetIntracellularAreaStimulus(), and GetStimulus().

double AbstractCardiacCellInterface::GetIntracellularAreaStimulus ( double  time  ) 

Get the value of the intracellular stimulus. This will always be in units of uA/cm^2.

Parameters:
time the time at which to evaluate the stimulus

Definition at line 85 of file AbstractCardiacCellInterface.cpp.

References GetIntracellularStimulus(), HeartConfig::GetSurfaceAreaToVolumeRatio(), HeartConfig::Instance(), and mIsUsedInTissue.

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.

Parameters:
tissue true if cell is in a tissue

Definition at line 100 of file AbstractCardiacCellInterface.cpp.

References mIsUsedInTissue.

void AbstractCardiacCellInterface::UseCellMLDefaultStimulus (  )  [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 105 of file AbstractCardiacCellInterface.cpp.

References EXCEPTION, and mHasDefaultStimulusFromCellML.

bool AbstractCardiacCellInterface::HasCellMLDefaultStimulus (  ) 

Returns:
Whether the cell was generated from a CellML file with stimulus metadata.

Definition at line 111 of file AbstractCardiacCellInterface.cpp.

References mHasDefaultStimulusFromCellML.

virtual void AbstractCardiacCellInterface::VerifyStateVariables (  )  [inline, virtual]

Empty method which can be over-ridden in concrete cell class which should go through the current state vector and go range checking on the values (eg check that concentrations are positive and gating variables are between zero and one). This method is called in the ComputeExceptVoltage method, unless NDEBUG has been defined.

Reimplemented in AbstractCardiacCell.

Definition at line 234 of file AbstractCardiacCellInterface.hpp.

Referenced by AbstractCvodeCell::Solve().

virtual AbstractLookupTableCollection* AbstractCardiacCellInterface::GetLookupTableCollection (  )  [inline, virtual]

If this cell uses lookup tables, returns the singleton object containing the tables for this cell's concrete class. Otherwise, returns NULL.

Must be implemented by subclasses iff they support lookup tables.

Definition at line 245 of file AbstractCardiacCellInterface.hpp.

boost::shared_ptr< AbstractStimulusFunction > AbstractCardiacCellInterface::GetStimulusFunction (  ) 

Returns:
The Intracellular stimulus function pointer

Definition at line 116 of file AbstractCardiacCellInterface.cpp.

References mpIntracellularStimulus.

const boost::shared_ptr< AbstractStimulusFunction > AbstractCardiacCellInterface::GetStimulusFunction (  )  const

For boost archiving use only (that's why the 'consts' are required)

Returns:
The Intracellular stimulus function pointer

Definition at line 122 of file AbstractCardiacCellInterface.cpp.

References mpIntracellularStimulus.

const boost::shared_ptr< AbstractIvpOdeSolver > AbstractCardiacCellInterface::GetSolver (  )  const

For boost archiving use only (that's why the 'consts' are required)

Returns:
pointer to the ODE solver being used

Definition at line 127 of file AbstractCardiacCellInterface.cpp.

References mpOdeSolver.

template<class Archive>
void AbstractCardiacCellInterface::serialize ( Archive &  archive,
const unsigned int  version 
) [inline, private]

This is needed to register the base-derived relationship between this class and AbstractCardiacCell; however, serialization of our members is done by AbstractCardiacCell in order to maintain archive backwards compatibility more easily, since we never archive AbstractCvodeCell s.

Parameters:
archive 
version 

Reimplemented in AbstractBackwardEulerCardiacCell< SIZE >, AbstractCardiacCell, FakeBathCell, and CML_noble_varghese_kohl_noble_1998_basic_with_sac.

Definition at line 312 of file AbstractCardiacCellInterface.hpp.


Friends And Related Function Documentation

friend class boost::serialization::access [friend]


Member Data Documentation

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 278 of file AbstractCardiacCellInterface.hpp.

Referenced by FitzHughNagumo1961OdeSystem::GetIIonic(), AbstractCvodeCell::GetVoltage(), AbstractCardiacCell::GetVoltage(), GetVoltageIndex(), AbstractCvodeCell::SetVoltage(), and AbstractCardiacCell::SetVoltage().

The intracellular stimulus current.

Definition at line 284 of file AbstractCardiacCellInterface.hpp.

Referenced by GetIntracellularStimulus(), GetStimulusFunction(), and SetIntracellularStimulusFunction().

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 291 of file AbstractCardiacCellInterface.hpp.

Referenced by AbstractCardiacCell::ComputeExceptVoltage(), CML_noble_varghese_kohl_noble_1998_basic_with_sac::EvaluateYDerivatives(), FitzHughNagumo1961OdeSystem::EvaluateYDerivatives(), AbstractCardiacCell::serialize(), and AbstractCvodeCell::SetVoltageDerivativeToZero().

Whether this cell exists in a tissue, or is an isolated cell.

Definition at line 294 of file AbstractCardiacCellInterface.hpp.

Referenced by GetIntracellularAreaStimulus(), AbstractCardiacCell::serialize(), and SetUsedInTissueSimulation().

Whether this cell has a default stimulus specified by CellML metadata

Definition at line 297 of file AbstractCardiacCellInterface.hpp.

Referenced by HasCellMLDefaultStimulus(), AbstractCardiacCell::serialize(), and UseCellMLDefaultStimulus().


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

Generated on Mon Apr 18 11:35:47 2011 for Chaste by  doxygen 1.5.5