Magnetic stimulation for non-homogeneous biological structures

Background Magnetic stimulation has gained relatively wide application in studying nervous system structures. This technology has the advantage of reduced excitation of sensory nerve endings, and hence results in quasi-painless action. It has become clinically accepted modality for brain stimulation. However, theoretical and practical solutions for assessment of induced current distribution need more detailed and accurate consideration. Some possible analyses are proposed for distribution of the current induced from excitation current contours of different shape and disposition. Relatively non-difficult solutions are shown, applicable for two- and three-dimensional analysis. Methods The boundary conditions for field analysis by the internal Dirichlet problem are introduced, based on the vector potential field excited by external current coils. The feedback from the induced eddy currents is neglected. Finite element modeling is applied for obtaining the electromagnetic fields distribution in a non-homogeneous domain. Results The distributions were obtained in a non-homogeneous structure comprised of homogeneous layers. A tendency was found of the induced currents to follow paths in lower resistivity layers, deviating from the expected theoretical course for a homogeneous domain. Current density concentrations occur at the boundary between layers, suggesting the possibility for focusing on, or predicting of, a zone of stimulation. Conclusion The theoretical basis and simplified approach for generation of 3D FEM networks for magnetic stimulation analysis are presented, applicable in non-homogeneous and non-linear media. The inconveniences of introducing external excitation currents are avoided. Thus, the possibilities are improved for analysis of distributions induced by time-varying currents from contours of various geometry and position with respect to the medium.


Background
The theory of electromagnetic fields can be considered as quite well developed. Nevertheless, some problems could be pointed-out in connection with its application to stimulation of biological structures: -distribution and spread of the time-varying magnetic field in non-linear and inhomogeneous domains of complex geometry; -the interaction of electromagnetic pulses with biophysical and physiological processes; -stationary or non-stationary processes of random and stray factors type, e.g. changes in the stimulated structure; -technical and technological aspects, e.g. coils shape and disposition with respect to the stimulated structures.
There may be other practical and theoretical problems connected with magnetic stimulation. However, in different medical applications, various approximate or heuristic solutions were found [1].
Published studies on electromagnetic field stimulation most often are based on distribution and propagation in homogeneous media. This approach could be correct for stationary magnetic field analysis, as biological tissues relative magnetic permeability is m r approximately equal to 1. In cases of time-varying magnetic field, feedback from the induced currents may be present. However, when propagation takes place in a weakly conductive medium, the intensity of the feedback-induced magnetic field would be several orders of magnitude lower than the external field. The induced field can then be ignored, and the medium considered to be homogeneous with respect to the external magnetic field.
On the other hand, the condition of continuity of the current lines requires taking into account all regions with different specific conductivities. An exception here is the case where the induced currents propagate tangentially to the layers in a non-homogeneous medium.
In some studies, the presence of non-homogeneities has been taken into account. Transcranial magnetic stimulation, where the scalp, skull, and brain represent a layered non-homogeneous structure, is a typical example [2,3]. A theoretical approach, applying a scalar magnetic potential in addition to the vector potential has been proposed as an alternative to finite-element modeling (FEM), but has considered only a simplified basic example [4]. FEM has been used for assessment of magnetically-induced currents for heart stimulation, with inclusion of the non-homogeneities of a canine thorax [5]. Field magnitudes for efficient myocardial stimulation were computed. Magnetic flux density in the upper arm, taking into account the skin and muscle layers and the bone, have been assessed for different coil designs [6]. A study of stimulation in the spinal region [7] has shown a tendency for current density concentration about the nerve roots. However, acceptable accuracy in evaluation of the induced currents values could only be attained when these currents were assessed on the basis of an appropriately obtained magnetic field distribution. A tendency of induced current concentrations at the interfaces between structures of strongly differing conductivity could be observed [3,7].
The present work investigates some possibilities for modeling and theoretical analysis of induced current distribution during magnetic stimulation, with excitation current contours of different shape and disposition. Relatively non-difficult solution approaches are proposed, applicable for two-dimensional and three-dimensional problems of field distribution in non-homogeneous media.

Initial conditions
Electrical current contours are positioned in the air above a partially homogeneous conductive medium that consists of homogeneous layers (Fig. 1), simulating part of the human body. The mean values of the electromagnetic parameters m (magnetic permeability) and g (specific conductivity) of the corresponding layers are assumed to be known. Each layer is taken as linear and isotropic. Effects of random factors are not considered.
The electrical currents in the several contours i k (t); k = 1,2,...,n are generated by an electromagnetic stimulator and create an alternating magnetic field. The distribution of the magnetic field in the examined region is described by the Maxwell equations in the form: where, besides the field vectors and , the electromagnetic potentials and V also take part. The last term in Eq. (1) represents the independent external excitation currents i k (t), generated by a current source.
It is well known that Eqs. (1) and (2)  Eq. (3) is valid for an arbitrary medium. Generally, the material constants in case of non-homogeneity should be considered as functions of the spacial co-ordinates and in anisotropic media -as tensor quantities. In cases of nonlinear media, the material constants are functions of the field vectors.
For homogeneous linear and isotropic medium, using the Lorentz gauge [9] gives: and the non-homogeneous diffusion equation is obtained: It is known that in harmonic mode, this equation can be reduced to the Helmholtz equation: rotH .
Equations (5) and (6) are also applicable for partially homogeneous media, as depicted in Fig. 1. In the internal part of the non-homogeneous region W, where , Eq. (5) and (6)  The solution for stationary mode and a homogeneous medium with respect to m, is: For current contours of relatively small cross section and for alternate current, a known formula is obtained: For the case shown in Fig. 1:

Figure 3
Magnetic vector-potential distribution, with potentials applied at the XOZ nodes. The excitation current I = 1000 A, f = 10 kHz flows in horizontally (a) or vertically (b) positioned coils. ; .  1) and for the external space. The non-homogeneous medium W, with specific conductivities g 1 ¹ g 2 ¹ ... ¹ g p ¹ ... ¹ g m , can be of layered structure, or include internal sub-spaces W p . It is evident, that the feedback from the eddy currents induced in the non-homogeneous medium W can be neglected, compared to the direct action of the external excitation.

Induced currents in the medium
The induced field and the eddy currents in the W-medium are obtained by the procedure described below.
The values of the vector-potential on the external surface of W are computed from Eq. (10) in all nodes necessary for the solution of the internal Dirichlet problem.
Using the homogeneous equation for the diffusion, relating to the magnetic vector-potential: or, for harmonic mode, the Helmholtz equation:

Figure 5
Distribution of the modulus of the magnetic field of the induced currents, excited by the horizontally (a) and vertically (b) positioned coils. Section and viewing angle -as in Fig. 4.

Figure 6
Distribution of the modulus of the current density vector in the homogeneous (a) and non-homogeneous (b) domains from a horizontally positioned coil. Section and viewing angle -as in Fig. 4,5.
the internal Dirichlet problem is solved. The boundary conditions for the vector-potential are defined by Eq. (10) for the boundary surface S W of the space W, subjected to electromagnetic stimulation. The feedback from the induced field was ignored, as discussed above. The respective conditions are valid for refraction of the field vector lines at the boundaries of two homogeneous regions [8,9].
The solution of Eq. (10), having in view the boundary conditions, does not present theoretical or computational difficulties. The problem of field distribution in the nonhomogeneous space W is solved by the finite element method (FEM) or by the boundary element method (BEM). Using appropriate software, the field and the currents induced by the external field in a weakly-conductive medium are obtained, together with the magnetic field of the induced currents.

Methods
In order to verify the applicability of the approach described in Section 2 above, a 3D spacial model consisting of five transverse layers was constructed, as a very approximate simulation of a biological object. The initial model consisted of 19,218 nodes and 104,693 eight-node brick elements. Two versions of the model were considered, al-lowing assignment of different resistivities to the respective layers, as follows: a) homogeneous medium W a of specific resistivity r = 20 Wm (Fig. 2a); the value was chosen as a rough average of very high and very low conductivities (soft tissue and blood; and bone and lung air, respectively). b) non-homogeneous medium W b containing isotrope layers of specific resistivities (Fig. 2b), chosen to correspond to body regions [10,11], as follows:r 1 = 20 Wmskin, r 2 ,r 4 = 7 Wm -skeletal muscle, r 3 = 160 Wm -bone, r 5 = 20 Wm -'buffer' zone needed to allow setting of adequate boundary conditions; for this test example, the layers were made of equal thickness, except for the buffer zone.
The external magnetic field was excited by a squareshaped coil. Two examples were considered initially (Fig.  2), with coils placed symmetrically on the upper surface of the 3D model: i) coil in horizontal position (XOZ plane), and ii) coil in vertical position (YOZ plane). The excitation current was arbitrarily chosen to be I = 1000 A rms , with a frequency f =10 kHz. The calculated values of the magnetic-vector-potentials were applied on the upper surface nodes of the model (1011 nodes), as magnetic boundary loads. Zero boundary conditions were adopted for the other surfaces. Harmonic magnetic analysis was performed to obtain the eddy currents distribution in the medium.

Figure 7
Distribution of the modulus of the current density vector in the homogeneous (a) and non-homogeneous (b) domains from a vertically positioned coil. Section and viewing angle -as in Fig. 4,5,6.

Results
The magnetic vector-potential distributions, induced by the horizontal and vertical current contours on the upper XOZ surface of the 3D model, are shown respectively in Fig. 3a,3b. The excited magnetic fields, viewed in a central XOY section and on the XOZ plane, are presented in Fig.  4a,4b for the corresponding current coil dispositions. Due to the homogeneous character of the layered structure with respect to the magnetic permeability m, and the neglected induced currents feedback, the external magnetic field distribution is not influenced by the non-homogeneities of the domain with respect to the specific conductivity g.
The magnetic fields excited by the induced currents in a non-homogeneous medium, from a horizontal and a vertical coil dispositions are shown in Figs. 5a,5b. The distribution of the eddy currents induced by a horizontal current contour in the homogeneous and non-homogeneous domains are shown in Fig. 6a,6b, respectively. The same distributions for the vertical coil are presented in Fig.  7a,7b.
It was found that the current lines tend to form loops in the low resistivity layers, in addition to running in opposition to the excitation current, as would be in an homogeneous medium. In order to better observe the induced current lines across layers of different resistivities, the non-homogeneity is modeled by vertical layers and the

Figure 8
Current density vector distribution obtained from a horizontally positioned (XOZ plane) coil in a vertically-sliced non-homogeneous volume.
coil position is on the upper surface of the medium, in the XOZ plane (Fig. 8). The effect of more complicated current contours can also be studied. The distribution of the induced vector-potential by the fan-like assembly (slinky coils, [12]), viewed from the upper XOZ surface, is shown in Fig. 9. The corresponding eddy currents distributions in the homogeneous and non-homogeneous media are shown in Fig. 10a,10b.

Discussion
The approach proposed in Section 2.2 considerably reduces the problems connected with generation of 3D FEM networks. However, this approach preserves the power for analysis in non-homogeneous and, in non-linear media. The inconveniences of introducing external excitation currents are avoided. Thus, the possibilities are improved for analysis of field distributions induced by different time-varying currents, with contours of various geometry and position with respect to the stimulated object.
A systematic error is introduced due to the ignored feedback from the eddy currents. However, this error is negligible, considering the low intensity of the eddy currents in the weakly conductive medium. These currents have their own magnetic field, which is several orders of magnitude lower than the external magnetic field from the stimulation coils, as can be seen in comparing the distributions in Fig. 5 and Fig. 4.
There exist theoretical possibilities by using an iteration procedure [e.g. [13]], for taking into account the feedback. Such a technique has been used in the development and practical application of numerical methods for field anal-

Figure 9
Magnetic vector-potential distribution on the domain surface (XOZ), excited by rectangular coils in a fan-like assembly (slinky coils).
ysis. However, in the case of weakly conductive media, the iteration just stops at the first cycle.
The different field structures in homogeneous and nonhomogeneous media are evident, as well as the relative difference in the current density amplitude values. The latter reach about 15% for horizontal and vertical contours. However, the difference is much more pronounced in some non-homogeneous regions (Figs. 6 and 7). It can be seen that the non-homogeneous structure has an essential influence on the induced currents distribution and on the differences in the current density vector.
An interesting event can be observed in relation to the distribution of the current density. The current lines follow low resistivity paths and in addition to crossing from layer to layer, flowing in opposition to the excitation current, they also rotate inside the low resistivity layer. In order to better illustrate this specific distribution, an example with vertical layers is shown in Fig. 8. Current lines forming loops in the vertical layers can be observed.
Analysis of the results shows, that the method and procedure adopted, introducing boundary conditions in the non-homogeneous domain, is workable. This method could also be applied for other coil shapes and domain structures, without any essential difficulties. Two examples with slinky coils ( Fig. 9 and 10) demonstrate the versatility of this method.
The fields resulting from the external coils and the fields created by the induced currents differ by many orders of magnitude (six in our examples), due to the relatively low conductivities of the medium and to the moderately high frequencies. Therefore, ignoring the feedback effects with respect to the introduced boundary condition in solving the internal Dirichlet problem, is justified.
The boundary conditions, related to the components A x ,A y ,A z of the magnetic vector-potential, can be considered invariant with respect to the co-ordinates.
Another point of interest is that a concentration of induced current densities can be observed at the interface (boundary) between structures of differing conductivities. Hence, a kind of concentration or selective stimulation (focusing) could be predicted or looked for, in structures with relatively well-defined layers of different resistivity (Figs. 6b, 7b and 10b). The concentration in question can be explained by the appearance of surface charges on these boundaries.
It should be reminded here, that the results were obtained for sinusoidal excitation currents. The application of the method for other types of time-varying currents could be addressed by adequate modeling software. The same was valid for studying field distributions in weakly conducting non-linear media.

Figure 10
Current density vectors distribution, excited by rectangular coils in a fan-like assembly in the homogeneous (a) and non-homogeneous (b) domains. Section and viewing angle as in Figs. 4,5,6,7.

Conclusion
The results obtained in this work do not provide data that would change the theory or practice of magnetic stimulation. However, our findings can be used to enhance concepts related to electromagnetic field distribution in nonhomogeneous structures.
For magnetic stimulation of biological targets embedded in electrically non-homogeneous media, it is helpful to decide upon external system geometry and parameters by considering an accurately obtained distribution of the induced current fields.
The complex interaction between external magnetic field and physiological processes should not deter the search for improved methods for more accurate assessment of the induced currents.

Authors' contributions
SP was responsible for the theoretical basis of the investigation. VK solved the problems of adequate FEM modeling and carried-out the generation and presentation of various field and current distributions. ID was responsible for setting and formulation of the problem and for the biomedical engineering aspect of the study. All authors read and approved the final manuscript.