LifeV
ElectroIonicModel Class Referenceabstract

#include <ElectroIonicModel.hpp>

+ Inheritance diagram for ElectroIonicModel:
+ Collaboration diagram for ElectroIonicModel:

Public Member Functions

virtual matrix_Type getJac (const vector_Type &v, Real h=1.0e-8)
 This methods computes the Jacobian numerically. More...
 
virtual std::vector< std::vector< Real > > getJac (const std::vector< Real > &v, Real h=1.0e-8)
 This methods computes the Jacobian numerically. More...
 
virtual void computeGatingRhs (const std::vector< vectorPtr_Type > &v, std::vector< vectorPtr_Type > &rhs)
 This methods computes the right hand side of the gating variables in the 3D case. More...
 
virtual void computeNonGatingRhs (const std::vector< vectorPtr_Type > &v, std::vector< vectorPtr_Type > &rhs)
 This methods computes the right hand side of the state variables that are not gating variables in the 3D case. More...
 
virtual void computeGatingVariablesWithRushLarsen (std::vector< vectorPtr_Type > &v, const Real dt)
 Compute the new value of the gating variables in 3D with the Rush Larsen method specified in the 0D version of the ionic model. More...
 
virtual void computeRhs (const std::vector< vectorPtr_Type > &v, std::vector< vectorPtr_Type > &rhs)
 Compute the right hand side of the ionic model in 3D. More...
 
virtual void computePotentialRhsICI (const std::vector< vectorPtr_Type > &v, std::vector< vectorPtr_Type > &rhs, matrix_Type &massMatrix)
 Compute the right hand side of the voltage equation linearly interpolating the ionic currents. More...
 
virtual void computePotentialRhsSVI (const std::vector< vectorPtr_Type > &v, std::vector< vectorPtr_Type > &rhs, FESpace< mesh_Type, MapEpetra > &uFESpace)
 Compute the right hand side of the voltage equation using SVI. More...
 
virtual void computePotentialRhsSVI (const std::vector< vectorPtr_Type > &v, std::vector< vectorPtr_Type > &rhs, FESpace< mesh_Type, MapEpetra > &uFESpace, const QuadratureRule &qr)
 Compute the right hand side of the voltage equation using SVI specifying the quadrature rule. More...
 
virtual void initialize (std::vector< Real > &v)
 Initialize the ionic model with a given vector of state variable (0D version) More...
 
virtual void initialize (std::vector< vectorPtr_Type > &v)
 Initialize the ionic model in 3D with a given vector of state variable vector pointers. More...
 
virtual void computeGatingRhs (const std::vector< Real > &v, std::vector< Real > &rhs)=0
 A new ionic model must have the following methods ///. More...
 
virtual Real computeLocalPotentialRhs (const std::vector< Real > &v)=0
 This methods contains the actual evaluation of the rhs of the voltage equation only (0D version) More...
 
virtual void computeRhs (const std::vector< Real > &v, std::vector< Real > &rhs)=0
 This methods contains the actual evaluation of the rhs of all state variablesin the model (0D version) More...
 
virtual void showMe ()=0
 This methods shows the parameters of the ionic model. More...
 
virtual void computeGatingVariablesWithRushLarsen (std::vector< Real > &, const Real)
 This methods contains the actual evaluation of the rhs of the voltage equation only (0D version) More...
 
virtual void computeNonGatingRhs (const std::vector< Real > &v, std::vector< Real > &rhs)
 In the case this method is improperly used, it should use this default implementation. More...
 

Protected Attributes

short int M_numberOfEquations
 
short int M_numberOfGatingVariables
 
std::vector< RealM_restingConditions
 
Real M_membraneCapacitance
 
Real M_appliedCurrent
 
vectorPtr_Type M_appliedCurrentPtr
 
function_Type M_pacingProtocol
 

Type definitions

typedef VectorEpetra vector_Type
 
typedef std::shared_ptr< VectorEpetravectorPtr_Type
 
typedef std::shared_ptr< VectorElementalelvecPtr_Type
 
typedef RegionMesh< LinearTetramesh_Type
 
typedef MatrixEpetra< Realmatrix_Type
 
typedef std::shared_ptr< matrix_TypematrixPtr_Type
 
typedef FESpace< mesh_Type, MapEpetrafeSpace_Type
 
typedef std::shared_ptr< feSpace_TypefeSpacePtr_Type
 
typedef std::function< Real(const Real &t, const Real &x, const Real &y, const Real &z, const ID &i) > function_Type
 
typedef FactorySingleton< Factory< ElectroIonicModel, std::string > > IonicModelFactory
 

Constructors & Destructor

 ElectroIonicModel ()
 Empty Constructor. More...
 
 ElectroIonicModel (int n)
 Constructor. More...
 
 ElectroIonicModel (int n, int g)
 Constructor. More...
 
 ElectroIonicModel (const ElectroIonicModel &Ionic)
 Copy Constructor. More...
 
virtual ~ElectroIonicModel ()
 Destructor. More...
 

Methods

virtual void setup (Teuchos::ParameterList &parameterList)
 setup the parameters of the ionic model from an xml file More...
 
const short int Size () const
 returns the number of equations of the ionic model More...
 
const short int numberOfGatingVariables () const
 returns the number of gating variables in the ionic model More...
 
const Real membraneCapacitance () const
 returns the value of the membrane capacitance in the model More...
 
const Real appliedCurrent () const
 returns the value of the applied current in the model/point More...
 
vectorPtr_Type appliedCurrentPtr ()
 returns the pointer to the applied current FE vector in the 3D case More...
 
const std::vector< RealrestingConditions () const
 returns the vector with the resting values of the variables in the ionic model More...
 
const function_Type pacaingProtocol () const
 returns the function describing the pacing protocol for the ionic model More...
 
void setMembraneCapacitance (const Real p)
 set the membrane capacitance in the ionic model More...
 
void setAppliedCurrent (const Real p)
 set the applied current in the ionic model/point More...
 
void setAppliedCurrentPtr (const vectorPtr_Type p)
 set the pointer to the applied current in the 3D ionic model More...
 
void setAppliedCurrent (const vector_Type &p)
 set the pointer to the applied current in the 3D ionic model More...
 
void setAppliedCurrentFromFunction (function_Type &f, feSpacePtr_Type feSpacePtr, Real time=0.0)
 Interpolate the function f on the FE space feSpacePtr at time time. More...
 
void setAppliedCurrentFromElectroStimulus (ElectroStimulus &stimulus, feSpacePtr_Type feSpacePtr, Real time=0.0)
 Interpolate the function of the electro stimulus. More...
 
void setPacingProtocol (function_Type pacingProtocol)
 Set the pacing protocol as boost function. More...
 
void setRestingCondtions (Real value, int j)
 Set component of the resting conditions. More...
 
void setRestingCondtions (std::vector< Real > &restingConditions)
 Set resting conditions. More...
 
void addAppliedCurrent (Real &rhs)
 Simple wrapper to add the applied current. More...
 
void addAppliedCurrent (std::vector< Real > &rhs)
 Simple wrapper to add the applied current. More...
 

Overloads

ElectroIonicModeloperator= (const ElectroIonicModel &Ionic)
 Assignment operator. More...
 

Detailed Description

Definition at line 96 of file ElectroIonicModel.hpp.

Member Typedef Documentation

◆ vector_Type

Definition at line 103 of file ElectroIonicModel.hpp.

◆ vectorPtr_Type

typedef std::shared_ptr<VectorEpetra> vectorPtr_Type

Definition at line 105 of file ElectroIonicModel.hpp.

◆ elvecPtr_Type

typedef std::shared_ptr<VectorElemental> elvecPtr_Type

Definition at line 107 of file ElectroIonicModel.hpp.

◆ mesh_Type

Definition at line 109 of file ElectroIonicModel.hpp.

◆ matrix_Type

Definition at line 111 of file ElectroIonicModel.hpp.

◆ matrixPtr_Type

typedef std::shared_ptr<matrix_Type> matrixPtr_Type

Definition at line 113 of file ElectroIonicModel.hpp.

◆ feSpace_Type

Definition at line 115 of file ElectroIonicModel.hpp.

◆ feSpacePtr_Type

typedef std::shared_ptr<feSpace_Type> feSpacePtr_Type

Definition at line 117 of file ElectroIonicModel.hpp.

◆ function_Type

typedef std::function< Real (const Real& t, const Real& x, const Real& y, const Real& z, const ID& i) > function_Type

Definition at line 123 of file ElectroIonicModel.hpp.

◆ IonicModelFactory

Definition at line 125 of file ElectroIonicModel.hpp.

Constructor & Destructor Documentation

◆ ElectroIonicModel() [1/4]

Empty Constructor.

Constructors.

Definition at line 47 of file ElectroIonicModel.cpp.

◆ ElectroIonicModel() [2/4]

ElectroIonicModel ( int  n)

Constructor.

Parameters
nnumber of equations in the ionic model

Definition at line 58 of file ElectroIonicModel.cpp.

◆ ElectroIonicModel() [3/4]

ElectroIonicModel ( int  n,
int  g 
)

Constructor.

If the number of gating variables is unknown use the method ElectroIonicModel ( int n );

Parameters
nnumber of equations in the ionic model
gnumber of gating variables in the ionic model

Definition at line 69 of file ElectroIonicModel.cpp.

◆ ElectroIonicModel() [4/4]

Copy Constructor.

Parameters
Ionican ionic model

Definition at line 80 of file ElectroIonicModel.cpp.

◆ ~ElectroIonicModel()

virtual ~ElectroIonicModel ( )
inlinevirtual

Destructor.

Definition at line 160 of file ElectroIonicModel.hpp.

Member Function Documentation

◆ setup()

virtual void setup ( Teuchos::ParameterList &  parameterList)
inlinevirtual

setup the parameters of the ionic model from an xml file

Parameters
parameterListTeuchos parameter list

Reimplemented in IonicFitzHughNagumo, and IonicAlievPanfilov.

Definition at line 171 of file ElectroIonicModel.hpp.

◆ Size()

const short int Size ( ) const
inline

returns the number of equations of the ionic model

Parameters

Definition at line 177 of file ElectroIonicModel.hpp.

◆ numberOfGatingVariables()

const short int numberOfGatingVariables ( ) const
inline

returns the number of gating variables in the ionic model

Parameters

Definition at line 186 of file ElectroIonicModel.hpp.

◆ membraneCapacitance()

const Real membraneCapacitance ( ) const
inline

returns the value of the membrane capacitance in the model

Parameters

Definition at line 195 of file ElectroIonicModel.hpp.

◆ appliedCurrent()

const Real appliedCurrent ( ) const
inline

returns the value of the applied current in the model/point

Parameters

Definition at line 204 of file ElectroIonicModel.hpp.

◆ appliedCurrentPtr()

vectorPtr_Type appliedCurrentPtr ( )
inline

returns the pointer to the applied current FE vector in the 3D case

Parameters

Definition at line 213 of file ElectroIonicModel.hpp.

◆ restingConditions()

const std::vector<Real> restingConditions ( ) const
inline

returns the vector with the resting values of the variables in the ionic model

Parameters

Definition at line 222 of file ElectroIonicModel.hpp.

◆ pacaingProtocol()

const function_Type pacaingProtocol ( ) const
inline

returns the function describing the pacing protocol for the ionic model

Parameters

Definition at line 231 of file ElectroIonicModel.hpp.

◆ setMembraneCapacitance()

void setMembraneCapacitance ( const Real  p)
inline

set the membrane capacitance in the ionic model

Parameters
pmembrane capacitance

Definition at line 240 of file ElectroIonicModel.hpp.

◆ setAppliedCurrent() [1/2]

void setAppliedCurrent ( const Real  p)
inline

set the applied current in the ionic model/point

Parameters
papplied current magnitude

Definition at line 249 of file ElectroIonicModel.hpp.

◆ setAppliedCurrentPtr()

void setAppliedCurrentPtr ( const vectorPtr_Type  p)
inline

set the pointer to the applied current in the 3D ionic model

Parameters
papplied current magnitude

Definition at line 258 of file ElectroIonicModel.hpp.

◆ setAppliedCurrent() [2/2]

void setAppliedCurrent ( const vector_Type p)
inline

set the pointer to the applied current in the 3D ionic model

Parameters
papplied current magnitude

Definition at line 267 of file ElectroIonicModel.hpp.

◆ setAppliedCurrentFromFunction()

void setAppliedCurrentFromFunction ( function_Type f,
feSpacePtr_Type  feSpacePtr,
Real  time = 0.0 
)
inline

Interpolate the function f on the FE space feSpacePtr at time time.

This method is a wrapper of the interpolation method from the FESpace class

Parameters
fboost function f(t,x,y,x,ID)
feSpacePtrpointer to the finite element space on which we interpolate
timetime at which we evaluate the boost function f

Definition at line 281 of file ElectroIonicModel.hpp.

◆ setAppliedCurrentFromElectroStimulus()

void setAppliedCurrentFromElectroStimulus ( ElectroStimulus stimulus,
feSpacePtr_Type  feSpacePtr,
Real  time = 0.0 
)
inline

Interpolate the function of the electro stimulus.

This method is a wrapper of the interpolation method from the FESpace class

Parameters
stimuluspacing protocol defined as ElectroStimulus boost function f(t,x,y,x,ID)
feSpacePtrpointer to the finite element space on which we interpolate
timetime at which we evaluate the boost function f

Definition at line 299 of file ElectroIonicModel.hpp.

◆ setPacingProtocol()

void setPacingProtocol ( function_Type  pacingProtocol)
inline

Set the pacing protocol as boost function.

Parameters
pacingProtocolboost function defining the pacing protocol

Definition at line 314 of file ElectroIonicModel.hpp.

◆ setRestingCondtions() [1/2]

void setRestingCondtions ( Real  value,
int  j 
)
inline

Set component of the resting conditions.

Parameters
valuevalue of the resting condition
jnumber of the component to be changed

Definition at line 324 of file ElectroIonicModel.hpp.

◆ setRestingCondtions() [2/2]

void setRestingCondtions ( std::vector< Real > &  restingConditions)
inline

Set resting conditions.

Parameters
restingConditionvalues of the resting condition

Definition at line 333 of file ElectroIonicModel.hpp.

◆ addAppliedCurrent() [1/2]

void addAppliedCurrent ( Real rhs)
inline

Simple wrapper to add the applied current.

Parameters
rhsright hand side of the voltage equation

Definition at line 342 of file ElectroIonicModel.hpp.

◆ addAppliedCurrent() [2/2]

void addAppliedCurrent ( std::vector< Real > &  rhs)
inline

Simple wrapper to add the applied current.

Parameters
rhsright hand side of the ionic model (in a point)

Definition at line 351 of file ElectroIonicModel.hpp.

◆ getJac() [1/2]

MatrixEpetra< Real > getJac ( const vector_Type v,
Real  h = 1.0e-8 
)
virtual

This methods computes the Jacobian numerically.

Parameters
vvector of pointers to the state variables vectors
hdifferentiation step

Reimplemented in IonicFitzHughNagumo.

Definition at line 142 of file ElectroIonicModel.cpp.

◆ getJac() [2/2]

std::vector< std::vector< Real > > getJac ( const std::vector< Real > &  v,
Real  h = 1.0e-8 
)
virtual

This methods computes the Jacobian numerically.

Parameters
vvector of the state variables
hdifferentiation step

Reimplemented in IonicFitzHughNagumo.

Definition at line 113 of file ElectroIonicModel.cpp.

+ Here is the caller graph for this function:

◆ computeGatingRhs() [1/2]

void computeGatingRhs ( const std::vector< vectorPtr_Type > &  v,
std::vector< vectorPtr_Type > &  rhs 
)
virtual

This methods computes the right hand side of the gating variables in the 3D case.

In 3D the gating variables are still treated nodewise (that is, as a local 0D system) It computes the right hand side of all state variables except for the voltage equation.

Parameters
vvector of pointers to the state variables vectors
rhsvector of pointers to the right hand side vectors of each variable

Definition at line 188 of file ElectroIonicModel.cpp.

◆ computeNonGatingRhs() [1/2]

void computeNonGatingRhs ( const std::vector< vectorPtr_Type > &  v,
std::vector< vectorPtr_Type > &  rhs 
)
virtual

This methods computes the right hand side of the state variables that are not gating variables in the 3D case.

This method is to be used together with Rush-Larsen time advancing scheme. In 3D the non gating variables are still treated nodewise (that is, as a local 0D system). It computes the right hand side of all non gating state variables except for the voltage equation.

Parameters
vvector of pointers to the state variables vectors
rhsvector of pointers to the right hand side vectors of each variable

Definition at line 230 of file ElectroIonicModel.cpp.

◆ computeGatingVariablesWithRushLarsen() [1/2]

void computeGatingVariablesWithRushLarsen ( std::vector< vectorPtr_Type > &  v,
const Real  dt 
)
virtual

Compute the new value of the gating variables in 3D with the Rush Larsen method specified in the 0D version of the ionic model.

This method wraps the 0D model to be used in 3D

Parameters
vvector of pointers to the state variables vectors
dttime step

Definition at line 471 of file ElectroIonicModel.cpp.

◆ computeRhs() [1/2]

void computeRhs ( const std::vector< vectorPtr_Type > &  v,
std::vector< vectorPtr_Type > &  rhs 
)
virtual

Compute the right hand side of the ionic model in 3D.

This method wraps the 0D model to be used in 3D

Parameters
vvector of pointers to the state variables vectors
rhsvector of right hand side state variables

Definition at line 274 of file ElectroIonicModel.cpp.

◆ computePotentialRhsICI()

void computePotentialRhsICI ( const std::vector< vectorPtr_Type > &  v,
std::vector< vectorPtr_Type > &  rhs,
matrix_Type massMatrix 
)
virtual

Compute the right hand side of the voltage equation linearly interpolating the ionic currents.

Parameters
vvector of pointers to the state variables vectors
rhsvector of right hand side state variables
massMatrixmass matrix of the monodomain system (may be lumped)

Definition at line 317 of file ElectroIonicModel.cpp.

◆ computePotentialRhsSVI() [1/2]

void computePotentialRhsSVI ( const std::vector< vectorPtr_Type > &  v,
std::vector< vectorPtr_Type > &  rhs,
FESpace< mesh_Type, MapEpetra > &  uFESpace 
)
virtual

Compute the right hand side of the voltage equation using SVI.

Parameters
vvector of pointers to the state variables vectors
rhsvector of right hand side state variables
uFESpacefinite element space of the voltage

Filling local elvec_u with potential values in the nodes

Reimplemented in IonicMinimalModel.

Definition at line 356 of file ElectroIonicModel.cpp.

◆ computePotentialRhsSVI() [2/2]

void computePotentialRhsSVI ( const std::vector< vectorPtr_Type > &  v,
std::vector< vectorPtr_Type > &  rhs,
FESpace< mesh_Type, MapEpetra > &  uFESpace,
const QuadratureRule qr 
)
virtual

Compute the right hand side of the voltage equation using SVI specifying the quadrature rule.

Parameters
vvector of pointers to the state variables vectors
rhsvector of right hand side state variables
uFESpacefinite element space of the voltage
qrquadrature rule for the integration

Definition at line 461 of file ElectroIonicModel.cpp.

◆ initialize() [1/2]

void initialize ( std::vector< Real > &  v)
virtual

Initialize the ionic model with a given vector of state variable (0D version)

Parameters
vvector of state variables initial conditions

Definition at line 501 of file ElectroIonicModel.cpp.

◆ initialize() [2/2]

void initialize ( std::vector< vectorPtr_Type > &  v)
virtual

Initialize the ionic model in 3D with a given vector of state variable vector pointers.

Parameters
vvector of state variables initial conditions

Definition at line 510 of file ElectroIonicModel.cpp.

◆ operator=()

ElectroIonicModel & operator= ( const ElectroIonicModel Ionic)

Assignment operator.

Methods.

Parameters
Ionicionic model

Definition at line 97 of file ElectroIonicModel.cpp.

◆ computeGatingRhs() [2/2]

virtual void computeGatingRhs ( const std::vector< Real > &  v,
std::vector< Real > &  rhs 
)
pure virtual

A new ionic model must have the following methods ///.

This methods contains the actual evaluation of the rhs of all state variables except for the voltage equation (0D version)

Parameters
vvector of state variables including the voltage (with n elements)
rhsvector of right hand side state variables excluding the voltage (with n-1 elements)

Implemented in IonicTenTusscher06, IonicFox, IonicLuoRudyI, IonicMinimalModel, IonicNoblePurkinje, IonicGoldbeter, IonicFitzHughNagumo, IonicHodgkinHuxley, IonicMitchellSchaeffer, and IonicAlievPanfilov.

◆ computeLocalPotentialRhs()

virtual Real computeLocalPotentialRhs ( const std::vector< Real > &  v)
pure virtual

This methods contains the actual evaluation of the rhs of the voltage equation only (0D version)

Parameters
vvector of state variables including the voltage (with n elements)

Implemented in IonicTenTusscher06, IonicFox, IonicLuoRudyI, IonicMinimalModel, IonicNoblePurkinje, IonicGoldbeter, IonicHodgkinHuxley, IonicFitzHughNagumo, IonicMitchellSchaeffer, and IonicAlievPanfilov.

◆ computeRhs() [2/2]

virtual void computeRhs ( const std::vector< Real > &  v,
std::vector< Real > &  rhs 
)
pure virtual

This methods contains the actual evaluation of the rhs of all state variablesin the model (0D version)

Although this method can just wrap the computeGatingRhs ( const std::vector<Real>& v, std::vector<Real>& rhs ) and the computeLocalPotentialRhs ( const std::vector<Real>& v ) methods, for efficiency it may be better to duplicate the code.

Parameters
vvector of state variables including the voltage (with n elements)
rhsvector of right hand side state variables including the voltage (with n elements)

Implemented in IonicTenTusscher06, IonicFox, IonicLuoRudyI, IonicMinimalModel, IonicNoblePurkinje, IonicGoldbeter, IonicHodgkinHuxley, IonicFitzHughNagumo, IonicMitchellSchaeffer, and IonicAlievPanfilov.

◆ showMe()

virtual void showMe ( )
pure virtual

◆ computeGatingVariablesWithRushLarsen() [2/2]

virtual void computeGatingVariablesWithRushLarsen ( std::vector< Real > &  ,
const Real   
)
inlinevirtual

This methods contains the actual evaluation of the rhs of the voltage equation only (0D version)

Overload this method in order to solve the ionic model with the Rush-Larsen method in the monodomain solver

Parameters
vvector of state variables
dttime step of the simulation

Reimplemented in IonicTenTusscher06, IonicFox, IonicLuoRudyI, IonicMinimalModel, IonicNoblePurkinje, and IonicHodgkinHuxley.

Definition at line 517 of file ElectroIonicModel.hpp.

◆ computeNonGatingRhs() [2/2]

virtual void computeNonGatingRhs ( const std::vector< Real > &  v,
std::vector< Real > &  rhs 
)
inlinevirtual

In the case this method is improperly used, it should use this default implementation.

This method should be used together with the Rush-Larsen time advancing scheme

Parameters
vvector of state variables including the voltage (with n elements)
rhsvector of right hand side state variables excluding the voltage and the gating variables (with n-g-1 elements)

Reimplemented in IonicTenTusscher06, IonicFox, and IonicLuoRudyI.

Definition at line 528 of file ElectroIonicModel.hpp.

Field Documentation

◆ M_numberOfEquations

short int M_numberOfEquations
protected

Definition at line 534 of file ElectroIonicModel.hpp.

◆ M_numberOfGatingVariables

short int M_numberOfGatingVariables
protected

Definition at line 544 of file ElectroIonicModel.hpp.

◆ M_restingConditions

std::vector<Real> M_restingConditions
protected

Definition at line 547 of file ElectroIonicModel.hpp.

◆ M_membraneCapacitance

Real M_membraneCapacitance
protected

Definition at line 550 of file ElectroIonicModel.hpp.

◆ M_appliedCurrent

Real M_appliedCurrent
protected

Definition at line 553 of file ElectroIonicModel.hpp.

◆ M_appliedCurrentPtr

vectorPtr_Type M_appliedCurrentPtr
protected

Definition at line 556 of file ElectroIonicModel.hpp.

◆ M_pacingProtocol

function_Type M_pacingProtocol
protected

Definition at line 559 of file ElectroIonicModel.hpp.


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