Dosimetry of Various Human Bodies Exposed to Microwave Broadband Electromagnetic Pulses

In this paper, human exposures to ultra-wideband (UWB) electromagnetic (EM) pulses in the microwave region are assessed using a frequency-dependent FDTD scheme previously proposed by the authors. Complex permittivity functions of all biological tissues used in the numerical analyses are accurately expressed by the four-term Cole–Cole model. In our method, we apply the fast inverse Laplace transform to determine the time-domain impulse response, utilize the Prony method to find the Z-domain representation, and extract residues and poles for use in the FDTD formulation. Update equations for the electric field are then derived via the Z-transformation. Firstly, we perform reflection and transmission analyses of a multilayer composed of six different biological tissues and then confirm the validity of the proposed method by comparison with analytical results. Finally, numerical dosimetry of various human bodies exposed to EM pulses from the front in the microwave frequency range is performed, and the specific energy absorption is evaluated and compared with that prescribed in international guidelines.

Introduction

In 2002, the Federal Communications Commission of the United States of America issued a ruling allowing the use of ultra-wideband (UWB) electromagnetic (EM) pulses in the frequency range between 3.1 and 10.6 GHz (1). Since then, numerous applications of UWB pulses have emerged as a result of the regulation, such as vital signal detection and locating moving human bodies, as well as ground-penetrating radar, remote sensing, non-destructive inspection, and so forth (2). Owing to many advantages of UWB pulses such as low power consumption and immunity to multipath of EM propagations and EM interferences, the widespread use of EM pulses is expected to continue. Some of these applications use UWB pulses in the vicinity of human bodies, such as wireless capsule endoscopy using broadband EM pulses for on/in-body communications (3, 4). These applications have led to public concern about the effect on the broadband EM pulses to human body.

Meanwhile, biological effects due to EM pulses have been numerically and experimentally investigated. The effects include the microwave hearing effect, microwave heating effect, and electroporation (5–9). Consequently, international organizations have prescribed exposure limits for the temporal peak of specific energy absorption (SA) in published guidelines to prevent adverse effects, particularly of microwave hearing, which is considered an acute biological effect (10, 11). The International Commission on Non-ionizing Radio Protection (ICNIRP) provided an SA limit of 2 mJ/kg in an arbitrary 10 g-averaged tissue for a single pulse illumination (10), while the SA limit is up to 576 J/kg for continuous exposure of 6 min in the regulation defined by the Institute of Electrical and Electronics Engineers (IEEE) (11). Recently, the ICNIRP has revised the Guidelines based on the relevant scientific knowledge and published them in 2020. In the guidelines, it is mentioned that there is no evidence that microwave hearing in any realistic exposure scenarios causes adverse health consequence and the microwave hearing effect is not considered in the guideline but it is still mandate to consider the heating effect from the pulse exposures (12). Although many biological effects due to EM pulses have been experimentally confirmed, there are few studies providing detailed exposure levels or showing the distribution of SA inside a human body.

To derive SA inside the human body, calculation of the interactions between EM pulses and biological bodies is necessary. In the earliest studies, most of the biological targets were objects having simple shapes such as a multilayer or a dielectric sphere, inside which SA or the induced electric field was derived theoretically (13, 14). However, there has been no detailed dosimetric information of the detailed human body exposed to EM pulses due to difficulties in the calculation of SA or the induced electric field inside biological bodies. These problems are mainly attributed to the frequency dependence of the dielectric properties of biological tissues, which are expressed by the four-term Cole-Cole model (15).

To perform numerical dosimetry of EM pulses, we need to consider the frequency dependence of the permittivity and conductivity of biological tissues over a broad frequency range. Many frequency-dependent finite-difference time-domain [(FD)²TD] approaches have been proposed, such as recursive convolution method (16), piecewise linear recursive convolution method (17), trapezoidal recursive convolution method (18), auxiliary differential equation method (19–21), and Z-transform method (22). However, these approaches have only been applied to materials having complex permittivity expressed by relatively simple models such as the Debye and Lorentz models. These approaches are not applicable to the Cole–Cole function, due to difficulties in finding the exact time-domain solution of a fractional-order differential equation. Nevertheless, many attempts have been made to address this problem, including those using the Riemann–Liouville theory to find the time-domain solution of the model (23, 24). Recently, our research group has proposed an FDTD formulation for analyses of arbitrary frequency-dependent materials via the use of the fast inverse Laplace transform (FILT) and the Prony method (25). The proposed method has also been extended to three-dimensional analyses of UWB antennas in the vicinity of the human body and the dosimetry of EM pulses incident to a human head (26, 27).

In this study, we extend our numerical models to whole-body human models which are exposed to broadband EM pulses. To the best of our knowledge, this is the first ever numerical dosimetry of EM pulses, and detailed information of SA inside the human body is provided for compliance with the SA limit prescribed in international guidelines. The proposed method involves two steps; firstly, we apply the FILT to permittivity functions in the complex frequency domain to transform them into time-domain impulse responses. Then, the Prony method is used to extract the model parameters and to determine expressions for the infinite impulse response (IIR) expressions in the Z-domain. The update equations of electric field are derived via the Z-transform.

This paper is outlined as follows. The proposed (FD)²TD formulation and the calculation of the update coefficients for the electric field are described in section (FD)²TD Formulation Using FILT and Prony Method. The validity of the method in calculating SA and internal electric field strength (IEFS) inside a multilayer model of biological tissues is demonstrated via comparison with the theoretical results in section Transmission Characteristics of EM Pulses Into Biological Bodies. Numerical dosimetry of anatomically detailed human body models exposed to UWB EM pulses is performed and physical quantities such as SA and IEFS are quantitatively derived and compared with those prescribed in the guidelines in section Transmission Characteristics of EM Pulses Into Biological Bodies. Finally, conclusions are drawn in section Conclusion.

(FD)²TD Formulation Using FILT and Prony Method

Methodology

In this study, all media used in the numerical analyses are biological tissues having complex relative permittivity expressed by the four-term Cole–Cole function as

$\begin{array}{ll}{\epsilon }_m(\omega )={\epsilon }_{\infty }+\frac{\sigma }{j\omega {\epsilon }_0}+{\displaystyle \sum _{q=1}^4}\frac{\Delta {\chi }_q}{1+{(j\omega {\tau }_q)}^{1-{\alpha }_q}},& (1)\end{array}$

where, ω, ε_∞, and σ are the angular frequency [rad/s], relative permittivity and conductivity [S/m] of a biological medium at infinite frequency, respectively. ε₀ is the free-space permittivity and Δχ_q represents the change in relative permittivity in the qth relaxation term. τ_q and α_q are the relaxation time and a parameter determining the broadness of the qth term, respectively. All parameters in (1) can be found in Gabriel's database of dielectric properties for biological tissues (15). Although Gabriel's permittivity data are de facto, it is noteworthy that different Cole–Cole parameters may be derived, depending on the method used in fitting the measurement data of the dielectric properties. By limiting the frequency range to between 1 MHz and 20 GHz, the number of Cole–Cole terms may be reduced from four terms to two terms while providing the best fit to the measurement data (28). The average deviations from the measurement results over a frequency range between 1 MHz and 20 GHz are higher than 20% for both relative permittivity and loss factor of some biological tissues but they are shown to be <15% for the two Cole-Cole terms used in Kensuke et al. (28).

Since the (FD)²TD formulation for biological bodies has been described in the literature (26), we hereby only show the numerical procedures of the proposed method and the derived update equation for electric fields. Firstly, we apply the FILT to electric susceptibility represented by the Cole–Cole model, with the relaxation time normalized by the time step interval used in the FDTD simulations and obtain impulse responses for each susceptibility term in the time domain. Then, we use the Prony method to transform the time-domain impulse response into that in the Z-domain. The permittivity is now expressed in the Z-domain as

$\begin{array}{ll}{\epsilon }_m(z)={\epsilon }_{\infty }+\frac{1+z^{-1}}{1-z^{-1}}\frac{\sigma \Delta t}{2{\epsilon }_0}+{\displaystyle \sum _{q=1}^4}{\displaystyle \sum _{l=1}^{L_q}}\frac{A_l^{(q)}}{1-p_l^{(q)}{z}^{-1}},& (2)\end{array}$

where Δt and L_q are the time step interval and number of poles for the qth Cole–Cole term, respectively. $A_l^{(q)}$ and $p_l^{(q)}$ are residues and poles, respectively. Note that the second term of the right hand side is obtained by applying the bilinear approximation, i.e., jω ≈ s = 2(1 − z⁻¹)/ (1 + z⁻¹)/ (Δt). Since the nested summation in the third term of the right-hand side of (2) can be merged into a single summation, (2) can simply be expressed as

$\begin{array}{ll}{\epsilon }_m(z)={\epsilon }_{\infty }+\frac{1+z^{-1}}{1-z^{-1}}\frac{\sigma \Delta t}{2{\epsilon }_0}+{\displaystyle \sum _{k=1}^N}\frac{A_k}{1-p_k{z}^{-1}},& (3)\end{array}$

where N is the total number of Debye terms, i.e., N = L₁ + L₂ + L₃ + L₄. Procedures for determining N, A_k, and p_k using the Prony method will be described in subsection Calculation of Specific Energy Absorption and Internal Electric Field Strength and can also be found in the literature (25). Substituting (3) into the discrete constitutional relation of Maxwell's equations, we obtain the update equation for the electric field as

$\begin{array}{ll}{E}^n=\frac{1}{L_0}\left[\frac{D^n}{{\epsilon }_0}-\frac{\sigma \Delta t}{2{\epsilon }_0}{E}^{n-1}-I^{n-1}-{\displaystyle \sum _{k=1}^N}{p}_k{P}_k^{n-1}\right],& (4)\end{array}$

where,

$\begin{array}{llll}{L}_0& =& {\epsilon }_{\infty }+\frac{\sigma \Delta t}{2{\epsilon }_0}+{\displaystyle \sum _{k=1}^N}{A}_k& (5)\end{array}$

$\begin{array}{llll}{I}^n& =& I^{n-1}+\frac{\sigma \Delta t}{2{\epsilon }_0}(E^n+E^{n-1})& (6)\end{array}$

$\begin{array}{llll}{P}_k^n& =& p_k{P}_k^{n-1}+A_k{E}^n,& (7)\end{array}$

Iⁿ and ${\text{P}}_k^n$ are the auxiliary variables which are initialized by setting I⁰ = 0 and ${\text{P}}_k^0$ = 0, respectively. Equations (4), (6), and (7) are the update equations for the electric field, the auxiliary field used for considering the conductivity term, and the auxiliary field used for considering the Cole–Cole terms, respectively. The update equations for the electric flux density and magnetic field can be obtained by applying the central difference to Maxwell's equations similarly to those in conventional FDTD procedures (20).

The main advantages of the proposed method are that we can avoid the formulation of fractional-order differential equations by using the FILT and the Prony method and it is straightforward to implement the proposed method into the conventional FDTD code.

Calculation of Specific Energy Absorption and Internal Electric Field Strength

The exposure level inside biological bodies illuminated by broadband EM pulses can be evaluated using SA, which has been used as a metric in the guidelines. In the FDTD simulations, SA can be calculated using the following equation:

$\begin{array}{l}{\text{SA}}^{n-\frac{1}{2}}=\Delta t{\displaystyle \sum _{m=1}^n{\left(\frac{E\left(t\right)}{\rho }·\frac{\partial D\left(t\right)}{\partial t}\right)|}_{t=\left(m-\frac{1}{2}\right)\Delta t}},\\ =\frac{1}{2\rho }{\displaystyle \sum _{m=1}^n\left(E^m+E^{m-1}\right)·\left(D^m-D^{m-1}\right)},& (8)\end{array}$

where ρ is the density of the biological tissue. The electric flux density in Equation (7) at the nth time step (Dⁿ) is updated using the magnetic field at the (n+1/2)^th time step (H^n+1/2) and Eⁿ is updated using (4).

Since a broadband pulse is utilized in our simulations, we can also obtain numerical solutions of the electric and magnetic fields at each frequency component within a single run. The electric field at a frequency is determined via Fourier transform of the waveform obtained at an observation location as

$\begin{array}{ll}E\left(\omega \right)={\displaystyle {\int }_0^{T}E\left(t\right)e^{-j\omega t}dt}={\displaystyle \sum _{n=0}^{N^T}{E}^n{e}^{-j\omega n\Delta t}}\Delta t,& (9)\end{array}$

where N^T is the total number of time steps. After the electric field at each frequency is obtained, the specific absorption rate (SAR) is then calculated as follows:

$\begin{array}{ll}SAR\left(\omega \right)=\frac{\sigma {\left|E\left(\omega \right)\right|}^2}{2\rho }.& (10)\end{array}$

Calculation of Update Coefficients for Electric Field

The procedures for determining coefficients A_k and p_k in the update equation for the electric field are described as follows. First, we transform relative permittivity represented in the frequency domain into that in the complex frequency domain by replacing jω with the complex frequency s and apply the FILT to find impulse response of the permittivity in the time domain. Then, the Prony method is used to extract the residues A_k and poles p_k from the expression for the IIR in the Z-domain.

As an example, we apply the FILT and the Prony method to the permittivity functions of biological tissues “Fat” and “Gray Matter.” Each Cole–Cole parameter is taken from the Gabriel's database and listed in Table 1. The time step interval Δt used to normalize the relaxation time in the Cole–Cole function before applying the FILT is set to 1.668 ps.

TABLE 1

Table 1. Cole–Cole parameters for “Fat” and “Gray Matter” from Gabriel's database.

Table 2 shows the update coefficients A_k and p_k used in the FDTD calculations for “Fat” and “Gray Matter.” These values are directly obtained from the Prony method and the number of coefficients for each Cole–Cole term is truncated when the ratio of {|A_k| /max|A_k|} for k = 1, 2, …, N is less than a tolerance value of 10⁻³. The update coefficients A_k and p_k physically correspond to the initial amplitude and the decreasing ratio of the time-domain impulse response, respectively. From Table 2, the total numbers of coefficients are 15 and 29 for “Fat” and “Gray Matter,” respectively.

TABLE 2

Table 2. Update coefficients for “Fat” and “Gray Matter.”

To demonstrate the validity of the update coefficients, we calculate the reflection coefficients from each biological medium by one-dimensional FDTD simulation using the model shown as Figure 4 in Chakarothai et al. (26) and compare their values with those obtained from the EM theory. The analysis model is half filled with biological tissues and truncated with perfectly matched layers in order to absorb the outgoing wave. Numerical results using a time step interval of 1.668 ps and a resolution of 0.5 mm are shown in Figure 1. The reflection coefficients are analytically calculated using $\Gamma =|1-\sqrt{{\epsilon }_m}|/|1+\sqrt{{\epsilon }_m}|,$ where ε_m is the complex relative permittivity expressed by Equation (1).

FIGURE 1

Figure 1. Reflection coefficients of “Fat” and “Gray Matter”.

As shown in Figure 1, the reflection coefficients of “Fat” and “Gray Matter” obtained via numerical simulations are within 2% of those obtained by the analytical method over a broad frequency range between 50 MHz and 10 GHz, demonstrating the validity of the update coefficients and our numerical approach. Update coefficients for the other type of biological tissues can also be calculated straightforwardly using the procedures described above. Note that when we change the time step interval, we also need to recalculate the update coefficients.

Transmission Characteristics of EM Pulses Into Biological Bodies

Multilayer Model

Figure 2 shows a multilayer model mimicking a human head, which comprises six biological tissues, similar to those used in the literature (29, 30). Table 3 indicates the thicknesses of biological tissues used in the analysis model and the number of the update coefficients for each biological tissue. The total size of the multilayer model is 180 mm. These coefficients are obtained by applying the FILT and the Prony method with a time step interval of Δt = 1.668 ps. The resolution and the total number of cells used in our analysis model are 0.5 mm and 5,000, respectively. CPMLs with eight layers are utilized on both sides of the analysis domain to absorb the outgoing EM waves. The total number of time steps is 100,000. The incident electric field is given by a Gaussian pulse expressed as

$\begin{array}{ll}{E}^{inc}(t)=exp\left(-{\left(\frac{t-T_0}{a_0}\right)}^2\right)u(t),& (11)\end{array}$

where T₀ = 0.250 ns, a₀ = 0.0633 ns, and u(t) is the unit step function, which is applied from the air region on the left side as shown in Figure 2. The applied Gaussian pulse contains frequency components from dc to approximately 9.3 GHz, where the power of the pulse decreases 1,000-fold from its maximum value. Number of sampling points used in fast Fourier transform to obtain the reflection coefficients and transmission characteristics is 120,000. Zero padding is used after 100,000 sampling data.

FIGURE 2

Figure 2. Multilayer model consisting of six biological tissues.

TABLE 3

Table 3. Tissue density and thickness of each tissue layer.

Figure 3 indicates the reflection coefficient and the transmission of the multilayer model as a percentage obtained by the FDTD and analytical methods from 50 MHz to 10 GHz. The transmission, indicating the power transmitting into a biological tissue layer, is calculated using the following equation:

$\begin{array}{ll}{P}^t(\omega )=(1-|\Gamma {|}^2)\times 100\left[\%\right],& (12)\end{array}$

where Γ is the reflection coefficient. It is shown that the numerical and analytical results are in good agreement, demonstrating the validity of our proposed FDTD method again. Note also that the broadband results are numerically obtained in the time domain by a single FDTD run and transformed to those in the frequency domain via the fast Fourier transform. From the results, it can be seen that the reflection coefficient decreases with increasing frequency and reaches a minimum value of 0.38 at ~2 GHz, while the transmission exhibits peak at 2 GHz. From Figure 3, more than 80% of the incident power penetrates into the multilayer model at the maximum transmission frequency of around 2 GHz. Note that the maximum transmission frequency depends on the thicknesses of the biological tissues in the model; thus, using a different model will yield different results from those shown in this study. Next, the transmission characteristics of the EM pulse are obtained from the ratio between the receiving power at an observation point inside a biological tissue layer and the power penetrating into the model from the leftmost boundary of the skin layer.

FIGURE 3

Figure 3. Reflection coefficient and transmission of multilayer mimicking a human head model.

Figure 4 indicates the transmission characteristics at the center of each layer in the left half of the analysis model as a percentage. From the figure, we observe a peak of the transmitted power in the skin at around 2.5 GHz. This peak shifts to a higher frequency with decreasing maximum value when the power penetrates into the subsequent layer, which is the “Bone” layer in this case. It is also shown that when the observation point is located inside the CSF layer, the transmission characteristics of this multilayer model are almost flat in the range between 300 and 800 MHz, having a maximum at ~500 MHz. In addition, when the frequency is larger than 1 GHz, most of the power is absorbed at the superficial layers before reaching the CSF layer and, therefore, the transmitted power monotonically decreases with increasing frequency in this region, except between 5 and 8 GHz, where the transmissions in the “Fat” layer is greater than that in the “Skin” layer. This may be due to a small loss in the “Bone,” compared to that of “Skin,” and multiple reflections occurring between the Skin–Fat and Fat–Bone boundaries that create a local maximum. Note that only a small proportion of the power reaches the other side of the analysis region. For example, as shown in Figure 4, the transmission power that reaches the CSF layer on the right side of the analysis model is <2% of the total transmission power and is concentrated at lower frequencies of below 100 MHz.

FIGURE 4

Figure 4. Transmission characteristics at observation points inside biological tissue layers.

Dosimetry of Various Human Bodies Exposed to EM Pulses

As discussed in the previous subsection, the transmission into the multilayer mimicking a human head model shows a high value in the frequency range between 300 and 800 MHz at the deep tissues such as “Dura” and “CSF” and decreases above a frequency of 1 GHz due to the superficial absorption of EM energy. Therefore, our target for numerical dosimetry is an EM pulse having broad frequency components below 1 GHz. As shown in Figure 5, numerical human models of anatomical adult male (TARO), adult female (HANAKO), 7-, 5-, and 3-year child models, which were developed by National Institute of Information and Communications Technology, Japan, are chosen as our targets for numerical dosimetry (31, 32). They contain 51 different biological tissues and have a spatial resolution of 2 mm. Figure 5 also shows exposure situations with an ungrounded model (free-space model) and a grounded model standing on a ground plane made of the perfectly electric conductor (PEC). The height and weight for each numerical human model are also indicated in Figure 5. In accordance with the Courant condition, the time step interval is determined as 3.85 ps. The time step interval used here is different from that used in subsection Multilayer Model; thus, we need to recalculate the update coefficients for all biological tissues used in the numerical simulations.

FIGURE 5

Figure 5. Numerical human models used in simulations. Shading color indicates the distance from the front of analysis domain. Darker for smaller distance.

The update coefficients are derived by applying the FILT and the Prony method as described above. The time step interval used to normalize the relaxation time in the Cole–Cole function is then set to 3.85 ps, the same as that used in the numerical simulations. After calculating the update coefficients, we validate them by performing one-dimensional simulations and computing the reflection coefficients. We have found that the reflection coefficients of all biological media having a permittivity function characterized by the Cole–Cole model match those obtained from the EM theory with a small difference of <2% in the frequency range between 10 MHz and 2 GHz, demonstrating the validity of the determined update coefficients used in FDTD simulations (33). The required number of terms of the update coefficients is different for each tissue, with the maximum number of N = 50 for “Infiltrated Cancellous Bone Barrow.” The reflection coefficients show values ranging from 0.70 to 0.90 at 10 MHz, which decrease to 0.40–0.80 at 2 GHz. The biological tissue with the lowest and highest reflection coefficients are “Non-infiltrated Bone Marrow” and “Cerebro-Spinal Fluid (CSF),” respectively. This is attributed to the fact that CSF has a higher conductivity than other tissues.

An incident electric field polarized in the z-direction, having the same Gaussian waveform with T₀ = 0.385 ns and a₀ = 0.146 ns impinges on each model from the front. The power of the incident EM pulse decreases by half at ~1.3 GHz. An EM plane wave has a polarization axis parallel with the human body axis. The amplitude of the incident electric field is 1 V/m. The calculation of SA can be carried out until the pulse strength decreased to almost zero at a specific time step. The calculation time and memory usage for each model are tabulated in Table 4. The calculation time also includes the computation time for an on-the-fly Fourier transform to obtain the electric field distributions inside the human models for use in deriving the whole-body average SAR at 24 different frequencies from 10 MHz to 1 GHz within a single FDTD run. Numerical simulations involving human models are carried out on a single calculation node (Intel Xeon E5-2680v4 @ 2.4 GHz, 256 GB memory) of TSUBAME3.0 supercomputer at Tokyo Institute of Technology, Japan, with 28 parallel threads. The total number of time steps is 10,000 steps. For our cases, total electric energy inside the human model was < -50 dB after 10,000 step, compared to its maximum during numerical simulations. The total size of the analysis region including perfectly matched layers (PMLs) are also indicated in the table.

TABLE 4

Table 4. Calculation time and memory usage for each model.

Figure 6 shows the SAR distribution inside the TARO model illuminated by a planewave with a 1 V/m electric field at 200 MHz. The results are obtained by using the proposed (FD)²TD approach and using constant dielectric properties for comparison. Both distributions are well-matched, demonstrating the validity of the (FD)²TD method. Figures 7A,B indicate the whole-body-average specific absorption rate (WBA-SAR) of various numerical human models for the ungrounded and grounded conditions, respectively. Note that our numerical results are obtained by a single run of FDTD computation for each model. The WBA-SAR at each frequency is calculated using an on-the-fly Fourier transform during the FDTD run. It is shown that the results obtained by the proposed method are in good agreement with those derived in the literature, which are also shown as a red dashed line in the figure (34), again demonstrating the validity of our proposed method. Note that the results in Figure 7 are normalized by the limits of the incident power density prescribed by the ICNIRP guidelines at each frequency. As indicated in Table 4, our proposed method requires almost three time longer than the conventional one using constant permittivity at a frequency, however, the conventional method can only provide WBA-SAR at only one frequency per FDTD computation and number of time steps required to obtain the converged solution is different for each frequency (35). Meanwhile, the proposed method requires only a single FDTD run to determine the solutions at 24 analysis frequencies. Although the proposed method requires more memory usage, in consequence, it is shown to be very computationally efficient. As number of frequency points increases, efficiency of the proposed method increases at the expense of using more memory. From Figure 7, the SAR peak is found at around 70 MHz for the ungrounded TARO model, which corresponds to the whole-body resonance and the peak frequency increases for the shorter HANAKO, 7-, 5-, and 3-y-child models. For grounded conditions, the SAR peaks occur at 40 and 70 MHz for the adult and 3-y-child models, respectively. These results are in good agreement with those indicated in the literature (32, 34).

FIGURE 6

Figure 6. SAR distribution at 500 MHz calculated using the proposed FDTD method and that using the conventional FDTD method with constant dielectric properties.

FIGURE 7

Figure 7. Whole-body average SAR over a frequency range between 10 MHz to 1 GHz, calculated using the proposed FDTD method for the ungrounded and grounded cases. (A) Ungrounded cases. (B) Grounded cases.

Figures 8, 9 show SA distributions on the human surface for ungrounded and grounded cases, respectively. For the ungrounded cases, SA peaks appear at the ankles, wrists (hand), forearm, and neck; SA increases by 53% (adult male) to 74% (3-y-child) at the ankle, whereas SA at the other parts remain almost the same in the ungrounded cases. Figure 10 shows the layer maximum SA at each height of various human bodies when exposed to an EM pulse of 1 V/m, polarized in the z-direction. The peak SA and the SA at the ankle for various human models, which are normalized to the ICNIRP-prescribed power density limit of 2 W/m², are summarized in Table 5 (10). Note that the power density limit is constant for 10–400 MHz and increases with respect to the frequency over 400 MHz. Therefore, there is no prescribed limit value of the power density for a wideband pulse such as one used in our study. However, we hereby use 2 W/m² for the normalization of SA in Table 5 as it should provide conservative evaluations. The maximum SA appears at the hands for the adult male, 7-, 5-, and 3-y-child models for the ungrounded case, while it appears at the neck for the adult female model. This may be attributed to the proportion of biological tissues in the male and female human models is different since the female model contains more fat in each body part. Note that 7-, 5-, and 3-y-child models are proportionally morphed using a morphing algorithm (32). The maximum SA among five human models is 0.437 pJ/kg. Note that this SA value also depends on the waveform (which is the Gaussian pulse in our study). To reach a dose of 2 mJ/kg, as prescribed in the ICNIRP guidelines, we need to increase the field strength from 1 V/m to more than 83 kV/m or 9.14 MW/m², which does not seem realistic in real life. Note that SA is obtained from the value at a voxel and it should be smaller for an average over 10 g tissues. Hence, the SA values shown in Table 5 assume a worst-case scenario. For compliance with the IEEE standards, the repetition rate of an incident pulse having a field strength of 87 kV/m must be <800 Hz or 800 pulses per second in order not to exceed the SA limit. In the grounded case, the maximum SA occurs at the ankle for the adult female, 5-, and 3-y-child while they are found at the same location (hands) in the ungrounded case for the adult male and 7-y-child. It is found that energy absorption at the ankle increases when the human body is grounded, while that at the other parts remains almost unchanged for all the models used in this study. Table 6 shows the peak 1 g-averaged and 10 g-averaged SAs normalized by a power density limit of 2 W/m². It is shown that the peak 1 g-averaged and 10 g-averaged SAs are higher in the grounded cases than those in the ungrounded cases for all models. The increases in the SAs are significant in a smaller model, e.g., for the grounded 3-y-child model, the peak 1 g-averaged and 10 g-averaged SAs are ~1.67 and 1.76 than those of the ungrounded cases. These results provide the first ever demonstration that the SA distribution due to broadband EM pulse illumination can be quantitatively evaluated in detail and compared with the SA limit prescribed in international guidelines or standards. Further detailed exposure levels for different incident angles, different pulse shapes, and postures will be investigated in the future.

FIGURE 8

Figure 8. SA distribution for one pulse illumination for various ungrounded human models.

FIGURE 9

Figure 9. SA distribution for one pulse illumination for various grounded human models.

FIGURE 10

Figure 10. Maximum SA at each height for various human models for 1 V/m electric field strength.

TABLE 5

Table 5. Maximum specific energy absorption normalized by an ICNIRP-prescribed power density limit of 2 W/m² for general exposures to various human model.

TABLE 6

Table 6. Peak 1 g-averaged and 10 g-averaged SA normalized by an ICNIRP-prescribed power density limit of 2 W/m² for general exposures to various human model.

Conclusion

We have performed numerical dosimetry on human bodies illuminated by an EM pulse from the front by using the (FD)²TD method, previously proposed by the authors. The method fully considers broadband characteristics of the complex relative permittivity of the biological media used in the analysis model via the application of the FILT and the Prony method. Firstly, we demonstrated the validity of the update coefficients, i.e., the residues and poles of the expression for the IIR in the z-domain, by comparing the numerical reflection coefficients with those derived from the EM theory. It was clarified that the numerical results within 2% of those obtained theoretically over a broad frequency range from 50 MHz to 10 GHz, demonstrating the validity of the proposed approach. It was also found that the transmission characteristics of the EM pulse into the CSF layer of a multilayer mimicking a human head are almost flat over a frequency range between 300 and 800 MHz and that the transmission decreases with increasing EM traveling distance from the skin boundary due to higher energy absorption at superficial biological tissues such as “Skin” and “Fat” when the frequency is higher than 1 GHz. Therefore, most of the pulse energy that penetrates into the biological body has a frequency below 1 GHz. Then, numerical dosimetry of various human models exposed to an EM pulse having a frequency component of up to ~1.3 GHz was performed. The whole-body average SAR at 24 frequencies was determined by a single run of broadband FDTD simulations. The results matched those published in the literature, demonstrating the validity and availability of the proposed FDTD method. Then the SA distribution of each numerical model was determined, and it was found that the maximum SA occurs at the hands and neck for the ungrounded model, while they appear at the hands and ankle when the model is grounded. The maximum SA value was 0.290 pJ/kg for an incident electric field strength of 1 V/m or 0.437 nJ/kg for an incident power density of 2 W/m². It has been shown for the first time that, by using our proposed FDTD approach with the FILT and the Prony method, we can obtain quantitatively detailed information on SA that can be compared with the limits prescribed in international guidelines or standards.

Data Availability Statement

The original contributions presented in the study are included in the article, further inquiries can be directed to the corresponding author.

Author Contributions

JC designed and performed all the simulations, analyzed the data, and wrote the paper. KW and KF gave advices and revised the paper. All authors have read and agreed to the published version of the manuscript.

Funding

This study was financially supported by a JSPS Grant-in-Aid for Scientific Research (JP18K18376), Japan, parts of the numerical calculations were carried out on the TSUBAME3.0 supercomputer at Tokyo Institute of Technology.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

1. FCC. Revision of Part 15 of the Commission's Rules Regarding Ultra-Wideband Transmission Systems (2002).

Google Scholar

2. Allen B, Dohler M, Okon E, Malik W, Brown A, Edwards D. Ultra-Wideband Antennas and Propagation for Communications, Radar and Imaging. West Sussex: John Wiley & Sons, Ltd. (2007).

Google Scholar

3. Thotahewa KMS, Redout JM, Yuce MR. Propagation, power absorption, and temperature analysis of UWB wireless capsule endoscopy devices operating in the human body. IEEE Trans Microwave Theory Tech. (2015) 63:3823–33. doi: 10.1109/TMTT.2015.2482492

CrossRef Full Text | Google Scholar

4. Hall PS. Antennas and Propagation for Body-Centric Wireless Communications. Norwood, MA: Artech House, Inc. (2006).

Google Scholar

5. Gandhi OP, Furse CM. Currents induced in the human body for exposure to ultrawideband electromagnetic pulses. IEEE Trans Electromagn Compat. (1997) 39:174–80. doi: 10.1109/15.584941

CrossRef Full Text | Google Scholar

6. Chou CK, Guy AW. Microwave-induced auditory responses in guinea pigs: relationship of threshold and microwave-pulse duration. Radio Sci. (1979) 14:193–7. doi: 10.1029/RS014i06Sp00193

CrossRef Full Text | Google Scholar

7. Weaver JC, Smith KC, Esser AT, Son RS, Gowrishankar TR. A brief overview of electroporation pulse strength-duration space: a region where additional intracellular effects are expected. Bioelectrochemistry. (2012) 87:236–43. doi: 10.1016/j.bioelechem.2012.02.007

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Lin JC. The microwave auditory effect. IEEE J Electromagn RF Microwaves Med Biol. (2021). doi: 10.1109/JERM.2021.3062826

CrossRef Full Text | Google Scholar

9. Weil CM. Absorption characteristics of multilayered sphere models exposed to UHF/microwave radiation. IEEE Trans Biomed Eng. (1975) 22:468–76. doi: 10.1109/TBME.1975.324467

PubMed Abstract | CrossRef Full Text | Google Scholar

10. ICNIRP. Guidelines for limiting exposure to time-varying electric, magnetic, and electromagnetic fields (up to 300 GHz). International Commission on Non-Ionizing Radiation Protection. Health Phys. (1998) 74:494–522.

Google Scholar

11. IEEE-C95.1. IEEE Standard for Safety Levels With Respect to Human Exposure to Radio Frequency Electromagnetic Fields, 3 kHz to 300 GHz. IEEE Std C951-2005 (Revision of IEEE Std C951-1991) (2006). p. 1–238.

Google Scholar

12. ICNIRP. Guidelines for limiting exposure to electromagnetic fields (100 kHz to 300 GHz). Health Phys. (2020) 118:483–524. doi: 10.1097/HP.0000000000001210

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Lin JC, Chuan-Lin W, Lam CK. Transmission of electromagnetic pulse into the head. Proc IEEE. (1975) 63:1726–7. doi: 10.1109/PROC.1975.10043

CrossRef Full Text | Google Scholar

14. Lin JC. Interaction of electromagnetic transient radiation with biological materials. IEEE Trans Electromagn Compat. (1975) 17:93–7.

Google Scholar

15. Gabriel S, Lau RW, Gabriel C. The dielectric properties of biological tissues: III. Parametric models for the dielectric spectrum of tissues. Phys Med Biol. (1996) 41:2271. doi: 10.1088/0031-9155/41/11/003

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Luebbers R, Hunsberger FP, Kunz KS, Standler RB, Schneider M. A frequency-dependent finite-difference time-domain formulation for dispersive materials. IEEE Trans Electromagn Compat. (1990) 32:222–7. doi: 10.1109/15.57116

CrossRef Full Text | Google Scholar

17. Kelley DF, Luebbers RJ. Piecewise linear recursive convolution for dispersive media using FDTD. IEEE Trans Antennas Propag. (1996) 44:792–7. doi: 10.1109/8.509882

CrossRef Full Text | Google Scholar

18. Shibayama J, Ando R, Nomura A, Yamauchi J, Nakano H. Simple trapezoidal recursive convolution technique for the frequency-dependent FDTD analysis of a Drude–Lorentz model. IEEE Photon Technol Lett. (2009) 21:100–2. doi: 10.1109/LPT.2008.2009003

CrossRef Full Text | Google Scholar

19. Goorjian PM, Taflove A. Direct time integration of Maxwell's equations in nonlinear dispersive media for propagation and scattering of femtosecond electromagnetic solitons. Opt Lett. (1992) 17:180–2. doi: 10.1364/OL.17.000180

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Taflove A, Hagness SC. Computational Electrodynamics: The Finite-Difference Time-Domain Method. Artech House (2005).

Google Scholar

21. Kashiwa T, Fukai I. A treatment by the FD-TD method of the dispersive characteristics associated with electronic polarization. Microwave Opt Technol Lett. (1990) 3:203–5. doi: 10.1002/mop.4650030606

CrossRef Full Text | Google Scholar

22. Sullivan DM. Frequency-dependent FDTD methods using Z transforms. IEEE Trans Antennas Propag. (1992) 40:1223–30. doi: 10.1109/8.182455

CrossRef Full Text | Google Scholar

23. Rekanos IT Papadopoulos TG. FDTD modeling of wave propagation in Cole-Cole media with multiple relaxation times. IEEE Antennas Wireless Propag Lett. (2010) 9:67–9. doi: 10.1109/LAWP.2010.2043410

CrossRef Full Text | Google Scholar

24. Mescia L, Bia P, Caratelli D. Fractional derivative based FDTD modeling of transient wave propagation in Havriliak-Negami media. IEEE Trans Microwave Theory Tech. (2014) 62:1920–9. doi: 10.1109/TMTT.2014.2327202

CrossRef Full Text | Google Scholar

25. Chakarothai J. Novel FDTD scheme for analysis of frequency-dependent medium using fast inverse Laplace transform and Prony's method. IEEE Trans Antennas Propag. (2019) 67:6076–89. doi: 10.1109/TAP.2018.2878077

CrossRef Full Text | Google Scholar

26. Chakarothai J, Watanabe S, Wake K. Numerical dosimetry of electromagnetic pulse exposures Using FDTD Method. IEEE Trans Antennas Propag. (2018) 66:5397–408. doi: 10.1109/TAP.2018.2862344

CrossRef Full Text | Google Scholar

27. Chakarothai J, Wake K, Watanabe S, Chen Q, Sawaya K. Frequency dependent FDTD method for ultra-wideband electromagnetic analyses. IEICE Trans Electr. (2019) J102-C:102–13.

Google Scholar

28. Kensuke S, Kanako W, Soichi W. Development of best fit Cole-Cole parameters for measurement data from biological tissues and organs between 1 MHz and 20 GHz. Radio Sci. (2014) 49:459–72. doi: 10.1002/2013RS005345

CrossRef Full Text | Google Scholar

29. Abdalla A, Teoh A. A multi layered model of human head irradiated by electromagnetic plane wave of 100MHz–300 GHz. Int J Sci Res. (2005) 15:1–7.

30. Sabbah AI, Dib N, Al-Nimr M. Evaluation of specific absorption rate and temperature elevation in a multi-layered human head model exposed to radio frequency radiation using the finite-difference time domain method. Microwaves Antennas Propag IET. (2011) 5:1073–80. doi: 10.1049/iet-map.2010.0172

CrossRef Full Text | Google Scholar

31. Nagaoka T, Watanabe S, Sakurai K, Kunieda E, Watanabe S, Taki M, et al. Development of realistic high-resolution whole-body voxel models of Japanese adult males and females of average height and weight, and application of models to radio-frequency electromagnetic-field dosimetry. Phys Med Biol. (2004) 49:1–15. doi: 10.1088/0031-9155/49/1/001

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Nagaoka T, Kunieda E, Watanabe S. Proportion-corrected scaled voxel models for Japanese children and their application to the numerical dosimetry of specific absorption rate for frequencies from 30 MHz to 3 GHz. Phys Med Biol. (2008) 53:6695. doi: 10.1088/0031-9155/53/23/004

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Chakarothai J, Fujii K, editors. A Unified approach for treatment of frequencydependent materials in FDTD method. In: 2019 International Symposium on Antennas and Propagation (ISAP). Xi'an (2019).

Google Scholar

34. Wang J, Fujiwara O, Kodera S, Watanabe S. FDTD calculation of whole-body average SAR in adult and child models for frequencies from 30 MHz to 3 GHz. Phys Med Biol. (2006) 51:4119. doi: 10.1088/0031-9155/51/17/001

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Chakarothai J, Wake K, Watanabe S. Convergence of a Single-Frequency FDTD Solution in Numerical Dosimetry. IEEE Trans Microwave Theory Tech. (2016) 64:707–14. doi: 10.1109/TMTT.2016.2518661

CrossRef Full Text | Google Scholar