1. cbcbeat package¶
1.1. Subpackages¶
- 1.1.1. cbcbeat.cellmodels package
- 1.1.1.1. Submodules
- 1.1.1.2. cbcbeat.cellmodels.beeler_reuter_1977 module
- 1.1.1.3. cbcbeat.cellmodels.cardiaccellmodel module
- 1.1.1.4. cbcbeat.cellmodels.fenton_karma_1998_BR_altered module
- 1.1.1.5. cbcbeat.cellmodels.fenton_karma_1998_MLR1_altered module
- 1.1.1.6. cbcbeat.cellmodels.fitzhughnagumo module
- 1.1.1.7. cbcbeat.cellmodels.fitzhughnagumo_manual module
- 1.1.1.8. cbcbeat.cellmodels.grandi_pasqualini_bers_2010 module
- 1.1.1.9. cbcbeat.cellmodels.nocellmodel module
- 1.1.1.10. cbcbeat.cellmodels.rogers_mcculloch_manual module
- 1.1.1.11. cbcbeat.cellmodels.tentusscher_2004_mcell module
- 1.1.1.12. cbcbeat.cellmodels.tentusscher_2004_mcell_cont module
- 1.1.1.13. cbcbeat.cellmodels.tentusscher_2004_mcell_disc module
- 1.1.1.14. cbcbeat.cellmodels.tentusscher_panfilov_2006_M_cell module
- 1.1.1.15. cbcbeat.cellmodels.tentusscher_panfilov_2006_epi_cell module
- 1.1.1.16. Module contents
1.2. Submodules¶
1.3. cbcbeat.bidomainsolver module¶
These solvers solve the (pure) bidomain equations on the form: find the transmembrane potential \(v = v(x, t)\) and the extracellular potential \(u = u(x, t)\) such that
where the subscript \(t\) denotes the time derivative; \(G_x\) denotes a weighted gradient: \(G_x = M_x \mathrm{grad}(v)\) for \(x \in \{i, e\}\), where \(M_i\) and \(M_e\) are the intracellular and extracellular cardiac conductivity tensors, respectively; \(I_s\) and \(I_a\) are prescribed input. In addition, initial conditions are given for \(v\):
Finally, boundary conditions must be prescribed. For now, this solver assumes pure homogeneous Neumann boundary conditions for \(v\) and \(u\) and enforces the additional average value zero constraint for u.
-
class
cbcbeat.bidomainsolver.
BasicBidomainSolver
(mesh, time, M_i, M_e, I_s=None, I_a=None, v_=None, params=None)[source]¶ Bases:
object
This solver is based on a theta-scheme discretization in time and CG_1 x CG_1 (x R) elements in space.
Note
For the sake of simplicity and consistency with other solver objects, this solver operates on its solution fields (as state variables) directly internally. More precisely, solve (and step) calls will act by updating the internal solution fields. It implies that initial conditions can be set (and are intended to be set) by modifying the solution fields prior to simulation.
- Arguments
- mesh (
dolfin.Mesh
) - The spatial domain (mesh)
- time (
dolfin.Constant
or None) - A constant holding the current time. If None is given, time is created for you, initialized to zero.
- M_i (
ufl.Expr
) - The intracellular conductivity tensor (as an UFL expression)
- M_e (
ufl.Expr
) - The extracellular conductivity tensor (as an UFL expression)
- I_s (
dict
, optional) - A typically time-dependent external stimulus given as a dict,
with domain markers as the key and a
dolfin.Expression
as values. NB: it is assumed that the time dependence of I_s is encoded via the ‘time’ Constant. - I_a (
dolfin.Expression
, optional) - A (typically time-dependent) external applied current
- v_ (
ufl.Expr
, optional) - Initial condition for v. A new
dolfin.Function
will be created if none is given. - params (
dolfin.Parameters
, optional) - Solver parameters
- mesh (
-
static
default_parameters
()[source]¶ Initialize and return a set of default parameters
- Returns
- A set of parameters (
dolfin.Parameters
)
To inspect all the default parameters, do:
info(BasicBidomainSolver.default_parameters(), True)
-
solution_fields
()[source]¶ Return tuple of previous and current solution objects.
Modifying these will modify the solution objects of the solver and thus provides a way for setting initial conditions for instance.
- Returns
- (previous v, current vur) (
tuple
ofdolfin.Function
)
-
solve
(interval, dt=None)[source]¶ Solve the discretization on a given time interval (t0, t1) with a given timestep dt and return generator for a tuple of the interval and the current solution.
- Arguments
- interval (
tuple
) - The time interval for the solve given by (t0, t1)
- dt (int, optional)
- The timestep for the solve. Defaults to length of interval
- interval (
- Returns
- (timestep, solution_fields) via (
genexpr
)
Example of usage:
# Create generator solutions = solver.solve((0.0, 1.0), 0.1) # Iterate over generator (computes solutions as you go) for (interval, solution_fields) in solutions: (t0, t1) = interval v_, vur = solution_fields # do something with the solutions
-
step
(interval)[source]¶ Solve on the given time interval (t0, t1).
- Arguments
- interval (
tuple
) - The time interval (t0, t1) for the step
- interval (
- Invariants
- Assuming that v_ is in the correct state for t0, gives self.vur in correct state at t1.
-
time
¶ The internal time of the solver.
-
class
cbcbeat.bidomainsolver.
BidomainSolver
(mesh, time, M_i, M_e, I_s=None, I_a=None, v_=None, params=None)[source]¶ Bases:
cbcbeat.bidomainsolver.BasicBidomainSolver
This solver is based on a theta-scheme discretization in time and CG_1 x CG_1 (x R) elements in space.
Note
For the sake of simplicity and consistency with other solver objects, this solver operates on its solution fields (as state variables) directly internally. More precisely, solve (and step) calls will act by updating the internal solution fields. It implies that initial conditions can be set (and are intended to be set) by modifying the solution fields prior to simulation.
- Arguments
- mesh (
dolfin.Mesh
) - The spatial domain (mesh)
- time (
dolfin.Constant
or None) - A constant holding the current time. If None is given, time is created for you, initialized to zero.
- M_i (
ufl.Expr
) - The intracellular conductivity tensor (as an UFL expression)
- M_e (
ufl.Expr
) - The extracellular conductivity tensor (as an UFL expression)
- I_s (
dict
, optional) - A typically time-dependent external stimulus given as a dict,
with domain markers as the key and a
dolfin.Expression
as values. NB: it is assumed that the time dependence of I_s is encoded via the ‘time’ Constant. - I_a (
dolfin.Expression
, optional) - A (typically time-dependent) external applied current
- v_ (
ufl.Expr
, optional) - Initial condition for v. A new
dolfin.Function
will be created if none is given. - params (
dolfin.Parameters
, optional) - Solver parameters
- mesh (
-
static
default_parameters
()[source]¶ Initialize and return a set of default parameters
- Returns
- A set of parameters (
dolfin.Parameters
)
To inspect all the default parameters, do:
info(BidomainSolver.default_parameters(), True)
-
linear_solver
¶ The linear solver (
dolfin.LUSolver
ordolfin.PETScKrylovSolver
).
-
nullspace
¶
1.4. cbcbeat.cardiacmodels module¶
This module contains a container class for cardiac models:
CardiacModel
. This class
should be instantiated for setting up specific cardiac simulation
scenarios.
-
class
cbcbeat.cardiacmodels.
CardiacModel
(domain, time, M_i, M_e, cell_models, stimulus=None, applied_current=None)[source]¶ Bases:
object
A container class for cardiac models. Objects of this class represent a specific cardiac simulation set-up and should provide
- A computational domain
- A cardiac cell model
- Intra-cellular and extra-cellular conductivities
- Various forms of stimulus (optional).
This container class is designed for use with the splitting solvers (
cbcbeat.splittingsolver
), see their documentation for more information on how the attributes are interpreted in that context.- Arguments
- domain (
dolfin.Mesh
) - the computational domain in space
- time (
dolfin.Constant
or None ) - A constant holding the current time.
- M_i (
ufl.Expr
) - the intra-cellular conductivity as an ufl Expression
- M_e (
ufl.Expr
) - the extra-cellular conductivity as an ufl Expression
- cell_models (
CardiacCellModel
) - a cell model or a dict with cell models associated with a cell model domain
- stimulus (
dict
, optional) - A typically time-dependent external stimulus given as a dict,
with domain markers as the key and a
dolfin.Expression
as values. NB: it is assumed that the time dependence of I_s is encoded via the ‘time’ Constant. - applied_current (
ufl.Expr
, optional) - an applied current as an ufl Expression
- domain (
1.5. cbcbeat.cellsolver module¶
This module contains solvers for (subclasses of) CardiacCellModel.
-
class
cbcbeat.cellsolver.
BasicSingleCellSolver
(model, time, params=None)[source]¶ Bases:
cbcbeat.cellsolver.BasicCardiacODESolver
A basic, non-optimised solver for systems of ODEs typically encountered in cardiac applications of the form: find a scalar field \(v = v(t)\) and a vector field \(s = s(t)\)
\[ \begin{align}\begin{aligned}v_t = - I_{ion}(v, s) + I_s\\s_t = F(v, s)\end{aligned}\end{align} \]where \(I_{ion}\) and \(F\) are given non-linear functions, \(I_s\) is some prescribed stimulus. If \(I_s\) depends on time, it is assumed that \(I_s\) is a
dolfin.Expression
with parameter ‘t’.Use this solver if you just want to test the results from a cardiac cell model without any spatial mesh dependence.
Here, this nonlinear ODE system is solved via a theta-scheme. By default theta=0.5, which corresponds to a Crank-Nicolson scheme. This can be changed by modifying the solver parameters.
Note
For the sake of simplicity and consistency with other solver objects, this solver operates on its solution fields (as state variables) directly internally. More precisely, solve (and step) calls will act by updating the internal solution fields. It implies that initial conditions can be set (and are intended to be set) by modifying the solution fields prior to simulation.
- Arguments
- model (
CardiacCellModel
) - A cardiac cell model
- time (
Constant
or None) - A constant holding the current time.
- params (
dolfin.Parameters
, optional) - Solver parameters
- model (
-
class
cbcbeat.cellsolver.
BasicCardiacODESolver
(mesh, time, model, I_s=None, params=None)[source]¶ Bases:
object
A basic, non-optimised solver for systems of ODEs typically encountered in cardiac applications of the form: find a scalar field \(v = v(x, t)\) and a vector field \(s = s(x, t)\)
\[ \begin{align}\begin{aligned}v_t = - I_{ion}(v, s) + I_s\\s_t = F(v, s)\end{aligned}\end{align} \]where \(I_{ion}\) and \(F\) are given non-linear functions, and \(I_s\) is some prescribed stimulus.
Here, this nonlinear ODE system is solved via a theta-scheme. By default theta=0.5, which corresponds to a Crank-Nicolson scheme. This can be changed by modifying the solver parameters.
Note
For the sake of simplicity and consistency with other solver objects, this solver operates on its solution fields (as state variables) directly internally. More precisely, solve (and step) calls will act by updating the internal solution fields. It implies that initial conditions can be set (and are intended to be set) by modifying the solution fields prior to simulation.
- Arguments
- mesh (
dolfin.Mesh
) - The spatial domain (mesh)
- time (
dolfin.Constant
or None) - A constant holding the current time. If None is given, time is created for you, initialized to zero.
- model (
cbcbeat.CardiacCellModel
) - A representation of the cardiac cell model(s)
- I_s (optional) A typically time-dependent external stimulus
- given as a
dolfin.GenericFunction
or a Markerwise. NB: it is assumed that the time dependence of I_s is encoded via the ‘time’ Constant. - params (
dolfin.Parameters
, optional) - Solver parameters
- mesh (
-
static
default_parameters
()[source]¶ Initialize and return a set of default parameters
- Returns
- A set of parameters (
dolfin.Parameters
)
-
solution_fields
()[source]¶ Return tuple of previous and current solution objects.
Modifying these will modify the solution objects of the solver and thus provides a way for setting initial conditions for instance.
- Returns
- (previous vs, current vs) (
tuple
ofdolfin.Function
)
-
solve
(interval, dt=None)[source]¶ Solve the problem given by the model on a given time interval (t0, t1) with a given timestep dt and return generator for a tuple of the interval and the current vs solution.
- Arguments
- interval (
tuple
) - The time interval for the solve given by (t0, t1)
- dt (int, optional)
- The timestep for the solve. Defaults to length of interval
- interval (
- Returns
- (timestep, current vs) via (
genexpr
)
Example of usage:
# Create generator solutions = solver.solve((0.0, 1.0), 0.1) # Iterate over generator (computes solutions as you go) for (interval, vs) in solutions: # do something with the solutions
-
step
(interval)[source]¶ Solve on the given time step (t0, t1).
End users are recommended to use solve instead.
- Arguments
- interval (
tuple
) - The time interval (t0, t1) for the step
- interval (
-
time
¶ The internal time of the solver.
-
class
cbcbeat.cellsolver.
CardiacODESolver
(mesh, time, model, I_s=None, params=None)[source]¶ Bases:
object
An optimised solver for systems of ODEs typically encountered in cardiac applications of the form: find a scalar field \(v = v(x, t)\) and a vector field \(s = s(x, t)\)
\[ \begin{align}\begin{aligned}v_t = - I_{ion}(v, s) + I_s\\s_t = F(v, s)\end{aligned}\end{align} \]where \(I_{ion}\) and \(F\) are given non-linear functions, and \(I_s\) is some prescribed stimulus.
Note
For the sake of simplicity and consistency with other solver objects, this solver operates on its solution fields (as state variables) directly internally. More precisely, solve (and step) calls will act by updating the internal solution fields. It implies that initial conditions can be set (and are intended to be set) by modifying the solution fields prior to simulation.
- Arguments
- mesh (
dolfin.Mesh
) - The spatial mesh (mesh)
- time (
dolfin.Constant
or None) - A constant holding the current time. If None is given, time is created for you, initialized to zero.
- model (
cbcbeat.CardiacCellModel
) - A representation of the cardiac cell model(s)
- I_s (
dolfin.Expression
, optional) - A typically time-dependent external stimulus. NB: it is assumed that the time dependence of I_s is encoded via the ‘time’ Constant.
- params (
dolfin.Parameters
, optional) - Solver parameters
- mesh (
-
static
default_parameters
()[source]¶ Initialize and return a set of default parameters
- Returns
- A set of parameters (
dolfin.Parameters
)
-
solution_fields
()[source]¶ Return current solution object.
Modifying this will modify the solution object of the solver and thus provides a way for setting initial conditions for instance.
- Returns
- (previous vs_, current vs) (
dolfin.Function
)
-
solve
(interval, dt=None)[source]¶ Solve the problem given by the model on a given time interval (t0, t1) with a given timestep dt and return generator for a tuple of the interval and the current vs solution.
- Arguments
- interval (
tuple
) - The time interval for the solve given by (t0, t1)
- dt (int, optional)
- The timestep for the solve. Defaults to length of interval
- interval (
- Returns
- (timestep, current vs) via (
genexpr
)
Example of usage:
# Create generator solutions = solver.solve((0.0, 1.0), 0.1) # Iterate over generator (computes solutions as you go) for (interval, vs) in solutions: # do something with the solutions
1.6. cbcbeat.dolfinimport module¶
This module handles all dolfin import in cbcbeat. Here dolfin and dolfin_adjoint gets imported. If dolfin_adjoint is not present it will not be imported.
1.7. cbcbeat.gossplittingsolver module¶
This module contains a CellSolver that uses JIT compiled Gotran models together with GOSS (General ODE System Solver), which can be interfaced by the GossSplittingSolver
-
class
cbcbeat.gossplittingsolver.
GOSSplittingSolver
(model, params=None)[source]¶ -
static
default_parameters
()[source]¶ Initialize and return a set of default parameters for the splitting solver
- Returns
- A set of parameters (
dolfin.Parameters
)
To inspect all the default parameters, do:
info(SplittingSolver.default_parameters(), True)
-
merge
(solution)[source]¶ Extract membrane potential from solutions from the PDE solve and put it into membrane potential used for the ODE solve.
- Arguments
- solution (
dolfin.Function
) - Function holding the combined result
- solution (
-
solution_fields
()[source]¶ Return tuple of previous and current solution objects.
Modifying these will modify the solution objects of the solver and thus provides a way for setting initial conditions for instance.
- Returns
- (current v, current vur) (
tuple
ofdolfin.Function
)
-
solve
(interval, dt)[source]¶ Solve the problem given by the model on a given time interval (t0, t1) with a given timestep dt and return generator for a tuple of the time step and the solution fields.
- Arguments
- interval (
tuple
) - The time interval for the solve given by (t0, t1)
- dt (int)
- The timestep for the solve
- interval (
- Returns
- (timestep, solution_fields) via (
genexpr
)
Example of usage:
# Create generator solutions = solver.solve((0.0, 1.0), 0.1) # Iterate over generator (computes solutions as you go) for ((t0, t1), (v, vur)) in solutions: # do something with the solutions
-
step
(interval)[source]¶ Solve the problem given by the model on a given time interval (t0, t1) with timestep given by the interval length.
- Arguments
- interval (
tuple
) - The time interval for the solve given by (t0, t1)
- interval (
- Invariants
- Given self._vs in a correct state at t0, provide v and s (in self.vs) and u (in self.vur) in a correct state at t1. (Note that self.vur[0] == self.vs[0] only if theta = 1.0.)
-
static
1.8. cbcbeat.gotran2cellmodel module¶
1.9. cbcbeat.gotran2dolfin module¶
1.10. cbcbeat.markerwisefield module¶
-
class
cbcbeat.markerwisefield.
Markerwise
(objects, keys, markers)[source]¶ Bases:
object
A container class representing an object defined by a number of objects combined with a mesh function defining mesh domains and a map between the these.
- Arguments
- objects (tuple)
- the different objects
- keys (tuple of ints)
- a map from the objects to the domains marked in markers
- markers (
dolfin.MeshFunction
) - a mesh function mapping which domains the mesh consist of
Example of usage
Given (g0, g1), (2, 5) and markers, let
g = g0 on domains marked by 2 in markers g = g1 on domains marked by 5 in markersletting:
g = Markerwise((g0, g1), (2, 5), markers)
1.11. cbcbeat.monodomainsolver module¶
These solvers solve the (pure) monodomain equations on the form: find the transmembrane potential \(v = v(x, t)\) such that
where the subscript \(t\) denotes the time derivative; \(G_i\) denotes a weighted gradient: \(G_i = M_i \mathrm{grad}(v)\) for, where \(M_i\) is the intracellular cardiac conductivity tensor; \(I_s\) ise prescribed input. In addition, initial conditions are given for \(v\):
Finally, boundary conditions must be prescribed. For now, this solver assumes pure homogeneous Neumann boundary conditions for \(v\).
-
class
cbcbeat.monodomainsolver.
BasicMonodomainSolver
(mesh, time, M_i, I_s=None, v_=None, params=None)[source]¶ Bases:
object
This solver is based on a theta-scheme discretization in time and CG_1 elements in space.
Note
For the sake of simplicity and consistency with other solver objects, this solver operates on its solution fields (as state variables) directly internally. More precisely, solve (and step) calls will act by updating the internal solution fields. It implies that initial conditions can be set (and are intended to be set) by modifying the solution fields prior to simulation.
- Arguments
- mesh (
dolfin.Mesh
) - The spatial domain (mesh)
- time (
dolfin.Constant
or None) - A constant holding the current time. If None is given, time is created for you, initialized to zero.
- M_i (
ufl.Expr
) - The intracellular conductivity tensor (as an UFL expression)
- I_s (
dict
, optional) - A typically time-dependent external stimulus given as a dict,
with domain markers as the key and a
dolfin.Expression
as values. NB: it is assumed that the time dependence of I_s is encoded via the ‘time’ Constant. - v_ (
ufl.Expr
, optional) - Initial condition for v. A new
dolfin.Function
will be created if none is given. - params (
dolfin.Parameters
, optional) - Solver parameters
- mesh (
-
static
default_parameters
()[source]¶ Initialize and return a set of default parameters
- Returns
- A set of parameters (
dolfin.Parameters
)
To inspect all the default parameters, do:
info(BasicMonodomainSolver.default_parameters(), True)
-
solution_fields
()[source]¶ Return tuple of previous and current solution objects.
Modifying these will modify the solution objects of the solver and thus provides a way for setting initial conditions for instance.
- Returns
- (previous v, current v) (
tuple
ofdolfin.Function
)
-
solve
(interval, dt=None)[source]¶ Solve the discretization on a given time interval (t0, t1) with a given timestep dt and return generator for a tuple of the interval and the current solution.
- Arguments
- interval (
tuple
) - The time interval for the solve given by (t0, t1)
- dt (int, optional)
- The timestep for the solve. Defaults to length of interval
- interval (
- Returns
- (timestep, solution_field) via (
genexpr
)
Example of usage:
# Create generator solutions = solver.solve((0.0, 1.0), 0.1) # Iterate over generator (computes solutions as you go) for (interval, solution_fields) in solutions: (t0, t1) = interval v_, v = solution_fields # do something with the solutions
-
step
(interval)[source]¶ Solve on the given time interval (t0, t1).
- Arguments
- interval (
tuple
) - The time interval (t0, t1) for the step
- interval (
- Invariants
- Assuming that v_ is in the correct state for t0, gives self.v in correct state at t1.
-
time
¶ The internal time of the solver.
-
class
cbcbeat.monodomainsolver.
MonodomainSolver
(mesh, time, M_i, I_s=None, v_=None, params=None)[source]¶ Bases:
cbcbeat.monodomainsolver.BasicMonodomainSolver
This solver is based on a theta-scheme discretization in time and CG_1 elements in space.
Note
For the sake of simplicity and consistency with other solver objects, this solver operates on its solution fields (as state variables) directly internally. More precisely, solve (and step) calls will act by updating the internal solution fields. It implies that initial conditions can be set (and are intended to be set) by modifying the solution fields prior to simulation.
- Arguments
- mesh (
dolfin.Mesh
) - The spatial domain (mesh)
- time (
dolfin.Constant
or None) - A constant holding the current time. If None is given, time is created for you, initialized to zero.
- M_i (
ufl.Expr
) - The intracellular conductivity tensor (as an UFL expression)
- I_s (
dict
, optional) - A typically time-dependent external stimulus given as a dict,
with domain markers as the key and a
dolfin.Expression
as values. NB: it is assumed that the time dependence of I_s is encoded via the ‘time’ Constant. - v_ (
ufl.Expr
, optional) - Initial condition for v. A new
dolfin.Function
will be created if none is given. - params (
dolfin.Parameters
, optional) - Solver parameters
- mesh (
-
static
default_parameters
()[source]¶ Initialize and return a set of default parameters
- Returns
- A set of parameters (
dolfin.Parameters
)
To inspect all the default parameters, do:
info(MonodomainSolver.default_parameters(), True)
-
linear_solver
¶ The linear solver (
dolfin.LUSolver
ordolfin.KrylovSolver
).
1.12. cbcbeat.splittingsolver module¶
This module contains splitting solvers for CardiacModel objects. In particular, the classes
- SplittingSolver
- BasicSplittingSolver
These solvers solve the bidomain (or monodomain) equations on the form: find the transmembrane potential \(v = v(x, t)\) in mV, the extracellular potential \(u = u(x, t)\) in mV, and any additional state variables \(s = s(x, t)\) such that
where
- the subscript \(t\) denotes the time derivative,
- \(M_i\) and \(M_e\) are conductivity tensors (in mm^2/ms)
- \(I_s\) is prescribed input current (in mV/ms)
- \(I_a\) is prescribed input current (in mV/ms)
- \(I_{ion}\) and \(F\) are typically specified by a cell model
- Note that M_i and M_e can be viewed as scaled by \(\chi*C_m\) where
- \(\chi\) is the surface-to volume ratio of cells (in 1/mm) ,
- \(C_m\) is the specific membrane capacitance (in mu F/(mm^2) ),
In addition, initial conditions are given for \(v\) and \(s\):
Finally, boundary conditions must be prescribed. These solvers assume pure Neumann boundary conditions for \(v\) and \(u\) and enforce the additional average value zero constraint for u.
The solvers take as input a
CardiacModel
providing the
required input specification of the problem. In particular, the
applied current \(I_a\) is extracted from the
applied_current
attribute, while the stimulus \(I_s\) is extracted from the
stimulus
attribute.
It should be possible to use the solvers interchangably. However, note that the BasicSplittingSolver is not optimised and should be used for testing or debugging purposes primarily.
-
class
cbcbeat.splittingsolver.
SplittingSolver
(model, params=None)[source]¶ Bases:
cbcbeat.splittingsolver.BasicSplittingSolver
An optimised solver for the bidomain equations based on the operator splitting scheme described in Sundnes et al 2006, p. 78 ff.
The solver computes as solutions:
- “vs” (
dolfin.Function
) representing the solution for the transmembrane potential and any additional state variables, and - “vur” (
dolfin.Function
) representing the transmembrane potential in combination with the extracellular potential and an additional Lagrange multiplier.
The algorithm can be controlled by a number of parameters. In particular, the splitting algorithm can be controlled by the parameter “theta”: “theta” set to 1.0 corresponds to a (1st order) Godunov splitting while “theta” set to 0.5 to a (2nd order) Strang splitting.
- Arguments
- model (
cbcbeat.cardiacmodels.CardiacModel
) - a CardiacModel object describing the simulation set-up
- params (
dolfin.Parameters
, optional) - a Parameters object controlling solver parameters
- model (
Example of usage:
from cbcbeat import * # Describe the cardiac model mesh = UnitSquareMesh(100, 100) time = Constant(0.0) cell_model = FitzHughNagumoManual() stimulus = Expression("10*t*x[0]", t=time, degree=1) cm = CardiacModel(mesh, time, 1.0, 1.0, cell_model, stimulus) # Extract default solver parameters ps = SplittingSolver.default_parameters() # Examine the default parameters info(ps, True) # Update parameter dictionary ps["theta"] = 1.0 # Use first order splitting ps["CardiacODESolver"]["scheme"] = "GRL1" # Use Generalized Rush-Larsen scheme ps["pde_solver"] = "monodomain" # Use monodomain equations as the PDE model ps["MonodomainSolver"]["linear_solver_type"] = "direct" # Use direct linear solver of the PDEs ps["MonodomainSolver"]["theta"] = 1.0 # Use backward Euler for temporal discretization for the PDEs solver = SplittingSolver(cm, params=ps) # Extract the solution fields and set the initial conditions (vs_, vs, vur) = solver.solution_fields() vs_.assign(cell_model.initial_conditions()) # Solve dt = 0.1 T = 1.0 interval = (0.0, T) for (timestep, fields) in solver.solve(interval, dt): (vs_, vs, vur) = fields # Do something with the solutions
- Assumptions
- The cardiac conductivities do not vary in time
- “vs” (
-
class
cbcbeat.splittingsolver.
BasicSplittingSolver
(model, params=None)[source]¶ A non-optimised solver for the bidomain equations based on the operator splitting scheme described in Sundnes et al 2006, p. 78 ff.
The solver computes as solutions:
- “vs” (
dolfin.Function
) representing the solution for the transmembrane potential and any additional state variables, and - “vur” (
dolfin.Function
) representing the transmembrane potential in combination with the extracellular potential and an additional Lagrange multiplier.
The algorithm can be controlled by a number of parameters. In particular, the splitting algorithm can be controlled by the parameter “theta”: “theta” set to 1.0 corresponds to a (1st order) Godunov splitting while “theta” set to 0.5 to a (2nd order) Strang splitting.
This solver has not been optimised for computational efficiency and should therefore primarily be used for debugging purposes. For an equivalent, but more efficient, solver, see
cbcbeat.splittingsolver.SplittingSolver
.- Arguments
- model (
cbcbeat.cardiacmodels.CardiacModel
) - a CardiacModel object describing the simulation set-up
- params (
dolfin.Parameters
, optional) - a Parameters object controlling solver parameters
- model (
- Assumptions
- The cardiac conductivities do not vary in time
-
static
default_parameters
()[source]¶ Initialize and return a set of default parameters for the splitting solver
- Returns
- A set of parameters (
dolfin.Parameters
)
To inspect all the default parameters, do:
info(BasicSplittingSolver.default_parameters(), True)
-
merge
(solution)[source]¶ Combine solutions from the PDE solve and the ODE solve to form a single mixed function.
- Arguments
- solution (
dolfin.Function
) - Function holding the combined result
- solution (
-
solution_fields
()[source]¶ Return tuple of previous and current solution objects.
Modifying these will modify the solution objects of the solver and thus provides a way for setting initial conditions for instance.
- Returns
- (previous vs, current vs, current vur) (
tuple
ofdolfin.Function
)
-
solve
(interval, dt)[source]¶ Solve the problem given by the model on a given time interval (t0, t1) with a given timestep dt and return generator for a tuple of the time step and the solution fields.
- Arguments
- interval (
tuple
) - The time interval for the solve given by (t0, t1)
- dt (int, list of tuples of floats)
- The timestep for the solve. A list of tuples of floats can also be passed. Each tuple should contain two floats where the first includes the start time and the second the dt.
- interval (
- Returns
- (timestep, solution_fields) via (
genexpr
)
Example of usage:
# Create generator dts = [(0., 0.1), (1.0, 0.05), (2.0, 0.1)] solutions = solver.solve((0.0, 1.0), dts) # Iterate over generator (computes solutions as you go) for ((t0, t1), (vs_, vs, vur)) in solutions: # do something with the solutions
-
step
(interval)[source]¶ Solve the problem given by the model on a given time interval (t0, t1) with timestep given by the interval length.
- Arguments
- interval (
tuple
) - The time interval for the solve given by (t0, t1)
- interval (
- Invariants
- Given self._vs in a correct state at t0, provide v and s (in self.vs) and u (in self.vur) in a correct state at t1. (Note that self.vur[0] == self.vs[0] only if theta = 1.0.)
- “vs” (
1.13. cbcbeat.utils module¶
This module provides various utilities for internal use.
-
cbcbeat.utils.
state_space
(domain, d, family=None, k=1)[source]¶ Return function space for the state variables.
- Arguments
- domain (
dolfin.Mesh
) - The computational domain
- d (int)
- The number of states
- family (string, optional)
- The finite element family, defaults to “CG” if None is given.
- k (int, optional)
- The finite element degree, defaults to 1
- domain (
- Returns
- a function space (
dolfin.FunctionSpace
)
-
cbcbeat.utils.
end_of_time
(T, t0, t1, dt)[source]¶ Return True if the interval (t0, t1) is the last before the end time T, otherwise False.
-
cbcbeat.utils.
convergence_rate
(hs, errors)[source]¶ Compute and return rates of convergence \(r_i\) such that
\[errors = C hs^r\]
-
class
cbcbeat.utils.
Projecter
(V, params=None)[source]¶ Bases:
object
Customized class for repeated projection.
- Arguments
- V (
dolfin.FunctionSpace
) - The function space to project into
- solver_type (string, optional)
- “iterative” (default) or “direct”
- V (
- Example of usage::
- my_project = Projecter(V, solver_type=”direct”) u = Function(V) f = Function(W) my_project(f, u)
1.14. Module contents¶
The cbcbeat Python module is a problem and solver collection for cardiac electrophysiology models.
To import the module, type:
from cbcbeat import *