LifeV
|
#include <ElectroIonicModel.hpp>
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< Real > | M_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< VectorEpetra > | vectorPtr_Type |
typedef std::shared_ptr< VectorElemental > | elvecPtr_Type |
typedef RegionMesh< LinearTetra > | mesh_Type |
typedef MatrixEpetra< Real > | matrix_Type |
typedef std::shared_ptr< matrix_Type > | matrixPtr_Type |
typedef FESpace< mesh_Type, MapEpetra > | feSpace_Type |
typedef std::shared_ptr< feSpace_Type > | feSpacePtr_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 ¶meterList) |
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< Real > | restingConditions () 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 | |
ElectroIonicModel & | operator= (const ElectroIonicModel &Ionic) |
Assignment operator. More... | |
Definition at line 96 of file ElectroIonicModel.hpp.
typedef VectorEpetra vector_Type |
Definition at line 103 of file ElectroIonicModel.hpp.
typedef std::shared_ptr<VectorEpetra> vectorPtr_Type |
Definition at line 105 of file ElectroIonicModel.hpp.
typedef std::shared_ptr<VectorElemental> elvecPtr_Type |
Definition at line 107 of file ElectroIonicModel.hpp.
typedef RegionMesh<LinearTetra> mesh_Type |
Definition at line 109 of file ElectroIonicModel.hpp.
typedef MatrixEpetra<Real> matrix_Type |
Definition at line 111 of file ElectroIonicModel.hpp.
typedef std::shared_ptr<matrix_Type> matrixPtr_Type |
Definition at line 113 of file ElectroIonicModel.hpp.
typedef FESpace<mesh_Type, MapEpetra> feSpace_Type |
Definition at line 115 of file ElectroIonicModel.hpp.
typedef std::shared_ptr<feSpace_Type> feSpacePtr_Type |
Definition at line 117 of file ElectroIonicModel.hpp.
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.
typedef FactorySingleton<Factory<ElectroIonicModel, std::string> > IonicModelFactory |
Definition at line 125 of file ElectroIonicModel.hpp.
ElectroIonicModel | ( | int | n | ) |
Constructor.
n | number of equations in the ionic model |
Definition at line 58 of file ElectroIonicModel.cpp.
ElectroIonicModel | ( | int | n, |
int | g | ||
) |
Constructor.
If the number of gating variables is unknown use the method ElectroIonicModel ( int n );
n | number of equations in the ionic model |
g | number of gating variables in the ionic model |
Definition at line 69 of file ElectroIonicModel.cpp.
ElectroIonicModel | ( | const ElectroIonicModel & | Ionic | ) |
Copy Constructor.
Ionic | an ionic model |
Definition at line 80 of file ElectroIonicModel.cpp.
|
inlinevirtual |
Destructor.
Definition at line 160 of file ElectroIonicModel.hpp.
|
inlinevirtual |
setup the parameters of the ionic model from an xml file
parameterList | Teuchos parameter list |
Reimplemented in IonicFitzHughNagumo, and IonicAlievPanfilov.
Definition at line 171 of file ElectroIonicModel.hpp.
|
inline |
returns the number of equations of the ionic model
Definition at line 177 of file ElectroIonicModel.hpp.
|
inline |
returns the number of gating variables in the ionic model
Definition at line 186 of file ElectroIonicModel.hpp.
|
inline |
returns the value of the membrane capacitance in the model
Definition at line 195 of file ElectroIonicModel.hpp.
|
inline |
returns the value of the applied current in the model/point
Definition at line 204 of file ElectroIonicModel.hpp.
|
inline |
returns the pointer to the applied current FE vector in the 3D case
Definition at line 213 of file ElectroIonicModel.hpp.
|
inline |
returns the vector with the resting values of the variables in the ionic model
Definition at line 222 of file ElectroIonicModel.hpp.
|
inline |
returns the function describing the pacing protocol for the ionic model
Definition at line 231 of file ElectroIonicModel.hpp.
|
inline |
set the membrane capacitance in the ionic model
p | membrane capacitance |
Definition at line 240 of file ElectroIonicModel.hpp.
|
inline |
set the applied current in the ionic model/point
p | applied current magnitude |
Definition at line 249 of file ElectroIonicModel.hpp.
|
inline |
set the pointer to the applied current in the 3D ionic model
p | applied current magnitude |
Definition at line 258 of file ElectroIonicModel.hpp.
|
inline |
set the pointer to the applied current in the 3D ionic model
p | applied current magnitude |
Definition at line 267 of file ElectroIonicModel.hpp.
|
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
f | boost function f(t,x,y,x,ID) |
feSpacePtr | pointer to the finite element space on which we interpolate |
time | time at which we evaluate the boost function f |
Definition at line 281 of file ElectroIonicModel.hpp.
|
inline |
Interpolate the function of the electro stimulus.
This method is a wrapper of the interpolation method from the FESpace class
stimulus | pacing protocol defined as ElectroStimulus boost function f(t,x,y,x,ID) |
feSpacePtr | pointer to the finite element space on which we interpolate |
time | time at which we evaluate the boost function f |
Definition at line 299 of file ElectroIonicModel.hpp.
|
inline |
Set the pacing protocol as boost function.
pacingProtocol | boost function defining the pacing protocol |
Definition at line 314 of file ElectroIonicModel.hpp.
|
inline |
Set component of the resting conditions.
value | value of the resting condition |
j | number of the component to be changed |
Definition at line 324 of file ElectroIonicModel.hpp.
|
inline |
Set resting conditions.
restingCondition | values of the resting condition |
Definition at line 333 of file ElectroIonicModel.hpp.
|
inline |
Simple wrapper to add the applied current.
rhs | right hand side of the voltage equation |
Definition at line 342 of file ElectroIonicModel.hpp.
|
inline |
Simple wrapper to add the applied current.
rhs | right hand side of the ionic model (in a point) |
Definition at line 351 of file ElectroIonicModel.hpp.
|
virtual |
This methods computes the Jacobian numerically.
v | vector of pointers to the state variables vectors |
h | differentiation step |
Reimplemented in IonicFitzHughNagumo.
Definition at line 142 of file ElectroIonicModel.cpp.
|
virtual |
This methods computes the Jacobian numerically.
v | vector of the state variables |
h | differentiation step |
Reimplemented in IonicFitzHughNagumo.
Definition at line 113 of file ElectroIonicModel.cpp.
|
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.
v | vector of pointers to the state variables vectors |
rhs | vector of pointers to the right hand side vectors of each variable |
Definition at line 188 of file ElectroIonicModel.cpp.
|
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.
v | vector of pointers to the state variables vectors |
rhs | vector of pointers to the right hand side vectors of each variable |
Definition at line 230 of file ElectroIonicModel.cpp.
|
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
v | vector of pointers to the state variables vectors |
dt | time step |
Definition at line 471 of file ElectroIonicModel.cpp.
|
virtual |
Compute the right hand side of the ionic model in 3D.
This method wraps the 0D model to be used in 3D
v | vector of pointers to the state variables vectors |
rhs | vector of right hand side state variables |
Definition at line 274 of file ElectroIonicModel.cpp.
|
virtual |
Compute the right hand side of the voltage equation linearly interpolating the ionic currents.
v | vector of pointers to the state variables vectors |
rhs | vector of right hand side state variables |
massMatrix | mass matrix of the monodomain system (may be lumped) |
Definition at line 317 of file ElectroIonicModel.cpp.
|
virtual |
Compute the right hand side of the voltage equation using SVI.
v | vector of pointers to the state variables vectors |
rhs | vector of right hand side state variables |
uFESpace | finite 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.
|
virtual |
Compute the right hand side of the voltage equation using SVI specifying the quadrature rule.
v | vector of pointers to the state variables vectors |
rhs | vector of right hand side state variables |
uFESpace | finite element space of the voltage |
qr | quadrature rule for the integration |
Definition at line 461 of file ElectroIonicModel.cpp.
|
virtual |
Initialize the ionic model with a given vector of state variable (0D version)
v | vector of state variables initial conditions |
Definition at line 501 of file ElectroIonicModel.cpp.
|
virtual |
Initialize the ionic model in 3D with a given vector of state variable vector pointers.
v | vector of state variables initial conditions |
Definition at line 510 of file ElectroIonicModel.cpp.
ElectroIonicModel & operator= | ( | const ElectroIonicModel & | Ionic | ) |
Assignment operator.
Methods.
Ionic | ionic model |
Definition at line 97 of file ElectroIonicModel.cpp.
|
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)
v | vector of state variables including the voltage (with n elements) |
rhs | vector 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.
This methods contains the actual evaluation of the rhs of the voltage equation only (0D version)
v | vector of state variables including the voltage (with n elements) |
Implemented in IonicTenTusscher06, IonicFox, IonicLuoRudyI, IonicMinimalModel, IonicNoblePurkinje, IonicGoldbeter, IonicHodgkinHuxley, IonicFitzHughNagumo, IonicMitchellSchaeffer, and IonicAlievPanfilov.
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.
v | vector of state variables including the voltage (with n elements) |
rhs | vector 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.
|
pure virtual |
This methods shows the parameters of the ionic model.
Implemented in IonicTenTusscher06, IonicFox, IonicLuoRudyI, IonicMinimalModel, IonicNoblePurkinje, IonicGoldbeter, IonicHodgkinHuxley, IonicFitzHughNagumo, IonicAlievPanfilov, and IonicMitchellSchaeffer.
|
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
v | vector of state variables |
dt | time step of the simulation |
Reimplemented in IonicTenTusscher06, IonicFox, IonicLuoRudyI, IonicMinimalModel, IonicNoblePurkinje, and IonicHodgkinHuxley.
Definition at line 517 of file ElectroIonicModel.hpp.
|
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
v | vector of state variables including the voltage (with n elements) |
rhs | vector 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.
|
protected |
Definition at line 534 of file ElectroIonicModel.hpp.
|
protected |
Definition at line 544 of file ElectroIonicModel.hpp.
|
protected |
Definition at line 547 of file ElectroIonicModel.hpp.
|
protected |
Definition at line 550 of file ElectroIonicModel.hpp.
|
protected |
Definition at line 553 of file ElectroIonicModel.hpp.
|
protected |
Definition at line 556 of file ElectroIonicModel.hpp.
|
protected |
Definition at line 559 of file ElectroIonicModel.hpp.