Non-Dimensional Analysis of Thermal Effect on Skin Exposure to an Electromagnetic Beam

Hongyun Wang¹, Wesley A. Burgei², Hong Zhou^3*

¹Department of Applied Mathematics, University of California, Santa Cruz, CA, USA.
²U.S. Department of Defense, Joint Intermediate Force Capabilities Office, Quantico, VA, USA.
³Department of Applied Mathematics, Naval Postgraduate School, Monterey, CA, USA.

DOI: 10.4236/ajor.2020.105011                

Abstract

We consider the problem of inducing withdrawal reflex on a test subject by exposing the subject’s skin to an electromagnetic beam. Heat-sensitive nociceptors in the skin are activated wherever the temperature is above the activation temperature. Withdrawal reflex occurs when the activated volume reaches a threshold. We non-dimensionalize the problem to write the temperature as the product of a parameter-free function of non-dimensional variables and a function of beam parameters. This formulation allows studying beam parameters without knowing skin material parameters. We examine the effects of spot size, total power and distribution type of the electromagnetic beam on 3 quantities at reflex: 1) the time to reflex, 2) the maximum temperature increase, and 3) the total energy consumption. We find that the flat-top beam is the best, with the lowest energy consumption and the smallest maximum temperature increase. The Super-Gaussian beam is only slightly inferior to the flat-top. The Gaussian beam has by far the worst performance among these three.

Keywords

Heat-Sensitive Nociceptor, Heat-Induced Withdrawal Reflex, Radius and Distribution Type of Electromagnetic Beam, Non-Dimensional Analysis

1. Introduction

Electromagnetic beams have been employed in many applications which include airport security screening, medical imaging and non-lethal crowd control weapons [1] [2]. The interaction of the electromagnetic fields with the human body is dependent upon the power density, beam spot size, frequency and duration. Particularly, when a high-powered microwave beam reaches a human subject, it quickly produces an intolerable heating sensation in the skin and compels the subject to withdraw from the beam [3].

In this paper we carry out non-dimensional analysis of thermal effect on skin caused by exposure to a stationary electromagnetic beam. We focus on how the beam radius, total power and distribution type of the electromagnetic beam affect the time to withdrawal reflex, the maximum temperature increase, and the total energy consumption. The theoretical results obtained in our study will help guide us to design safe, effective and efficient electromagnetic heating experiments.

2. Mathematical Formulation

We consider the situation where a skin area of the test subject is exposed to an electromagnetic beam. Let $T\left(z\mathrm{,}r\mathrm{,}t\right)$ denote the temperature of skin as a function of 3-D spatial coordinates $\left(z\mathrm{,}r\right)$ and time t. Here z is the depth from the skin surface and $r$ is the 2-D coordinates on the skin surface. We assume 1) the electromagnetic beam is perpendicular to the skin surface (i.e., incident angle = 0); 2) before the exposure, the 3-D skin has a uniform initial temperature $T_{\text{base}}$ (baseline temperature); and 3) heat conduction is included only in the depth direction. Under these assumptions, which correspond to case B in our previous study [4], the temperature distribution $T\left(z\mathrm{,}r\mathrm{,}t\right)$ can be solved analytically.

At each surface location $r$, the temperature along the depth is governed by

$\{\begin{array}{l}{\rho }_{\text{m}}{C}_p\frac{\partial T}{\partial t}=K\frac{{\partial }^{2}T}{\partial z^2}+P\left(r\right)\mu exp\left(-\mu z\right)\hfill \\ {\frac{\partial T}{\partial z}|}_{z=0}=\mathrm{0,}T\left(z\mathrm{,0}\right)=T_{\text{base}}\hfill \end{array}$

where

· ${\rho }_{\text{m}}$ is the mass density of the skin,

· $C_p$ is the specific heat capacity of the skin,

· K is the heat conductivity of the skin,

· $\mu$ is the absorption coefficient of the skin at the beam frequency, and

· $P\left(r\right)$ is the beam power density absorbed into the skin at location $r$.

In the ADT CHEETEH model [5], the absorbed power density $P\left(r\right)$ is modeled as a given fraction of the incident power density at the skin surface:

$P\left(r\right)=\left(1-\gamma \right)P_{\text{incident}}(\; r\; )$

where $\gamma$ fraction of the incident power is reflected by the skin.

The absorbed electromagnetic energy increases the temperature in skin. Heat-sensitive nociceptors are activated wherever the local temperature is above the nociceptor activation temperature $T_{\text{act}}$. When the combined signal from nociceptors is over a certain threshold, withdrawal reflex occurs and the test subject moves away from the beam. Our study is aimed at understanding the process of beam power absorbed at the skin surface to withdrawal reflex. We first introduce characteristic scales for the depth, time and temperature, and then use them to carry out non-dimensionalization.

· Length scale in the depth direction and time scale:

$z_{\text{s}}\equiv \frac{1}{\mu }\mathrm{,}{t}_{\text{s}}\equiv \frac{{\rho }_{\text{m}}{C}_p}{K{\mu }^2}$

· Non-dimensional depth and time:

$z_{\text{nd}}\equiv \frac{z}{z_{\text{s}}}=\mu z\mathrm{,}{t}_{\text{nd}}\equiv \frac{t}{t_{\text{s}}}=t\frac{K{\mu }^2}{{\rho }_{\text{m}}{C}_p}$

· Temperature scale and power density scale:

$\Delta T\equiv T_{\text{act}}-T_{\text{base}}\mathrm{,}{P}_{\text{s}}\equiv K\mu \Delta T$

· Non-dimensional temperature as a function of $\left(z_{\text{nd}}\mathrm{,}{t}_{\text{nd}}\right)$ :

$T_{\text{nd}}\left(z_{\text{nd}}\mathrm{,}{t}_{\text{nd}}\right)\equiv \frac{T\left(z\mathrm{,}t\right)-T_{\text{base}}}{\Delta T}$

· Non-dimensional power density (absorbed into the skin) at surface location $r$ :

$P_{\text{nd,1}}\left(r\right)\equiv \frac{P\left(r\right)}{P_{\text{s}}}=\frac{P\left(r\right)}{K\mu \Delta T}$

Note that the function value of $P_{\text{nd,1}}\left(r\right)$ is dimensionless but the independent variable $r$ is not, which will be non-dimensionalized shortly. At each surface location $r$, the non-dimensional temperature is governed by

$\{\begin{array}{l}\frac{\partial T_{\text{nd}}}{\partial t_{\text{nd}}}=\frac{{\partial }^2{T}_{\text{nd}}}{\partial z_{\text{nd}}^2}+P_{\text{nd,1}}\left(r\right)exp\left(-z_{\text{nd}}\right)\hfill \\ {\frac{\partial T_{\text{nd}}}{\partial z_{\text{nd}}}|}_{z_{\text{nd}}=0}=\mathrm{0,}T\left(z_{\text{nd}}\mathrm{,0}\right)=0\hfill \end{array}$

We adopt the same general assumption used in the ADT CHEETEH model [5] that the occurrence of withdrawal reflex is solely determined by the number of heat-sensitive nociceptors activated at the given time. When the nociceptor density in the skin is uniform, withdrawal reflex is solely determined by the volume of activated region at the given time, which is highly influenced by the beam spot size.

To study the effect of beam spot size in a proper mathematical setting, we need to introduce a length scale for non-dimensionalizing the 2-D surface coordinates $r$. Let $v_{\text{c}}$ be the critical threshold on activated volume for triggering withdrawal reflex. In a general formulation, we allow the realized value of threshold $v_{\text{c}}$ in individual tests to be a random variable, fluctuating from one test to another to model the effect of biovariability and environmental uncertainty [6]. To keep the formulation simple, we assume that all other parameters are deterministic and that the effect of uncertainty from all sources is summarized in the randomness of threshold $v_{\text{c}}$. By definition, the reflex time $t_{\text{ref}}$ is the time it takes for the activated volume to reach threshold $v_{\text{c}}$. The uncertainty in threshold $v_{\text{c}}$ leads to fluctuations in the observed reflex time $t_{\text{ref}}$ in individual tests.

Let $E\left\{v_{\text{c}}\right\}$ denote the average of threshold $v_{\text{c}}$, which provides us a well-defined physical volume. For surface coordinates $r$, we introduce the characteristic radius $r_{\text{s}}$ in the equation below, based on the physical volume $E\left\{v_{\text{c}}\right\}$.

$\left(\pi r_{\text{s}}^2\right)z_{\text{s}}=E\left\{v_{\text{c}}\right\}$

Geometrically, a cylinder of base radius $r_{\text{s}}$ and height $z_{\text{s}}$ has volume $E\left\{v_{\text{c}}\right\}$. We use $r_{\text{s}}$ as the length scale to non-dimensionalize surface coordinates $r$.

· Length scale for surface coordinates and volume scale:

$r_{\text{s}}\equiv \sqrt{\frac{E\left\{v_{\text{c}}\right\}}{\pi z_{\text{s}}}}\mathrm{,}{v}_{\text{s}}\equiv r_{\text{s}}^2{z}_{\text{s}}=\frac{E\left\{v_{\text{c}}\right\}}{\pi }$

· Non-dimensional surface coordinates:

$r_{\text{nd}}\equiv \frac{r}{r_{\text{s}}}$

· Non-dimensional volume and non-dimensional volume threshold:

$V_{\text{nd}}\equiv \frac{V\pi }{E\left\{v_{\text{c}}\right\}}\mathrm{,}{v}_{\text{c,nd}}\equiv \frac{v_{\text{c}}\pi }{E\left\{v_{\text{c}}\right\}}$

· Non-dimensional power density as a function of $r_{\text{nd}}$ :

$P_{\text{nd}}\left(r_{\text{nd}}\right)\equiv P_{\text{nd,1}}\left(r_s{r}_{\text{nd}}\right)$

Note that the volume scale is $v_{\text{s}}\equiv \frac{E\left\{v_{\text{c}}\right\}}{\pi }$. Function $P_{\text{nd}}(\cdot )$ is the same as $P_{\text{nd,1}}(\cdot )$ except that it has non-dimensional $r_{\text{nd}}$ as the independent variable. After non-dimensionalization, parameters satisfy

$E\left\{v_{\text{c,nd}}\right\}=\pi \mathrm{,}{T}_{\text{base,nd}}=\mathrm{0,}{T}_{\text{act,nd}}=1$

The non-dimensional temperature has the analytical expression

$\{\begin{array}{l}T\left(z\mathrm{,}r\mathrm{,}t\right)=P\left(r\right)H\left(z\mathrm{,}t\right)\\ H\left(z\mathrm{,}t\right)={\displaystyle {\int }_0^t}G\left(z\mathrm{,}s\right)\text{d}s\\ G\left(z\mathrm{,}t\right)=\frac{1}{2}\text{erfc}\left(\frac{2t-z}{\sqrt{4t}}\right){\text{e}}^{t-z}+\frac{1}{2}\text{erfc}\left(\frac{2t+z}{\sqrt{4t}}\right){\text{e}}^{t+z}\end{array}$ (1)

Here, for conciseness, we have dropped the subscript “nd” and used the simple notations for all non-dimensional quantities. For example, $P\left(r\right)$ in (1) means $P_{\text{nd}}\left(r_{\text{nd}}\right)$. Based on the non-dimensional temperature, we calculate the non-dimensional activated volume

$\begin{array}{c}{V}_{\text{act}}\left(t\right)=\text{Volume}\left\{\left(z\mathrm{,}r\right)|T\left(z\mathrm{,}r\mathrm{,}t\right)\ge T_{\text{act}}\right\}\mathrm{,}{T}_{\text{act}}=1\\ =\text{Volume}\left\{\left(z\mathrm{,}r\right)|P\left(r\right)H\left(z\mathrm{,}t\right)\ge 1\right\}\end{array}$

The non-dimensional reflex time is governed by $V_{\text{act}}\left(t_{\text{ref}}\right)=v_{\text{c}}$ :

$\text{Volume}\left\{\left(z\mathrm{,}r\right)|P\left(r\right)H\left(z\mathrm{,}{t}_{\text{ref}}\right)\ge 1\right\}=v_{\text{c}}$ (2)

where non-dimensional volume threshold $v_{\text{c}}$ satisfies $E\left\{v_{\text{c}}\right\}=\pi$. Notice that in (1), function $H\left(z\mathrm{,}t\right)$ is parameter-free. In this non-dimensional formulation, $T\left(z\mathrm{,}r\mathrm{,}t\right)$ is completely determined by $P\left(r\right)$, there are no other parameters. The non-dimensional reflex time $t_{\text{ref}}$ is influenced by two factors in the non-dimensional beam power density $P\left(r\right)$ :

1) The non-dimensional total power (absorbed into the skin) $P_{\text{tot}}\equiv {\displaystyle \int P\left(r\right)\text{d}r}$.

2) The spatial distribution of beam power. We will examine the effect of distribution type (Gaussian, super-Gaussian, or flat-top) and beam spot size.

One convenient measure of beam spot size is the half width at half maximum (HWHM) which is half of the full width at half maximum (FWHM) for all symmetric beam power distributions. Let $r_{\text{hw}}$ denote the HWHM. We study three symmetric distributions in the form of $P\left(r\right)=P\left(r\right)$ where $r\equiv \Vert r\Vert$.

Gaussian:

$P\left(r\right)=P\left(0\right)exp\left[-\left(ln2\right){\left(\frac{r}{r_{\text{hw}}}\right)}^2\right]$

Super-Gaussian ( $n=10$ ):

$P\left(r\right)=P\left(0\right)exp\left[-\left(ln2\right){\left(\frac{r}{r_{\text{hw}}}\right)}^n\right]$

Flat-top:

$P\left(r\right)=P\left(0\right)\theta \left(r_{\text{hw}}-r\right)\mathrm{,}\theta \left(s\right)\equiv \{\begin{array}{ll}0,\hfill & s<0\hfill \\ 1,\hfill & s\ge 0\hfill \end{array}$

where $P\left(0\right)$ is the beam center power density (absorbed into the skin) and $\theta \left(s\right)$ is the Heaviside step function. Given a distribution type (Gaussian, super-Gaussian or flat-top), the expression above allows us to describe the beam power distribution using $P\left(0\right)$ and $r_{\text{hw}}$. As we will see in Figure 1, at a given total power (absorbed into the skin), $r_{\text{hw}}$ does not provide a good description of how concentrated the beam power is, when comparing different types of distribution. To measure the effective area of the beam, we use the effective mode radius $r_{\text{eff}}$ defined as [7]

$r_{\text{eff}}\equiv \sqrt{\frac{1}{\pi }\cdot \frac{{\left({\displaystyle \int P\left(r\right)2\pi r\text{d}r}\right)}^2}{{\displaystyle \int P{\left(r\right)}^{2}2\pi r\text{d}r}}}\mathrm{<mo>\; =\; </mo>}\frac{{\displaystyle \int P\left(r\right)2r\text{d}r}}{\sqrt{{\displaystyle \int P{\left(r\right)}^{2}2r\text{d}r}}}$ (3)

For conciseness, we shall refer to $r_{\text{eff}}$ simply as the beam radius. We like to specify the beam power distribution in terms of the total power $P_{\text{tot}}\equiv {\displaystyle \int P\left(r\right)2\pi r\text{d}r}$ and the radius $r_{\text{eff}}$. For that purpose, we express $P_{\text{tot}}$ and $r_{\text{eff}}$ in terms of $P\left(0\right)$ and $r_{\text{hw}}$, for each of the three distribution types above. We first derive a general integral formula

${\displaystyle {\int }_0^{+\infty }}exp\left(-\alpha x^n\right)2x\text{d}x=\frac{2}{n}{\alpha }^{-\frac{2}{n}}{\displaystyle {\int }_0^{+\infty }}exp\left(-s\right)s^{\frac{2}{n}-1}\text{d}s=\frac{2}{n}{\alpha }^{-\frac{2}{n}}\Gamma (\; 2\; n\; )$

Applying this integral formula to the three distributions, we obtain

Gaussian:

$P_{\text{tot}}=P\left(0\right){\displaystyle \int }exp\left[-\left(ln2\right){\left(\frac{r}{r_{\text{hw}}}\right)}^2\right]2\pi r\text{d}r=P\left(0\right)\pi r_{\text{hw}}^2\cdot \frac{1}{ln2}$

$r_{\text{eff}}=r_{\text{hw}}\cdot \sqrt{\frac{2}{ln2}}$

Super-Gaussian ( $n=10$ ):

$P_{\text{tot}}=P\left(0\right){\displaystyle \int }exp\left[-\left(ln2\right){\left(\frac{r}{r_{\text{hw}}}\right)}^n\right]2\pi r\text{d}r=P\left(0\right)\pi r_{\text{hw}}^2\cdot \frac{\frac{2}{n}\Gamma \left(\frac{2}{n}\right)}{{\left(ln2\right)}^{\frac{2}{n}}}$

$r_{\text{eff}}=r_{\text{hw}}\cdot \sqrt{\frac{2}{n}\Gamma \left(\frac{2}{n}\right)}{\left(\frac{2}{ln2}\right)}^{\frac{1}{n}}$

Flat-top:

$P_{\text{tot}}={\displaystyle \int P\left(r\right)2\pi r\text{d}r}=P\left(0\right)\pi r_{\text{hw}}^2$

$r_{\text{eff}}=r_{\text{hw}}$

We write the three distributions in terms of $P_{\text{tot}}$ and $r_{\text{eff}}$

Gaussian:

$P\left(r\right)=\frac{2P_{\text{tot}}}{\pi r_{\text{eff}}^2}exp\left[-2{\left(\frac{r}{r_{\text{eff}}}\right)}^2\right]$ (4)

Super-Gaussian:

$P\left(r\right)=\frac{2^{\frac{2}{n}}{P}_{\text{tot}}}{\pi r_{\text{eff}}^2}exp\left[-2{\left(\frac{2}{n}\Gamma \left(\frac{2}{n}\right)\right)}^{\frac{n}{2}}{\left(\frac{r}{r_{\text{eff}}}\right)}^n\right]$ (5)

Flat-top:

$P\left(r\right)=\frac{P_{\text{tot}}}{\pi r_{\text{eff}}^2}\theta \left(r_{\text{eff}}-r\right)$ (6)

Figure 1 compares the 3 beam power distributions at total power $P_{\text{tot}}=\pi$. In the left panel, the HWHM is fixed at $r_{\text{hw}}=1$ for the three distributions. It is clear that $r_{\text{hw}}$ does not give a robust measure of how spread out the beam power is when applied to different distribution types. The right panel plots the three distributions with the (effective mode) radius fixed at $r_{\text{eff}}=1$. The results suggest that $r_{\text{eff}}$ is more robust when applied to different distribution types. Thus, $r_{\text{eff}}$ offers a better measure of how spread out the beam power is. Below, we shall use $r_{\text{eff}}$ to describe the beam spot size.

The reflex time $t_{\text{ref}}$ is determined in Equation (2) with beam power density $P\left(r\right)$ given by one of (4), (5) or (6), specified in terms of the total power $P_{\text{tot}}$ and the beam radius $r_{\text{eff}}$. With the mathematical formulation established, we investigate the effect of distribution type and parameters $\left(P_{\text{tot}}\mathrm{,}{r}_{\text{eff}}\right)$ on the reflex time.

Before we end this section, we clarify the relation between the non-dimensional and the physical versions of power density $P\left(0\right)$, total power $P_{\text{tot}}$, and total energy $E_{\text{tot}}\equiv t_{\text{ref}}{P}_{\text{tot}}$. This clarification is necessary because we have two different length scales in the non-dimensionalization, and there is more than one way to non-dimensionalize these quantities. All scales used in non-dimensionalization are summarized in Table 1.

Figure 1. Comparison of 3 beam power distributions at $P_{\text{tot}}=\pi$. (a) The HWHM is fixed at $r_{\text{hw}}=1$ ; (b) The (effective mode) radius is fixed at $r_{\text{eff}}=1$.

Table 1. Scales used in the non-dimensionalization.

The physical and the non-dimensional versions of $\left(P\left(0\right)\mathrm{,}{P}_{\text{tot}}\mathrm{,}{E}_{\text{tot}}\right)$ are related by

$P_{\text{nd}}\left(0\right)=\frac{P_{\text{phy}}\left(0\right)}{K\mu \Delta T}$

$P_{\text{tot,nd}}=\pi \frac{P_{\text{tot,phy}}}{K{\mu }^2\Delta TE\left\{v_{\text{c}}\right\}}$

$E_{\text{tot,nd}}=\pi \frac{E_{\text{tot,phy}}}{{\rho }_{\text{m}}{C}_p\Delta TE\left\{v_{\text{c}}\right\}}$

3. Results of Non-Dimensional Analysis

We carry out non-dimensional analysis of the deterministic case, $v_{\text{c}}\equiv \pi$. In this section, all quantities are non-dimensional even though they are not denoted so explicitly.

3.1. Effect of Beam Radius Given Total Power and Distribution Type

In this sub-section, we use formulation (2) to examine the effect of beam radius $r_{\text{eff}}$ on the reflex time, when the distribution type of beam power and the total power $P_{\text{tot}}$ are both fixed. The two panels of Figure 2 show plots of radius vs reflex time, respectively, for $P_{\text{tot}}=8$ and $P_{\text{tot}}=2$.

For each value of total power, the reflex time attains a minimum at a certain value of radius, which we shall call the optimal radius and denote it by $r_{\text{opt}}$.

$r_{\text{opt}}\equiv \arg \underset{r_{\text{eff}}}{\min }{t}_{\text{ref}}\left(r_{\text{eff}}\right)$

Given the total power, with a beam radius below or above the optimal radius $r_{\text{opt}}$, it takes longer time for the beam to induce withdrawal reflex on the test

Figure 2. Beam radius vs reflex time for the three distribution types. (a) Total beam power: $P_{\text{tot}}=8$ ; (b) $P_{\text{tot}}=2$.

subject. This is true for all three distribution types. Results in Figure 2 indicate that the optimal radius $r_{\text{opt}}$ varies with the total power.

3.2. Effect of Total Power on Minimization of Reflex Time with Respect to Beam Radius

In Figure 3, we examine how the optimal radius $r_{\text{opt}}$ and the corresponding minimum reflex time $t_{\text{ref}}\left(r_{\text{opt}}\right)$ vary with the total beam power.

In the left panel of Figure 3, the optimal beam radius increases with $P_{\text{tot}}$, and approaches 1 (non-dimensional) at large $P_{\text{tot}}$. In the right panel, $min{t}_{\text{ref}}\left(r_{\text{eff}}\right)$ is approximately inversely proportional to $P_{\text{tot}}$. We fit $t_{\text{ref}}=\frac{c}{P_{\text{tot}}}$ to the results for each distribution and we obtain

Gaussian: $\underset{r_{\text{eff}}}{min}{t}_{\text{ref}}\left(r_{\text{eff}}\right)\approx \frac{10.8}{P_{\text{tot}}}$

Super-Gaussian: $\underset{r_{\text{eff}}}{min}{t}_{\text{ref}}\left(r_{\text{eff}}\right)\approx \frac{8.2}{P_{\text{tot}}}$

Flat-top: $\underset{r_{\text{eff}}}{min}{t}_{\text{ref}}\left(r_{\text{eff}}\right)\approx \frac{7.4}{P_{\text{tot}}}$

The relation $t_{\text{ref}}\propto 1/P_{\text{tot}}$ is true only when the beam radius set to $r_{\text{opt}}\left(P_{\text{tot}}\right)$ at each individual value of $P_{\text{tot}}$ to minimize $t_{\text{ref}}$. The optimal radius varies with $P_{\text{tot}}$, as demonstrated in the left panel of Figure 3. When the beam radius is set to a fixed value for all $P_{\text{tot}}$, the relation $t_{\text{ref}}\propto 1/P_{\text{tot}}$ is no longer valid. This is demonstrated in Figure 4.

3.3. Effect of Total Beam Power and Distribution Type on Energy Consumption by the Time of Reflex

Recall that $P_{\text{tot}}$ is the total power absorbed into the skin. Let $E_{\text{tot}}$ denote the

Figure 3. Minimization of reflex time with respect to beam radius. (a) Total power vs optimal radius; (b) Total power vs reflex time at optimal radius.

Figure 4. Total beam power $P_{\text{tot}}$ vs reflex time $t_{\text{ref}}$, at fixed radius $r_{\text{eff}}=1$.

total energy absorbed into the skin by the time of reflex. It has the expression

$E_{\text{tot}}=t_{\text{ref}}{P}_{\text{tot}}$

The total energy absorbed is proportional to the total energy incident at the skin, which, in turn, may be viewed as being proportional to the total energy consumed at the beam generator. In this simplified model, minimizing the energy consumption at beam generator is equivalent to minimizing $E_{\text{tot}}$. At each value of $P_{\text{tot}}$, total energy $E_{\text{tot}}$ is minimized when $t_{\text{ref}}$ attains its minimum at the optimal beam radius $r_{\text{opt}}\left(P_{\text{tot}}\right)$. Let $E_{\text{tot,opt}}\left(P_{\text{tot}}\right)$ denote the minimum of total energy with respect to beam radius $r_{\text{eff}}$ at the given $P_{\text{tot}}$.

$E_{\text{tot,opt}}\left(P_{\text{tot}}\right)\equiv \underset{r_{\text{eff}}}{min}{t}_{\text{ref}}\left(r_{\text{eff}}\right)P_{\text{tot}}=t_{\text{ref}}\left(r_{\text{opt}}\left(P_{\text{tot}}\right)\right)P_{\text{tot}}$ (7)

If $min{t}_{\text{ref}}\left(r_{\text{eff}}\right)$ is exactly proportional to $1/P_{\text{tot}}$, then $E_{\text{tot,opt}}\left(P_{\text{tot}}\right)$ in (7) would be a constant, independent of $P_{\text{tot}}$. Figure 5 examines the detailed behavior of $E_{\text{tot,opt}}$ vs $P_{\text{tot}}$.

The total energy at the optimal beam radius varies slightly when the total beam power increases from 1/8 to 32 (non-dimensional). The differences among the three distribution types are more prominent than the variations over the range of $P_{\text{tot}}$ within each type. For each of the three distribution types, $E_{\text{tot,opt}}$ attains a minimum around $P_{\text{tot}}=3$ (non-dimensional). More precisely, values of the optimal total beam power $P_{\text{tot,opt}}$, the corresponding optimal beam radius ${r_{\text{opt}}|}_{P_{\text{tot,opt}}}$, and the corresponding minimum energy consumption ${E_{\text{tot,opt}}|}_{P_{\text{tot,opt}}}$ are listed in Table 2 for the three distribution types.

Table 2 tells us that in order to minimize the energy consumption for inducing withdrawal reflex, we should use the flat-top distribution with total beam power = 2.62 (non-dimensional) and beam radius = 0.742 (nom-dimensional). The corresponding energy consumption is 7.17 (non-dimensional). The scales used in non-dimensionalization are summarized in Table 1. The physical value

Figure 5. Total beam power vs energy consumption at the optimal beam radius.

Table 2. Optimal values for minimizing energy consumption.

of the minimum energy consumption for the flat-top beam has the expression

$\underset{\text{Flat-top}}{min}{E}_{\text{phy}}=\frac{7.17}{\pi }{\rho }_{\text{m}}{C}_p\Delta TE\left\{v_{\text{c}}\right\}=2.28{\rho }_{\text{m}}{C}_p\Delta TE\left\{v_{\text{c}}\right\}$ (8)

(8) gives a tight lower bound on the energy needed to induce withdrawal reflex. No matter which of the three distribution types we adopt, no matter what total beam power and beam radius we select, the energy consumption, measured as the total energy absorbed into the skin by the reflex time, cannot be lower than the amount on the right hand side of (8). Note that ${\rho }_{\text{m}}{C}_p\Delta TE\left\{v_{\text{c}}\right\}$ is the amount of energy needed to heat up a 3-D region of volume $E\left\{v_{\text{c}}\right\}$ by $\Delta T$ degree if 1) all energy is used solely to heat the confined region, 2) the heating within the region is uniform, and 3) no energy goes into the surrounding outside the region. In reality, this idealized heating is impossible. The heating within the region won’t be uniform and it is inevitable that some electromagnetic energy will end up heating the surrounding outside the confined region. The coefficient 2.28 in (8) reflects this intrinsic deficiency in utilizing the beam energy, and it represents the best that can be achieved.

If the flat-top beam is not available, we should use the super Gaussian distribution with total beam power = 2.73 and beam radius = 0.746. The corresponding energy consumption is 7.99, not much higher than 7.17 for the flat-top beam. Among the three distribution types, the Gaussian beam has the highest energy consumption (see Table 2). As we will see in the next sub-section, high energy consumption is not the only drawback of the Gaussian beam.

3.4. Minimizing Reflex Time for the Gaussian Beam Leads to a Very Large Surface Temperature Increase

Since the electromagnetic heating decays exponentially with the depth, at any time the largest temperature increase occurs on the skin surface and at the beam center. As a function of time, the temperature at beam center attains its maximum at withdrawal reflex since the electromagnetic heating is discontinued (moved away from the test subject) upon withdrawal reflex. In this sub-section, we focus on this maximum surface temperature increase. The physical temperature increase has the expression

${\left(T\left(\mathrm{0,}{t}_{\text{ref}}\right)-T_{\text{base}}\right)}_{\text{phy}}=T_{\text{nd}}\left(t_{\text{ref,nd}}\right)\Delta T$

where $\Delta T\equiv T_{\text{act}}-T_{\text{base}}$ is the temperature scale (see Table 1). In the case of temperature scale $\Delta T=8˚\text{C}$, for example, a non-dimensional temperature of $T_{\text{nd}}=5$ corresponds to a physical temperature increase of 40˚C. Large temperature increases lead to burn injuries in tests, which we need to avoid by all means.

We examine the non-dimensional surface temperature $T_{\text{nd}}$ when the beam is set to its optimal radius for minimizing the reflex time. Figure 6 plots total beam power vs maximum $T_{\text{nd}}$ at the optimal beam radius $r_{\text{opt}}$. The flat-top distribution has the smallest surface temperature, around 2; the super-Gaussian distribution has a slightly higher surface temperature, around 2.5. In contrast, the Gaussian distribution has the highest surface temperature, well above 5, much

Figure 6. Total beam power vs maximum temperature when the beam radius is set to the optimal radius $r_{\text{opt}}$, for each of the three distribution types.

larger than those of flat-top and super-Gaussian. In the case of temperature scale $\Delta T=8˚\text{C}$, the Gaussian beam with the optimal radius $r_{\text{opt}}$ poses a significant risk of burn injury in tests. For this reason, the Gaussian distribution should be avoided if a better alternative is available in applications. When neither the flat-top nor the super-Gaussian distribution is available, the Gaussian beam should not be set to its optimal radius. Rather, the radius of the Gaussian beam should be selected to ensure $max{T}_{\text{nd}}\le T_{\text{tolerance}}$ where $T_{\text{tolerance}}$ is a prescribed temperature tolerance.

In Figure 7, we study the performance of the Gaussian beam with beam radius set to a fixed multiple of optimal radius ( $\beta r_{\text{opt}}$ ). In particular, we want to investigate how the maximum temperature and the energy consumption change in response to an increase in the radius multiplier ( $\beta$ ). Intuitively, a larger beam radius will produce a more smooth heating and thus result in a smaller maximum temperature at withdrawal reflex. The total energy consumption at beam radius $\beta r_{\text{opt}}$ has the expression

$E_{\text{tot},\beta }\left(P_{\text{tot}}\right)=t_{\text{ref}}\left(\beta r_{\text{opt}}\left(P_{\text{tot}}\right)\right)P_{\text{tot}}$ (9)

The left panel of Figure 7 plots the total power vs the maximum temperature with the beam radius set to $\beta r_{\text{opt}}$, for several fixed values of $\beta \ge 1$. When the radius of the Gaussian beam is increased from $r_{\text{opt}}$ to $1.6r_{\text{opt}}$, the maximum $T_{\text{nd}}$ (non-dimensional) drops from above 5 to below 3, comparable to those of the super-Gaussian beam at its optimal radius (see Figure 6). The right panel of Figure 7 shows the total power vs the energy consumption $E_{\text{tot,}\beta }$, described in (9). When the beam radius is increased to bring down the maximum temperature, the beam heats up a larger surrounding volume and leads to a larger energy consumption. With beam radius $1.6r_{\text{opt}}$, the energy consumption of the Gaussian

Figure 7. Behavior of the Gaussian beam with beam radius $\beta r_{\text{opt}}$ for several values of $\beta \ge 1$. (a) $P_{\text{tot}}$ vs maximum $T_{\text{nd}}$ ; (b) $P_{\text{tot}}$ vs energy consumption $E_{\text{tot,}\beta }$.

beam is around 14 (non-dimensional). In contrast, the super-Gaussian beam ( $n=10$ ) at its optimal radius has a similar maximum temperature but has an energy consumption significantly lower than 14, in the range of 8 - 8.5 (see Figure 5).

4. Concluding Remarks

In this paper, we consider the problem of inducing withdrawal reflex on a test subject by exposing the subject to an electromagnetic beam of a millimeter wavelength. During the exposure to the electromagnetic beam, heat-sensitive nociceptors in the skin are activated wherever the local temperature is above the nociceptor activation temperature. Withdrawal reflex occurs when the combined signal from all nociceptors reaches a certain threshold. In the case of uniform nociceptor density in the skin, withdrawal reflex occurs when the activated volume reaches the volume threshold.

For an electromagnetic beam of given frequency, its thermal effect on the skin is specified by 1) total beam power absorbed into the skin, 2) beam effective model radius, and 3) power distribution type. While the skin material properties and the volume threshold may not be known, we expect these skin parameters to be independent of the beam parameters, which we can vary in designing tests as inout variables. Each test has several key outcome quantities that we like to monitor and/or optimize. These quantities are: a) the time to withdrawal reflex, b) the maximum temperature, and c) the energy consumption by the reflex time. The main task of this study is to analyze the effect of beam parameters on these outcome quantities in a proper setting that is independent of skin material properties.

In the non-dimensionalization process, the length scale in the depth direction is naturally defined by the absorption coefficient ( $\mu$ ), which describes the spatial decay of electromagnetic heating in the depth direction. We introduce a separate length scale for skin surface coordinates based on the volume threshold for withdrawal reflex. The temperature evolution of the 3-D skin is governed by the electromagnetic heating and the heat conduction in the skin depth direction. The heat conduction perpendicular to the depth direction is neglected because the length scale for surface coordinates is much larger than the depth scale. In addition to the two length scales, we introduce scales for the time, temperature increase, and power density. Using these scales, we carry out non-dimensionalization. The non-dimensional temperature is expressed as the product of the non-dimensional beam power density and a parameter-free function. This form of analytical solution provides a suitable mathematical framework for investigating the effect of beam parameters without the interference of skin parameters.

We examine three distribution types for beam power density: Gaussian, super-Gaussian ( $n=10$ ), and flat-top. For each distribution type, at each value of total beam power, the reflex time varies with the beam radius and attains its minimum at a certain radius, called the optimal radius. The optimal radius increases with the total power and approaches 1 at large total power. The minimum reflex time is approximately inversely proportional to the total power. We use the total energy absorbed into the skin by the reflex time as a measure of the energy consumption at the beam generator. The energy consumption with the optimal beam radius varies slightly over the range of total beam power, attaining a minimum around total power = 3.0 (non-dimensional). The inter-distribution-type differences of energy consumption are much more prominent than the variations within each type. Among the three distribution types, the flat-top beam has the lowest energy consumption, attaining a minimum of 7.17 (non-dimensional) with respect to total beam power and beam radius. This minimum energy consumption gives a tight lower bound on the energy needed to induce withdrawal reflex. No matter which of the three distribution types we adopt, no matter what total beam power and beam radius we select, the energy consumption cannot be lower than 7.17 (non-dimensional).

Among the three distribution types, the Gaussian beam has the worst performance in two aspects: a) the Gaussian beam has the highest energy consumption even with the optimal radius; b) the Gaussian beam, when set to its optimal radius, produces a very large temperature increase at withdrawal reflex, which may lead to burn injury in tests. The temperature increase can be attenuated by selecting a beam radius larger than the optimal radius. But doing so will further increase the energy consumption. In summary, for the purpose of minimizing energy consumption and minimizing burn injury risk, the Gaussian beam should be avoided. If the Gaussian beam is the only option available to bring down its maximum temperature to match that of the super-Gaussian beam ( $n=10$ ), the radius of the Gaussian beam should be set to about 1.6 times the optimal radius.

Acknowledgements and Disclaimer

The authors acknowledge the Joint Intermediate Force Capabilities Office of U.S. Department of Defense and the Naval Postgraduate School for supporting this work. The views expressed in this document are those of the authors and do not reflect the official policy or position of the Department of Defense or the U.S. Government.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Topfer, F. and Oberhammer, J. (2015) Millimeter-Wave Tissue Diagnosis: The Most Promising Fields for Medical Applications. IEEE Microwave Magazine, 16, 97-113. https://doi.org/10.1109/MMM.2015.2394020
[2]	Zhadobov, M., Chahat, N., Sauleau, R., Le Quement, C. and Le Drean, Y. (2011) Millimeter-Wave Interactions with the Human Body: State of Knowledge and Recent Advances. International Journal of Microwave and Wireless Technologies, 3, 237-247. https://doi.org/10.1017/S1759078711000122
[3]	Walters, T.J., Blick, D.W., Johnson, L.R., Adair, E.R. and Foster, K.R. (2000) Heating and Pain Sensation Produced in Human Skin by Millimeter Waves: Comparison to a Simple Thermal Model. Health Physics, 78, 259-267. https://doi.org/10.1097/00004032-200003000-00003
[4]	Wang, H., Burgei, W.A. and Zhou, H. (2020) A Concise Model and Analysis for Heat-Induced Withdrawal Reflex Caused by Millimeter Wave Radiation. American Journal of Operations Research, 10, 31-81. https://doi.org/10.4236/ajor.2020.102004
[5]	Cazares, S.M., Snyder, J.A., Belanich, J., Biddle, J.C., Buytendyk, A.M., Teng, S.H.M. and O’Connor, K. (2019) Active Denial Technology Computational Human Effects End-to-End Hypermodel for Effectiveness (ADT CHEETEH-E). Human Factors and Mechanical Engineering for Defense and Safety, 3, Article No. 13. https://doi.org/10.1007/s41314-019-0023-7
[6]	Wang, H., Burgei, W.A. and Zhou, H. (2018) Dose-Injury Relation as a Model for Uncertainty Propagation from Input Dose to Target Dose. American Journal of Operations Research, 8, 360-385. https://doi.org/10.4236/ajor.2018.85021
[7]	Paschotta, R. (2008) Effective Mode Area. In: The Encyclopedia of Laser Physics and Technology, Wiley-VCH, Weinheim.