 Research
 Open Access
OpenMEEG: opensource software for quasistatic bioelectromagnetics
 Alexandre Gramfort^{1, 2}Email author,
 Théodore Papadopoulo^{1},
 Emmanuel Olivi^{1} and
 Maureen Clerc^{1}Email author
https://doi.org/10.1186/1475925X945
© Gramfort et al; licensee BioMed Central Ltd. 2010
 Received: 4 May 2010
 Accepted: 6 September 2010
 Published: 6 September 2010
Abstract
Background
Interpreting and controlling bioelectromagnetic phenomena require realistic physiological models and accurate numerical solvers. A semirealistic model often used in practise is the piecewise constant conductivity model, for which only the interfaces have to be meshed. This simplified model makes it possible to use Boundary Element Methods. Unfortunately, most Boundary Element solutions are confronted with accuracy issues when the conductivity ratio between neighboring tissues is high, as for instance the scalp/skull conductivity ratio in electroencephalography. To overcome this difficulty, we proposed a new method called the symmetric BEM, which is implemented in the OpenMEEG software. The aim of this paper is to present OpenMEEG, both from the theoretical and the practical point of view, and to compare its performances with other competing software packages.
Methods
We have run a benchmark study in the field of electro and magnetoencephalography, in order to compare the accuracy of OpenMEEG with other freely distributed forward solvers. We considered spherical models, for which analytical solutions exist, and we designed randomized meshes to assess the variability of the accuracy. Two measures were used to characterize the accuracy. the Relative Difference Measure and the Magnitude ratio. The comparisons were run, either with a constant number of mesh nodes, or a constant number of unknowns across methods. Computing times were also compared.
Results
We observed more pronounced differences in accuracy in electroencephalography than in magnetoencephalography. The methods could be classified in three categories: the linear collocation methods, that run very fast but with low accuracy, the linear collocation methods with isolated skull approach for which the accuracy is improved, and OpenMEEG that clearly outperforms the others. As far as speed is concerned, OpenMEEG is on par with the other methods for a constant number of unknowns, and is hence faster for a prescribed accuracy level.
Conclusions
This study clearly shows that OpenMEEG represents the state of the art for forward computations. Moreover, our software development strategies have made it handy to use and to integrate with other packages. The bioelectromagnetic research community should therefore be able to benefit from OpenMEEG with a limited development effort.
Keywords
 Boundary Element Method
 Electrical Impedance Tomography
 Functional Electrical Stimulation
 Head Model
 Forward Problem
Introduction
Many devices used in the clinical or the cognitive science domain perform electromagnetic measurements, or stimulation, on the human body. Devices measuring electric fields include electroencephalograms (EEG), cardiograms (ECG), myograms (EMG), while magnetoencephalograms (MEG) or magnetocardiograms (MCG) measure magnetic fields. Among stimulating devices, transcranial magnetic stimulation (TMS) uses magnetic coils to stimulate brain regions, while functional electric stimulation (FES) and electrical impedance tomography (EIT) impose an electric current or potential through contact electrodes.
To interpret and control the bioelectromagnetic phenomena involved with these devices, realistic physiological modeling is required, in terms of geometry and conductivity [1]. Accurate numerical solutions of the governing equations must be computed: obtaining the best accuracy possible for a given computational model is one of the goals of OpenMEEG, the opensource software package introduced in this article.
written here in the case where j = 0 on the boundary.
A forward problem consists of simulating V and/or B when σ, J^{p}, and boundary current j are known. The forward electro magnetostatics problem is wellposed, in contrast to the far more difficult, illposed inverse problem of estimating σ, or J^{p}, from partial boundary data. Still, obtaining an accurate solution for the forward problem is far from trivial. This paper presents a software package, OpenMEEG, that makes available to the community recent research efforts to improve the accuracy of forward solvers.
Forward solutions provide the relationship between the quantities of interest and the measurements. To obtain a good description of this relationship, one must model the sources, the conductivity, and the sensors.
 1.
if the conductivity field can be described using simple geometries (with axilinear, planar, cylindrical, spherical or ellipsoidal symmetry), analytical methods can be derived, and fast algorithms have been proposed that converge to the analytical solutions for EEG [2, 3] and MEG [4];
 2.
for piecewise constant conductivity fields as in Figure 1, Boundary Element Methods (BEMs) can be applied, resulting in a simplified description of the geometry only on the boundaries [5];
 3.
general nonhomogeneous and anisotropic conductivity fields are handled by volumic methods; Finite Element Methods and Finite Difference Methods belong to this category.
This paper deals with Boundary Element Methods, whose advantage over volumic methods is to use an economic representation of the conductivity field (a few conductivity parameters  one per tissue, and a few triangular meshes to represent the interfaces). Until recently, all Boundary Element Methods used in bioelectromagnetics were based on a Green representation theorem, involving operators called single and doublelayer potentials. In OpenMEEG the approach is to consider an extended version of the Green representation theorem [6], involving more operators, and leading to a new BEM formulation, called the Symmetric Boundary Element Method [7].
The structure of the paper is as follows: the applications targeted by OpenMEEG are first presented, then some mathematical aspects of the Symmetric BEM are explained (details are presented in an appendix). Then, to motivate the use of OpenMEEG, a comparison study with four other solvers in the context of EEG forward modeling is presented. The accuracy of the solvers is tested on multiple random sphere models. The accuracy of OpenMEEG is also assessed by numerical experiments for MEG forward modeling. Finally a section provides technical details on the OpenMEEG software package while giving sample code for using OpenMEEG via Python or via the Fieldtrip Matlab toolbox. A more complete presentation of the usage of OpenMEEG via a command line interface is available in Additional file 1.
Target applications of OpenMEEG
The main purpose of OpenMEEG is to solve the forward problems arising in magneto and electroencephalography. OpenMEEG is therefore primarily developed in the field of brain research, but its development could directly benefit to other problems dealing with electromagnetic biosignals (see for example [8] in the context of electrocardiography). A geometrical model of nested meshes representing tissue interfaces must be provided to OpenMEEG, which does not perform any segmentation or meshing. OpenMEEG handles nested regions, and could in principle be extended to more general, disjoint regions [9]. The main output of a forward problem is a leadfield, i.e., the linear application relating sources at specific locations to sensor measurements. Although the principal target of OpenMEEG is magneto and electroencephalography, other types of bioelectromagnetic problems have also been handled with OpenMEEG: we hereforth describe its current scope.
Electroencephalography (EEG) is concerned with the variations of electric potential on the scalp, due to sources within the brain. At frequencies of interest, the quasistatic regime is valid, resulting in the electrostatics relation (1). The air surrounding the scalp is supposed nonconductive, hence the normal current vanishes on the scalp: j = 0 in boundary condition (2). The sources within the brain are represented by dipoles: in equation (1), the sources are ${J}^{\text{p}}=\overrightarrow{q}{\delta}_{p}$ where p is a dipole position and $\overrightarrow{q}$ the associated dipolar moment. OpenMEEG computes the electric potential and the normal current on each interface between two homogeneous tissues, due to electric sources within the brain. EEG sensors are electrodes, modeled in OpenMEEG as discrete positions on the scalp at which the potential can be measured (infinite impedance assumption). On these sensors, OpenMEEG computes the EEG leadfield, representing the linear relationship between source amplitude (for fixed position and orientation) and sensor values.
Magnetoencephalography (MEG) is concerned with magnetic fields produced by sources within the brain, which, like for EEG, are modeled as dipoles [10]. OpenMEEG makes use of the Biot and Savart relation (3) to compute the magnetic field, and hence requires the electric potential to be computed beforehand on all the interfaces [11]. Magnetometers or gradiometers can be modeled in OpenMEEG. Magnetometers are defined by their position and the direction of the field they measure. Gradiometers and more general sensors are handled by providing to OpenMEEG the positions, the orientations and the weights of integration points. For example with axial gradiometers present in CTF MEG systems, the forward field for a sensor is obtained by subtracting the two leadfields computed at the locations of two nearby magnetometers.
Electrical Impedance Tomography (EIT) infers characteristics of a conductive domain, by analyzing the potential resulting from the application of a current on the boundary. This method has been applied to calibrate the conductivity of EEG head models [12–14]. OpenMEEG computes the forward problem associated to EIT: given a conductive domain Ω defined by the interfaces between homogeneous regions, and their conductivity, and for a prescribed normal current j, OpenMEEG computes the potential V and the normal current σ∂_{ n } V on each interface by solving (1) for J^{p} = 0, and with boundary condition (2). The electrodes are modeled with a P0 approximation over the triangles of the scalp; the injected current is assumed constant over a triangle.
Intracranial electric potentials are measured in certain clinical settings, either on the surface (ElectroCorticography) or within the brain (intracranial EEG, or stereotaxic EEG). OpenMEEG is able to compute leadfields for such intracranial electrodes, in realistic head models.
Functional Electrical Stimulation (FES) provides a way to restore movement of paralyzed body regions by activating the efferent somatic axons. For this a current is applied to a nerve sheath, using specially designed electrodes. The precise location and time course of the applied current must be optimized in order to achieve the best selectivity, and to minimize the current intensity for a desired outcome. Optimizing the stimulation parameter in a realistic nerve model requires a forward model for FES. OpenMEEG provides such a tool, by combination of the concepts of EIT and of internal potential simulation (see section on the applications targeted by OpenMEEG) [15, 16].
Cortical Mapping is an inverse problem that aims to recover the potential and the normal current occurring on the surface of the cortex (i.e., under the skull bone), given EEG measurements on electrodes [17]. A particularly elegant solution to this problem has been proposed with the symmetric BEM [18], making it possible e.g. to solve further source localization problems.
Methods: implementation
OpenMEEG uses a Galerkin Boundary Element formulation, that jointly considers the electric potential and the normal current on each interface as unknowns of the problem. A P1 (piecewise linear) approximation is used for the electric potential, whereas the normal current $\sigma \frac{\partial V}{\partial n}$ is approximated with P0 (piecewise constant) boundary elements. The most intricate part of the implementation concerns the assembly of the system matrix, requiring singular kernel integration over triangles. Those are double integrals. The inner integrals are computed with analytical schemes [19], whereas the outer integrals are computed with a 7point Gauss quadrature scheme [20]. An adaptive integration scheme recursively subdivides the triangles until the required precision is achieved. This adaptive integration has an influence on the accuracy, as will be exposed in the benchmark results further on.
Since the electric potential can only be computed up to a constant, the system matrix is deflated to make it invertible. In practice, it consists of imposing the constraint that the integral of the potential over the external layer is 0.
The magnetic field is computed by using the Biot and Savart equation with the Galerkin piecewise linear approximation for the potential [11].
As stated above, given a set of dipole positions and orientations, a set of sensors and a head model defining homogeneous conductivity regions, the M/EEG forward problem produces as output the leadfield. OpenMEEG computes such lead fields with the following procedure (see the tutorial in Additional file 1 for a detailed graphical representation of this flowchart):

assemble the system matrix involving boundary integral operators on the discretized surfaces;

for a specified set of dipole positions and orientations, assemble a discretized, sourcerelated term;

solve the linear system relating the two above matrices (providing V and $\sigma \frac{\partial V}{\partial n}$ on each interface)

interpolate the scalp potential (for EEG) or apply the discretized Biot and Savart relation (for MEG).
Mathematical details, as well as the practical usage of OpenMEEG can be found in the Appendix.
Methods: benchmark comparison study
In order to motivate the use of the OpenMEEG software, we have conducted a set of numerical experiments that compare the accuracy and the robustness of OpenMEEG with alternative M/EEG forward solvers. The comparison in the context of EEG forward modeling is run with the 4 alternative BEM solvers (see [21] for a comparison with FEM). The precision of the MEG forward solver is demonstrated using known analytical properties of the magnetic field when considering sphere models and with two other solvers.
Publicly available M/EEG forward solvers
Several software projects have the ability to solve the M/EEG forward problems. MNE, BrainStorm, EEGLAB (via the NFT Toolbox), Fieldtrip, Simbio, OpenMEEG and SPM, which shares with Fieldtrip the same M/EEG forward solvers.
Fieldtrip and SPM offer two implementations of the BEM. The first one, called Dipoli, was written by Oostendorp [22] and is not open source (only binary files for UNIX systems are available), while the second one, called BEMCP, is opensource and was written by Christophe Phillips during his PhD [23]. The Dipoli solver implements a linear collocation method with Isolated Skull Approach (ISA) [24], whereas BEMCP implements a simple linear collocation method. The Dipoli implementation details can be found in [22]. Another implementation with linear collocation (with and without ISA) can be found in the Matlab toolbox called Helsinki BEM [25]. The MNE package written in pure C code by Hämäläinen also offers a linear collocation (with ISA) implementation of the BEM. The Brainstorm toolbox in its latest version uses only sphere models for both EEG and MEG. Recently the EEGLAB toolbox also has provided a package called NFT that is based on a BEM called METU [26]. Finally, the Simbio forward solver also consists of a BEM with linear collocation and ISA.
Commercial software packages are not listed here. However, to our knowledge, commercial products that provide a BEM solver for forward modeling implement a linear collocation method with ISA. This is for example the case of ASA [27]. Note also that beyond the field of brain research, alternative BEM implementation exist. The SCIRun/BioPSE project contains for example a BEM implementation based on [19] that can solve the forward problem of electrocardiography [28].
Accuracy measures
The accuracy of forward solvers can be assessed for simple geometries such as nested spheres, by comparison with analytical results. The precision of a forward solution is tested with two measures. the RDM (Relative Difference Measure) and the MAG (Magnitude ratio) [29].
In both of these expression, the norm is the discrete ℓ^{2} norm over the set of sensor measurements. The closer to 0 (resp. to 1) the RDM (resp. the MAG), the better it is.
Geometrical models
Results
Results: accuracy of the electric potential simulations
The implementations tested for EEG are: OpenMEEG with and without adaptive numerical integration (abbr. OM and OMNA), Simbio (abbr. SB), Helsinki BEM with and without ISA (abbr. HBI and HB), Dipoli (abbr. DP) and BEMCP (abbr. CP). The METU solver was also tested but we were unable to obtain the precisions advertised in [26], so it was decided not to include it in the comparison.
From these simulations it can be observed that:

HB and CP, that implement a simple linear collocation method, have similar results and are clearly the less precise solvers.

HBI, SB and DP, that implement a linear collocation method with ISA, have very similar results. SB and HBI are however slightly more accurate than DP.

OpenMEEG provides the most precise solutions even when no adaptive integration is used. The adaptive integration significantly improves the results, particularly when the meshes are coarsely sampled (42 and 162 vertices per layer).
Results: accuracy of the magnetic field simulations
We next present results for MEG forward modeling. MEG manufacturers propose 3 kinds of sensors (magnetometers, axial gradiometers, and planar gradiometers), all of which are oriented radially with respect to the helmet.
Results: computation speed
The OpenMEEG Software Package
All the tests performed in this paper were made using the version 2.0 of OpenMEEG.
Licence
OpenMEEG is distributed under the French opensource license CeCILLB. It is intended to give users the freedom to modify and redistribute the software. It is therefore compatible with popular opensource licenses such as the GPL and BSD licenses. The CeCILLB license imposes to anybody distributing a software incorporating OpenMEEG the obligation to give credits (by citing the appropriate publications), in order for all contributions to be properly identified and acknowledged. The references to be acknowledged are [7], and the present article.
Source code
OpenMEEG is implemented in C/C++ with limited external dependencies. It uses the Intel MKL libraries on Windows and ATLAS (BLAS/LAPACK) on Unix systems for fast and accurate linear algebra. A modified version of the MATIO library http://sourceforge.net/projects/matio has been integrated in OpenMEEG for Matlab compatible IOs for vectors and matrices. The source code of OpenMEEG is hosted on the INRIA GForge platform and is accessible to an anonymous user via a public version control system. OpenMEEG binaries and source code are both available from http://openmeeg.gforge.inria.fr.
Multiplatform
OpenMEEG is available as precompiled binaries for GNULinux systems, Mac OS and Windows. OpenMEEG's build and packaging system is based on CMake/CPack http://www.cmake.org allowing easy development and deployment on all architectures.
Parallel processing
Testing
Deployment on multiple architectures with heterogenous hardware and software environments requires testing procedures to assess the stability of the solutions provided by compiled binaries. This testing procedure is run through the CMake/CTest testing software. OpenMEEG test suite guarantees the integrity of the results obtained by MEG, EEG and EIT forward solvers. A part of the tests consists of running OpenMEEG on a 3layer spherical geometry like the one presented in Figure 2. Outputs are then compared to analytical results to test that the accuracy is not degraded by a modification of the code.
Integration
Considerable efforts have been made to facilitate the use of OpenMEEG by the M/EEG community. OpenMEEG can be simply invoked via a command line interface, or via higher levels languages. OpenMEEG can also be called from Python or via the Fieldtrip Toolbox, where it has been fully integrated in the M/EEG forward modeling routines.
OpenMEEG demo script in Python
Python code 

Import openmeeg as om 
# Load data 
condFile = 'om_demo.cond' 
geomFile = 'om_demo.geom' 
dipoleFile = 'cortex.dip' 
squidsFile = 'meg_squids.txt' 
electrodesFile = 'eeg_electrodes.txt' 
geom = om.Geometry() 
geom.read(geomFile,condFile) 
dipoles = om.Matrix() 
dipoles.load(dipoleFile) 
squids = om.Sensors() 
squids.load(squidsFile) 
electrodes = om.Matrix() 
electrodes.load(electrodesFile) 
# Compute forward problem 
gaussOrder = 3; # Integration order 
hm = om.HeadMat(geom,gaussOrder) 
hminv = hm.inverse() 
dsm = om.DipSourceMat (geom, dipoles, gaussOrder) 
ds2 mm = om.DipSource2MEGMat (dipoles, squids). 
h2 mm = om.Head2MEGMat (geom, squids) 
h2em = om.Head2EEGMat (geom, electrodes). 
gain_meg = om.GainMEG (hminv, dsm, h2 mm, ds2mm) 
gain_eeg = om.GainEEG (hminv, dsm, h2em) 
OpenMEEG demo script in Matlab
Matlab code 

%% Structure for BEM volume conduction model 
%% Each layer mesh is indexed by k 
% vol.bnd(k).pnt. : vertices for mesh "k" 
% vol.bnd (k).tri : triangles for mesh "k" 
% %% Set the conductivities of each domain 
% vol .cond : conductivities 
%% EEG electrodes 
% sens.pnt : locations of electrodes 
%% Positions of the dipoles 
% pos :locations of dipoles 
%% Compute the BEM 
% choose BEM method (OpenMEEG, BEMCP or Dipoli) 
cfg.method = 'openmeeg'; 
% Compute the BEM matrix 
vol = ft_prepare_bemmodel(cfg, vol); 
cfg. vol = vol; 
cfg. grid. pos = pos; 
cfg. elec = sens; 
% Compute leadfield 
% with no orientation constraint 
lf_openmeeg = ft_prepare_leadfield(cfg); 
Benchmark
The benchmark presented in the previous section was run within the Fieldtrip environment and its MATLAB source code can be obtained for noncommercial use from the authors. Since SPM uses the same code as Fieldtrip for forward modeling, SPM can now benefit from integration of OpenMEEG. Moreover, the binary format used by OpenMEEG is that of Matlab, by use of the opensource MATIO library.
A sample dataset for M/EEG forward modeling can be downloaded from http://openmeeg.gforge.inria.fr. The sample dataset is provided with scripts that can be run to compute MEG and EEG leadfields on a realistic 3layer model.
Documentation
A tutorial for OpenMEEG is available on the web site, and it is briefly summarised in Additional file 1.The tutorial describes the lowlevel interface and details the different steps to be followed when computing M/EEG lead fields (Figure 9). Developper documentation can be generated via doxygen http://www.stack.nl/~dimitri/doxygen/.
Conclusion
In this paper, the OpenMEEG software project has been detailed, from the mathematical grounds of the symmetric BEM to more practical aspects.
The relevance of the OpenMEEG solver for quasistatic bioelectromagnetics has been demonstrated by a benchmark incorporating many alternative solvers in the context of M/EEG forward modeling. According to the results of this simulation study, OpenMEEG outperforms all the alternative solvers tested. By providing stateoftheart solutions for both EEG and MEG forward problems, OpenMEEG enables the combined use of these two complementary modalities.
It should be mentioned that OpenMEEG is being used for many problems in the field of quasistatic bioelectromagnetics, including Electrical Impedance Tomography, Intracranial electric potentials, Functional Electrical Stimulation and Cortical Mapping. This wide range of application domains, as well as its integration into highlevel languages make OpenMEEG unique and particularly valuable for basic and clinical research purposes.
Appendix
The symmetric BEM
In early numerical experiments to compare a Boundary Element and a Finite Element Method (FEM) for forward electroencephalography, we found a superior accuracy of the FEM [21]. This triggered a quest to improve the precision of Boundary Element Methods and led us to study the extended Green representation theorem [6]. We proposed a common formalism for the integral formulations of the forward EEG problem, and derived three different Boundary Element Methods within the same framework [7]. In this section we recall the mathematical background of Boundary Element Methods, and present both the doublelayer BEM, which is the most widespread method, and the symmetric BEM, which is a new formulation.
Green Representation
A fundamental result in potential theory shows that a harmonic function (i.e., such that Δu = 0) is uniquely determined within a domain Ω from its value on the boundary ∂Ω (Dirichlet condition), or the value of its normal derivative (Neumann condition). The Green Representation Theorem gives an explicit representation of a piecewiseharmonic function as a combination of boundary integrals of its jumps and the jumps of its normal derivative across interfaces. Before stating this theorem, some notation must be defined.

The restriction of a function f to a surface S _{ j }is indicated by f _{ sj }.

The functions ${f}_{{S}_{j}}^{}$ and ${f}_{{S}_{j}}^{+}$ represent the interior and exterior limits of f on S _{ j }:$\text{for}r\in {S}_{j},{f}_{{S}_{j}}^{\pm}(r)=\underset{\alpha \to {0}^{\pm}}{\mathrm{lim}}f(r+\alpha n).$

The jump of a function f across S _{ j }is denoted by:${[f]}_{{S}_{j}}={f}_{{S}_{j}}^{}{f}_{{S}_{j}}^{+},$

∂_{ n } V = n·∇V denotes the partial derivative of V in the direction of a unit vector n,

The function $G(r)=\frac{1}{4\pi \Vert r\Vert}$ is the fundamental solution of the Laplacian in ℝ^{3}, such that ΔG = δ _{0}.
The notation is simplified by introducing two integral operators, which map a scalar function f on ∂Ω to another scalar function on ∂Ω: the "doublelayer" operator $(\mathcal{D}f)(r)={\displaystyle \underset{\partial \Omega}{\int}{\partial}_{{n}^{\prime}}G(r{r}^{\prime})}f({r}^{\prime})\text{d}s({r}^{\prime})$, and the "singlelayer" operators $(\mathcal{S}f)(r)={\displaystyle \underset{\partial \Omega}{\int}G(r{r}^{\prime})}f({r}^{\prime})\text{d}s({r}^{\prime})$. For an operator $\mathcal{D}$, its restriction ${\mathcal{D}}_{ij}$ maps a function of S _{ j } to a function of S _{ i } .
where ${\mathcal{D}}^{\ast}$ is the adjoint of the operator $\mathcal{D}$. The Geselowitz formula exploits only the first boundary integral representation equation (5), while it is possible to exploit both (5) and (6). Thus three Boundary Element Methods can be derived within a unified setting: a BEM involving only singlelayer potentials, a BEM involving only doublelayer potentials, and a symmetric BEM combining single and doublelayer potentials [7]. We concentrate hereforth on the doublelayer and on the symmetric BEMs.
The doublelayer BEM
The above formula was established by Geselowitz [5], and was the only one used to model electroencephalography or electrocardiography, until recently, when [7] showed the diversity of BEMs that can be derived. This classical BEM is called a doublelayer BEM because it only involves the doublelayer operator $\mathcal{D}$.
The symmetric BEM
The originality of the symmetric Boundary Element Method is to consider a different piecewise harmonic function for each domain: the function u Ω _{ i } equal to $V\frac{{v}_{{\Omega}_{i}}}{{\sigma}_{i}}$ within Ω _{ i } and to $\frac{{v}_{{\Omega}_{i}}}{{\sigma}_{i}}$ outside of Ω _{ i } . This ${u}_{{\Omega}_{i}}$ is indeed harmonic in ℝ^{3}\∂Ω _{ i } , and the representation equations (5) and (6) can be applied, leading to a system of integral equations involving two types of unknowns: the potential V _{ i } and the normal current (σ∂_{ n } V) _{ i } on each interface.
where ${\hat{\sigma}}_{ij}$ (resp.${\stackrel{\u2323}{\sigma}}_{ij}$) is defined as σ _{ i } + σ _{ j } (resp. ${\sigma}_{i}^{1}+{\sigma}_{j}^{1}$) and where b _{1} and c _{1} are the coefficients of the P0 (resp. P1) boundary element decomposition of the source term ${\partial}_{n}{v}_{{\Omega}_{1}}$ (resp. ${\sigma}_{1}^{1}{v}_{{\Omega}_{1}}$).
Blocks N _{ ij } and D _{ ij } map a potential V _{ j } on S _{ j } to a function defined on S _{ i } . Block S _{ ij } maps a normal current p _{ j } on S _{ j } to a function defined on S _{ i } . The resulting matrix is blockdiagonal, and symmetric, hence the name "symmetric BEM".
The magnetic field is computed from the electric field and the primary source distribution using the Biot and Savart equation (3), as proposed by Ferguson, Zhang and Stroink [11].
In summary, the symmetric BEM introduces an additional unknown into the problem. the normal current, and uses an additional set of representation equations linking the normal current and the potential. The symmetric BEM departs from the doublelayer BEM in several ways:

the normal current to each surface is explicitely modeled;

only the surfaces which bound a common compartment have an interaction (whence the blocks of zeros in the matrix);

only the surfaces which bound a compartment containing sources have a source term;

the matrix is symmetric;

the matrix is larger for a given head model.
Declarations
Acknowledgements
The authors acknowledge support from the ANR grant ViMAGINE ANR08BLAN025002, and for EO from a PhD grant from the Regional Council of Provence Alpes Cote d'Azur. Help is gratefully acknowledged from C. Micheli and R. Oostenveld regarding the integration of OpenMEEG within Fieldtrip, and from A. Anwander, M. Dannhauer and C. Wolters regarding the use of Simbio.
Authors’ Affiliations
References
 Soleimani M, Shipley R, Smith N, Mitchell C: Medical imaging and physiological modelling. linking physics and biology. BioMedical Engineering OnLine 2009, 8: 1. 10.1186/1475925X81View ArticleGoogle Scholar
 de Munck J, Peters M: A fast method to compute the potential in the multisphere model. IEEE Trans on Biomed Engin 1993, 40(11):1163–1174.View ArticleGoogle Scholar
 Kariotou F: Electroencephalography in ellipsoidal geometry. J Math Anal Appl 2004, 290: 324–42. 10.1016/j.jmaa.2003.09.066MathSciNetView ArticleGoogle Scholar
 Sarvas J: Basic mathematical and electromagnetic concepts of the biomagnetic inverse problem. Phys Med Biol 1987, 32: 11–22. 10.1088/00319155/32/1/004View ArticleGoogle Scholar
 Geselowitz DB: On bioelectric potentials in an homogeneous volume conductor. Biophysics Journal 1967, 7: 1–11. 10.1016/S00063495(67)865718View ArticleGoogle Scholar
 Nédélec JC: Acoustic and Electromagnetic Equations. Springer Verlag; 2001.View ArticleGoogle Scholar
 Kybic J, Clerc M, Abboud T, Faugeras O, Keriven R, Papadopoulo T: A Common Formalism for the Integral Formulations of the Forward EEG Problem. IEEE Transactions on Medical Imaging 2005, 24: 12–28. 10.1109/TMI.2004.837363View ArticleGoogle Scholar
 Voth EJ: The inverse problem of electrocardiography. industrial solutions and simulations. IJBEM 2005., 7(2):Google Scholar
 Kybic J, Clerc M, Faugeras O, Keriven R, Papadopoulo T: Generalized head models for MEG/EEG. boundary element method beyond nested volumes. Physics in Medicine and Biology 2006, 51: 1333–1346. 10.1088/00319155/51/5/021View ArticleGoogle Scholar
 Hämäläinen M, Hari R, IImoniemi RJ, Knuutila J, Lounasmaa OV: Magnetoencephalographytheory, instrumentation, and applications to noninvasive studies of the working human brain. Reviews of Modern Physics 1993, 65(2):413–497. 10.1103/RevModPhys.65.413View ArticleGoogle Scholar
 Ferguson AS, Zhang X, Stroink G: A Complete Linear Discretization for Calculating the Magnetic Field Using the Boundary Element Method. IEEE Trans Biomed Eng 1994, 41(5):455–459. 10.1109/10.293220View ArticleGoogle Scholar
 Goncalves S, de Munck J, Verbunt J, Heethaar R, Lopes da Silva F: In vivo measurement of the brain and skull resistivities using an EITbased method and the combined analysis of SEF/SEP data. IEEE Transactions on Biomedical Engineering 2003, 50(9):1124–1128. 10.1109/TBME.2003.816072View ArticleGoogle Scholar
 Clerc M, Badier JM, Adde G, Kybic J, Papadopoulo T: Boundary Element formulation for Electrical Impedance Tomography. ESAIM. Proceedings, EDP Sciences 2005, 14: 63–71.MathSciNetGoogle Scholar
 Clerc M, Adde G, Kybic J, Papadopoulo T, Badier JM: In vivo conductivity estimation with symmetric boundary elements. In International Journal of Bioelectromagnetism Edited by: Malmivuo J. 2005, 7: 307–310.Google Scholar
 Jacquir S, Fruitet J, Guiraud D, Clerc M: Computation of the electrical potential inside the nerve induced by an electrical stimulus. EMBC 2007: Proceedings of the 29th Annual International Conference of the IEEE Engineering in Medicine and Biology Society 2007.Google Scholar
 Laforêt J, Guiraud D, Clerc M: A toolchain to simulate and investigate selective stimulation strategies for FES. EMBC 2009: Proceedings of the 30th Annual International Conference of the IEEE Engineering in Medicine and Biology Society 2009.Google Scholar
 He B, Wang Y, Wu D: Estimating Cortical Potentials from Scalp EEGs in a Realistically Shaped Inhomogeneous Head Model by Means of the Boundary Element Method. IEEE Transactions on Biomedical Engineering 1999, 46(10):1264–1268. 10.1109/10.790505View ArticleGoogle Scholar
 Clerc M, Kybic J: Cortical mapping by LaplaceCauchy transmission using a boundary element method. Inverse Problems 2007, 23: 2589–2601. 10.1088/02665611/23/6/020MathSciNetView ArticleGoogle Scholar
 de Munck JC: A linear discretization of the Volume Conductor Boundary Integral Equation using Analytically Integrated Elements. IEEE Trans Biomed Eng 1992, 39(9):986–990. 10.1109/10.256433View ArticleGoogle Scholar
 Bonnet M: Boundary Integral Equations Methods for Solids and Fluids. John Wiley and Sons; 1999.Google Scholar
 Clerc M, Dervieux A, Faugeras O, Keriven R, Kybic J, Papadopoulo T: Comparison of BEM and FEM methods for the E/MEG problem. Proceedings of BIOMAG Conference 2002.Google Scholar
 Oostendorp T, van Oosterom A: Source parameter estimation in inhomogeneous volume conductors of arbitrary shape. IEEE Trans Bioned Engin 1989, 36: 382–391. 10.1109/10.19859View ArticleGoogle Scholar
 Phillips C: Source estimation in EEG. PhD thesis. Université de Liège; 2000.Google Scholar
 Hämäläinen MS, Sarvas J: Realistic Conductivity Geometry Model of the Human Head for Interpretation of Neuromagnetic Data. IEEE Trans Biomed Eng 1989, 36(2):165–171. 10.1109/10.16463View ArticleGoogle Scholar
 Stenroos M, Mäntynen V, Nenonen J: A Matlab library for solving quasistatic volume conduction problems using the boundary element method. Comput Methods Prog Biomed 2007, 88(3):256–263. 10.1016/j.cmpb.2007.09.004View ArticleGoogle Scholar
 AkalinAcar Z, Gençer NG: An advanced boundary element method (BEM) implementation for the forward problem of electromagnetic source imaging. Physics in medicine and biology 2004, 49(21):5011–28. 10.1088/00319155/49/21/012View ArticleGoogle Scholar
 Zanow F, Knosche T: ASAAdvanced Source Analysis of Continuous and EventRelated EEG/MEG Signals. Brain Topography 2004, 16(4):287–290. 10.1023/B:BRAT.0000032867.41555.d0View ArticleGoogle Scholar
 Parker S, Johnson C: SCIRun: A Scientific Programming Environment for Computational Steering. Proceedings of the IEEE/ACM Supercomputing Conference, IEEE 1995, 52.View ArticleGoogle Scholar
 Meijs JWH, Weier OW, Peters MJ, van Oosterom A: On the numerical accuracy of the boundary element method. IEEE Trans Biomed Eng 1989, 36: 1038–1049. 10.1109/10.40805View ArticleGoogle Scholar
 Van Uitert R, Johnson C, Zhukov L: Influence of head tissue conductivity in forward and inverse magnetoencephalographic simulations using realistic head models. IEEE Trans on Biomed Engin 2004, 51(12):2129–2137. 10.1109/TBME.2004.836490View ArticleGoogle Scholar
 Nolte G: The magnetic lead field theorem in the quasistatic approximation and its use for magnetoencephalography forward calculation in realistic volume conductors. Physics in Medicine and Biology 2003, 48: 3637–3652. 10.1088/00319155/48/22/002View ArticleGoogle Scholar
Copyright
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.