>Formulation of a Forward Problem of Electrocardiography using the Boundary Element Method

Rob MacLeod

In this discussion, the forward and inverse problems of electrocardiography have been solved in terms of epicardial and body surface potentials. The method of the forward solution is the boundary element method applied to a realistic homogeneous human torso the result is a forward transfer coefficient matrix ZBH. Direct inversion of ZBH is not possible since the problem is ill-posed and any resulting inverse solution would be unstable. Therefore, we have applied the technique of regularization to compute a constrained inverse transfer coefficient matrix, ZHB, with which epicardial potentials can be safely predicted from body surface potential maps. This section outlines the formulation of the problem and describes the steps involved in the generation of the forward solution matrix.

Analytical derivation

The starting point of the derivation is Green's second identity for two scalar, piecewise continuous functions, f and g, which states that

$\displaystyle \int_{S}^{}$(f$\displaystyle \nabla$g - g$\displaystyle \nabla$f ) . d$\displaystyle \vec{A}\,$ = $\displaystyle \int_{V}^{}$(f$\displaystyle \nabla^{2}_{}$g - g$\displaystyle \nabla^{2}_{}$f )dV, (1)

where, the Green's volume V is bounded by surface S whose outward directed element is denoted d$ \vec{A}\,$. We can define the function f as 1/r, where r is the distance from the source point q to an arbitrary field point p, called the observation point, and the function g as the scalar electric potential, $ \phi$. The area element, d$ \vec{A}\,$ is defined as $ \vec{n}\,$ dA, where $ \vec{n}\,$ is the outward normal to the surface and dAis the scalar element of area of the integrating surface. With these substitutions, the Green's identity now reads

$\displaystyle \int_{S}^{}$($\displaystyle {\frac{1}{r}}$$\displaystyle \nabla$$\displaystyle \phi$ - $\displaystyle \phi$$\displaystyle \nabla$$\displaystyle {\frac{1}{r}}$) . d$\displaystyle \vec{A}\,$ = $\displaystyle \int_{V}^{}$($\displaystyle {\frac{1}{r}}$$\displaystyle \nabla^{2}_{}$$\displaystyle \phi$ - $\displaystyle \phi$$\displaystyle \nabla^{2}_{}$$\displaystyle {\frac{1}{r}}$) dV. (2)

If there are no sources within the volume, V, the scalar potential $ \phi$satisfies the Laplace equation $ \nabla^{2}_{}$$ \phi$ = 0, and the first term of the volume integral disappears. The second term of the volume integral includes the expression $ \nabla^{2}_{}$r-1, which equals 0 for all r $ \neq$ 0. In the immediate vicinity of r = 0,$ \phi$ can be considered constant so that the volume integral becomes

- $\displaystyle \int_{V}^{}$$\displaystyle \phi$$\displaystyle \nabla^{2}_{}$$\displaystyle {\frac{1}{r}}$ dV = $\displaystyle \phi$(ps)C, (3)

where ps represents the location of a singularity at which observation point and integration volume meet (r = 0) and

C = - $\displaystyle \int_{V_{s}}^{}$$\displaystyle \nabla^{2}_{}$$\displaystyle {\frac{1}{r}}$ dV, (4)

where Vs is the reduced volume which just includes the singularity. It can be shown [1] that in general the integrand may be written as

$\displaystyle \nabla^{2}_{}$$\displaystyle {\frac{1}{r}}$ dV = - 4$\displaystyle \pi$$\displaystyle \delta$($\displaystyle \vec{r}\,$ - $\displaystyle \vec{r^\prime}\,$), (5)

where $ \vec{r}\,$ and $ \vec{r^{\prime}}\,$ are vectors from the origin to the observation point and any point p within the volume of integration, respectively. The integral can now be evaluated as

- $\displaystyle \int_{V_{s}}^{}$$\displaystyle \nabla^{2}_{}$$\displaystyle {\frac{1}{r}}$ dV = $\displaystyle \int_{V_{s}}^{}$4$\displaystyle \pi$$\displaystyle \delta$($\displaystyle \vec{r}\,$ - $\displaystyle \vec{r^\prime}\,$dV = $\displaystyle \left\{\vphantom{ \begin{array}{cl}
0 & \mbox{($p$\space outside...
... )} \\
4\pi & \mbox{($p$\space inside ${\rm V_{s}}$ ).}
\end{array} }\right.$$\displaystyle \begin{array}{cl}
0 & \mbox{($p$\space outside ${\rm V_{s}}$ )} \\
4\pi & \mbox{($p$\space inside ${\rm V_{s}}$ ).}
\end{array}$ $\displaystyle \left.\vphantom{ \begin{array}{cl}
0 & \mbox{($p$\space outside ...
... )} \\
4\pi & \mbox{($p$\space inside ${\rm V_{s}}$ ).}
\end{array} }\right.$ (6)

It is possible to choose an observation point anywhere within the Green's volume to satisfy equation 2. The standard approach is to select the observation point very close to the boundary of, but still within, the Green's volume. Physically, we can assume that the potential at such a point is virtually identical to that on the real surface in question, while mathematically, the points remain inside the volume and the result -4$ \pi$ can be used for the integral in equation 6. This approach was taken by both Barr et al. [2] and Messinger-Rapport and Rudy [3] and from equation 2 yields the following equation:

4$\displaystyle \pi$$\displaystyle \phi$(p) = $\displaystyle \int_{S}^{}$($\displaystyle {\frac{1}{r}}$$\displaystyle \nabla$$\displaystyle \phi$ - $\displaystyle \phi$$\displaystyle \nabla$$\displaystyle {\frac{1}{r}}$) . d$\displaystyle \vec{A}\,$. (7)

A second strategy is possible, however, in which the observation point is located not just inside, but on the surface. It can be shown [4] that the integral in equation 3 can be evaluated at a point on the surface with the result -2$ \pi$. This results in a slightly different form of equation 2 than that given in equation 7:

2$\displaystyle \pi$$\displaystyle \phi$(p) = $\displaystyle \int_{S^{-}}^{}$($\displaystyle {\frac{1}{r}}$$\displaystyle \nabla$$\displaystyle \phi$ - $\displaystyle \phi$$\displaystyle \nabla$$\displaystyle {\frac{1}{r}}$) . $\displaystyle \vec{dA}\,$, (8)

the obvious difference being an additive factor of 2$ \pi$$ \phi$(p). S- in this equation is the surface area excluding the singularity, that is, the integral is proper-valued and must be evaluated excluding the singularity at r = 0. To maintain consistency between the two equations 8 and 7, if the observation point is taken close to but still inside the surface of integration, then the volume integral of $ \nabla^{2}_{}$(1/r)dV yields the result -4$ \pi$, as in equation 6, and the surface integration is also performed over the entire surface; there is no singularity. If, on the other hand, the observation point is placed on the surface of integration, the volume integral yields the value -2$ \pi$ and the surface integral must be calculated by excluding the singularity, that is, as a proper-valued integral. Conceptually, placing the observation point on the surface is simpler and more accurate, and therefore, this approach is used in further discussion.

Figure 1 depicts the geometry to which we now apply Green's second identity. In this particular case, the Green's volume is the region enclosed by the epicardial and body surfaces, SH and SB, which together comprise S. We proceed from equation 8 by separating the surface which bounds the Green's volume V into the heart surface, SH, and the body surface, SB. The sense of the normal to the heart surface is redirected into the Green's volume so that by rearranging the integrals, we now have

2$\displaystyle \pi$$\displaystyle \phi$(p) = $\displaystyle \int_{S_{H}}^{}$$\displaystyle \phi$ d$\displaystyle \Omega$ - $\displaystyle \int_{S_{B}}^{}$$\displaystyle \phi$ d$\displaystyle \Omega$ - $\displaystyle \int_{S_{H}}^{}$$\displaystyle {\frac{\nabla \phi}{r}}$ . d$\displaystyle \vec{A}\,$ + $\displaystyle \int_{S_{B}}^{}$$\displaystyle {\frac{\nabla \phi}{r}}$ . d$\displaystyle \vec{A}\,$. (9)

We have introduced d$ \Omega$, the incremental solid angle, for the $ \nabla$(1/r) . d$ \vec{A}\,$

Derivation of the coefficient matrices

The general approach to finding solutions to integral equations like equation 9 is to write one equation for each of a number of points on both of the surfaces and solve these equations simultaneously. This is known as the collocation method in numerical mathematics and provides a means of reducing an integration over an arbitrary smooth surface to a sum of (somewhat simpler) integrals, each of which can be evaluated separately [5]. The particular application of the collocation method to this problem originates with work of Barr, Ramsey and Spach [2] and their notation will be used throughout.

We start by rewriting equation 9 for two observation points placed on the body and epicardial surface, respectively, as

$\displaystyle \phi_{B}^{i}$ - $\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{H}}^{}$$\displaystyle \phi_{H}^{}$ d$\displaystyle \Omega_{BH}^{i}$ + $\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{B}}^{}$$\displaystyle \phi_{B}^{}$ d$\displaystyle \Omega_{BB}^{i}$ + $\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{H}}^{}$$\displaystyle {\frac{\nabla \phi_{H}}{r^{i}}}$ . d$\displaystyle \vec{A}\,$ = 0 (10)


$\displaystyle \phi_{H}^{i}$ - $\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{H}}^{}$$\displaystyle \phi_{H}^{}$ d$\displaystyle \Omega_{HH}^{i}$ + $\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{B}}^{}$$\displaystyle \phi_{B}^{}$ d$\displaystyle \Omega_{HB}^{i}$ + $\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{H}}^{}$$\displaystyle {\frac{\nabla \phi_{H}}{r^{i}}}$ . d$\displaystyle \vec{A}\,$ = 0. (11)

Here we have defined the differential solid angle d$ \Omega_{PQ}^{i}$ as that subtended by an elemental area of the integration surface Q at the ithobservation point on surface P. Likewise, $ \phi^{i}_{P}$ is the potential at the ith location on the surface P.

For NB points defined on the body surface and NH points on the epicardium, we can write equations 10 and 11, NB and NH times, respectively (collocation method). Let us write discretized expressions for each of the resulting terms in equations 10 and 11; these expressions will be examined in detail below. From equation 10:

$\displaystyle \phi_{B}^{i}$ + $\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{B}}^{}$$\displaystyle \phi_{B}^{}$ d$\displaystyle \Omega_{BB}^{i}$ = $\displaystyle \sum_{j=1}^{N_{B}}$pijBB$\displaystyle \Phi^{j}_{B}$, (12)

- $\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{H}}^{}$$\displaystyle \phi_{H}^{}$d$\displaystyle \Omega_{BH}^{i}$ = $\displaystyle \sum_{j=1}^{N_{H}}$pijBH$\displaystyle \Phi^{j}_{H}$, (13)

$\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{H}}^{}$$\displaystyle {\frac{\nabla \phi_{H}}{r^{i}}}$ . d$\displaystyle \vec{A}_{H}^{}$ = $\displaystyle \sum_{j=1}^{N_{H}}$gijBH$\displaystyle \Gamma^{j}_{H}$, (14)

and from equation  11:

$\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{B}}^{}$$\displaystyle \phi_{B}^{}$ d$\displaystyle \Omega_{HB}^{i}$ = $\displaystyle \sum_{j=1}^{N_{B}}$pijHB$\displaystyle \Phi^{j}_{B}$, (15)

$\displaystyle \phi_{H}^{i}$ - $\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{H}}^{}$$\displaystyle \phi_{H}^{}$d$\displaystyle \Omega_{HH}^{i}$ = $\displaystyle \sum_{j=1}^{N_{H}}$pijHH$\displaystyle \Phi^{j}_{H}$, (16)


$\displaystyle {\frac{1}{2\pi}}$$\displaystyle \int_{S_{H}}^{}$$\displaystyle {\frac{\nabla \phi_{H}}{r^{i}}}$ . d$\displaystyle \vec{A}_{H}^{}$ = $\displaystyle \sum_{j=1}^{N_{H}}$gijHH$\displaystyle \Gamma^{j}_{H}$. (17)

The argument of each of the summations can be separated into the product of a potential ( $ \Phi_{B}^{j}$ or $ \Phi_{H}^{j}$) or the gradient of a potential ( $ \Gamma_{H}^{j}$) at a specific point j on either one of the surfaces and a second factor (the p's and g's) based entirely on the geometry of the torso and the heart. $ \Phi^{j}_{B}$ and $ \Phi^{j}_{H}$ are the potentials at node j on the body and heart surfaces, respectively; $ \Gamma^{j}_{H}$ is the normal component of the potential gradient for point j on the heart surface (the corresponding quantity is zero on the body surface). In general the gPQij term links the value of the potential gradient ( $ \Gamma^{j}_{}$) at point j on surface P to the observation point i on surface Q while pPQij is the geometrical coefficient which weights the contribution of the potential at node j of surface Q to the potential at observation point i on surface P. The first subscript and superscript of each p or g term indicate the observation point, the second subscript and superscript, the element of the surface of integration; thus, for example, pijHB is the coefficient for the observation point i on the epicardial surface and point j on the surface of integration, the body surface.

The Forward Transform Matrix

Now by inserting the appropriate right-hand sides of equations  1217, we get the discretized equivalent equations to 10 and 11,

piBB$\displaystyle \Phi_{B}^{}$ + piBH$\displaystyle \Phi_{H}^{}$ + giBH$\displaystyle \Gamma_{H}^{}$ = 0 (18)


piHB$\displaystyle \Phi_{B}^{}$ + piHH$\displaystyle \Phi_{H}^{}$ + giHH$\displaystyle \Gamma_{H}^{}$ = 0, (19)

where the summations over j are implicit in each term and i refers to a specific observation point on either the heart (equation 19) or the body (equation 18) surface. The p and g terms are the row vectors which express the geometrical contribution of each point on the surface of integration to the potential at the observation point i. If we write equation 18 for each point on the body surface and equation 19 for each point on the heart surface, two sets of equations result, which in matrix notation can be written as:

PBB$\displaystyle \Phi_{B}^{}$ + PBH$\displaystyle \Phi_{H}^{}$ + GBH$\displaystyle \Gamma_{H}^{}$ = 0 (20)


PHB$\displaystyle \Phi_{B}^{}$ + PHH$\displaystyle \Phi_{H}^{}$ + GHH$\displaystyle \Gamma_{H}^{}$ = 0. (21)

The P's and G's are the matrices formed by collocating all the elements of the associated pi and gi row vectors, one row for each observation point. For the PHB matrix, for example, each row contains NH elements from one pHBi vector, and there are NH rows, each representing a different value of i. Here again, the first subscript represents the surface containing the observation points, the second subscript the surface of integration. PHH and GHH are square matrices of size NH x NH, PBH and GBH are sized NB x NH, PBB is another square matrix of size NB x NB, and PHBis sized NH x NB.

By solving the equation 21 for $ \Gamma_{H}^{}$ and substituting the result into equation 20 we remove the need for explicit knowledge of the potential gradients. This leads, after sorting of variables, to

(PBB - GBHGHH-1PHB)$\displaystyle \Phi_{B}^{}$ = (GBHGHH-1PHH - PBH)$\displaystyle \Phi_{H}^{}$, (22)

which can be rewritten as

$\displaystyle \Phi_{B}^{}$ = ZBH$\displaystyle \Phi_{H}^{}$, (23)

with ZBH defined as


Equations 23 and 24 define the solution to the forward problem in the desired form; ZBH is the transfer coefficient matrix which directly relates epicardial potentials to body surface potentials. It remains to examine each term of equation 22 and develop accurate computational methods for evaluating the elements of all the matrices involved. Once this has been achieved, we may convert epicardial potentials into body surface potentials and, thus, solve the forward problem.


J.D. Jackson.
Classical Electrodynamics.
John Wiley & Sons, New York, 1975.
Second edition.

R.C. Barr, M. Ramsey, and M.S. Spach.
Relating epicardial to body surface potential distributions by means of transfer coefficients based on geometry measurements.
IEEE Trans Biomed. Eng., BME-24:1-11, 1977.

Y. Rudy and B.J. Messinger-Rapport.
The inverse solution in electrocardiography: Solutions in terms of epicardial potentials.
Crit. Rev. Biomed. Eng., 16:215-268, 1988.

A.C.L. Barnard, I.M. Duck, and M.S. Lynn.
The application of electromagnetic theory to electrocardiology: I Derivation of the integral equation.
Biophys. J., 7:443-462, 1967.

C.W. Steele.
Numerical Computation of Electrical and Magnetic Fields.
Van Nostrand Reinhold, New York, 1987.

About this document ...

Solutions to the Forward and Problem of Electrocardiography using the Boundary Element Method

This document was generated using the LaTeX2HTML translator Version 98.1p1 release (March 2nd, 1998)

Copyright © 1993, 1994, 1995, 1996, 1997, Nikos Drakos, Computer Based Learning Unit, University of Leeds.

The command line arguments were:
latex2html -no_math -html_version 3.2,math -split 3 -link 1 -no_navigation bem-form.

The translation was initiated by Rob MacLeod on 1999-06-10


Rob MacLeod