CardioWave | Dev / DiscreteMathematics browse

The computer needs discrete quantities to work with and thus we cannot directly apply computational methods to the Bidomain Equations. We first must grid or mesh the domain and then run some discretization method on the "chopped up" domain. The simulation system does not directly impose any restrictions on the gridding or discretization method. Regular, irregular, and warped regular grids may all be used in the gridding stage. Finite differences, finite elements, or box (finite volume) methods may then be used on any of the grids. The result should be a matrix representation of the elliptic operators, one for each domain. The simulation system will then read in these matrices and begin the simulation. Note, of course, that if the matrices are generated improperly, the CardioWave codes have no real way of knowing this other than the computed results running off to infinity. There are also a number of "hidden" links between the modules and the matrix. In particular, the Conjugate Gradient is only guaranteed to work on symmetric matrices - but again, there is no easy way for the simulation code to know if the matrix is symmetric (other than to crank through the actual work of looking at A(i,j) and A(j,i) for all i and j). Similarly, if a matrix is very poorly conditioned, it may require 500 or more iterations to solve. However, the "maximum iteration" parameter for the solver is settable by the user at run time and thus could easily be set too small and cause apparent non-convergence.

Some simple tests could perhaps be performed by various modules to insure that the matrix or matrices meet certain criteria. However, such testing may be costly especially if the matrix is changing frequently as is the case for adaptive time-steppers and implicit methods.

Taking a look at the discretized form of the Monodomain Equations:

The discretization process produces a matrix representation of the elliptic operator, Mi, and brings out a term to scale for the different Voronoi areas or volumes, A. Note that this A term scales both parts of the transmembrane current (capacitive and ionic). The lower case vm is now the discrete approximation to the continuous system, hence the partial derivatives now become total derivatives. We have also rewritten the membrane functions as Fq() and Fv() since they will be computed together, in one function call, to improve computational efficiency. Note that we have also switched the signs on the functions. For the sake of simplicity, we will always assume the monodomain to have ve=0 everywhere, thus vm=vi and hence the subscripts on the matrix M and stimulus S are considered to be intracellular.

For the Bidomain Equations:

We see similar changes. Both elliptic operators are represented by their matrix approximations, Mi and Me. Again, the Voronoi area or volume term is used to scale the transmembrane current pieces. The lower case variables, vm, ve, and vi are used to designate the discrete approximations to Vm, Phi-e, and Phi-i, respectively.