1. cbcbeat package

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

\[ \begin{align}\begin{aligned}v_t - \mathrm{div} ( G_i v + G_i u) = I_s\\\mathrm{div} (G_i v + (G_i + G_e) u) = I_a\end{aligned}\end{align} \]

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\):

\[v(x, 0) = v_0\]

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
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 of dolfin.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
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
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
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 or dolfin.PETScKrylovSolver).

nullspace
step(interval)[source]

Solve on the given time step (t0, t1).

Arguments
interval (tuple)
The time interval (t0, t1) for the step
Invariants
Assuming that v_ is in the correct state for t0, gives self.vur in correct state at t1.
variational_forms(k_n)[source]

Create the variational forms corresponding to the given discretization of the given system of equations.

Arguments
k_n (ufl.Expr or float)
The time step
Returns
(lhs, rhs) (tuple of ufl.Form)

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
applied_current()[source]

An applied current: used as a source in the elliptic bidomain equation

cell_models()[source]

Return the cell models

conductivities()[source]

Return the intracellular and extracellular conductivities as a tuple of UFL Expressions.

Returns (M_i, M_e) (tuple of ufl.Expr)

domain()[source]

The spatial domain (dolfin.Mesh).

extracellular_conductivity()[source]

The intracellular conductivity (ufl.Expr).

intracellular_conductivity()[source]

The intracellular conductivity (ufl.Expr).

stimulus()[source]

A stimulus: used as a source in the parabolic bidomain equation

time()[source]

The current time (dolfin.Constant or None).

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
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
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 of 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
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
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
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
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
class cbcbeat.cellsolver.SingleCellSolver(model, time, params=None)[source]

Bases: cbcbeat.cellsolver.CardiacODESolver

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_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 of dolfin.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
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)
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.)

1.8. cbcbeat.gotran2cellmodel module

1.9. cbcbeat.gotran2dolfin module

class cbcbeat.gotran2dolfin.DOLFINCodeGenerator(code_params=None)[source]

Bases: PythonCodeGenerator

Class for generating a DOLFIN compatible declarations of an ODE from a gotran file

static default_parameters()[source]
function_code(comp, indent=0, default_arguments=None, include_signature=True)[source]
init_parameters_code(ode, indent=0)[source]

Generate code for setting parameters

init_states_code(ode, indent=0)[source]

Generate code for setting initial condition

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 markers

letting:

g = Markerwise((g0, g1), (2, 5), markers)
keys()[source]

The keys or domain numbers

markers()[source]

The markers

values()[source]

The objects

cbcbeat.markerwisefield.handle_markerwise(g, classtype)[source]
cbcbeat.markerwisefield.rhs_with_markerwise_field(g, mesh, v)[source]

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

\[v_t - \mathrm{div} ( G_i v) = I_s\]

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\):

\[v(x, 0) = v_0\]

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
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 of dolfin.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
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
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
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 or dolfin.KrylovSolver).

step(interval)[source]

Solve on the given time step (t0, t1).

Arguments
interval (tuple)
The time interval (t0, t1) for the step
Invariants
Assuming that v_ is in the correct state for t0, gives self.v in correct state at t1.
variational_forms(k_n)[source]

Create the variational forms corresponding to the given discretization of the given system of equations.

Arguments
k_n (ufl.Expr or float)
The time step
Returns
(lhs, rhs, prec) (tuple of ufl.Form)

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

\[ \begin{align}\begin{aligned}v_t - \mathrm{div} (M_i \mathrm{grad} v + M_i \mathrm{grad} u) = - I_{ion}(v, s) + I_s\\ \mathrm{div} (M_i \mathrm{grad} v + (M_i + M_e) \mathrm{grad} u) = I_a\\s_t = F(v, s)\end{aligned}\end{align} \]

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\):

\[ \begin{align}\begin{aligned}v(x, 0) = v_0\\s(x, 0) = s_0\end{aligned}\end{align} \]

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

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
static default_parameters()[source]

Initialize and return a set of default parameters for the splitting solver

Returns
The set of default parameters (dolfin.Parameters)

Example of usage:

info(SplittingSolver.default_parameters(), True)
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
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_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 of dolfin.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.
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)
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.)

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
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”
Example of usage::
my_project = Projecter(V, solver_type=”direct”) u = Function(V) f = Function(W) my_project(f, u)
static default_parameters()[source]

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 *