Electromagnetic Wave Propagation Fundamental Field Equations

#### **Chapter 1**

## Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical Structures with Impedance Boundary Condition

*Shiva Hayati Raad and Zahra Atlasbaf*

### **Abstract**

The integral equation (IE) method is one of the efficient approaches for solving electromagnetic problems, where dyadic Green's function (DGF) plays an important role as the Kernel of the integrals. In general, a layered medium with planar, cylindrical, or spherical geometry can be used to model different biomedical media such as human skin, body, or head. Therefore, in this chapter, different approaches for the derivation of Green's function for these structures will be introduced. Due to the recent great interest in two-dimensional (2D) materials, the chapter will also discuss the generalization of the technique to the same structures with interfaces made of isotropic and anisotropic surface impedances. To this end, general formulas for the dyadic Green's function of the aforementioned structures are extracted based on the scattering superposition method by considering field and source points in the arbitrary locations. Apparently, by setting the surface conductivity of the interfaces equal to zero, the formulations will turn into the associated problem with dielectric boundaries. This section will also aid in the design of various biomedical devices such as sensors, cloaks, and spectrometers, with improved functionality. Finally, the Purcell factor of a dipole emitter in the presence of the layered structures will be discussed as another biomedical application of the formulation.

**Keywords:** biomedical, dyadic Green's function, integral equation, 2D materials, sensors, Purcell factor

#### **1. Introduction**

Planarly, cylindrically and spherically layered media have been widely used to model the human skin, body, or head. In particular, a rectangular slab is proposed to model and analyze skin temperature distribution [1]. Moreover, in the multilayered skin model, three stacked layers are exploited to simulate the performance of the epidermis, dermis, and sub-cutis parts of the skin [2]. In other research, a planarly layered medium has been proposed as the simplified human body model by considering the impact of skin, fat, and muscle in the electromagnetic performance [3]. Cylindrical-shaped equivalent phantom of the skin is alternatively used to characterize the interactions between an antenna and the human body [4]. For a more precise investigation, multilayered cylinders are proposed to model a biological system with different tissues [5]. Considering spherical geometrics, the interactions of a

three-layered spherical human head model with a finite-length dipole and a cellular phone helical antenna are investigated using Green's function [6, 7]. The research is extended to the six-layered model, including skin, fat, bone, dura, CSF, and brain layers of the head [8]. This model has also been used to explore the impact of the shapes and positions of coils in the MRI system [9]. For the analysis of the aforementioned structures, integral equation (IE) based methods have been widely used. Dyadic Green's function, which is a powerful method to calculate the electromagnetic response for different excitation sources, plays an essential role in the IE method [10–12].

Different mathematical approaches can be utilized to calculate the dyadic Green's function. The dyadic Green's function in a homogeneous environment can be expressed in terms of vector functions *M*, *N*, and *L*. The Green's function is singular when the field and observation points coincide. In this case, Green's function can be written as the main part plus the portion proportional to the impulse function. This type of decomposition of the Green's function is not unique and depends on the shape of the volume that separates the environment [13]. In general, to obtain the dyadic Green's function in a layered medium, the Green's function of the homogeneous environment can be used. Then, the effect of inhomogeneity can be considered by adding reflected and transmitted waves to the dyadic Green's function [14]. In a similar approach, called the scattering superposition method, the scattering Green's function in each layer is expressed in terms of vector wave functions with the unknown coefficients that are obtained by applying the boundary conditions [15]. Moreover, to calculate the Green's function of a medium by impedance method, instead of each dielectric layer a transmission line, and instead of each metallic layer a current source can be considered in the equivalent model in the rectangular, cylindrical, and spherical coordinates [15–17].

Considering medical diagnostics and treatment, planarly, cylindrically and spherically layered media are engineered mainly with plasmonic materials. In this regard, the interaction of an environment with a nanometer-scale dipole emitter is of interest in different biomedical fields. For example, the optical activity of the proteins can be investigated using Green's tensor approach. Furthermore, a single excited molecule in the vicinity of a metallic structure can be potentially used in the sensors because of behaving as a resonant filter. Moreover, the point source is an appropriate model for the concentrated light sources which are used in medical applications [18]. Following these trends, different plasmonic structures have been presented theoretically and realized in the real environment. For example, the effect of the size of the gold nanoparticles on the decay rate and the energy transfer of dipole emitters is investigated using Green's tensor formulation and compared with the results obtained from the fabrication [19]. The same analysis is carried out for a metallic cylinder coated with a dielectric layer and the resulted integrals are solved numerically [20]. Also, Fermat's golden rule is used to connect the imaginary part of the Green's function to the radiation impedance of the dipole antenna adjacent to a medium, and the Purcell factor is extracted for the planar meta-surface, which is then experimentally characterized in the microwave frequency band [21].

Recently, graphene's plasmons are proposed as the low loss and reconfigurable alternative to the plasmons of the noble metals. There are two approaches to use graphene in analytical and numerical methods. In the first one, the graphene boundary is modeled with a 2D surface, characterized by its surface conductivity, whereas in the second one, the graphene layer is replaced by a very thin dielectric [22]. Although the latter can be analyzed using the available formulas previously presented for the multilayered dielectric structures in the literature [13, 23], graphene analysis using surface conductivity model has many advantages. First, in the dielectric model, each graphene interface adds an extra layer to the structure. Therefore, when the number of graphene layers is large, the problem becomes very complicated. Secondly,

#### *Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

in the dielectric model, it is necessary to compute the special functions of the cylindrical and spherical coordinates with complex arguments, which requires the implementation of specific algorithms for their effective calculation [24], while in the impedance boundary condition method, the surface conductivity of graphene appears as a coefficient for special functions. Third, in the dielectric model, the thickness of the graphene layer is considered to be about 0.335 nm, which is often very small compared to other geometrical parameters and makes the convergence of the analytical functions slow [25]. Also, if the goal is using the dielectric model in numerical methods, it is necessary to use a dense mesh for the equivalent dielectric of the graphene that is not optimal in terms of time and memory [26].

Applying the graphene surface conductivity model for the derivation of the Green's function has been considered in recent years. For example, for the graphene sheet under electric bias, the dyadic Green's function is derived by the Hertzian potential and plane-wave expansion methods and the corresponding integrals are solved with the saddle point method [27, 28]. Also, the method is expanded for the analysis of graphene with tensor surface conductivity which can be used to analyze graphene with magnetic bias or spatial dispersion. Romberg's integration procedure is proposed for the numerical solution of the resulting integrals [29]. In another research, the analysis of the electric dipole in the proximity of the parallel plate waveguide with graphene walls is studied by extracting the Green's function and calculating the spontaneous emission. It is found that symmetric and asymmetric plasmonic modes lead to a sharp increase in this parameter [30]. As another instance, a point source is taken into account in the vicinity of the infinite cylinder with graphene cover and it is observed that by changing the distance of the source from the cylinder as well as changing the chemical potential of the graphene layer, the Fano resonances can be controlled [31]. In this chapter, the dyadic Green's function of various planar, cylindrical, and spherical geometries with impedance boundary conditions will be calculated using the scattering superposition method. Specifically, we have focused on the graphene material due to its wide range of applications. Apparently, another 2D material can be considered by replacing the graphene surface conductivity with the surface conductivity of the desired material. The presented formulas can be potentially sued to design various biomedical devices. Moreover, by approaching the surface impedance to zero, these structures can be potentially used to investigate the interaction of the human body with different electromagnetic sources. For the Green's function calculation of complex media in different coordinates using vector wave functions, the reader is referred to [32–35].

#### **2. Surface conductivity of graphene material under different conditions**

Graphene is a two-dimensional material made of carbon atoms and can be considered in the solution of Maxwell's equations using surface conductivity boundary condition [23]. Depending on the geometrical and optical conditions, graphene surface conductivity can be isotropic or anisotropic. The purpose of this section is to provide an overview of the graphene surface conductivity under different conditions (electric bias, magnetic bias, and spatial dispersion) and for different geometries (continuous or patterned sheets).

It should be noted that graphene material is mainly synthesized through four methods, including 1) mechanical exfoliation of highly ordered pyrolytic graphite (HOPG), 2) the epitaxial growth of graphene on silicon carbide (SiC), 3) the reduction of graphene oxide, and 4) chemical vapor deposition (CVD) technique. The comparison of these methods in terms of quality and the fabricated area is provided in **Table 1** [36].


**Table 1.**

*Comparison of the quality and area of synthesized graphene using different techniques [36].*

In this regard, a monolayer graphene film with metallic electrodes is transferred to a high impedance surface (HIS) which is realized by the printed circuit board (PCB) technology at microwave frequencies. The measurement is conducted in a small microwave chamber [37]. In another research in the same spectrum, an infrared laser is used to etch the CVD grown graphene sheet with predefined periodicity and later investigate the absorption of the designed structure inside a rectangular waveguide [38]. The electron beam lithography is another approach used for pattering the graphene sheet for enhanced light matter interaction [39]. Moreover, a transparent graphene millimeter wave absorber constructed by multiple transfer-etch processing is characterized by reflectometery technique at 140 GHz [40]. Also, CVD-grown graphene is used to enhance the sensitivity of the surface-enhanced Raman spectroscopy (SERS)-based chemical sensor. The measurement is done using a Raman spectrometer at the laser wavelength of 785 nm (red) [41]. In another sensor chip, DNA is hybridized to the graphene-based substrate under UV light with the wavelength of 260 nm. In this sensor, the atomic force microscopy (AFM) is used to ensure the continuity and uniformity of the synthesized graphene, and Raman characterization is used to investigate its quality and the number of layers [42].

#### **2.1 Graphene material under electric bias**

When a graphene sheet is under electric bias, its surface conductivity is isotropic and can be approximately calculated by using Kubo's formulas as [43]:

$$\sigma\_{\rm intra} = \frac{2ie^2k\_BT}{\hbar^2\pi(\alpha + i/\pi)}\ln\left[2\cosh\left(\frac{\mu\_c}{2k\_BT}\right)\right] \tag{1}$$

$$\sigma\_{\text{inter}} = \frac{e^2}{4\hbar} \times \left[ \frac{1}{2} + \frac{1}{\pi} \arctan\left(\frac{\hbar\omega - 2\mu\_c}{2k\_B T}\right) - \frac{i}{2\pi} \ln\left(\frac{(\hbar\omega + 2\mu\_c)^2}{(\hbar\omega - 2\mu\_c)^2 + (2k\_B T)^2}\right) \right] \tag{2}$$

In the above equations,*T* is the temperature, *μ<sup>c</sup>* is the chemical potential of graphene, ℏ is the reduced Planck's constant, *KB* is the Boltzmann's constant, and ω is the angular frequency. At low-THz frequencies, the inter-band contribution of the surface conductivity can be neglected. Also, the graphene layer can be modeled as a dielectric with very low thickness *δ* and equivalent dielectric constant *ε* ¼ <sup>1</sup> � *<sup>i</sup> <sup>σ</sup> ωε*0*<sup>δ</sup>* [22]. Note that in the above equations, graphene is assumed to be in the linear region, otherwise, other terms proportional to 1*=ω*<sup>3</sup> and 1*=ω*<sup>4</sup> should be added respectively for the frequency range of ℏ*ω* <2*μ<sup>c</sup>* and ℏ*ω*≥2*μ<sup>c</sup>* [44].

To increase the light-matter interaction of graphene material, the nano-pattering method has been proposed [45]. The surface conductivity of the periodic graphene

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

**Figure 1.**

*(a) Graphene square patches with the periodicity of* D*<sup>0</sup> and the air gap distance of* g *and (b) densely packed graphene strips with the width* W *and the periodicity* L *[46].*

elements under the electric bias is also isotropic. For the square patches shown in **Figure 1(a)**, by considering the periodicity of *D*<sup>0</sup> and the air gap distance of *g*, closed-form surface impedance (inverse of the surface conductivity) is as [47]:

$$Z\_s = \frac{D\_0}{\sigma\_\circ (D\_0 - \text{g})} + i \frac{\pi}{2o\epsilon \epsilon\_0 \left(\frac{c\_2 + 1}{2}\right) D\_0 \ln\left[\csc\left(\frac{\pi \text{g}}{2D\_0}\right)\right]} \tag{3}$$

The effective surface conductivity of periodic graphene elements with arbitrary pattern can be measured or extracted using the parameter retrieval method through full-wave simulation [48]. It should be noted that for computing the electric field required for each considered chemical potential, approximate equations can be derived as [49]:

$$\mu\_c(eV) = \begin{cases} \lambda\_1 E\_0^{\lambda\_2}, E\_0(V/nm) \ge 0 \\ -\lambda\_1 E\_0^{\lambda\_2}, E\_0(V/nm) < 0 \end{cases} \tag{4}$$

where, *λ*<sup>1</sup> ¼ 0*:*3677 and *λ*<sup>2</sup> ¼ 0*:*5010. For the chemical potentials in the range of [�1,1] eV, the required bias fields are in the order of several volts per nanometer, which can be implemented practically [29].

The surface conductivity of densely packed graphene strips, as illustrated in **Figure 1(b)**, is anisotropic and can be approximated in the form of a diagonal tensor using the effective medium formulation as [50]:

$$\begin{aligned} \sigma\_{\infty} &= \frac{W \sigma \sigma\_C}{L \sigma\_C + W \sigma} \\ \sigma\_{\mathcal{W}} &= \sigma \frac{W}{L} \\ \sigma\_{\text{xy}} &= \sigma\_{\text{yx}} = \mathbf{0} \end{aligned} \tag{5}$$

In the above equations, *W* and *L* are the width and periodicity of the strips, respectively. Also, *σ* is the surface conductivity of graphene under electric bias, and *σ<sup>c</sup>* is the static conductivity of the surface. Two very important properties of this environment are the existence of near-zero surface conductivity and hyperbolic dispersion region with the potential applications in 2D lens structures and the spontaneous emission enhancement of the dipole emitters [50, 51].

#### **2.2 Graphene sheet with spatial dispersion effects**

When the graphene sheet is placed on a substrate with a high dielectric constant, its surface conductivity is a tensor in which the elements depend on the wave propagation constant in the structure. This is called the spatial dispersion effect and the associated formulae for calculating surface conductivity are [29]:

$$
\sigma\_{\infty} = \sigma + a \frac{d^2}{d\mathfrak{x}^2} + \beta \frac{d^2}{d\mathfrak{y}^2} \tag{6}
$$

$$
\sigma\_{\mathcal{Y}} = \sigma + \beta \frac{d^2}{d\mathbf{x}^2} + \alpha \frac{d^2}{d\mathbf{y}^2} \tag{7}
$$

$$
\sigma\_{\text{xy}} = \sigma\_{\text{yx}} = 2\beta \frac{d^2}{d\text{xdy}} \tag{8}
$$

The parameters *σ*, *α*, and *β* are extracted for an unbiased sheet ð Þ *μ<sup>c</sup>* ¼ 0 thorough perturbation theory [29]. The resulting equations are valid for the electrically biased sheet ð Þ *μ<sup>c</sup>* 6¼ 0 , as well [52].

#### **2.3 Graphene sheet under magnetic bias**

When a graphene sheet is under magnetic bias, its surface conductivity is also a tensor. The diagonal elements of this tensor are equal and the off-diagonal elements are opposite in sign defined as [53].

$$\sigma\_{\text{xx}}(o\omega, B\_0) = \sigma\_0 \frac{1 - io\sigma}{\left(o\iota\_\iota\pi\right)^2 + \left(1 - io\sigma\right)^2} \tag{9}$$

$$\sigma\_{\rm yx}(o\sigma, B\_0) = \sigma\_0 \frac{o\iota\_\mathfrak{t}\mathfrak{r}}{\left(o\iota\_\mathfrak{t}\mathfrak{r}\right)^2 + \left(\mathfrak{1} - i\iota\mathfrak{r}\right)^2} \tag{10}$$

where *<sup>σ</sup>*<sup>0</sup> <sup>¼</sup> <sup>2</sup>*e*2*<sup>τ</sup> <sup>π</sup>*ℏ<sup>2</sup> *kBT* ln 2 cosh *<sup>μ</sup><sup>c</sup>* 2*kBT* � � and *<sup>ω</sup><sup>c</sup>* <sup>¼</sup> *eB*0*ν*<sup>2</sup> *F μc* . The approximate formulas for the calculation of the surface impedance of the square graphene elements under magnetic bias, shown in **Figure 1(a)**, are as follows [54]:

$$\overline{\overline{Z}}\_p = F\_\text{G} \overline{\overline{Z}}\_\text{g} + \frac{i}{2a} \sqrt{\frac{\mu\_0}{\varepsilon\_0 \varepsilon\_\text{eff}}} \begin{bmatrix} 1 & 0\\ 0 & 1 \end{bmatrix} \tag{11}$$

where *FG* <sup>¼</sup> <sup>0</sup>*:*6½ � *<sup>D</sup>*0*=*ð Þ *<sup>D</sup>*<sup>0</sup> � *<sup>g</sup>* <sup>3</sup> <sup>þ</sup> <sup>0</sup>*:*4 and *<sup>α</sup>* ¼ � *<sup>k</sup>*0*D*<sup>0</sup> ffiffiffiffi *εeff* p *<sup>π</sup>* ln sin *<sup>π</sup><sup>g</sup>* 2*D*<sup>0</sup> h i � � . Similar to the electrically biased patterned elements with arbitrary shapes, the parameter retrieval method can be used under the applied magnetic bias [55]. Given the discussion of the above three sections, it is observed that assuming the surface conductivity of graphene as:

$$
\overline{\overline{\sigma}} = \begin{bmatrix}
\sigma\_{\text{xx}} & \sigma\_{\text{xy}} \\
\sigma\_{\text{yx}} & \sigma\_{\text{yy}}
\end{bmatrix} \tag{12}
$$

All items expressed above can be extracted as a special case.

#### **3. Analysis of graphene-based structures using dyadic Green's function**

Dyadic Green's functions for the planarly, cylindrically, and spherically layered structures with graphene interfaces will be derived in this section. To this end, the

boundary conditions of the continuity and discontinuity of tangential electric and magnetic fields are respectively satisfied regarding the considered surface conductivity model for the graphene.

#### **3.1 Graphene-based planarly layered media**

This section aims to obtain dyadic Green's functions for planar structures with graphene boundaries. This problem can be solved either in the rectangular or cylindrical coordinates. In the first sub-section, a graphene sheet with the tensor surface conductivity boundary condition (TSCBC) is considered and its dyadic Green's function is calculated in the rectangular coordinates. In this case, the anisotropy of the surface impedance causes the coupling of the transverse electric (TE) and transverse magnetic (TM) fields. As a result of coupling, the number of unknown coefficients in the expansion of the dyadic Green's function is increased concerning the electrically biased sheet. In the second part of this section, a graphene-dielectric stack with an arbitrary number of layers is investigated considering the electric bias for the graphene sheets. This problem is solved in the cylindrical coordinates to simplify the calculation of the resulted Sommerfeld integrals.

#### *3.1.1 Graphene sheet with the tensor surface conductivity boundary condition*

The purpose of this section is to obtain the dyadic Green's function of a graphene sheet with the tensor surface conductivity in the interface of half-spaces, as shown in **Figure 2(a)**. The constitutive parameters of the top and bottom regions are considered as (*ε*1, *μ*1) and (*ε*2, *μ*2), respectively. Without losing the generality of the problem, the source is assumed to be in the first environment and the graphene boundary is considered in *z* = 0 interface. Dyadic Green's function of this structure will be calculated using the scattering superposition method. For this purpose, the Green's function in each region of the problem is written in the form of the Green's function in the absence and presence of structure. Thus [23]:

$$
\overline{\overline{\mathbf{G}}}\_{\varepsilon}^{(11)} = \overline{\overline{\mathbf{G}}}\_{\varepsilon 0} \left( \overline{\overline{\mathbf{R}}}, \overline{\mathbf{R}'} \right) + \overline{\overline{\mathbf{G}}}\_{\varepsilon \varepsilon}^{(11)} \left( \overline{\overline{\mathbf{R}}}, \overline{\mathbf{R}'} \right) \quad z > 0 \tag{13}
$$

$$
\overline{\overline{\mathbf{G}}}\_{e}^{(21)} = \overline{\overline{\mathbf{G}}}\_{e}^{(21)} \left( \overline{\mathbf{R}}, \overline{\mathbf{R}'} \right) \quad z < 0 \tag{14}
$$

where *G* ð Þ 11 *es* and *G* ð Þ 21 *es* are respectively the scattering Green's function in regions 1 and 2 and they are expanded in terms of *M* and *N* vector wave functions.

**Figure 2.**

*(a) Graphene sheet with the tensor surface conductivity boundary condition at the interface of half-spaces [46] and (b) graphene-dielectric stack.*

Also, *G<sup>e</sup>*<sup>0</sup> is the free-space Green's function which can be computed using the *Gm* method as:

$$\overline{\overline{\mathbf{G}}}\_{e0}\left(\overline{\mathbf{R}}, \overline{\mathbf{R}'}\right) = -\frac{1}{\left(k\_1\right)^2} \hat{\mathbf{z}} \hat{\mathbf{z}} \,\delta\left(\overline{\mathbf{R}}, \overline{\mathbf{R}'}\right) + \int\_{-\infty}^{\infty} \int\_{-\infty}^{\infty} dk\_x dk\_y \,\mathbf{C}^{(1)}$$

$$\times \left\{ \overline{\mathbf{M}}(-h\_1) \overline{\mathbf{M}}'(h\_1) + \overline{\mathbf{N}}(-h\_1) \overline{\mathbf{N}}'(h\_1) \right\} \,\,\mathbf{z} < \mathbf{z}' \tag{15}$$

where *<sup>C</sup>*ð Þ<sup>1</sup> <sup>¼</sup> *<sup>i</sup>* 8*π*2*h*<sup>1</sup> *k*<sup>2</sup> *<sup>x</sup>*þ*k*<sup>2</sup> ð Þ*<sup>y</sup>* . Using the scalar wave function *ψ k* � � <sup>¼</sup> exp *ikxx* <sup>þ</sup> *ikyy* <sup>þ</sup> *ikzz* � �, it can be readily found that:

$$\overline{\mathbf{M}}(\overline{k}) = i \left( k\_{\gamma} \hat{\mathbf{x}} - k\_{\mathbf{x}} \hat{\mathbf{y}} \right) \boldsymbol{\upmu} \left( \overline{k} \right) \tag{16}$$

$$\overline{\mathbf{N}}(\overline{\mathbf{k}}) = \left( -\frac{k\_x}{k\_j} \left( k\_x \hat{\mathbf{x}} + k\_\gamma \hat{\mathbf{y}} \right) + \frac{k\_\times^2 + k\_\gamma^2}{k\_j} \hat{\mathbf{z}} \right) \boldsymbol{\upmu}(\overline{\mathbf{k}}) \tag{17}$$

The parameters *kx*, *ky*, and *kz* are the wavenumbers in *x*, *y*, and *z* directions, respectively, and *kj* shows the wavenumber for *j* = 1, 2. These parameters are not independent and are related to each other via *kz* ¼ � ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi *k*2 *<sup>j</sup>* � *<sup>k</sup>*<sup>2</sup> *<sup>x</sup>* � *<sup>k</sup>*<sup>2</sup> *y* q ¼ � ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi *k*2 *<sup>j</sup>* � *<sup>k</sup>*<sup>2</sup> *ρ* q . The vector wave function *M* represents the electric field of TE modes and the vector wave function *N* shows the electric field of the TM waves. In the structure under consideration, the anisotropy of surface conductivity leads to the coupling of TE and TM fields. Therefore [56]:

$$
\begin{split}
\overline{\overline{\mathbf{G}}}\_{\text{es}}^{(11)}\left(\overline{\mathbf{R}},\overline{\mathbf{R}}'\right) &= \int^{\stackrel{\text{es}}{}} \int^{\text{es}} dk\_{\text{x}} dk\_{\text{y}} \mathbf{C}^{(1)} \times \left[ \left( a\_{1} \overline{\mathbf{M}}(h\_{1}) + a\_{1}' \overline{\mathbf{N}}(h\_{1}) \right) \overline{\mathbf{M}}'(h\_{1}) \right. \\ &\left. + \left( b\_{1} \overline{\mathbf{N}}(h\_{1}) + b\_{1}' \overline{\mathbf{M}}(h\_{1}) \right) \overline{\mathbf{N}}'(h\_{1}) \right] \end{split} \tag{18}
$$

$$\begin{split} \overline{\overline{\mathbf{G}}}\_{\alpha}^{(21)} \left( \overline{\mathbf{R}}, \overline{\mathbf{R}}' \right) &= \int\_{-\infty}^{\infty} \int dk\_{x} dk\_{y} \mathbf{C}^{(1)} \\ &\times \left[ \left( a\_{2} \overline{\mathbf{M}}(-h\_{2}) + a\_{2}' \overline{\mathbf{N}}(-h\_{2}) \right) \overline{\mathbf{M}}'(h\_{1}) + \left( b\_{2} \overline{\mathbf{N}}(-h\_{2}) + b\_{2}' \overline{\mathbf{M}}(-h\_{2}) \right) \overline{\mathbf{N}}'(h\_{1}) \right] \end{split} \tag{19}$$

The unknown coefficients *a*1, *a*<sup>0</sup> 1*, b*1*, b*<sup>0</sup> <sup>1</sup>*, a*2, *a*<sup>0</sup> 2, *b*2*,* and *b*<sup>0</sup> <sup>2</sup> will be obtained by applying the boundary conditions. Using the self and mutual orthogonality of the vector wave functions, the above equations can be divided into two systems of equations, each with four unknown coefficients. The boundary conditions on the electric and magnetic Green's functions respectively state that:

$$
\hat{\mathbf{z}} \times \left( \overline{\overline{\mathbf{G}}}\_{\varepsilon}^{(11)} - \overline{\overline{\mathbf{G}}}\_{\varepsilon}^{(21)} \right) = \mathbf{0} \tag{20}
$$

$$\hat{\mathbf{z}} \times \left( \frac{\nabla \times \overline{\overline{\mathbf{G}}}\_{\varepsilon}^{(11)}}{i o \mu\_1} - \frac{\nabla \times \overline{\overline{\mathbf{G}}}\_{\varepsilon}^{(21)}}{i o \mu\_2} \right) = \overline{\overline{\sigma}} \,\overline{\overline{\mathbf{G}}}\_{\varepsilon}^{(21)} \tag{21}$$

After applying the above boundary conditions and removing the coupling effect from the tangential components of the electric field, and by defining, *A* ¼ *σxxkx* þ *σxyky*, *B* ¼ *σyxkx* þ *σyyky*, *C* ¼ *σxxky* � *σxykx*, and, *D* ¼ *σyykx* � *σyxky*,

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

the unknown coefficients of the TE waves, *<sup>a</sup>*<sup>1</sup> <sup>¼</sup> <sup>Δ</sup>ð Þ<sup>1</sup> *TE* <sup>Δ</sup>*TE* and *a*<sup>0</sup> <sup>1</sup> <sup>¼</sup> <sup>Δ</sup>ð Þ<sup>2</sup> *TE* <sup>Δ</sup>*TE*, can be obtained as [46]:

$$\Delta\_{\rm TE} = ik\_{\rho}^2 P^+ Q^+ + P^+ \left(Ak\_\mathbf{x} + Bk\_\mathbf{y}\right) i\rho \frac{h\_\mathbf{1}}{k\_\mathbf{1}} + Q^+ i\rho \left(\mathbf{C}k\_\mathbf{y} + Dk\_\mathbf{x}\right) + \alpha^2 i \frac{h\_\mathbf{1}}{k\_\mathbf{1}} (\mathbf{BC} + \mathbf{AD}) \tag{22}$$

$$\Delta\_{\rm TE}^{(1)} = ik\_{\rho}^2 Q^+ P^- + \left(Ak\_x + Bk\_\gamma\right) P^- i\nu \frac{h\_1}{k\_1} - ai\left(Ck\_\gamma + Dk\_x\right) Q^+ - \alpha^2 i \frac{h\_1}{k\_1} (BC + AD) \tag{23}$$

$$
\Delta\_{\rm TE}^{(2)} = \alpha \frac{2h\_1}{\mu\_1} \left( Dk\_\mathcal{y} - Ck\_\mathcal{x} \right) \tag{24}
$$

where, *<sup>P</sup>*<sup>þ</sup> <sup>¼</sup> *<sup>h</sup>*<sup>1</sup> *μ*1 <sup>þ</sup> *<sup>h</sup>*<sup>2</sup> *μ*2 , *<sup>Q</sup>*<sup>þ</sup> <sup>¼</sup> *<sup>k</sup>*<sup>1</sup> *μ*1 <sup>þ</sup> *<sup>k</sup>*<sup>2</sup> <sup>2</sup>*h*<sup>1</sup> *k*1*h*2*μ*<sup>2</sup> , *<sup>P</sup>*� <sup>¼</sup> *<sup>h</sup>*<sup>1</sup> *μ*1 � *<sup>h</sup>*<sup>2</sup> *μ*2 , and *<sup>Q</sup>*� <sup>¼</sup> *<sup>k</sup>*<sup>1</sup> *μ*1 � *<sup>k</sup>*<sup>2</sup> <sup>2</sup>*h*<sup>1</sup> *k*1*h*2*μ*<sup>2</sup> . By defining the TM waves expansion coefficients as *<sup>b</sup>*<sup>1</sup> <sup>¼</sup> <sup>Δ</sup>ð Þ<sup>1</sup> *TM* <sup>Δ</sup>*TM* and *b*<sup>0</sup> <sup>1</sup> <sup>¼</sup> <sup>Δ</sup>ð Þ<sup>2</sup> *TM* <sup>Δ</sup>*TM*, it can be shown:

$$
\Delta\_{\rm TM} = ik\_{\rho}^{2}Q^{+}P^{+} + i\nu Q^{+} \left(\mathrm{Ck}\_{\gamma} + Dk\_{\mathrm{x}}\right) + i\nu \frac{h\_{1}}{k\_{1}}P^{+} \left(Ak\_{\mathrm{x}} + Bk\_{\gamma}\right) + i\nu^{2} \frac{h\_{1}}{k\_{1}}(AD + BC) \tag{25}
$$

$$
\Delta\_{\rm TM}^{(1)} = -ik\_{\rho}^{2}Q^{-}P^{+} + i\nu \frac{h\_{1}}{k\_{1}}P^{+} + \left(Ak\_{\mathrm{x}} + Bk\_{\mathrm{y}}\right) - i\nu \left(Dk\_{\mathrm{x}} + Ck\_{\gamma}\right)Q^{-}
$$

$$
+ i\nu \frac{h\_{1}}{k\_{1}}(AD + BC) \tag{26}
$$

$$
\Delta\_{TM}^{(2)} = -2a \frac{h\_1}{\mu\_1} \left( Ak\_\gamma - Bk\_x \right) \tag{27}
$$

Other unknown coefficients can be obtained using decoupling equations in [46]. To validate the obtained coefficients, the structure of **Figure 1(a)** consisting of square patches with *D*<sup>0</sup> = 5 μm, *g* = 0.5 μm, *μ<sup>c</sup>* = 0.5 eV and τ = 0.5 ps with plane wave illumination is considered under electric and magnetic biases. Since Green's function coefficients are the same as reflection and transmission coefficients of the plane wave, the results of Green's function are compared with the results of the circuit model as [54]:

$$\text{LS}\_{\text{31}} = \frac{\text{2}(\text{2} + \eta\_0 \sigma\_d)}{\text{(2} + \eta\_0 \sigma\_d)^2 + \left(\eta\_0 \sigma\_0\right)^2} \tag{28}$$

$$\mathcal{S}\_{\mathsf{41}} = \frac{2\eta\_0 \sigma\_0}{\left(2 + \eta\_0 \sigma\_d\right)^2 + \left(\eta\_0 \sigma\_0\right)^2} \tag{29}$$

**Figure 3** shows the magnitude of the transmission coefficient by considering the electric bias for the graphene layer. The results of the two methods are identical,

**Figure 3.**

*The magnitude of the transmission coefficient for the graphene nano-patch with the parameters* D*<sup>0</sup> = 5 μm,* g *= 0.5 μm,* μ<sup>c</sup> *= 0.5 eV, and τ = 0.5 ps [46].*

#### **Figure 4.**

*(a) The magnitude and (b) phase of the transmission coefficient for the graphene nano-patch with* D*<sup>0</sup> = 5 μm,* g *= 0.5 μm,* μ<sup>c</sup> *= 0.5 eV, ps. τ =0.5, and* B*<sup>0</sup> = 0.5 T [46].*

and because of the absence of electromagnetic coupling under electric bias, the transmission coefficient due to mutual coupling is zero. In **Figure 4** the same results are illustrated for the applied magnetic bias of 0.5 Tesla. There is good agreement in the magnitude and phase of the transmission coefficient in the abovementioned two methods. Also, by finding the poles of the coefficients, the electromagnetic wave propagation constants for the electrically and magnetically biased graphene sheets can be obtained which are in full agreement with [28, 57], respectively. The correctness of the extracted coefficients confirms the validity of the dyadic Green's function formulation.

#### *3.1.2 Graphene-dielectric stack*

Dyadic Green's function for an *N*-layer dielectric environment has been previously formulated using the scattering superposition method [58]. In this section, the above equations are extended to the environment with the electrically biased graphene boundaries, as shown in **Figure 2(b)**. The graphene boundary can be either continuous or periodically patterned as discussed in **section 2**. To start the analysis, the layers are numbered by starting from the top layer, and an arbitrary field point *i* and source point *j* are assumed. The problem is solved in the cylindrical coordinates. Since the cylindrical wave functions are discussed in detail in the next section, they are not mentioned here. The dyadic Green's function can be expanded as [58]:

*G* ð Þ*ij es <sup>r</sup>*, *<sup>r</sup>*<sup>0</sup> ð Þ¼ *<sup>i</sup>* 4*π* ð<sup>∞</sup> �∞ *dh*X<sup>∞</sup> *n*¼0 <sup>2</sup> � *<sup>δ</sup>*<sup>0</sup> *n* � � *λh <sup>j</sup>* <sup>1</sup> � *<sup>δ</sup><sup>N</sup> i* � �*Mn<sup>λ</sup>* � ð Þ *hi* <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *j* � �*Aij MM*<sup>0</sup> *<sup>n</sup><sup>λ</sup>* �*h <sup>j</sup>* � � <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> j* � �*Bij MM*<sup>0</sup> *<sup>n</sup><sup>λ</sup> h <sup>j</sup>* � � h i<sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> i* � �*Nn<sup>λ</sup>*ð Þ *hi* <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *j* � �*Aij NN*<sup>0</sup> *<sup>n</sup><sup>λ</sup>* �*h <sup>j</sup>* � � <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> j* � �*Bij NN*<sup>0</sup> *<sup>n</sup><sup>λ</sup> h <sup>j</sup>* � � h i<sup>þ</sup> <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *i* � �*Mn<sup>η</sup>*ð Þ �*hi* <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *j* � �*Cij MM*<sup>0</sup> *<sup>n</sup><sup>λ</sup>* �*h <sup>j</sup>* � � <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> j* � �*Dij MM*<sup>0</sup> *<sup>n</sup><sup>λ</sup> h <sup>j</sup>* � � h i<sup>þ</sup> <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *i* � �*Nn<sup>η</sup>*ð Þ� �*hi* <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *j* � �*Cij NN*<sup>0</sup> *<sup>n</sup><sup>λ</sup>* �*h <sup>j</sup>* � � <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> j* � �*Dij NN*<sup>0</sup> *<sup>n</sup><sup>λ</sup> h <sup>j</sup>* � � h io (30)

In the above equations, *M* and *N* are vector wave functions in the cylindrical coordinate system, and *Aij <sup>M</sup>*,*<sup>N</sup>*, *Bij <sup>M</sup>*,*<sup>N</sup>*, *<sup>C</sup>ij <sup>M</sup>*,*<sup>N</sup>*, *Dij <sup>M</sup>*,*<sup>N</sup>* are unknown coefficients. The Kranoker delta function is used in the field expansion to generalize the formulation for the arbitrary locations of the field and source points. By applying the boundary conditions on the tangential components of the electric field at the interface of the arbitrary two layers, denoted by *i* and *i* + 1, it can be shown that:

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

$$\begin{aligned} \begin{bmatrix} \boldsymbol{A}\_{M}^{\boldsymbol{\tilde{y}}} \\ \boldsymbol{B}\_{M}^{\boldsymbol{\tilde{y}}} \\ \boldsymbol{O}\_{i} \boldsymbol{A}\_{N}^{\boldsymbol{\tilde{y}}} \\ \boldsymbol{O}\_{i} \boldsymbol{B}\_{N}^{\boldsymbol{\tilde{y}}} \end{bmatrix} e^{i\boldsymbol{h}\boldsymbol{x}\_{i}} + \begin{bmatrix} \boldsymbol{C}\_{M}^{\boldsymbol{\tilde{y}}} \\ \boldsymbol{D}\_{M}^{\boldsymbol{\tilde{y}}} + \boldsymbol{\delta}\_{j}^{\boldsymbol{\tilde{y}}} \\ \boldsymbol{O}\_{i} \boldsymbol{C}\_{N}^{\boldsymbol{\tilde{y}}} \\ \boldsymbol{O}\_{i} \left(\boldsymbol{D}\_{N}^{\boldsymbol{\tilde{y}}} + \boldsymbol{\delta}\_{j}^{\boldsymbol{\tilde{y}}}\right) \end{bmatrix} e^{-i\boldsymbol{h}\boldsymbol{x}\_{i}} = \begin{bmatrix} \left(\boldsymbol{A}\_{M}^{(i+1)\boldsymbol{\tilde{y}}} + \boldsymbol{\delta}\_{i+1}^{\boldsymbol{\tilde{y}}}\right) \\ \boldsymbol{B}\_{M}^{(i+1)\boldsymbol{\tilde{y}}} \\ \boldsymbol{O}\_{i+1} \left(\boldsymbol{A}\_{N}^{(i+1)\boldsymbol{\tilde{y}}} + \boldsymbol{\delta}\_{i+1}^{\boldsymbol{\tilde{y}}}\right) \\ \boldsymbol{O}\_{i+1} \boldsymbol{B}\_{N}^{(i+1)\boldsymbol{\tilde{y}}} \\ \boldsymbol{O}\_{i} \left(\boldsymbol{D}\_{M}^{(i+1)\boldsymbol{\tilde{y}}}\right) \\ \boldsymbol{D}\_{i+1} \boldsymbol{C}\_{N}^{(i+1)\boldsymbol{\tilde{y}}} \\ \boldsymbol{O}\_{i+1} \boldsymbol{D}\_{N}^{(i+1)\boldsymbol{\tilde{y}}} \end{bmatrix} e^{-i\boldsymbol{h}\_{i+1}\boldsymbol{x}\_{i}} \end{aligned} \tag{31}$$

where *Oi* <sup>¼</sup> *hi ki* , *Pi* <sup>¼</sup> *hi μi* , *Qi* <sup>¼</sup> *ki μi* . The boundary condition on the tangential components of the magnetic field yields:

$$
\begin{bmatrix} P\_{i+1} \left( A\_{M}^{(i+1)} + \delta\_{i+1}^{\hat{I}} \right) \\\\ P\_{i+1} B\_{M}^{(i+1)\hat{\}} \\\\ Q\_{i+1} \left( A\_{N}^{(i+1)\hat{\}} + \delta\_{i+1}^{\hat{I}} \right) \\\\ Q\_{i+1} B\_{N}^{(i+1)\hat{\}} \\\\ \end{bmatrix} e^{\hat{\mu}\_{i+1}\varepsilon\_{i}} - \begin{bmatrix} P\_{i+1} C\_{M}^{(i+1)\hat{\}} \\\\ P\_{i+1} D\_{M}^{(i+1)\hat{\}} \\\\ Q\_{i+1} C\_{N}^{(i+1)\hat{\}} \\\\ Q\_{i+1} D\_{N}^{(i+1)\hat{\}} \end{bmatrix} e^{-\hat{\mu}\_{i+1}\varepsilon\_{i}} - \begin{bmatrix} P\_{i} A\_{M}^{\hat{\mu}} \\\\ P\_{i} B\_{M}^{\hat{\mu}} \\\\ Q\_{i} A\_{N}^{\hat{\mu}} \\\\ Q\_{i} B\_{N}^{\hat{\mu}} \end{bmatrix} e^{\hat{\mu}\_{i+1}\varepsilon\_{i}} \\\\ + \begin{bmatrix} P\_{i} C\_{M}^{\hat{\mu}} \\\\ P\_{i} (D\_{M}^{\hat{\mu}} + \delta\_{i}^{\hat{\mu}}) \\\\ Q\_{i} (D\_{N}^{\hat{\mu}} + \delta\_{i}^{\hat{\mu}}) \\\\ P\_{i} B\_{N}^{\hat{\mu}} \end{bmatrix} e^{-\hat{\mu}\_{i+1}\varepsilon\_{i}} - \begin{bmatrix} A\_{M}^{\hat{\mu}} \\\\ Q\_{M}^{\hat{\mu}} \\\\ P\_{i} B\_{N}^{\hat{\mu}} \\\\ P\_{i} B\_{N}^{\hat{\mu}} \end{bmatrix} e^{\hat{\mu}\_{i+1}\varepsilon\_{i}} \\\\ \end{bmatrix} e^{-\hat{\mu}\_{i$$

By re-writing the coefficients as a matrix:

$$
\begin{bmatrix}
\mathbf{A}\_{M,N}^{(i+1)j} + \delta\_{i+1}^{j} & \mathbf{B}\_{M,N}^{(i+1)j} \\
\mathbf{C}\_{M,N}^{(i+1)j} & \mathbf{D}\_{M,N}^{(i+1)j} \\
\end{bmatrix} = \begin{bmatrix}
\mathbf{1} & \mathbf{R}\_{Fi}^{H,V} \\
\mathbf{T}\_{Fi}^{H,V} & \mathbf{T}\_{Fi}^{H,V} \\
\mathbf{R}\_{Pi}^{H,V} & \mathbf{1} \\
\mathbf{T}\_{Pi}^{H,V} & \mathbf{T}\_{Pi}^{H,V} \\
\end{bmatrix} \begin{bmatrix}
\mathbf{A}\_{M,N}^{ij} & \mathbf{B}\_{M,N}^{ij} \\
\mathbf{C}\_{M,N}^{ij} & \mathbf{D}\_{M,N}^{ij} + \delta\_i^{j} \\
\end{bmatrix} \tag{33}
$$

the outgoing and incoming reflection and transmission coefficients can be defined and used to extract the recursive relations. The coefficients for *M*<sup>0</sup> sources are:

$$R\_{Fi}^H = \frac{\mu\_i h\_{i+1} - \mu\_{i+1} h\_i + \mathbf{g}}{\mu\_i h\_{i+1} + \mu\_{i+1} h\_i + \mathbf{g}} \tag{34}$$

$$R\_{Pi}^{H} = \frac{\mu\_i h\_{i+1} - \mu\_{i+1} h\_i - \mathbf{g}}{\mu\_i h\_{i+1} + \mu\_{i+1} h\_i - \mathbf{g}} \tag{35}$$

$$T\_{Fi}^H = \frac{2\mu\_i h\_{i+1}}{\mu\_i h\_{i+1} + \mu\_{i+1} h\_i + \mathbf{g}} \tag{36}$$

*Electromagnetic Wave Propagation for Industry and Biomedical Applications*

$$T\_{Pi}^{H} = \frac{2\mu\_i h\_{i+1}}{\mu\_i h\_{i+1} + \mu\_{i+1} h\_i - \mathbf{g}} \tag{37}$$

The coefficients for *N*<sup>0</sup> sources are:

$$R\_{Fi}^V = \frac{\mu\_i h\_i k\_{i+1}^2 - k\_i^2 \mu\_{i+1} h\_{i+1} + g h\_i h\_{i+1}}{\mu\_i h\_i k\_{i+1}^2 + \mu\_{i+1} h\_{i+1} k\_i^2 + g h\_i h\_{i+1}} \tag{38}$$

$$R\_{Pi}^V = \frac{\mu\_i h\_i k\_{i+1}^2 - \mu\_{i+1} k\_i^2 h\_{i+1} - g h\_i h\_{i+1}}{\mu\_i h\_i k\_{i+1}^2 + \mu\_{i+1} k\_i^2 h\_{i+1} - g h\_i h\_{i+1}} \tag{39}$$

$$T\_{Pi}^V = \frac{2\mu\_{i+1}k\_ik\_{i+1}h\_{i+1}}{\mu\_ih\_ik\_{i+1}^2 + \mu\_{i+1}k\_i^2h\_{i+1} - gh\_ih\_{i+1}}\tag{40}$$

$$T\_{Fi}^V = \frac{2\mu\_i h\_{i+1} k\_i k\_{i+1}}{\mu\_i h\_i k\_{i+1}^2 + \mu\_{i+1} h\_{i+1} k\_i^2 + g h\_i h\_{i+1}} \tag{41}$$

The superscripts *H* and *V* respectively denote TE and TM sources. Also, subscripts *F* and *P* are used to show the incoming and outgoing waves, respectively. The procedure of extracting the unknown coefficients using (34)-(41) is discussed in [58]. To validate the proposed formulas, a parallel plate waveguide with graphene walls is considered. To extract the characteristic equation using the proposed formulations, it is necessary to force the denominator of the coefficients equal to zero. For this three-layer medium:

$$T^{(1)} = T^2 . T^1 = \cdot \begin{bmatrix} \frac{1}{T\_{F2}^H} e^{i(h\_2 - h\_1)d} & \frac{R\_{F2}^H}{T\_{F2}^H} e^{-i(h\_2 - h\_1)d} \\\\ \frac{R\_{F2}^H}{T\_{F2}^H} e^{i(h\_2 + h\_1)d} & \frac{1}{T\_{F2}^H} e^{-i(h\_2 - h\_1)d} \\\\ \frac{R\_{F2}^H}{T\_{F2}^H} e^{i(h\_2 + h\_1)d} & \frac{1}{T\_{F2}^H} e^{-i(h\_2 - h\_1)d} \end{bmatrix} \cdot \begin{bmatrix} \frac{1}{T\_{F1}^H} & \frac{R\_{F1}^H}{T\_{F1}^H} \\\\ \frac{R\_{F1}^H}{T\_{F1}^H} & \frac{1}{T\_{F1}^H} \end{bmatrix} \tag{42}$$

From which it can be deduced that *B*<sup>11</sup> *<sup>M</sup>*,*<sup>N</sup>* ¼ � *<sup>T</sup>*ð Þ<sup>1</sup> 12 *T*ð Þ<sup>1</sup> 11 . By setting *T*ð Þ<sup>1</sup> <sup>11</sup> equal to zero and assuming that the medium (1) and (3) are the same, and also defining *h* ¼ *iωμ*1*σ*, for the *H* coefficients it can be concluded that:

$$\frac{e^{ih\_2d}-1}{e^{ih\_2d}+1}h\_1 - (h\_2 + h) = 0 \implies h\_2 + h + j\tan\left(\frac{h\_2d}{2}\right)h\_1 = 0\tag{43}$$

This procedure is repeatable for *V* sources. Also, to calculate the reflection coefficient from a multilayer structure, it is necessary to consider the field and source points in region 1. In this case, the only non-zero coefficient in Green's function expansion is *B*<sup>11</sup> *<sup>M</sup>*,*<sup>N</sup>* coefficient representing the plane wave reflection coefficient from the multilayer structure.

#### **3.2 Graphene-based cylindrical structures**

In this section, the dyadic Green's function of a cylindrical structure with the tensor surface conductivity boundary condition will be extracted. Later, different examples of guiding and scattering problems are provided to investigate the validity *Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

of the formulation. In general, in cylindrical structures, TE and TM modes are coupled, which leads to the complexity of mathematical relations. Therefore, the generalization of the formulation to the multilayered cylinders is not considered here. Note that graphene sheets can be wrapped around cylindrical particles due to the presence of van der Waals force [59]. For this purpose, tape-assist transfer under micromanipulation and spin-coating methods are proposed [60].

#### *3.2.1 Dyadic Green's function for a cylinder with tensor surface conductivity boundary condition*

The dyadic Green's function of the single-layer cylinder with the tensor surface conductivity boundary condition, as considered in **Figure 5.**, will be extracted in the following. The interior region of the cylinder is made of dielectric material and its cover is considered as a full tensor surface conductivity. To solve the problem, the vector wave functions are defined as [23]:

$$\overline{\mathbf{M}}\_{\mu}(h) = e^{im\phi} e^{ihx} \left[ im \frac{Z\_m(\mu r)}{r} \hat{r} - \frac{\partial Z\_m(\mu r)}{\partial r} \hat{\Phi} \right] \tag{44}$$

$$\overline{\mathbf{N}}\_{\mu}(h) = \frac{1}{k\_{\circ}} e^{im\phi} e^{ihx} \left[ ih \frac{\partial Z\_{m}(\mu r)}{\partial r} \hat{r} - \frac{hm}{r} Z\_{m}(\mu r) \hat{\Phi} - \mu^{2} Z\_{m}(\mu r) \hat{\mathbf{z}} \right] \tag{45}$$

In the above equations, *Zm* is the cylindrical Bessel function in the inner layer and the cylindrical Hankel function in the outer layer, both with the orders of *m*. The wavenumber in the radial direction is *μ* and the wavenumber along the length is *h*. The free-space Green's function in the cylindrical coordinates is [23]:

$$\begin{split} \overline{\mathbf{G}}\_{\epsilon 0} \left( \overline{\mathbf{R}}, \overline{\mathbf{R}}' \right) &= -\frac{1}{k^2} \hat{r} \hat{r} \delta \left( \overline{\mathbf{R}} - \overline{\mathbf{R}}' \right) \\\\ &+ \frac{i}{8\pi} \int\_{-\infty}^{\infty} dh \sum\_{m=-\infty}^{\infty} \frac{1}{\eta^2} \times \left[ \overline{\mathbf{M}}\_{\eta}(h) \overline{\mathbf{M}}\_{\eta}'^{(1)}(-h) + \overline{\mathbf{N}}\_{\eta}(h) \overline{\mathbf{N}}\_{\eta}'^{(1)}(-h) \right] \quad R \,< \mathsf{R}' \end{split} \tag{46}$$

#### **Figure 5.**

*(a) A monolayer cylinder with a 2D cover with the tensor surface conductivity boundary condition and (b) its special cases constructed by densely packed strips and square nano-patches [61].*

In the cylindrical coordinate system, TE and TM modes are coupled and the expansion of the fields are [23].

$$\begin{split} \overline{\mathbf{G}}\_{ee}^{(11)} \left( \overline{\mathbf{R}}, \overline{\mathbf{R}}' \right) &= \frac{i}{8\pi} \int\_{-\infty}^{\infty} dh \sum\_{m=-\infty}^{\infty} \frac{1}{\eta^2} \times \\\\ \left\{ \left[ A\_{\eta} \overline{\mathbf{M}}\_{\eta}^{(1)}(h) + B\_{\eta} \overline{\mathbf{N}}\_{\eta}^{(1)}(h) \right] \overline{\mathbf{M}}\_{\eta}^{(1)}(-h) + \left[ C\_{\eta} \overline{\mathbf{N}}\_{\eta}^{(1)}(h) + D\_{\eta} \overline{\mathbf{M}}\_{\eta}^{(1)}(h) \right] \overline{\mathbf{N}}\_{\eta}^{(1)}(-h) \right\} \\\\ &\dots \end{split} \tag{47}$$

$$\begin{split} \overline{\overline{\mathbf{G}}}\_{\varepsilon\varepsilon}^{(21)} \left( \overline{\mathbf{R}}, \overline{\mathbf{R}}' \right) &= \frac{i}{8\pi} \int\_{-\infty}^{\infty} dh \sum\_{m=-\infty}^{\infty} \frac{1}{\eta^2} \times \\\\ \left\{ \left[ a\_{\xi} \overline{\mathbf{M}}\_{\xi}(h) + b\_{\xi} \overline{\mathbf{N}}\_{\xi}(h) \right] \overline{\mathbf{M}}\_{\eta}^{(1)}(-h) + \left[ c\_{\xi} \overline{\mathbf{N}}\_{\xi}(h) + d\_{\xi} \overline{\mathbf{M}}\_{\xi}(h) \right] \overline{\mathbf{N}}\_{\eta}^{(1)}(-h) \right\} \end{split} \tag{48}$$

where *Aη*, *Bη*, *Cη*, and *D<sup>η</sup>* are the unknown coefficients of the DGF expansion in region 1. Also, *aξ*, *bξ*, *cξ*, and *d<sup>ξ</sup>* are DGF expansion coefficients in the region 2. The boundary condition on the tangential components of the magnetic Green's function is given by [61]:

$$\hat{r} \times \left( \frac{\nabla \times \overline{\overline{\mathbf{G}}}\_{\varepsilon}^{(11)}}{\mu\_1} - \frac{\nabla \times \overline{\overline{\mathbf{G}}}\_{\varepsilon}^{(21)}}{\mu\_2} \right) = i\alpha \overline{\overline{\sigma}}\_{\phi\varepsilon} \cdot \overline{\overline{\mathbf{G}}}\_{\varepsilon}^{(21)} \tag{49}$$

By applying the above-mentioned boundary condition along with the boundary condition regarding the continuity of the electric Green's function to **(47)–(48),** the system of equations to determine the unknown coefficients can be obtained as:

$$\begin{bmatrix} -\frac{\partial H\_m^0(\eta\alpha)}{\partial a} & -\frac{1}{k\_1}\frac{hm}{a}H\_m^{(1)}(\eta a) & \frac{\partial J\_m(\xi a)}{\partial a} & \frac{1}{k\_2}\frac{hm}{a}J\_m(\xi a) \\\\ \mathbf{0} & -\frac{1}{k\_1}\eta^2H\_m^{(1)}(\eta a) & \mathbf{0} & \frac{1}{k\_2}\hat{\varepsilon}^2J\_m(\xi a) \\\\ -\frac{1}{\mu\_1}\frac{hm}{a}H\_m^{(1)}(\eta a) & -\frac{k\_1}{\mu\_1}\frac{H\_m^{(1)}(\eta a)}{\partial a} & \frac{1}{\mu\_2}\frac{hm}{a}J\_m(\xi a) + \sigma\_1 & \frac{k\_2}{\mu\_2}\frac{J\_m(\xi a)}{\partial a} + \sigma\_2 \\\\ \frac{1}{\mu\_1}\eta^2H\_m^{(1)}(\eta a) & \mathbf{0} & -\frac{1}{\mu\_2}\hat{\varepsilon}^2J\_m(\xi a) + \sigma\_3 & \sigma\_4 \\\\ \mathbf{0} & \mathbf{0} & \\\\ \mathbf{0} & \mathbf{1} & \mathbf{1}\_{m=1} \\\\ \frac{1}{\mu\_1}\frac{hm}{a}J\_m(\eta a) & & \\\\ -\frac{1}{\mu\_1}\eta^2J\_m(\eta a) & & \\\\ \end{bmatrix}$$

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

$$
\begin{bmatrix}
0 & \frac{1}{\mu\_1}\eta^2H\_m^{(1)}(\eta a) & \sigma\_3' & -\frac{1}{\mu\_2}\xi^2J\_m(\xi a) + \sigma\_4' \\
\end{bmatrix} \times \begin{bmatrix} \mathcal{C}\_\eta \\ \mathcal{I}\_\eta \\ \mathcal{C}\_\xi \\ \mathcal{I}\_\xi \\ \mathcal{I}\_\xi \end{bmatrix}
$$
 
$$
\begin{bmatrix}
\frac{1}{k\_1}\frac{hm}{a}J\_m(\eta a) \\
\frac{1}{k\_1}\eta^2J\_m(\eta a) \\
\mu\_1 \end{bmatrix}.
$$

(51)

where:

$$\sigma\_1 = i\alpha\sigma\_{x\phi}\frac{\partial f\_m(\xi a)}{\partial a}, \sigma\_2 = i\alpha\frac{\sigma\_{xx}}{k\_2}\xi^2 f\_m(\xi a) + i\alpha\sigma\_{z\phi}\frac{hm}{k\_2a}f\_m(\xi a)$$

$$\sigma\_3 = i\alpha\sigma\_{\phi\phi}\frac{\partial f\_m(\xi a)}{\partial a}, \sigma\_4 = i\alpha\sigma\_{\phi\phi}\frac{1}{k\_2}\frac{hm}{a}f\_m(\xi a) + i\alpha\sigma\_{\phi x}\frac{1}{k\_2}\xi^2 f\_m(\xi a)$$

$$\sigma'\_1 = i\alpha\sigma\_{xz}\frac{1}{k\_2}\xi^2 f\_m(\xi a) + i\alpha\sigma\_{z\phi}\frac{hm}{k\_2a}f\_m(\xi a), \sigma'\_2 = i\alpha\sigma\_{z\phi}\frac{\partial f\_m(\xi a)}{\partial a}$$

$$\sigma'\_3 = i\alpha\sigma\_{\phi\phi}\frac{1}{k\_2}\frac{hm}{a}f\_m(\xi a) + i\alpha\sigma\_{\phi x}\frac{1}{k\_2}\xi^2 f\_m(\xi a), \sigma'\_4 = i\alpha\sigma\_{\phi\phi}\frac{\partial f\_m(\xi a)}{a}$$

Nullifying the determinant of the matrices **(50)–(51)** for *m* = 0, each of the above matrices can be separated as the multiplication of:

$$\frac{\alpha \mathfrak{e}\_2}{\mathfrak{z}} \frac{J\_1(\mathfrak{xi}a)}{J\_0(\mathfrak{z}a)} - \frac{\alpha \mathfrak{e}\_1}{\eta} \frac{H\_1^{(1)}(\eta a)}{H\_0^{(1)}(\eta a)} + i \sigma\_d = \mathbf{0} \tag{52}$$

$$\frac{\xi}{\rho \alpha \mu\_2} \frac{J\_0(\xi a)}{J\_1(\xi a)} - \frac{\eta}{\alpha \mu\_1} \frac{H\_0^{(1)}(\eta a)}{H\_1^{(1)}(\eta a)} + i \sigma\_d = 0 \tag{53}$$

For the graphene shell under magnetic bias, both matrices in **(50)–(51)** result in the following equation for the propagation constant of hybrid TE and TM waves:

$$\left[\frac{\alpha \mathfrak{e}\_2}{\xi} \frac{H\_1^{(1)}(-\xi a)}{H\_0^{(1)}(-\xi a)} - \frac{\alpha \mathfrak{e}\_1}{\eta} \frac{J\_1(-\eta a)}{J\_0(-\eta a)} + i \sigma\_d \right] \times \left[\frac{\eta}{\alpha \mu\_1} \frac{J\_0(-\eta a)}{J\_1(-\eta a)} - \frac{\xi}{\alpha \mu\_2} \frac{H\_0^{(1)}(-\xi a)}{H\_1^{(1)}(-\xi a)} + i \sigma\_d \right] = -\sigma\_0^2 \tag{54}$$

Which is in agreement with [57]. Also, the total scattering cross-section (TSCS) of a graphene-coated cylinder with the parameters *a* = 50 μm, *ε*<sup>2</sup> = 2.4, *τ* = 1 ps, *μ<sup>c</sup>* = 0.25 eV under electric bias is calculated for both TE and TM polarizations in

**Figure 6**, and compared with the results of the CST2017 software package. As can be seen, both methods have resulted in the same results.

As another example, the densely packed graphene strips with the parameters *L* = 420 nm, *W* = 400 nm, *μ<sup>c</sup>* = 0.5 eV, and *τ* = 1 ps are considered around the dielectric cylinder as in **Figure 5.** It is assumed that the strips are wrapped around a hollow cylinder with the radius of *a* = 50 μm. **Figure 7** shows the TSCS of this structure for the magnetic biases with the strength in the range of 20-40 T. As observed, by increasing the magnetic bias, the resonant frequency of the surface plasmons blue shifts. The associated planar structure behaves as a hyperbolic metasurface [62]. Under locally flat consideration of the curvature, this structure can also be considered as a hyperbolic medium. In the cylindrical geometries, hyperbolic meta-surfaces can be obtained using graphene-dielectric stacks [63]. The advantage of this hyperbolic structure is its two-dimensional nature and reconfigurability. It is also demonstrated that covering the surface of nanotubes with the hyperbolic meta-surface increases the interaction of the light with dipole emitters [64].

Finally, as **Figure 5(b)** illustrates, graphene-based square patches around the cylinder are considered under magnetic bias. Geometrical and optical parameters are as follows: *τ* = 1 ps, *g* = 0.5 μm, and *D* = 0.5 μm and the TSCS is illustrated in **Figure 8.** for *B*<sup>0</sup> = 0 T and *B*<sup>0</sup> = 10 T. As can be seen, by changing the magnetic bias, the optical state changes from the maximum scattering to the minimum scattering. Such capability has recently been proposed by using a phase change material for switching between these two situations [65]. In this structure, the operating frequency can also be adjusted by changing the electric bias of graphene.

#### **Figure 6.**

*TSCS for an infinite cylinder covered with electrically biased graphene shell with* a *= 50 μm,* ε*<sup>2</sup> = 2.4,* τ *= 1 ps, and* μ<sup>c</sup> *= 0.25 eV considering (a) TE and (b) TM polarization for the incident wave [61].*

#### **Figure 7.**

*TSCS for an infinite length hollow cylinder (*a *= 50 μm) coated with densely packed graphene strips with* L *= 420 nm,* W *= 400 nm,* μ<sup>c</sup> *= 0.5 eV, and* τ *= 1 ps considering different magnetic biases [61].*

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

**Figure 8.**

*TSCS for an infinite cylinder coated with graphene patches with* D *= 0.5 μm,* g *= 0.5 μm, and* τ *= 0.5 ps under different electric and magnetic biases [61].*

It is essential to note that in the cylindrical structures, the spatial domain Green's function consists of an integral on the real axis. Due to the existence of poles in the integration path, the integration path is usually deformed into a triangular shape. Also, a common method for calculating the spatial domain Green's function is the generalized pencil of function (GPOF) method in which the Green's function is expanded in terms of the complex exponential functions [66]. The unknowns can be found via the algorithm provided in [67]. Also, as mentioned earlier, the dyadic Green's function of the graphene-based multilayered cylindrical structures will not be considered here. For calculating the scattering cross-section of graphene-based multilayer cylindrical structures, a simple approach based on the transfer matrix method (TMM) is proposed that will be suitable for establishing novel optical devices [68, 69].

#### **3.3 Graphene-based spherically layered medium**

In this section, a multi-layered spherical structure with the graphene boundaries is considered and its Green's function is extracted by assuming different locations for the source and observation points. The relationship between the Green's function expansion coefficients and the modified Mie-Lorentz coefficients is exhibited to discuss how to solve the scattering problems using the Green's function. Scattering analysis of graphene-based layered structures is of great importance in the design of novel optical devices [70]. Finally, the procedure for calculating the Purcell factor is considered as an important application. Instances of experimentally realized graphene-coated spherical particles can be found in [71, 72], where improved template method and hydrothermal method are proposed for the synthesis. Also, transmission electron microscopy (TEM) and field emission scanning electron microscopy (FE-SEM) are used for characterization.

#### *3.3.1 Dyadic Green's function of a graphene-based spherically layered structure*

Let us consider an *N*-layer spherical medium with the graphene boundaries as shown in **Figure 9.** The purpose of this section is to compute the dyadic Green's function of this structure with the assumption of arbitrary locations for the field and source points. For this purpose, the Green's function in each layer is expanded in terms of vector wave functions with unknown coefficients. These functions are calculated using the scalar wave function of [74]:

**Figure 9.** *Spherically layered medium with graphene boundaries (a) 2D and (b) 3D views* [73]*.*

$$\varphi\_{mn}(r,\theta,\phi) = z\_n(k\_p r) P\_n^m(\cos\theta) e^{im\phi} \tag{55}$$

where *zn*ð Þ*:* represents the spherical Bessel or Hankel functions of order *n* and *Pm <sup>n</sup>* ð Þ*:* is the associated Legendre function with degree *n* and order *m*. It can be readily shown that:

$$\overline{\mathbf{M}}\_{mn}(k\_p) = z\_n(k\_p r) \, e^{im\phi} \left[ \frac{im}{\sin \theta} P\_n^m(\cos \theta) \hat{\boldsymbol{\theta}} - \frac{d P\_n^m(\cos \theta)}{d \theta} \hat{\boldsymbol{\phi}} \right] \tag{56}$$

$$\begin{aligned} \overline{\mathbf{N}}\_{mn}(\mathbf{k}\_p) &= \frac{n(n+1)}{k\_p r} z\_n(\mathbf{k}\_p r) P\_n^m(\cos \theta) e^{im\phi} \hat{r} + \\\ \frac{1}{k\_p r} \frac{d \left[ r z\_n(\mathbf{k}\_p r) \right]}{dr} \left[ \frac{d P\_n^m(\cos \theta)}{d \theta} \hat{\theta} + \frac{im}{\sin \theta} P\_n^m(\cos \theta) \hat{\phi} \right] e^{im\phi} \end{aligned} \tag{57}$$

The above vector functions are self and mutually orthogonal. This feature will be used to decouple the equation when computing the unknown coefficients.

As mentioned earlier, in the scattering superposition method, the dyadic Green's function is written as the sum of the free-space and scattering Green's functions. The free-space Green's function is related to the source in an infinite homogeneous medium while the scattering Green's function is due to the source in the presence of the layered medium. The expansion of the Green's function for the spherical structure with *N* concentric layers, assuming that the source and field points are respectively located in the desired layers with the labels *p* and *q*, can be written as [23]:

$$
\overline{\hat{\mathbf{G}}}\_{\epsilon}^{(pq)}(\mathbf{r},\mathbf{r}') = \overline{\hat{\mathbf{G}}}\_{0\epsilon}^{(pq)}(\mathbf{r},\mathbf{r}')\delta\_{q}^{p} + \overline{\hat{\mathbf{G}}}\_{\epsilon\epsilon}^{(pq)}(\mathbf{r},\mathbf{r}') \tag{58}
$$

The free-space Green's function can be obtained using the residue theorem as [23]:

$$\begin{split} \overline{\overline{\mathbf{G}}}\_{0\epsilon}(r,r') &= -\frac{\hat{r}\hat{r}}{k\_q^2}\delta(r-r') + \frac{ik\_q}{4\pi} \sum\_{n=1}^{\infty} \sum\_{m=-n}^{n} \left(2-\delta\_m^0\right) \frac{2n+1}{n(n+1)} \times \\ &\frac{(n-m)!}{(n+m)!} \times \left\{ \begin{array}{c} \overline{\mathbf{M}}\_{mn}^{(1)}(k\_q)\overline{\mathbf{M}}\_{mn}^{'}(k\_q) + \overline{\mathbf{N}}\_{mn}^{(1)}(k\_q)\overline{\mathbf{N}}\_{mn}^{'}(k\_q) \end{array} \right. \\ &\left. \left(\overline{\mathbf{M}}\_{mn}(k\_q)\overline{\mathbf{M}}^{(1)}{}\_{mn}\left(k\_q\right) + \overline{\mathbf{N}}\_{mn}\left(k\_q\right)\overline{\mathbf{N}}^{(1)}{}\_{mn}\left(k\_q\right)r < r' \right. \end{split} \tag{59}$$

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

Moreover, scattering Green's function in each layer can be expanded as [75]:

*G*^ ð Þ *pq es <sup>r</sup>*, *<sup>r</sup>*<sup>0</sup> ð Þ¼ *ikq* 4*π* X∞ *n*¼1 X*n m*¼�*n* <sup>2</sup> � *<sup>δ</sup>*<sup>0</sup> *m* � � 2*n* þ 1 *n n*ð Þ þ 1 ð Þ *n* � *m* ! ð Þ *n* þ *m* ! � <sup>1</sup> � *<sup>δ</sup><sup>N</sup> p* � �*M*ð Þ<sup>1</sup> *mn* <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *q* � �*Apq <sup>H</sup> M*<sup>0</sup> *mn* <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> q* � �*Bpq <sup>H</sup> <sup>M</sup>*0ð Þ<sup>1</sup> *mn* <sup>n</sup> h i <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> p* � �*N*ð Þ<sup>1</sup> *mn* <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *q* � �*Apq <sup>V</sup> N*<sup>0</sup> *mn* <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> q* � �*Bpq <sup>V</sup> <sup>N</sup>*0ð Þ<sup>1</sup> *mn* h i <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *p* � �*Mmn* <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *q* � �*Cpq <sup>H</sup> M*<sup>0</sup> *mn* <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> q* � �*Dpq <sup>H</sup> <sup>M</sup>*0ð Þ<sup>1</sup> *mn* h i <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *p* � �*Nmn* <sup>1</sup> � *<sup>δ</sup>*<sup>1</sup> *q* � �*Cpq <sup>V</sup> N*<sup>0</sup> *mn* <sup>þ</sup> <sup>1</sup> � *<sup>δ</sup><sup>N</sup> q* � �*Dpq <sup>V</sup> <sup>N</sup>*0ð Þ<sup>1</sup> *mn* <sup>h</sup> io (60)

where *Apq <sup>H</sup>*,*<sup>V</sup>*, *Bpq <sup>H</sup>*,*<sup>V</sup>*, *Cpq <sup>H</sup>*,*<sup>V</sup>*, and *<sup>D</sup>pq <sup>H</sup>*,*<sup>V</sup>* are the unknown coefficients of the Green's function to be obtained. In the above formulation, the superscript (1) represents that the spherical Hankel functions of the first type are chosen to represent the spherical vector functions. For the other vector wave functions, the first-kind Bessel function should be selected. As observed, in the spherical structures the TE and TM waves are not coupled. Thus, the expansion of the fields only includes the interaction of the functions *M* and *M*<sup>0</sup> as well as the functions *N* and *N*<sup>0</sup> . Moreover, to solve the problem for the structures with the arbitrary number of layers, the Kronecker delta function is used in the expansions of the fields. In the middle layers of the spherical structure, the electric field is a linear combination of the Bessel and Hankel functions, whereas in the outermost layer, only the Hankel functions, and in the innermost layer, only the Bessel functions will exist. So we have used 1 � *<sup>δ</sup>*<sup>1</sup> *p* � � and 1 � *<sup>δ</sup><sup>N</sup> p* � � functions. To determine the delta functions related to the source

terms, free-space Green's function can be used. For *q* = 1 (external layer), the source functions related to the second criterion of the free-space Green's function is used. Also, for *q* = *N* (internal layer) the source functions related to the first criterion should be used. In the middle layers, a linear combination of different source functions must be used.

It should be noted that by using the addition theorem in Legendre functions, the internal series in the Green's function expansion can be eliminated [11]. To select the number of terms required for the convergence of the external series, the conditions of the problem must be considered. In the other words, in the structures whose electrical size is very much smaller than that of the wavelength, the term *n* = 1 is sufficient for the convergence [76]. Also, in the case where the distance between the source and observation points is large, the series can be truncated in the number *Nt* <sup>¼</sup> *<sup>x</sup>* <sup>þ</sup> <sup>3</sup> ffiffiffi *<sup>x</sup>* <sup>p</sup> <sup>þ</sup> 2 , where *<sup>x</sup>* <sup>=</sup> *<sup>k</sup>*0*R*1. Otherwise, the convergence of the series is weak, and a large number of terms in the range of 20 *k*0*R*<sup>1</sup> should be considered. In this case, the series acceleration techniques will be highly efficient in terms of computational efficiency [77, 78].

Boundary conditions on tangential components of electric and magnetic Green's functions in the interface of two adjacent layers are [73]:

$$
\hat{\mathbf{r}} \times \overline{\hat{\mathbf{G}}}^{(pq)}\_{\epsilon} = \hat{\mathbf{r}} \times \overline{\hat{\mathbf{G}}}^{(p+1)q}\_{\epsilon} \tag{61}
$$

$$\frac{1}{i\alpha\mu\_{p+1}}\hat{\mathbf{r}} \times \nabla \times \overline{\hat{\mathbf{G}}}\_{\varepsilon}^{[(p+1)q]} - \frac{1}{i\alpha\mu\_{p}}\hat{\mathbf{r}} \times \nabla \times \overline{\hat{\mathbf{G}}}\_{\varepsilon}^{[pq)} = -\sigma\_{(p+1)p}\hat{\mathbf{r}} \times \left(\hat{\mathbf{r}} \times \overline{\overline{\mathbf{G}}}\_{\varepsilon}^{[pq)}\right) \tag{62}$$

which respectively lead to the following equations:

$$\begin{aligned} \left[\begin{array}{l} \xi\_{n}^{pp}A\_{H}^{pq} \\ \xi\_{n}^{pp}B\_{H}^{q} \\ \delta\_{n}^{pp}A\_{V}^{pq} \\ \delta\_{n}^{pp}B\_{V}^{pq} \end{array}\right] + \left[\begin{array}{l} \Psi\_{n}^{pp}C\_{H}^{q} \\ \Psi\_{n}^{pp}\left(D\_{H}^{pq} + \delta\_{p}^{q}\right) \\ \partial\nu\_{n}^{pp}C\_{V}^{pq} \\ \partial\nu\_{n}^{pp}\left(D\_{V}^{pq} + \delta\_{p}^{q}\right) \end{array}\right] = \left[\begin{array}{l} \xi\_{n}^{(p+1)p}\left(A\_{H}^{(p+1)q} + \delta\_{p}^{q}\right) \\ \xi\_{n}^{(p+1)p}B\_{H}^{(p+1)q} \\ \delta\_{n}^{\xi\_{n}^{(p+1)p}}\left(A\_{V}^{(p+1)q} + \delta\_{p+1}^{q}\right) \\ \delta\_{n}^{\xi\_{n}^{(p+1)p}}C\_{V}^{(p+1)q} \\ \end{array}\right] \\ + \left[\begin{array}{l} \Psi\_{n}^{(p+1)p}C\_{H}^{(p+1)q} \\ \Psi\_{n}^{(p+1)p}D\_{H}^{(p+1)q} \\ \partial\nu\_{n}^{(p+1)p}C\_{V}^{(p+1)q} \\ \end{array}\right] \tag{63}$$

*kp*þ<sup>1</sup> *μ<sup>p</sup>*þ<sup>1</sup> *<sup>∂</sup>ξ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> <sup>A</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>H</sup>* þ *δ q p*þ1 � � *<sup>∂</sup>ξ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> <sup>B</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> H <sup>ξ</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> <sup>A</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>q</sup> <sup>V</sup>* þ *δ q p*þ1 � � *<sup>ξ</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> <sup>B</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>q</sup> V* 2 6 6 6 6 6 6 6 6 4 3 7 7 7 7 7 7 7 7 5 þ *kp*þ<sup>1</sup> *μ<sup>p</sup>*þ<sup>1</sup> *∂ψ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> C*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>q</sup> H ∂ψ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> D*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>q</sup> H ψ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> C*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>q</sup> V ψ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> D*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>q</sup> V* 2 6 6 6 6 6 6 4 3 7 7 7 7 7 7 5 � *kp μp ∂ξpp <sup>n</sup> Apq H ∂ξpp <sup>n</sup> Bpq H ξpp <sup>n</sup> <sup>A</sup>pq V ξpp <sup>n</sup> Bpq V* 2 6 6 6 6 6 4 3 7 7 7 7 7 5 � *kp μp ∂ψpp <sup>n</sup> Cpq H ∂ψpp <sup>n</sup> Dpq <sup>H</sup>* <sup>þ</sup> *<sup>δ</sup><sup>q</sup> p* � � *ψpp <sup>n</sup> Cpq V ψpp <sup>n</sup> Dpq <sup>V</sup>* <sup>þ</sup> *<sup>δ</sup><sup>q</sup> p* � � 2 6 6 6 6 6 6 6 4 3 7 7 7 7 7 7 7 5 ¼ *σ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> ξpp <sup>n</sup> Apq H ξpp <sup>n</sup> Bpq H ∂ξpp <sup>n</sup> Apq V ∂ξpp <sup>n</sup> Bpq V* 2 6 6 6 6 6 4 3 7 7 7 7 7 5 þ *σ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> ψpp <sup>n</sup> Cpq H ψpp <sup>n</sup> Dpq <sup>H</sup>* <sup>þ</sup> *<sup>δ</sup><sup>q</sup> p* � � *∂ψpp <sup>n</sup> Cpq V ∂ψpp <sup>n</sup> Dpq <sup>V</sup>* <sup>þ</sup> *<sup>δ</sup><sup>q</sup> p* � � 2 6 6 6 6 6 6 6 4 3 7 7 7 7 7 7 7 5 (64)

In the above equations, *ψpq <sup>n</sup>* ¼ *j <sup>n</sup> kpRq* � � , *ξpq <sup>n</sup>* <sup>¼</sup> *<sup>h</sup>*ð Þ<sup>1</sup> *<sup>n</sup> kpRq* � �, *∂ψpq <sup>n</sup>* ¼ 1*=ρ* d *ρ j <sup>n</sup>*ð Þ *<sup>ρ</sup>* � �� � *<sup>ρ</sup>*¼*kpRq* and *<sup>∂</sup> <sup>ξ</sup>pq <sup>n</sup>* <sup>¼</sup> <sup>1</sup>*=<sup>ρ</sup>* <sup>d</sup> *<sup>ρ</sup>h*ð Þ<sup>1</sup> *<sup>n</sup>* ð Þ *ρ* h i� � � *ρ*¼*kpRq* . To obtain recursive relations for unknown coefficients of Green's function, separating the above equations is necessary. Therefore [75]:

$$
\begin{bmatrix}
\boldsymbol{A}\_{H,V}^{(p+1)q} + \boldsymbol{\delta}\_{p+1}^{q} \\
\boldsymbol{B}\_{H,V}^{(p+1)q}
\end{bmatrix} = \frac{\mathbf{1}}{T\_{Fp}^{H,V}} \begin{bmatrix}
\boldsymbol{A}\_{H,V}^{pq} \\
\boldsymbol{B}\_{H,V}^{pq}
\end{bmatrix} + \frac{\boldsymbol{R}\_{Fp}^{H,V}}{T\_{Fp}^{H,V}} \begin{bmatrix}
\boldsymbol{C}\_{H,V}^{pq} \\
\boldsymbol{D}\_{H,V}^{pq} + \boldsymbol{\delta}\_{p}^{q}
\end{bmatrix} \tag{65}
$$

$$
\begin{bmatrix}
\mathbf{C}\_{H,V}^{(p+1)q} \\
\mathbf{D}\_{H,V}^{(p+1)q}
\end{bmatrix} = \frac{\mathbf{R}\_{p p}^{H,V}}{T\_{P p}^{H,V}} \begin{bmatrix}
\mathbf{A}\_{H,V}^{pq} \\
\mathbf{B}\_{H,V}^{pq}
\end{bmatrix} + \frac{\mathbf{1}}{T\_{P p}^{H,V}} \begin{bmatrix}
\mathbf{C}\_{H,V}^{pq} \\
\mathbf{D}\_{H,V}^{pq} + \delta\_p^q
\end{bmatrix} \tag{66}
$$

where the reflection and transmission coefficients for TE waves are computed as:

$$R\_{fp}^H = \frac{k\_{p+1}\mu\_p \partial \boldsymbol{\nu}\_n^{(p+1)p} \boldsymbol{\nu}\_n^{pp} - k\_p \mu\_{p+1} \partial \boldsymbol{\nu}\_n^{pp} \boldsymbol{\nu}\_n^{(p+1)p} + \mathbf{g} \boldsymbol{\nu}\_n^{pp} \boldsymbol{\nu}\_n^{(p+1)p}}{k\_{p+1}\mu\_p \partial \boldsymbol{\nu}\_n^{(p+1)p} \boldsymbol{\xi}\_n^{pp} - k\_p \mu\_{p+1} \partial \boldsymbol{\xi}\_n^{pp} \boldsymbol{\nu}\_n^{(p+1)p} + \mathbf{g} \boldsymbol{\xi}\_n^{pp} \boldsymbol{\nu}\_n^{(p+1)p}} \tag{67}$$

$$R\_{p\_p}^H = \frac{k\_{p+1}\mu\_p \partial \xi\_n^{(p+1)p} \xi\_n^{pp} - k\_p \mu\_{p+1} \partial \xi\_n^{pp} \xi\_n^{(p+1)p} + \mathbf{g}\_n \xi\_n^{pp} \xi\_n^{(p+1)p}}{k\_{p+1}\mu\_p \partial \xi\_n^{(p+1)p} \eta\_n^{pp} - k\_p \mu\_{p+1} \partial \eta\_n^{pp} \ \xi\_n^{(p+1)p} + \mathbf{g}\_n \nu\_n^{pp} \ \xi\_n^{(p+1)p}} \tag{68}$$

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

*T<sup>H</sup> Fp* ¼ *kp*þ1*μ<sup>p</sup> <sup>∂</sup>ψ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> <sup>ξ</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup>* � *<sup>ψ</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> <sup>∂</sup>ξ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup>* � � *kp*þ1*μq∂ψ*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup> <sup>ξ</sup>pp <sup>n</sup>* � *kpμp*þ<sup>1</sup>*∂ξpp <sup>n</sup> <sup>ψ</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup>* <sup>þ</sup> *<sup>g</sup> <sup>ξ</sup>pp <sup>n</sup> <sup>ψ</sup>*ð Þ *<sup>p</sup>*þ<sup>1</sup> *<sup>p</sup> <sup>n</sup>* (69)

$$T\_{p\_p}^H = \frac{k\_{p+1}\mu\_p \left(\partial \xi\_n^{(p+1)p} \,\mu\_n^{(p+1)p} - \partial \varphi\_n^{(p+1)p} \xi\_n^{(p+1)p}\right)}{k\_{p+1}\mu\_p \partial \xi\_n^{(p+1)p} \,\mu\_n^{pp} - k\_p \mu\_{p+1} \partial \varphi\_n^{pp} \xi\_n^{(p+1)p} + \mathbf{g} \,\varphi\_n^{pp} \xi\_n^{(p+1)p}} \tag{70}$$

Also, the reflection and transmission coefficients for the TM waves are:

$$R\_{Fp}^V = \frac{k\_{p+1}\mu\_p \nu\_n^{(p+1)p} \partial \nu\_n^{pp} - k\_p \mu\_{p+1} \nu\_n^{pp} \partial \nu\_n^{(p+1)p} + \mathbf{g} \, \partial \nu\_n^{pp} \partial \nu\_n^{(p+1)p}}{k\_{p+1}\mu\_p \nu\_n^{(p+1)p} \, \partial \xi\_n^{pp} - k\_p \mu\_{p+1} \xi\_n^{pp} \partial \nu\_n^{(p+1)p} + \mathbf{g} \, \partial \xi\_n^{pp} \partial \nu\_n^{(p+1)p}} \tag{71}$$

$$R\_{fp}^V = \frac{k\_{p+1}\mu\_p \wp\_n^{(p+1)p} \partial \wp\_n^{pp} - k\_p \mu\_{p+1} \wp\_n^{pp} \partial \wp\_n^{(p+1)p} + \mathbf{g} \, \partial \wp\_n^{pp} \partial \wp\_n^{(p+1)p}}{k\_{p+1}\mu\_p \wp\_n^{(p+1)p} \, \partial \xi\_n^{pp} - k\_p \mu\_{p+1} \xi\_n^{pp} \, \partial \wp\_n^{(p+1)p} + \mathbf{g} \, \partial \xi\_n^{pp} \, \partial \wp\_n^{(p+1)p}} \tag{72}$$

$$T\_{Fp}^V = \frac{k\_{p+1}\mu\_p \left(\nu\_n^{(p+1)p} \partial \xi\_n^{(p+1)p} - \xi\_n^{(p+1)p} \partial \nu\_n^{(p+1)p}\right)}{k\_{p+1}\mu\_p \nu\_n^{(p+1)p} \partial \xi\_n^{pp} - k\_p \mu\_{p+1} \xi\_n^{pp} \partial \nu\_n^{(p+1)p} + \mathfrak{g} \ \partial \xi\_n^{pp} \partial \nu\_n^{(p+1)p}} \tag{73}$$

$$T\_{Pp}^V = \frac{k\_{p+1}\mu\_p \left(\xi\_n^{(p+1)p} \partial \mathsf{w}\_n^{(p+1)p} - \mathsf{w}\_n^{(p+1)p} \partial \xi\_n^{(p+1)p}\right)}{k\_{p+1}\mu\_p \xi\_n^{(p+1)p} \partial \mathsf{w}\_n^{pp} - k\_p \mu\_{p+1} \mathsf{w}\_n^{pp} \partial \xi\_n^{(p+1)p} + \mathsf{g} \ \partial \mathsf{w}\_n^{pp} \partial \xi\_n^{(p+1)p}} \tag{74}$$

In the above formulas, the subscripts *F* and *P* represent the outgoing and incoming waves, respectively. The symbols *T<sup>H</sup>* ð Þ *<sup>P</sup>*,*<sup>F</sup> <sup>p</sup>* and *<sup>R</sup><sup>H</sup>* ð Þ *<sup>P</sup>*,*<sup>F</sup> <sup>p</sup>* express the transmission and reflection of TE waves (due to the presence of superscript *H*), whereas the expression *T<sup>V</sup>* ð Þ *<sup>P</sup>*,*<sup>F</sup> <sup>p</sup>* and *<sup>R</sup><sup>V</sup>* ð Þ *<sup>P</sup>*,*<sup>F</sup> <sup>p</sup>* express the transmission and reflection of the TM waves (due to the presence of superscript *V*). Using the matrix form:

$$\begin{bmatrix} A\_{H,V}^{(p+1)q} + \delta\_{p+1}^q & B\_{H,V}^{(p+1)q} \\ & C\_{H,V}^{(p+1)q} & D\_{H,V}^{(p+1)q} \end{bmatrix} = \begin{bmatrix} \mathbf{1} & R\_{fp}^{H,V} \\ T\_{fp}^{H,V} & T\_{fp}^{H,V} \\ R\_{fp}^{H,V} & \mathbf{1} \\ T\_{fp}^{H,V} & T\_{fp}^{H,V} \end{bmatrix} \begin{bmatrix} A\_{H,V}^{pq} & B\_{H,V}^{pq} \\ C\_{H,V}^{pq} & D\_{H,V}^{pq} + \delta\_{p}^q \end{bmatrix} \tag{75}$$

Recursive relations can use to start through: *ANq <sup>H</sup>*,*<sup>V</sup>* <sup>¼</sup> *<sup>B</sup>Nq <sup>H</sup>*,*<sup>V</sup>* <sup>¼</sup> *<sup>C</sup>*<sup>1</sup>*<sup>q</sup> <sup>H</sup>*,*<sup>V</sup>* <sup>¼</sup> *<sup>D</sup>*<sup>1</sup>*<sup>q</sup> <sup>H</sup>*,*<sup>V</sup>* ¼ 0. It should be noted that due to the difficulty of constructing multiple concentric graphene shells, one can consider each boundary as a dielectric or perfect electric conductor (PEC) in the optical design. For example, an optical absorber consisting of a metal-dielectric spherical resonator whose outermost layer is coated with graphene is shown to enhance the absorption of the resonator [79]. To design optical structures such as absorber and invisible cloaks using Green's function formulation, it is necessary to transfer the point source to infinity to resemble a plane wave. In this case, considering the expansion of Green's function, it can be observed that the only unknown coefficients of expansion are *Bpq <sup>H</sup>*,*<sup>V</sup>* and *<sup>D</sup>pq H*,*V* coefficients. Using the convolution integral it can be shown that [23]:

$$
\overline{\overline{\mathbf{G}}}\_{0\epsilon} J + \overline{\overline{\mathbf{G}}}\_{\epsilon\epsilon}^{(11)} J = \overline{\overline{\mathbf{G}}}\_{\epsilon\epsilon}^{(21)} J \tag{76}
$$

The source dependency is the same in all of the above Green's functions and can be simplified from both sides of the equation. It is observed that the equation is converted into the equation resulting from the Mie analysis of the spherically layered structures. By directly starting from the Mie-Lorentz theory, the same results can be obtained [80]. Moreover, due to the sub-wavelength nature of the localized graphene plasmons, the final formulas can be simplified using the polynomial approximation of the special functions [81]. These structures can also be used as the building blocks of optical meta-materials [82].

#### *3.3.2 Purcell factor and energy transfer between donor-acceptor emitters*

As mentioned earlier, one of the important applications of Green's function is studying the interaction of dipole emitters in the vicinity of nanostructures. For this purpose, a vertical dipole in the vicinity of the graphene-based spherical structure, which was introduced in the previous section, is considered and its Purcell factor is calculated using the Green's function. Assuming that the field point and observation points for the dipole with moment *d*<sup>⊥</sup> <sup>0</sup> ¼ *d*0^*r* are in the same location of*r*<sup>0</sup> ¼ Δ >*R*1, *θ*<sup>0</sup> ¼ 0, and *ϕ*<sup>0</sup> ¼ 0, by calculating the scattered field using the convolution integral and the use of [18]:

$$\frac{\Gamma\_{total}^{d0}}{\Gamma\_0} = \mathbf{1} + \frac{6\pi\varepsilon\_0}{k\_1^3 d\_0^2} \text{Im} \left[ \mathbf{d}\_0, E\_{sca}^{d\_0}(r') \right] \tag{77}$$

where symbol *Im* represents the imaginary part of the complex function. The decay rate can be calculated. For this purpose, using relationships:

$$P\_n^m(\cos \theta)|\_{\theta \to 0} \simeq \sin^m \theta = \begin{cases} 1 \ m = 0 \\ 0 \ m \neq 0 \end{cases} \tag{78}$$

$$\frac{dP\_n^{\mathfrak{m}}(\cos \theta)}{d\theta}|\_{\theta \to 0} \simeq m \sin^{(m-1)}\theta = \begin{cases} 1 & m = 1 \\ 0 & m \neq 1 \end{cases} \tag{79}$$

The scattered field can be calculated. Thus [73]:

$$\frac{\Gamma\_{total}^{\perp}}{\Gamma\_0} = 1 - \frac{3}{2} \sum\_{n} n(n+1)(2n+1)\Re\left[\left(\frac{z\_n(k\_1\Delta)}{k\_1\Delta}\right)^2 B\_V^{11}\right] \tag{80}$$

As can be seen, the above equation is in full agreement with [83] which has extracted the decay rate for a core-shell plasmonic sphere, whereas in this case, it is necessary to use the Mie coefficients of graphene-based structure, namely, *B*<sup>11</sup> *<sup>V</sup>* . Using the derived formulas, the positions of the dipole can be considered arbitrarily. The transferred energy between the donor-acceptor pairs can be calculated straightforwardly using the Green's function *G* as [84]:

$$\frac{\Gamma(\boldsymbol{\alpha})}{\Gamma\_0(\boldsymbol{\alpha})} = \frac{|\mathbf{d}\_A \cdot \mathbf{G}(\boldsymbol{r}\_A, \boldsymbol{r}\_B, \boldsymbol{\alpha}) . d\_B|^2}{\left|\mathbf{d}\_A . \mathbf{G}\_0(\boldsymbol{r}\_A, \boldsymbol{r}\_B, \boldsymbol{\alpha}) . d\_B\right|^2} \tag{81}$$

where the subscript 0 refers to the free-space parameters and *d*<sup>A</sup> and *d*<sup>B</sup> are respectively the dipole moments of the acceptor and donor.

#### **4. Conclusion**

In conclusion, dyadic Green's function extraction for planarly, cylindrically, and spherically layered medium based on scattering superposition method is a unified

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

approach to deal with a wide range of electromagnetic problems in the realm of biomedicine. Specifically, the interaction of the human skin, body, and head with the electromagnetic sources with arbitrary distributions can be studied. Moreover, by engineering the constitutive parameters of the layers, a variety of novel devices for medical diagnostics and treatment can be proposed. Plasmonic metals and 2D materials are two main categories of such materials. For the sake of efficient analytical analysis, the impedance boundary condition is satisfied in the case of 2D materials to be used in the design of compact devices.

### **Author details**

Shiva Hayati Raad\* and Zahra Atlasbaf\* Department of Electrical and Computer Engineering, Tarbiat Modares University, Tehran, Iran

\*Address all correspondence to: shiva.hayati@modares.ac.ir and atlasbaf@modares.ac.ir

© 2021 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

### **References**

[1] J. R. Monds and A. G. McDonald, "Determination of skin temperature distribution and heat flux during simulated fires using Green's functions over finite-length scales," *Applied Thermal Engineering*, vol. 50, no. 1, pp. 593–603, 2013.

[2] J. Schmitt, G. Zhou, E. Walker, and R. Wall, "Multilayer model of photon diffusion in skin," *JOSA A,* vol. 7, no. 11, pp. 2141–2153, 1990.

[3] A. Christ, T. Samaras, A. Klingenböck, and N. Kuster, "Characterization of the electromagnetic near-field absorption in layered biological tissue in the frequency range from 30 MHz to 6000 MHz," *Physics in Medicine Biology*, vol. 51, no. 19, p. 4951, 2006.

[4] N. Chahat, M. Zhadobov, L. Le Coq, S. I. Alekseev, and R. Sauleau, "Characterization of the interactions between a 60-GHz antenna and the human body in an off-body scenario," *IEEE Transactions on Antennas Propagation*, vol. 60, no. 12, pp. 5958– 5965, 2012.

[5] S. Caorsi, M. Pastorino, and M. Raffetto, "Analytic SAR computation in a multilayer elliptic cylinder for bioelectromagnetic applications," *Bioelectromagnetics: Journal of the Bioelectromagnetics Society,The Society for Physical Regulation in Biology Medicine,The European Bioelectromagnetics Association*, vol. 20, no. 6, pp. 365–371, 1999.

[6] K. S. Nikita, G. S. Stamatakos, N. K. Uzunoglu, and A. Karafotias, "Analysis of the interaction between a layered spherical human head model and a finite-length dipole," *IEEE Transactions on Microwave Theory Techniques,* vol. 48, no. 11, pp. 2003–2013, 2000.

[7] S. Koulouridis and K. S. Nikita, "Study of the coupling between human head and cellular phone helical antennas," *IEEE*

*Transactions on electromagnetic Compatibility*, vol. 46, no. 1, pp. 62-70, 2004.

[8] J. Kim and Y. Rahmat-Samii, "Implanted antennas inside a human body: Simulations, designs, and characterizations," *IEEE Transactions on microwave theory techniques*, vol. 52, no. 8, pp. 1934-1943, 2004.

[9] F. Liu and S. Crozier, "Electromagnetic fields inside a lossy, multilayered spherical head phantom excited by MRI coils: models and methods," Physics in Medicine Biology, vol. 49, no. 10, p. 1835, 2004.

[10] F. Keshmiri and C. Craeye, "Moment-method analysis of normalto-body antennas using a Green's function approach," *IEEE transactions on antennas propagation*, vol. 60, no. 9, pp. 4259–4270, 2012.

[11] J. S. Meiguni, M. Kamyab, and A. Hosseinbeig, "Theory and experiment of spherical aperture-coupled antennas," *IEEE Transactions on Antennas and Propagation,* vol. 61, no. 5, pp. 2397–2403, 2013.

[12] H. Khodabakhshi and A. Cheldavi, "Irradiation of a six-layered spherical model of human head in the near field of a half-wave dipole antenna," *IEEE transactions on microwave theory techniques*, vol. 58, no. 3, pp. 680-690, 2010.

[13] W. C. Chew, *Waves and fields in inhomogeneous media*. IEEE press, 1995.

[14] J. Sun, C.-F. Wang, J. L.-W. Li, and M.-S. Leong, "Mixed potential spatial domain Green's functions in fast computational form for cylindrically stratified media," *Progress In Electromagnetics Research,* vol. 45, pp. 181–199, 2004.

[15] T. Itoh, "Spectral domain immitance approach for dispersion characteristics

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

of generalized printed transmission lines," *IEEE transactions on Microwave Theory and Techniques,* vol. 28, no. 7, pp. 733–736, 1980.

[16] M. Thiel and A. Dreher, "Dyadic Green's function of multilayer cylindrical closed and sector-structures for waveguide, microstrip-antenna, and network analysis," *IEEE transactions on microwave theory and techniques,* vol. 50, no. 11, pp. 2576–2579, 2002.

[17] T. V. B. Giang, M. Thiel, and A. Dreher, "A unified approach to the analysis of radial waveguides, dielectric resonators, and microstrip antennas on spherical multilayer structures," *IEEE transactions on microwave theory and techniques,* vol. 53, no. 1, pp. 404–409, 2005.

[18] L. Novotny and B. Hecht, *Principles of nano-optics*. Cambridge university press, 2012.

[19] C. A. Marocico, X. Zhang, and A. L. Bradley, "A theoretical investigation of the influence of gold nanosphere size on the decay and energy transfer rates and efficiencies of quantum emitters," *The Journal of chemical physics,* vol. 144, no. 2, p. 024108, 2016.

[20] V. Karanikolas, C. A. Marocico, and A. L. Bradley, "Spontaneous emission and energy transfer rates near a coated metallic cylinder," *Physical Review A,* vol. 89, no. 6, p. 063817, 2014.

[21] K. Rustomji, R. Abdeddaim, C. M. de Sterke, B. Kuhlmey, and S. Enoch, "Measurement and simulation of the polarization-dependent Purcell factor in a microwave fishnet metamaterial," *Physical Review B,* vol. 95, no. 3, p. 035156, 2017.

[22] A. Vakil, "Transformation optics using graphene: one-atom-thick optical devices based on graphene," 2012.

[23] C.-T. Tai, *Dyadic Green functions in electromagnetic theory*. Institute of Electrical & Electronics Engineers (IEEE), 1994.

[24] L.-W. Cai, "On the computation of spherical Bessel functions of complex arguments," *Computer Physics Communications,* vol. 182, no. 3, pp. 663–668, 2011.

[25] D. Bhattacharya, B. Ghosh, P. K. Goswami, and K. Sarabandi, "Evaluation of Efficient Green's Functions for Spherically Stratified Media," *IEEE Transactions on Antennas and Propagation,* vol. 66, no. 3, pp. 1613–1618, 2018.

[26] Y. S. Cao, L. J. Jiang, and A. E. Ruehli, "An equivalent circuit model for graphene-based terahertz antenna using the PEEC method," *IEEE Transactions on Antennas and Propagation,* vol. 64, no. 4, pp. 1385–1393, 2016.

[27] A. Y. Nikitin, F. J. Garcia-Vidal, and L. Martin-Moreno, "Analytical expressions for the electromagnetic dyadic Green's function in graphene and thin layers," *IEEE Journal of Selected Topics in Quantum Electronics,* vol. 19, no. 3, pp. 4600611–4600611, 2012.

[28] G. W. Hanson, "Dyadic Green's functions and guided surface waves for a surface conductivity model of graphene," *Journal of Applied Physics,* vol. 103, no. 6, p. 064302, 2008.

[29] G. W. Hanson, "Dyadic Green's functions for an anisotropic, non-local model of biased graphene," *IEEE Transactions on antennas and propagation,* vol. 56, no. 3, pp. 747–757, 2008.

[30] M. Cuevas, "Surface plasmon enhancement of spontaneous emission in graphene waveguides," *Journal of Optics,* vol. 18, no. 10, p. 105003, 2016

[31] T. J. Arruda, R. Bachelard, J. Weiner, and P. W. Courteille, "Tunable Fano resonances in the decay rates of a pointlike emitter near a graphenecoated nanowire," *Physical Review B,* vol. 98, no. 24, p. 245419, 2018.

[32] E. L. Tan and S. Y. Tan, "A unified representation of the dyadic Green's

functions for planar, cylindrical and spherical multilayered biisotropic media," *Progress In Electromagnetics Research,* vol. 20, pp. 75–100, 1998

[33] E. L. Tan and S. Y. Tan, "Coordinate-independent dyadic formulation of the dispersion relation for bianisotropic media," *IEEE Transactions on Antennas Propagation,* vol. 47, no. 12, pp. 1820– 1824, 1999.

[34] E. L. Tan and S. Y. Tan, "Concise spectral formalism in the electromagnetics of bianisotropic media," *Progress In Electromagnetics Research,* vol. 25, pp. 309–331, 2000.

[35] E. Tan, "Vector wave function expansions of dyadic Green's functions for bianisotropic media," *IEE Proceedings-Microwaves, Antennas propagation,* vol. 149, no. 1, pp. 57–63, 2002.

[36] W.-B. Lu, H. Chen, and Z.-G. Liu, "A review of microwave devices based on CVD-grown graphene with experimental demonstration," *EPJ Applied Metamaterials*, vol. 6, p. 8, 2019.

[37] J. Zhang, Z. Liu, W. Lu, H. Chen, B. Wu, and Q. Liu, "A low profile tunable microwave absorber based on graphene sandwich structure and high impedance surface," *International Journal of RF Microwave Computer-Aided Engineering,* vol. 30, no. 2, p. e22022, 2020.

[38] D. Yi, X.-C. Wei, and Y.-L. Xu, "Tunable microwave absorber based on patterned graphene," *IEEE Transactions on Microwave Theory Techniques,* vol. 65, no. 8, pp. 2819–2826, 2017.

[39] Z. Fang *et al.*, "Active tunable absorption enhancement with graphene nanodisk arrays," vol. 14, no. 1, pp. 299– 304, 2014.

[40] B. Wu *et al.*, "Experimental demonstration of a transparent graphene millimetre wave absorber with 28% fractional bandwidth at 140 GHz," *Scientific reports,* vol. 4, p. 4130, 2014.

[41] C. Srichan *et al.*, "Highly-sensitive surface-enhanced Raman spectroscopy (SERS)-based chemical sensor using 3D graphene foam decorated with silver nanoparticles as SERS substrate," *Scientific reports,* vol. 6, p. 23733, 2016.

[42] J. B. Maurya and Y. K. Prajapati, "Experimental Demonstration of DNA Hybridization Using Graphene-Based Plasmonic Sensor Chip," *Journal of Lightwave Technology*, 2020.

[43] L. A. Falkovsky, "Optical properties of graphene and IV–VI semiconductors," *Physics-Uspekhi,* vol. 51, no. 9, p. 887, 2008.

[44] H. Nasari and M. Abrishamian, "Nonlinear terahertz frequency conversion via graphene microribbon array," *Nanotechnology,* vol. 27, no. 30, p. 305202, 2016.

[45] S. Thongrattanasiri, F. H. Koppens, and F. J. G. De Abajo, "Complete optical absorption in periodically patterned graphene," *Physical review letters,* vol. 108, no. 4, p. 047401, 2012.

[46] S. H. Raad, Z. Atlasbaf, M. Shahabadi, and J. Rashed-Mohassel, "Dyadic Green's Function for the Tensor Surface Conductivity Boundary Condition," *IEEE Transactions on Magnetics,* vol. 55, no. 11, pp. 1–7, 2019.

[47] Y. R. Padooru, A. B. Yakovlev, C. S. Kaipa, G. W. Hanson, F. Medina, and F. Mesa, "Dual capacitive-inductive nature of periodic graphene patches: Transmission characteristics at lowterahertz frequencies," *Physical Review B,* vol. 87, no. 11, p. 115401, 2013.

[48] A. Andryieuski and A. V. Lavrinenko, "Graphene metamaterials based tunable terahertz absorber: effective surface conductivity approach," *Optics express,* vol. 21, no. 7, pp. 9144–9155, 2013.

*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

[49] Y. Guo, T. Zhang, W.-Y. Yin, and X.- H. Wang, "Improved hybrid FDTD method for studying tunable graphene frequency-selective surfaces (GFSS) for THz-wave applications," *IEEE Transactions on Terahertz Science and Technology,* vol. 5, no. 3, pp. 358–367, 2015.

[50] J. S. Gomez-Diaz, M. Tymchenko, and A. Alu, "Hyperbolic plasmons and topological transitions over uniaxial metasurfaces," *Physical Review Letters,* vol. 114, no. 23, p. 233901, 2015.

[51] J. Gomez-Diaz and A. Alu, "Flatland optics with hyperbolic metasurfaces," *ACS Photonics,* vol. 3, no. 12, pp. 2211– 2224, 2016.

[52] A. Fallahi and J. Perruisseau-Carrier, "Design of tunable biperiodic graphene metasurfaces," *Physical Review B,* vol. 86, no. 19, p. 195408, 2012.

[53] D. L. Sounas and C. Caloz, "Gyrotropy and nonreciprocity of graphene for microwave applications," *IEEE Transactions on Microwave Theory and Techniques,* vol. 60, no. 4, pp. 901– 914, 2012.

[54] L. Lin, L.-S. Wu, W.-Y. Yin, and J.-F. Mao, "Modeling of magnetically biased graphene patch frequency selective surface (FSS)," in *2015 IEEE MTT-S InternationalMicrowaveWorkshop Series on AdvancedMaterials and Processes for RF and THz Applications (IMWS-AMP)*, 2015, pp. 1–3: IEEE.

[55] Y. Yermakov *et al.*, "Effective surface conductivity of optical hyperbolic metasurfaces: from far-field characterization to surface wave analysis," *Scientific reports,* vol. 8, no. 1, p. 14135, 2018.

[56] J. Komijani, J. Rashed-Mohassel, and A. Mirkamali, "Dyadic Green Functions for a Dielectric Layer on a PEMC Plane," *Progress In Electromagnetics Research,* vol. 6, pp. 9–22, 2009.

[57] J. Gómez-Díaz and J. Perruisseau-Carrier, "Propagation of hybrid transverse magnetic-transverse electric plasmons on magnetically biased graphene sheets," *Journal of Applied Physics,* vol. 112, no. 12, p. 124906, 2012.

[58] Raad, Shiva Hayati, Zahra Atlasbaf, and Mauro Cuevas. "Dyadic Green's function for the graphene-dielectric stack with arbitrary field and source points." arXiv preprint arXiv: 2104.03581 (2021).

[59] M. Cuevas, M. A. Riso, and R. A. Depine, "Complex frequencies and field distributions of localized surface plasmon modes in graphene-coated subwavelength wires," *Journal of Quantitative Spectroscopy and Radiative Transfer,* vol. 173, pp. 26–33, 2016.

[60] M. Hajati and Y. Hajati, "Plasmonic characteristics of two vertically coupled graphene-coated nanowires integrated with substrate," *Applied Optics*, vol. 56, no. 4, pp. 870–875, 2017.

[61] S. H. Raad and Z. Atlasbaf, "Dyadic analysis of a cylindrical wire consisting of a cover with fully-populated surface conductivity tensor," *Optics express,* vol. 27, no. 15, pp. 21214–21225, 2019.

[62] J. Gomez-Diaz, M. Tymchenko, and A. Alù, "Hyperbolic metasurfaces: surface plasmons, light-matter interactions, and physical implementation using graphene strips," *Optical Materials Express,* vol. 5, no. 10, pp. 2313–2329, 2015.

[63] K.-H. Kim, Y.-S. No, S. Chang, J.-H. Choi, and H.-G. Park, "Invisible hyperbolic metamaterial nanotube at visible frequency," *Scientific reports,* vol. 5, p. 16027, 2015.

[64] Y. Mazor and A. Alù, "Angularmomentum selectivity and asymmetry in highly confined wave propagation along sheath-helical metasurface tubes," *Physical Review B,* vol. 99, no. 15, p. 155425, 2019.

[65] Y. Huang, Y. Shen, C. Min, and G. Veronis, "Switching photonic nanostructures between cloaking and superscattering regimes using phasechange materials," *Optical Materials Express,* vol. 8, no. 6, pp. 1672–1685, 2018.

[66] Y. Hua and T. K. Sarkar, "Generalized pencil-of-function method for extracting poles of an EM system from its transient response," 1989.

[67] C. Tokgoz and G. Dural, "Closedform Green's functions for cylindrically stratified media," *IEEE Transactions on Microwave theory and Techniques,* vol. 48, no. 1, pp. 40–49, 2000.

[68] S. H. Raad, C. J. Zapata-Rodríguez, and Z. Atlasbaf, "Multi-frequency super-scattering from sub-wavelength graphene-coated nanotubes," *JOSA B,* vol. 36, no. 8, pp. 2292–2298, 2019.

[69] S. H. Raad, C. J. Zapata-Rodríguez, and Z. Atlasbaf, "Graphene-coated resonators with frequency-selective super-scattering and super-cloaking," *Journal of Physics D: Applied Physics,* vol. 52, no. 49, p. 495101, 2019.

[70] S. H. Raad, Z. Atlasbaf, and M. Cuevas, "Scattering from Multilayered Graphene-Based Cylindrical and Spherical Particles," in *Nanoplasmonics*: IntechOpen, 2020.

[71] X. Ma *et al.*, "Graphene oxide wrapped gold nanoparticles for intracellular Raman imaging and drug delivery," vol. 1, no. 47, pp. 6495–6500, 2013.

[72] D. Cai, L. Ding, S. Wang, Z. Li, M. Zhu, and H. Wang, "Facile synthesis of ultrathin-shell graphene hollow spheres for high-performance lithium-ion batteries," *Electrochimica Acta*, vol. 139, pp. 96-103, 2014.

[73] S. H. Raad, Z. Atlasbaf, C. J. Zapata-Rodríguez, M. Shahabadi, and J. Rashed-Mohassel, "Dyadic Green's function for the electrically biased graphene-based

multilayered spherical structures," *Journal of Quantitative Spectroscopy Radiative Transfer,* vol. 256, p. 107251, 2020.

[74] Z.-S. Wu, Z.-J. Li, H. Li, Q.-K. Yuan, and H.-Y. Li, "Off-axis Gaussian beam scattering by an anisotropic coated sphere," *IEEE Transactions on Antennas and Propagation,* vol. 59, no. 12, pp. 4740–4748, 2011.

[75] L.-W. Li, P.-S. Kooi, M.-S. Leong, and T.-S. Yee, "Electromagnetic dyadic Green's function in spherically multilayered media," *IEEE Transactions on Microwave Theory and Techniques,* vol. 42, no. 12, pp. 2302–2310, 1994.

[76] M. Farhat, C. Rockstuhl, and H. Bağcı, "A 3D tunable and multi-frequency graphene plasmonic cloak," *Optics express,* vol. 21, no. 10, pp. 12592–12603, 2013.

[77] F. TING, "Fast solutions of electromagnetic fields in layered media," 2008.

[78] J. D. Shumpert, "Modeling of Periodic Dielectric Structures (electomagnetic Crystals)," University of Michigan, 2001.

[79] M. Wan *et al.*, "Strong tunable absorption enhancement in graphene using dielectric-metal core-shell resonators," *Scientific reports,* vol. 7, no. 1, p. 32, 2017.

[80] S. H. Raad, Z. Atlasbaf, J. Rashed-Mohassel, and M. Shahabadi, "Scattering from graphene-based multilayered spherical structures," *IEEE Transactions on Nanotechnology,* vol. 18, pp. 1129–1136, 2019.

[81] S. H. Raad and Z. Atlasbaf, "Equivalent RLC ladder circuit for scattering by graphene-coated nanospheres," *IEEE Transactions on Nanotechnology,* vol. 18, pp. 212–219, 2019.

[82] S. H. Raad and Z. Atlasbaf, "Tunable optical meta-surface using graphene*Dyadic Green's Function for Multilayered Planar, Cylindrical, and Spherical… DOI: http://dx.doi.org/10.5772/intechopen.95834*

coated spherical nanoparticles," *AIP Advances,* vol. 9, no. 7, p. 075224, 2019.

[83] T. J. Arruda, A. S. Martinez, F. A. Pinheiro, R. Bachelard, S. Slama, and P. W. Courteille, "Fano resonances in plasmonic core-shell particles and the Purcell effect," in *Fano Resonances in Optics and Microwaves*: Springer, 2018, pp. 445–472.

[84] V. D. Karanikolas and A. L. Bradley, "Interactions Between Quantum Emitters in the Presence of Plasmonic Nanostructures," Trinity College Dublin, 2016.

#### **Chapter 2**

## A TLM Formulation Based on Fractional Derivatives for Dispersive Cole-Cole Media

*Mohammed Kanjaa, Otman El Mrabet and Mohsine Khalladi*

#### **Abstract**

An auxiliary differential equation (ADE) transmission line method (TLM) is proposed for broadband modeling of electromagnetic (EM) wave propagation in biological tissues with the Cole-Cole dispersion Model. The fractional derivative problem is surmounted by assuming a linear behavior of the polarization current when the time discretization is short enough. The polarization current density is approached using Lagrange extrapolation polynomial and the fractional derivation is obtained according to Riemann definition of a fractional *α*-order derivative. Reflection coefficients at an air/muscle and air/fat tissues interfaces simulated in a 1-D domain are found to be in good agreement with those obtained from the analytic model over a broad frequency range, demonstrating the validity of the proposed approach.

**Keywords:** computational electromagnetics, numerical methods, transmission line matrix, biological system modeling, Cole-Cole medium, fractional derivative

#### **1. Introduction**

In the last two decades there has been a growing interest in the interaction between biological tissues and electromagnetic field at microwave frequencies. New promising applications of this technology in biomedical engineering like microwave imaging [1, 2], minimal invasive cancer therapies as thermal ablations [3], ultrawide band temperature dependent dielectric spectroscopy [4] and EM dosimetry [5], rely heavily on an accurate mathematical model of the response of these tissues to an external electromagnetic field. Numerical resolution of the propagation problem within these tissues requires a robust mathematical model and the previous incorporation of the dielectric data that at all the working frequencies.

A first model of the time response in a time varying electric field of biological tissues was formulated by Debye [6] through a time decaying polarization current *j t*ðÞ¼ *<sup>∂</sup><sup>P</sup> <sup>∂</sup><sup>t</sup>* resulting from the time variation of the polarization vector *P t*ð Þ, the permittivity for this model is given in Eq. (1) and the corresponding argand diagram is depicted in **Figure 1**.

$$
\varepsilon\_r^\*(o) - \varepsilon\_{r\infty} = \varepsilon\_r(o) - \varepsilon\_{r\infty} - j\varepsilon\_r'(o) = \frac{\Delta\varepsilon\_p}{1 + jo\tau\_p} \tag{1}
$$

**Figure 1.** *Argand diagram for a unique pole Debye model obtained by representing Eq. (1) in the complexe plane.*

This model even if it fits the experimental results in liquids it loses its accuracy when applied over a large band of frequency or in the presence of more than one type of polar molecule. This non-Debye relaxation is attributed to the existence of different relaxation processes [7] each with its own relaxation time *τ* and its amplitude Δ*ε*. A more accurate approximation to the behavior of this category of dielectrics is given by the sum of the individual relaxation processes leading to a multipole model, the permittivity is given as the sum of the *p* terms corresponding to the *p* poles each with its own amplitude. In the case of the biological tissue there is a multiple contribution of each relaxation process resulting in a broadening of the relaxation zone. Cole-Cole [8] proposed an Argand diagram in which in which *ε<sup>r</sup>* the imaginary part of the complex permittivity is plotted as a function of its real part *εr*, following the empirical relation:

$$
\varepsilon\_r^\*\left(\boldsymbol{\alpha}\right) - \varepsilon\_{r\boldsymbol{\alpha}} = \varepsilon\_r(\boldsymbol{\alpha}) - \varepsilon\_{r\boldsymbol{\alpha}} - j\varepsilon\_r'(\boldsymbol{\alpha}) = \sum\_{p=1}^p \frac{\Delta\varepsilon\_p}{1 + \left(\boldsymbol{\,}j\boldsymbol{\sigma}\boldsymbol{\tau}\_p\right)^{a\_p}} \tag{2}
$$

where *ε<sup>r</sup>*<sup>∞</sup> is the optical relative permittivity, Δ*ε<sup>p</sup>* ¼ *εrs* � *ε<sup>r</sup>*<sup>∞</sup> is the amplitude of the *p* � *th* pole, *ε<sup>s</sup>* is the static relative permittivity, *ω* is the angular frequency , *τ<sup>p</sup>* is the relaxation time, *α<sup>p</sup>* the parameter that indicates the broadening of the dispersion for this pole, which in the Cole-Cole model must satisfy 0< *α<sup>p</sup>* < 1, in the case of *α<sup>p</sup>* ¼ 1 the model is simplified to a Debye dispersion problem and *j* ¼ ffiffiffiffi �p 1 (**Figure 2**).

The difficulty in the numerical implementation of such a model arises from the *α<sup>p</sup>* parameter which is a noninteger. The consequence is a fractional order differentiation in the time domain which is more challenging compared to the Debye time domain solution. To address this problem, authors in [9] used the Letnikov fractional derivative by introducing a Stirling asymptotic formula and a recursive relation for the polarization vector. The computational demands of the FDTD scheme were considerably reduced, but nevertheless it required the storage of a large number of previous values of the polarization vector. In [10], authors used the Z transform to formulate the frequency dependence between the electric flux density and the electric field, the fractional derivative was approximated by using a polynomial method. Rekanos et al. [11] used the Pade approximation, where the

*A TLM Formulation Based on Fractional Derivatives for Dispersive Cole-Cole Media DOI: http://dx.doi.org/10.5772/intechopen.96378*

**Figure 2.**

*Argand diagram for a simple pole Cole-Cole model obtained by representing Eq. (2) in the complexe plane.*

fractional power term is approximated using a Pade technique. All the latter approximations were made in the context of the FDTD method.

The TLM method, even if it is flexible and wide band and being a time domain method, cannot deal with the dispersive aspect of the Cole-Cole medium directly. In [12] an approach based on a convolution product between the susceptibility and the electric field and the temporal behavior is deduced by a DFT and a nonrecursive summation leading to a considerable computational cost and a nonnegligible error at high frequencies. The causality principle is used to justify the minor dependence of the recent susceptibility values on previous ones, but the problem of the fractional Differentegration wasn't addressed.

In this work an ADE-TLM algorithm to model the Cole-Cole dispersion. The polarization current density is approached using an extrapolation with Lagrange polynomial method and the fractional derivation is obtained using the Riemann definition of the *α*-order derivative. The auxiliary differential equation is used to establish the update equation of the polarization currents which are included later in the general structure of the SCN-TLM node.

#### **2. Formulation of the method**

In the Cole-Cole model the relationship between the polarization current related to the *pth* pole and the electric field is given by:

$$\mathbf{J}\_p(\boldsymbol{\omega}) = \varepsilon\_0 \Delta \varepsilon\_p \frac{j\boldsymbol{\alpha}}{\mathbf{1} + \left(\boldsymbol{j}\boldsymbol{\alpha}\tau\_p\right)^{a\_p}} \mathbf{E}(\boldsymbol{\omega}) \tag{3}$$

where *jωτ<sup>p</sup> <sup>α</sup><sup>p</sup>* is the power law function of frequency which in time domain results in a fractional derivative of order *α*:

$$\mathbf{J}\_p(t) + \tau^{a\_p} D^{a\_p} \mathbf{J}\_p(t) = \varepsilon\_0 \Delta \varepsilon\_p \frac{\partial \mathbf{E}(t)}{\partial t} \tag{4}$$

One could consider an analytic solution for this equation, but due to the fractional derivative this can only be done for particular values of *αp*, hence a numerical solution is necessary with a previous discretization of the temporal values of the vectors.

In order to obtain the discretized expression of **J***p*ð Þ*t* , the polarization current vector can be approached using a Lagrange interpolation in the time interval *ti* <*t*< *tj* with cardinal functions:

$$L\_i(t) = \prod\_{j=1; i \neq j}^{n} \frac{(t - t\_j)}{(t\_i - t\_j)} \tag{5}$$

the interpolated expression of **J***p*ð Þ*t* is then obtained:

$$\mathbf{J}\_p(t) = \sum\_{i=1}^n \mathbf{J}^i L\_i(t) \tag{6}$$

and by applying (5) in the time interval *n*Δ*t*<*t* <ð Þ *n* þ 1 Δ*t* we obtain the cardinal functions:

$$L\_n(t) = \frac{t - (n+1)\Delta t}{n\Delta t - (n+1)\Delta t} \qquad \qquad L\_{n+1}(t) = \frac{t - n\Delta t}{(n+1)\Delta t - n\Delta t} \tag{7}$$

By substituting (7) in (6) the polynomial extrapolation of the polarization current density between time steps *n*Δ*t* and ð Þ *n* þ 1 Δ*t* can be found in a straightforward manner as follow:

$$\mathbf{J}\_p(t) = \frac{t - (n+1)\Delta t}{n\Delta t - (n+1)\Delta t} \mathbf{J}\_p^n + \frac{t - n\Delta t}{(n+1)\Delta t - n\Delta t} \mathbf{J}\_p^{n+1} \tag{8}$$

which after simplification can be rewritten as:

$$\mathbf{J}\_p(t) = \frac{\mathbf{J}\_p^{n+1} - \mathbf{J}\_p^n}{\Delta t} t + (n+1)\mathbf{J}\_p^n - n\mathbf{J}\_p^{n+1} \tag{9}$$

Furthermore, by substituting (8) into (3) we derive the auxiliary differential equation with a fractional order derivative given by:

$$\frac{\mathbf{J}\_p^{n+1} + \mathbf{J}\_p^n}{2} + \boldsymbol{\tau}^{a\_p} \frac{\mathbf{J}\_p^{n+1} - \mathbf{J}\_p^n}{\Delta t} D^{a\_p} t + D^{a\_p} \left( (n+1) \mathbf{J}\_p^n - n \mathbf{J}\_p^{n+1} \right) = \varepsilon\_0 \Delta \varepsilon\_p \frac{\partial \mathbf{E}}{\partial t} \tag{10}$$

The generalization of the derivation operator to arbitrary non-integer orders, has been subject to intensive research from mathematicians [13, 14] and in electromagnetism [15]. One of the definitions in [16] for the fractional derivative, symbolized by the operator **D***<sup>α</sup> <sup>x</sup>* is given for the power series with fractional exponents by considering a function *h x*ð Þ defined as a power series *h x*ð Þ¼ <sup>P</sup>*<sup>n</sup> <sup>i</sup>*¼<sup>1</sup>*Ai*ð Þ *<sup>x</sup>* � *<sup>a</sup> <sup>ν</sup>*þ*i=<sup>n</sup>* for *<sup>ν</sup>*<sup>&</sup>gt; � 1 and *<sup>n</sup>* a positive integer, so each term ð Þ *<sup>x</sup>* � *<sup>a</sup> <sup>p</sup>* could have a noninteger exponent *p* ¼ *ν* þ *i=n*, then according to Reimann [15] a fractional derivative of order *α* for this term is given by:

$$\_aD\_x^a(\mathbf{x} - a)^p = \frac{d^a(\mathbf{x} - a)^p}{d(\mathbf{x} - a)^a} = \frac{\Gamma(p + 1)}{\Gamma(p - a + 1)}(\mathbf{x} - a)^{p - a} \quad for \quad \mathbf{x} > a \tag{11}$$

where in our case *α* is a noninteger number [16].

*A TLM Formulation Based on Fractional Derivatives for Dispersive Cole-Cole Media DOI: http://dx.doi.org/10.5772/intechopen.96378*

Hence by applying the fractional derivation operator given in (11) to (10) and for *p* ¼ 1 and *a* ¼ 0 we obtain the update equation for the polarization current:

$$\mathbf{J}\_p^{n+1} = -\frac{\phi\_p}{\nu\_p} \mathbf{J}\_p^n + \frac{\zeta\_p}{\nu\_p} \left( \mathbf{E}^{n+1} - \mathbf{E}^n \right) \tag{12}$$

where the update parameters *ψp*, *ϕ<sup>p</sup>* and *ζ<sup>p</sup>* are:

$$\varphi\_p = \frac{\mathbf{1}}{2} + \frac{\tau\_p^{a\_p}}{\Delta t} \frac{\Gamma(2)}{\Gamma(2 - a\_p)} t^{1 - a\_p} - n \frac{\tau\_p^{a\_p}}{\Gamma(1 - a\_p)} t^{-a\_p} \tag{13}$$

$$\phi\_p = \frac{1}{2} - \frac{\tau\_p^{a\_p}}{\Delta t} \frac{\Gamma(2)}{\Gamma(2 - a\_p)} t^{1 - a\_p} + (n + 1) \frac{\tau\_p^{a\_p}}{\Gamma(1 - a\_p)} t^{-a\_p} \tag{14}$$

$$
\zeta\_p = \frac{\varepsilon\_0 \Delta \varepsilon\_p}{\Delta t} \tag{15}
$$

To get the update equation of the electric field components we start from the Maxwell-Ampere equation, and by including the conductivity term:

$$\nabla \times \mathbf{H} = \varepsilon\_0 \varepsilon\_\infty \frac{\partial \mathbf{E}}{\partial t} + \sigma\_0 \mathbf{E} + \sum\_{p=1}^p \mathbf{J}\_p \tag{16}$$

where *<sup>σ</sup>*<sup>0</sup> is the static ionic conductivity, Eq. (16) formulated at time step *<sup>n</sup>* <sup>þ</sup> <sup>1</sup> 2 and by approaching **J** *<sup>n</sup>*þ<sup>1</sup> <sup>2</sup> by the average of its values at time steps *n* and *n* þ 1

$$\mathbf{J}\_p^{n+\frac{1}{2}} = \frac{\mathbf{J}\_p^{n+1} + \mathbf{J}\_p^n}{2} \tag{17}$$

Finally, the updated electric field is given by:

$$E\_{x,y,x}^{n+1} = E\_{x,y,x}^{n} \left( 1 - \frac{2\sigma\_0 \Delta t}{D} \right) - \sum\_{p} \frac{\left( 1 - \frac{\phi\_p}{\nu\_p} \right) \Delta t}{D} I\_{px,y,x}^{n} + \frac{2\Delta t}{D} (\nabla \times \mathbf{H})\_{x,y,x}^{n+\frac{1}{2}} \tag{18}$$

$$D = 2\varepsilon\_0 \varepsilon\_\infty + \sum\_p \frac{\zeta\_p \Delta t}{\nu\_p} + \sigma\_0 \Delta t \tag{19}$$

#### **3. The TLM Formalisme**

The TLM method is a numerical technique based on the discretization of the computational domain according to Huygens principle as an alternative to the Maxwell equations used in the FDTD method [17]. In this method the simulation domain is discretized in cells where a series of uniform transmission lines, parallel and series stubs and additional sources are used to take account of the real characteristics of the propagation through the medium as given by Maxwell equations. Therefore, instead of electric and magnetic field components the electromagnetic field is represented by voltage and current waves, propagating through the unit cell circuit referred to as symmetrical condensed node (SCN). The relationship between the electromagnetic field and the voltage and current waves at the time step *n*Δ*t* and at the center of the node are formulated as follows [18].

**Figure 3.**

*Structure of the TLM symmetrical condensed node, (a) unit volume of thedielectric material with dimensions Δx, Δy,and Δz, (b) equivalent SCN.*

As can be seen in **Figure 3(b)**, the SCN consists of interconnected transmission lines. This structure models a unit cell of the propagation medium as in **Figure 3(a)**. Each face of the cell corresponds to two ports orthogonal to each other and labeled from 1 to 12, these stubs model the propagation through free space and have therefore its characteristic impedance *Z*<sup>0</sup> ¼ ffiffiffiffiffiffiffiffiffi *ε*0*μ*<sup>0</sup> p . Three additional open circuit (13,14,15) stubs must be added to the node to account for the dispersive behavior of the medium. The SCN nodes are connected to each other to form the simulation domain given in **Figure 4** for a 1-D propagation.

The relationship between the electromagnetic field and the voltage and current waves at the time step *n*Δ*t* and at the center of the node are formulated as follows [18]:

$$E\_{\mathbf{x},y,x}^{n} = \frac{V\_{\mathbf{x},y,x}^{n}}{\Delta l} \qquad \qquad H\_{\mathbf{x},y,x}^{n} = \frac{I\_{\mathbf{x},y,x}^{n}}{\Delta l} \tag{20}$$

In this algorithm as in [11, 19] the dispersive properties of the medium are accounted for by adding voltage sources *sV* to each node, these sources constitute the counterparts of the terms in the FDTD formulation that have no equivalents in the TLM formalism. Therefore the voltage at the center of the node, as deduced from the conservation laws [20], becomes:

$$V\_{x,y,x}^{n+1} = V\_{x,y,x}^{n} + \frac{1}{4 + Y\_{ox,y,x}} \left( {}\_{t}V\_{x,y,x}^{n+1} + {}\_{t}V\_{x,y,x}^{n} \right) + \frac{4}{4 + Y\_{ox,y,x}} \frac{\Delta l \Delta t}{\varepsilon\_{0}} (\nabla \times \mathbf{H})\_{x,y,x}^{n+\frac{1}{2}} \tag{21}$$

where *Yocx*,*y*,*<sup>z</sup>* indicates the normalized admittance of the open circuit stub added to the node.

When expressed in terms of incident voltages on the corresponding stubs and accounting for the voltage sources injected in the stubs (16,17,18) to complete the model of the Cole-Cole medium:

$$V\_{x,y,x}^{n+1} = \frac{2}{4 + Y\_{\alpha x,y,x}} \left[ V\_{1,3,5}^{i} + V\_{2,4,6}^{i} + V\_{9,8,7}^{i} + V\_{12,11,10}^{i} + Y\_{\alpha;x,y,x} V\_{13,14,15}^{i} + \frac{1}{2} V\_{16,17,18}^{n+1} \right] \tag{22}$$

where the subscript *i* indicates the incident pulses on the indicated stubs.

*A TLM Formulation Based on Fractional Derivatives for Dispersive Cole-Cole Media DOI: http://dx.doi.org/10.5772/intechopen.96378*

**Figure 4.** *Structure of the 1-D SCN grid used in this work to model the simulation domain.*

In the TLM formalism the update equation issued from the ADE in the Eq. (12) becomes:

$$J\_{p\mathbf{x},y,x}^{n+1} = -\frac{\phi\_p}{\nu\_p} J\_{p\mathbf{x},y,x}^n + \frac{\zeta\_p}{\nu\_p \Delta l} \left( V\_{\mathbf{x},y,x}^{n+1} - V\_{\mathbf{x},y,x}^n \right) \tag{23}$$

and the electric field update equation is also formulated in the TLM formalism as:

$$V\_{x,y,x}^{n+1} = V\_{x,y,x}^{n} \left(1 - \frac{2\sigma\_{0}\Delta t}{D}\right) - \sum\_{p} \frac{\left(1 - \frac{\phi\_{p}}{\nu\_{p}}\right)\Delta t\Delta l}{D} J\_{p x,y,x}^{n} + \frac{2\Delta t\Delta l}{D} (\nabla \times \mathbf{H})\_{x,y,x}^{n+1} \tag{24}$$

The analogy between Eqs. (22) and (24) the expression of the normalized admittance of the stub added to the SCN node is obtained straightforward:

$$Y\_{\text{oc}} = 4\left(\frac{D}{2\varepsilon\_0} - 1\right) \tag{25}$$

and the update equation for voltage sources at the center of the node:

$$\mathbf{V}\_{s}V\_{\mathbf{x},\mathbf{y},\mathbf{z}}^{n+1} = -{}\_{s}V\_{\mathbf{x},\mathbf{y},\mathbf{z}}^{n} - \frac{4\sigma\_{0}\Delta t}{\varepsilon\_{0}}V\_{\mathbf{x},\mathbf{y},\mathbf{z}}^{n} - \sum\_{\mathbf{p}} \frac{2\left(\mathbf{1} - \frac{\phi\_{\mathbf{p}}}{\mathbf{y}\_{\mathbf{p}}}\right)\Delta t\Delta l}{\varepsilon\_{0}}I\_{\mathbf{p}\mathbf{x},\mathbf{y},\mathbf{z}}\tag{26}$$

that can be simplified to

$$\mathbf{C}\_{s}\mathbf{V}\_{\mathbf{x},\mathbf{y},\mathbf{z}}^{n+1} = -\,\_{s}\mathbf{V}\_{\mathbf{x},\mathbf{y},\mathbf{z}}^{n} + \mathbf{C}\_{1}\mathbf{V}\_{\mathbf{x},\mathbf{y},\mathbf{z}}^{n} + \sum\_{\mathbf{p}} \mathbf{C}\_{2\mathbf{p}}I\_{\mathbf{p}\mathbf{x},\mathbf{y},\mathbf{z}}\tag{27}$$

with update constants:

$$C\_1 = -\frac{4\sigma\_0 \Delta t}{\varepsilon\_0} \tag{28}$$

$$C\_{2p} = -\frac{2\left(1 - \frac{\phi\_p}{\nu\_p}\right)\Delta t \Delta l}{\varepsilon\_0} \tag{29}$$

A final step to establish the TLM formulation is to express the scattering process on the capacitive stubs:

$$\left[\boldsymbol{V}\_{\mathbf{1}3,\mathbf{1}4,\mathbf{15}}^{r}\right]^{n} = \left[\boldsymbol{V}\_{\mathbf{x},\mathbf{y},\mathbf{z}}\right]^{n} - \left[\boldsymbol{V}\_{\mathbf{13},\mathbf{14},\mathbf{15}}^{i}\right]^{n} \tag{30}$$

#### **4. Simulation and results**

In order to verify the new proposed algorithm, and for the sake of simplicity a one dimensional problem is simulated, where a Cole-Cole medium (fat or muscle) occupies the region *z*≥0, while the rest of space is air, 2 dimensional and 3 dimensional propagation and reflection models can be naturally deducted by ensuring the connexion between TLM cells in the *y* axis using ports 1,5,7 and 12, to account for the propagation along the *x* axis ports 3,6,10 and 11 must be connected to adjacent cells as in **Figure 5**.

The incident wave is a derivative Gaussian pulse given by *E* !*inc* ¼

ð Þ *t*�*t*<sup>0</sup> *<sup>τ</sup>*<sup>2</sup> *exp* �4*<sup>π</sup>* ð Þ *<sup>t</sup>*�*t*<sup>0</sup> <sup>2</sup> *τ*2 ^*ex* where *t*<sup>0</sup> ¼ 200Δ*t* and *τ* ¼ 220Δ*t* polarized in the *x* direction, and propagates along the *z* axis. In the TLM implementation, the considered 1- D simulation domain (see **Figure 6**) consists of a grid of 1000 cells interconnected on the z axis as in **Figure 4**, 500 of which were used to model the Cole-Cole medium

**Figure 5.**

*2-D ADE-TLM simulation of an incident and reflected pulse on the air/muscle, the tissue is modeled by a 4 pole Cole-Cole model, the interface is located at the TLM cell z* ¼ 100*, at time steps 160 (a) and 218 (b).*

#### **Figure 6.**

*2-D ADE-TLM simulation of a reflected and transmitted pulse on the air/muscle , the tissue is modeled by a 4 pole Cole-Cole model, the interface is located at the TLM cell z = 100 , at (a) time step 305, and (b) at time step 446.*

#### *A TLM Formulation Based on Fractional Derivatives for Dispersive Cole-Cole Media DOI: http://dx.doi.org/10.5772/intechopen.96378*

and the remaining cells were used to model the air. The cell dimension Δ*l* satisfies the stability condition Δ*t* ¼ Δ*l=*2*c*, the time step is set as Δ*t* ¼ 0*:*125*ps*. The simulation domain is truncated by an unsplit PML layer [21] of 10 cells at the beginning and the end.

The fourth order Cole-Cole parameters for fat tissue as well as for muscle tissue are listed in **Table 1** [9]. The values of the electric field are recorded at a point P located at the center of the SCN node 10 cells before the air-medium interface. **Figure 7(a)** shows the incident impulse at iteration 399 propagating in the air and **Figure 7(b)** at iteration 927 depicts the reflected and transmitted impulse on the air/muscle interface.


**Table 1.**

*Cole-Cole constants for tissues of muscle and fat.*

**Figure 7.**

*Propagation of an incident derivative Gaussian pulse through Cole-Cole medium model of human muscle (a) at iteration 399 (incident pulse) (b) at iteration 677 (reflected and transmitted pulse).*

**Figure 8.** *Reflection coefficient (dB) at normal incidence on an air/muscle interface.*

**Figure 9.** *Reflection coefficient (dB) at normal incidence on an air/fat interface.*

**Figure 10.**

*2-D ADE-TLM simulation of a reflected and transmitted pulse on the air/muscle, the tissue is modeled by a 4 pole Cole-Cole model, the interface is located at the TLM cell z* ¼ 100*, at time steps are 305 and 446.*

In **Figures 8** and **9**, the simulation results for both cases (fat and muscle) are compared to the theoretical reflection at the air/Cole-Cole medium interface which can be obtained by using the following equation [22]:

$$|R| = \frac{|\mathbf{1} - \sqrt{\varepsilon}\_r|}{|\mathbf{1} + \sqrt{\varepsilon}\_r|} \tag{31}$$

A good agreement over the whole frequency band is observed. The slight discrepancy between the numerical and theoretical results can be ascribed to the approximation made to the polarization current value when performing the Lagrange extrapolation on each iteration. The propagation and reflection through a 2-D TLM

#### **Figure 11.**

*Reflection coefficient (dB) at normal incidence on an air/adenoma (a), air/Normal tyroidian tissue interface (b).*

#### **Figure 12.** *Reflection coefficient (dB) at normal incidence on an air.*

grid modeling the incidence on an air/muscle interface is presented in **Figures 5** and **10** at different time steps, the interface is located at *z* ¼ 100 in a 10 by 200 cells grid, **Figure 5(a)** shows the Gaussian wave before the arrival on the interface while **Figure 10(a)** and **(b)** depicts it after reflection and transmission. A more precise simulation that targeted the Dielectric properties of pathological thyroid tissue types including adenoma, thyroiditis and the properties of the healthy thyroid based on the results of the experimental study presented in a previous published work in [2], and for the healthy thyroid tissue [23], in both works experimental results are used to extract the 2 poles Cole-Cole model parameters. Based upon this parameters we conducted a simulation for the 4 tissue types, the results of this simulation are presented in **Figure 11(a)** and **(b)** for the air/adenoma and the air/normal thyroidian interface respectively, **Figure 12(a)** and **(b)** depicts the reflection coefficient on an air/ thyroyditis affected thyroidian tissue and an air/tumor respectively also in this case a perfect agreement over the working frequency is observed.

#### **5. Conclusion**

Numerical methods are an essential part of the modeling process, new propagation media with anomalous relaxation properties impose innovative modeling methods. An effective ADE-TLM formulation way to model electromagnetic wave propagation in biological tissues using Cole-Cole dispersion model. The fractionalorder derivative in the Cole-Cole model is tackled by using a polynomial extrapolation of the polarization current. This approximation reduces considerably both the numerical cost of the simulation and its complexity since only two previous values of the *J* field are needed at each iteration. a Reimann derivation on the obtained polynomial extrapolation of the polarization current is then performed to solve the fractional order derivative. The presented results indicate the accuracy of our approach.

### **Conflict of interest**

The authors declare no conflict of interest.

### **Author details**

Mohammed Kanjaa\*†, Otman El Mrabet† and Mohsine Khalladi† University of Abdelmalek Essadi, Tetouan, Morocco

\*Address all correspondence to: mohammedkanjaa@uae.ac.ma

† These authors contributed equally.

© 2021 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

*A TLM Formulation Based on Fractional Derivatives for Dispersive Cole-Cole Media DOI: http://dx.doi.org/10.5772/intechopen.96378*

#### **References**

[1] Martellosio A, Pasian M, Bozzi M, Perregrini A, Mazzanti F, Svelto PE, et al. L: 0.5-50 ghz dielectric characterisation of breast cancer tissues. Electronics Letters. 2015;**51**(13):974-975

[2] Gavazzi S, Limone P, De Rosa G, Molinari F, and Vecchi G: Comparison of microwave dielectric properties of human normal, benign and malignant thyroid tissues obtained from surgeries: A preliminary study, Biomedical Physics Engineering Express, vol. 4, no. 4, p. 047 003, 2018.

[3] Ruvio G, Eaton-Evans J, Shahzad A, OHalloran M. Numerical evaluation of microwave thermal ablation to treat small adrenocortical masses. International Journal of RF and Microwave Computer-Aided Engineering. 2018;**28**(3):e21236

[4] Ley S, Schilling S, Fiser O, Vrba J, Sachs J, Helbig M. Ultra-wideband temperature dependent dielectric spectroscopy of porcine tissue and blood in the microwave frequency range. Sensors. 2019;**19**(7):1707

[5] Chakarothai J, Wake K, Watanabe S. Convergence of a single-frequency FDTD solution in numerical dosimetry. IEEE Transactions on Microwave Theory and Techniques. 2016;**64**(3): 707-714

[6] Debye P. Part i. dielectric constant. energy absorption in dielectrics with polar molecules. Transactions of the Faraday Society. 1934;**30**:679-684

[7] Cole KS, Cole RH. Dispersion and absorption in dielectrics i. alternating current characteristics. The Journal of Chemical Physics. 1941;**9**(4):341-351

[8] Rekanos IT, Yioultsis TV. Approximation of grunwald-letnikov fractional derivative for FDTD modeling of Cole-Cole media. IEEE

Transactions on magnetics. 2014;**50**(2): 181-184

[9] Guo B, Li J, Zmuda H. A new FDTD formulation for wave propagation in biological media with Cole-Cole model. IEEE Microwave and Wireless Components Letters. 2006;**16**(12): 633-635

[10] Rekanos IT, Papadopoulos TG. An auxiliary differential equation method for FDTD modeling of wave propagation in Cole-Cole dispersive media. IEEE Transactions on Antennas and Propagation. 2010;**58**(11): 3666-3674

[11] Barba I. Cabeceira A C L, Panizo M, and Represa J: Modelling dispersive dielectrics in TLM method, International Journal of Numerical Modelling: Electronic Networks. Devices and Fields. 2001;**14**(1):15-30

[12] Samko SG, Kilbas AA, Marichev OI, et al. Fractional integrals and derivatives. Gordon and Breach Science Publishers, Yverdon Yverdon-les-Bains, Switzerland, 1993. In: vol. 1993

[13] Techniques in the fractional calculus," in The Fractional Calculus, ser. Mathematics in Science and Engineering, K. B. Oldham and J. Spanier, Eds., vol. 111, Elsevier, 1974, pp. 133–160.

[14] N. Engheta: n the role of fractional calculus in electromagnetic theory," IEEE Antennas and Propagation Magazine, vol. 39, no. 4, pp. 35–46, Aug. 1997.

[15] Engheta N. On fractional calculus and fractional multipoles in electromagnetism. IEEE Transactions on Antennas and Propagation. Apr. 1996;**44**(4):554-566

[16] Wharmby AW, Bagley RL. The application of the fractional calculus model for dispersion and absorption in dielectrics in terahertz waves. International Journal of Engineering Science. 2015;**93**:1-12

[17] Jin H, Vahldieck R. Direct derivations of TLM symmetrical condensed node and hybrid symmetrical condensed node from maxwell's equations using centered differencing and averaging. IEEE Transactions on Microwave Theory and Techniques. Dec. 1994;**42**(12):2554-2561

[18] Christopoulos C, The Transmission-Line Modeling (TLM) Method in Electromagnetics. Morgan Claypool, 2006.

[19] Cabeceira ACL, Barba I, Grande A, Represa J. A 2D-TLM model for electromagnetic wave propagation in chiral media. Microwave and Optical Technology Letters. 2005;**46**(2):180-182

[20] Yaich MI, Khalladi M. The far-zone scattering calculation of frequencydependent materials objects using the tlm method. IEEE Transactions on Antennas and Propagation. Nov. 2002; **50**(11):1605-1608

[21] Yaich MI, Kanjaa M, Adraoui S, Mounirh K, Khalladi M. An unsplit formulation of the 3D-PML absorbing boundary conditions for TLM-method in time domain. In: 2018 6th International Conference on Multimedia Computing and Systems (ICMCS). May 2018. pp. 1-5

[22] Taflove A, Hagness SC. Computational electrodynamics: the finite-difference time-domain method. Artech house. 2005

[23] Gabriel S, Lau RWand Gabriel C1996b The dielectric properties of biological tissues: III. Parametric models for the dielectric spectrum of tissues Phys. Med. Biol. 41 227193

#### **Chapter 3**

## Averaged No-Regret Control for an Electromagnetic Wave Equation Depending upon a Parameter with Incomplete Initial Conditions

*Abdelhak Hafdallah and Mouna Abdelli*

#### **Abstract**

This chapter concerns the optimal control problem for an electromagnetic wave equation with a potential term depending on a real parameter and with missing initial conditions. By using both the average control notion introduced recently by E. Zuazua to control parameter depending systems and the no-regret method introduced for the optimal control of systems with missing data. The relaxation of averaged no-regret control by the averaged low-regret control sequence transforms the problem into a standard optimal control problem. We prove that the problem of average optimal control admits a unique averaged no-regret control that we characterize by means of optimality systems.

**Keywords:** optimal control, averaged no-regret control, electromagnetic wave equation, parameter depending equation, systems with missing data

#### **1. Introduction**

The research in the field of electromagnetism is set to become a vital factor in biomedical technologies. Those studies included several areas like the usage of electromagnetic waves for probing organs and advanced MRI techniques, microwave biosensors, non-invasive electromagnetic diagnostic tools, therapeutic applications of electromagnetic waves, radar technologies for biosensing, the adoption of electromagnetic waves in medical sensing, cancer detection using ultra-wideband signal, the interaction of electromagnetic waves with biological tissues and living systems, theoretical modeling of electromagnetic propagation through human body and tissues and imaging applications of electromagnetic.

Actually, the principal goal of the study is to control such electromagnetic waves to be compatible with some biomedical needs like X-rays in the framework of medical screening and wireless power transfer of electromagnetic waves through the human body [1] where we want to make waves closer to a desired distribution.

In this chapter, we consider a linear wave equation with a potential term *p x*ð Þ , *σ* supposed dependent on space variable *x* and real parameter *σ* ∈ ð Þ 0, 1 , this term generally comprises the dielectric permittivity of the medium which has different

properties and cannot be exactly presented, this is because of the difference or lack of knowledge of the physical properties of the material penetrated from the electromagnetic waves. The initial position and velocity are also supposed unknown.

In this study, we consider an optimal control problem for electromagnetic wave equation depending upon a parameter and with missing initial conditions. We use the method of no-regret control which was introduced firstly in statistics by Savage [2] and later by Lions [3, 4] where he used this concept in optimal control theory, and its related idea is "low-regret" control to apply it to control distributed systems of incomplete data which has the attention of many scholars [5–12], motivated by various applications in ecology, and economics as well [13]. Also, we use the notion of average control because our system depends upon a parameter, Zuazua was the first who introduced this new concept in [14].

The rest of this chapter is arranged as follows. Section 2, lists the definition of the problem we are studying. Section 3, is devoted to the study of the averaged noregret control and the averaged low-regret control for the electromagnetic wave equation. Ultimately, we prove the existence of a unique average low-regret control, and the characterization of the average optimal is given in Section 4. Finally, we make a conclusion in Section 5.

#### **2. Statement of the problem**

Consider a bounded open domain Ω with a smooth boundary ∂Ω. We set Σ ¼ ð Þ� 0, *T* ∂Ω and *Q* ¼ Ω � ð Þ 0, *T* . We introduce the following linear electromagnetic wave equation depending on a parameter

$$\begin{cases} \frac{\partial^2 y}{\partial t^2} - \Delta y + p(\mathbf{x}, \boldsymbol{\sigma}) y = \mathbf{0} & \text{in } Q \\\\ y = \begin{cases} v & \text{on } \Sigma\_0 \\ \mathbf{0} & \text{on } \Sigma \backslash \Sigma\_0 \end{cases} \\\\ y(\mathbf{x}, \mathbf{0}) = y\_0(\mathbf{x}); \frac{\partial y}{\partial t}(\mathbf{x}, \mathbf{0}) = y\_1(\mathbf{x}) \end{cases} \tag{1}$$

where *<sup>p</sup>* <sup>∈</sup>*L*<sup>∞</sup>ð Þ <sup>Ω</sup> is the potential term supposed dependent on a real parameter *σ* ∈ð Þ 0, 1 presents the dielectric permittivity and permeability of the medium and such that 0<sup>&</sup>lt; *<sup>α</sup>*<sup>1</sup> <sup>≤</sup> *p x*ð Þ , *<sup>σ</sup>* <sup>≤</sup>*α*<sup>2</sup> *<sup>a</sup>:e:*in <sup>Ω</sup>, *<sup>v</sup>* is a boundary control in *<sup>L</sup>*<sup>2</sup> ð Þ <sup>Σ</sup><sup>0</sup> , *<sup>y</sup>*<sup>0</sup> <sup>∈</sup> *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ Ω , *y*<sup>1</sup> ∈ *L*<sup>2</sup> ð Þ Ω are the initial position and velocity respectively, both supposed unknown. For all *<sup>σ</sup>* <sup>∈</sup>ð Þ 0, 1 , the wave Eq. (1) has a unique solution *y v*, *<sup>y</sup>*0, *<sup>y</sup>*1, *<sup>σ</sup>* � � in *<sup>C</sup>* ½ � 0, *<sup>T</sup>* ; *<sup>H</sup>*<sup>1</sup> ð Þ <sup>Ω</sup> � � <sup>∩</sup>*C*<sup>1</sup> ½ � 0, *<sup>T</sup>* ; *<sup>L</sup>*<sup>2</sup> ð Þ <sup>Ω</sup> � � [15].

Denote by *g* ¼ *y*0, *y*<sup>1</sup> � �∈ *H*<sup>1</sup> <sup>0</sup>ð Þ� <sup>Ω</sup> *<sup>L</sup>*<sup>2</sup> ð Þ Ω the initial data. We want to choose a control *u* independently of *σ* and *g* in a way such that the average state function *y* approaches a given observation *yd* ∈*L*<sup>2</sup> ð Þ *Q* . To achieve our goal, let' associate to (1) the following quadratic cost functional

$$J(\boldsymbol{v}, \mathbf{g}) = \left\| \int\_{0}^{1} \boldsymbol{\mathcal{y}}(\boldsymbol{v}, \mathbf{g}, \boldsymbol{\sigma}) d\boldsymbol{\sigma} - \boldsymbol{\mathcal{y}}\_{d} \right\|\_{L^{2}(Q)}^{2} + N \left\| \boldsymbol{v} \right\|\_{L^{2}(\Sigma\_{0})}^{2} \tag{2}$$

where *N* ∈ <sup>∗</sup> þ .

In this work, we aim to characterize the solution *u* of the optimal control problem with missing data given by

*Averaged No-Regret Control for an Electromagnetic Wave Equation Depending upon… DOI: http://dx.doi.org/10.5772/intechopen.95447*

$$\inf\_{v \in L^2(\Sigma\_0)} J(v, \mathbf{g}) \text{ subject to} (\mathbf{1}) \tag{3}$$

independently of *g* and *σ*.

#### **3. Averaged no-regret control and averaged low-regret control for the electromagnetic wave equation**

A classical method to obtain the optimality system is then to solve the minmax problem

$$\inf\_{v \in L^2(\Sigma\_0)} \left( \sup\_{\mathbf{g} \in H\_0^1(\Omega) \times L^2(\Omega)} (J(v, \mathbf{g})) \right),\tag{4}$$

but *J v*ð Þ , *g* is not upper bounded since sup*<sup>g</sup>* <sup>∈</sup> *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ� <sup>Ω</sup> *<sup>L</sup>*2ð Þ <sup>Ω</sup> ð Þ¼þ *J v*ð Þ , *<sup>g</sup>* <sup>∞</sup>. A natural idea of Lions [3] is to search for controls *v* such that

$$J(v, \mathbf{g}) - J(\mathbf{0}, \mathbf{g}) \le \mathbf{0}, \forall \mathbf{g} \in H\_0^1(\Omega) \times L^2(\Omega) \tag{5}$$

Those controls *v* are called averaged no-regret controls.

As in [16, 17], we introduce the averaged no-regret control defined by.

Definition 1 [1] We say that *v*∈*L*<sup>2</sup> ð Þ Σ<sup>0</sup> is an averaged no-regret control for (1) if *v* is a solution of

$$\inf\_{v, v \in L^2(\Sigma\_0)} \left( \sup\_{\mathbf{g} \in H\_0^1(\mathfrak{A}) \times L^2(\mathfrak{A})} (f(v, \mathbf{g}) - f(\mathbf{0}, \mathbf{g})) \right). \tag{6}$$

Let us start by giving the following important lemma. Lemma 1 For all *v*∈ *L*<sup>2</sup> ð Þ Σ<sup>0</sup> and *g* ∈ *G* we have

$$(J(\nu, \mathbf{g}) - J(\mathbf{0}, \mathbf{g})) \quad = \quad J(\nu, \mathbf{0}) - J(\mathbf{0}, \mathbf{0}) \tag{7}$$

$$-2\int\_{\Omega} \mathcal{y}\_0(\mathbf{x}) \int\_0^1 \frac{\partial \zeta}{\partial t}(\mathbf{x}, \mathbf{0}) d\sigma d\mathbf{x} + 2\int\_{\Omega} \mathcal{y}\_1(\mathbf{x}) \int\_0^1 \zeta(\mathbf{x}, \mathbf{0}) d\sigma d\mathbf{x} \tag{8}$$

where *ζ* is given by the following backward wave equation

$$\begin{cases} \frac{\partial^2 \zeta}{\partial t^2} - \Delta \zeta + p(\varkappa, \sigma) \zeta = \int\_0^1 \jmath(v, \mathbf{0}, \sigma) d\sigma & \text{in } Q \\ \zeta = \mathbf{0} & \text{on } \Sigma \\ \zeta(\varkappa, T) = \mathbf{0}, \frac{\partial \zeta}{\partial t}(\varkappa, T) = \mathbf{0} \end{cases} \tag{9}$$

which has a unique solution in *<sup>C</sup>* ½ � 0, *<sup>T</sup>* ; *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ <sup>Ω</sup> � � <sup>∩</sup>*C*<sup>1</sup> ½ � 0, *<sup>T</sup>* ; *<sup>L</sup>*<sup>2</sup> ð Þ <sup>Ω</sup> � � [15]. **Proof.** It's easy to check that for all ð Þ *<sup>v</sup>*, *<sup>g</sup>* <sup>∈</sup> *<sup>L</sup>*<sup>2</sup> ð Þ� Σ<sup>0</sup> *G*

$$f(v, \mathbf{g}) - f(0, \mathbf{g}) = f(v, \mathbf{0}) - f(0, \mathbf{0}) + 2 \int\_{Q} \left( \int\_{0}^{1} y(v, \mathbf{0}) d\sigma \right) \left( \int\_{0}^{1} y(0, \mathbf{g}) d\sigma \right) d\mathbf{x} dt. \tag{10}$$

Use (9) and apply Green formula to get

$$\int\_{Q} \left( \int\_{0}^{1} \mathbf{y}(\nu, \mathbf{0}) d\sigma \right) \left( \int\_{0}^{1} \mathbf{y}(\mathbf{0}, \mathbf{g}) d\sigma \right) dxdt = \int\_{0}^{T} \left[ \left( \frac{\partial^{2} \zeta}{\partial t^{2}} - \Delta \zeta + p(\mathbf{x}, \sigma) \zeta \right) \left( \int\_{0}^{1} \mathbf{y}(\mathbf{0}, \mathbf{g}) d\sigma \right) dxdt \right] \tag{11}$$

$$=-2\int\_{\mathfrak{U}}\mathcal{y}\_{0}(\boldsymbol{\varkappa})\int\_{0}^{1}\frac{\partial\boldsymbol{\zeta}}{\partial t}(\boldsymbol{v},\boldsymbol{0})d\boldsymbol{\sigma}d\boldsymbol{\varkappa}+2\int\_{\mathfrak{U}}\mathcal{y}\_{1}(\boldsymbol{\varkappa})\int\_{0}^{1}\boldsymbol{\zeta}(\boldsymbol{v},\boldsymbol{0})d\boldsymbol{\sigma}d\boldsymbol{\varkappa}.\tag{12}$$

■ The no-regret control seems to be hard to characterize (see [11]), for this. reason we relax the no-regret control problem by making some quadratic perturbation as follows.

Definition 2 [17] We say that *u<sup>γ</sup>* ∈*L*<sup>2</sup> ð Þ Σ<sup>0</sup> is an averaged low-regret control for (1) if *u<sup>γ</sup>* is a solution of

$$\inf\_{\nu,\nu \in L^2(\Sigma\_0)} \left( \sup\_{\mathbf{g} \in H\_0^1(\mathfrak{Q}) \times L^2(\mathfrak{Q})} (f(\nu,\mathbf{g}) - f(\mathbf{0},\mathbf{g}) - \gamma \left\| \mathcal{y}\_0 \right\|\_{H\_0^1(\mathfrak{Q})}^2 - \gamma \left\| \mathcal{y}\_1 \right\|\_{L^2(\mathfrak{Q})}^2 \right), \eta > 0. \tag{13}$$

Using (9) the problem (13) can be written as

$$\inf\_{\boldsymbol{v},\boldsymbol{\sigma}\in\boldsymbol{L}^{2}(\boldsymbol{\Sigma}\_{0})} \left( \begin{aligned} &f(\boldsymbol{v},\boldsymbol{0}) - f(\boldsymbol{0},\boldsymbol{0}) + \sup\_{\boldsymbol{\sigma}\in\boldsymbol{G}} 2(\int\_{\Omega} \boldsymbol{y}\_{1}(\boldsymbol{x}) \int\_{0}^{1} \boldsymbol{\zeta}(\boldsymbol{v},\boldsymbol{\sigma})(\boldsymbol{x},\boldsymbol{0}) d\boldsymbol{\sigma} d\boldsymbol{x} - \int\_{\Omega} \boldsymbol{y}\_{0}(\boldsymbol{x}) \int\_{0}^{1} \frac{d\boldsymbol{\zeta}}{\boldsymbol{\mathcal{U}}}(\boldsymbol{v},\boldsymbol{\sigma})(\boldsymbol{x},\boldsymbol{0}) d\boldsymbol{\sigma} d\boldsymbol{x} \\ & \quad - \boldsymbol{\eta} \|\boldsymbol{y}\_{0}\|\_{H^{1}\_{0}(\Omega)}^{2} - \boldsymbol{\eta} \|\boldsymbol{y}\_{1}\|\_{L^{2}(\Omega)}^{2}). \end{aligned} \right), \boldsymbol{\eta} > 0. \tag{14}$$

And thanks to Legendre transform (see [18, 19]), we have

$$\sup\_{\xi \in G} 2 \left( - \int\_{\Omega} \boldsymbol{y}\_0(\mathbf{x}) \int\_0^1 \frac{\partial \zeta}{\partial t}(\mathbf{v}, \mathbf{0}) d\sigma d\mathbf{x} + \int\_{\Omega} \boldsymbol{y}\_1(\mathbf{x}) \int\_0^1 \zeta(\mathbf{v}, \mathbf{0}) d\sigma d\mathbf{x} - \gamma \left\| \boldsymbol{y}\_0 \right\|\_{H^1\_0(\Omega)}^2 - \gamma \left\| \boldsymbol{y}\_1 \right\|\_{L^2(\Omega)}^2 \right), \tag{15}$$

$$=\frac{1}{\mathcal{Y}}\left\|\int\_{0}^{1}\frac{\partial\zeta(\boldsymbol{\nu},\boldsymbol{\sigma})}{\partial t}(\boldsymbol{\varkappa},\mathbf{0})d\boldsymbol{\sigma}\right\|\_{H\_{0}^{1}(\Omega)}^{2}+\frac{1}{\mathcal{Y}}\left\|\int\_{0}^{1}\zeta(\boldsymbol{\nu},\boldsymbol{\sigma})(\boldsymbol{\varkappa},\mathbf{0})d\boldsymbol{\sigma}\right\|\_{L^{2}(\Omega)}^{2}.\tag{16}$$

Then, the averaged low-regret control problem (9) is equivalent to the following classical optimal control problem

$$\inf\_{\nu \in L^2(\Sigma\_0)} J\_{\mathcal{I}}(\nu) \tag{17}$$

where

$$J\_{\mathcal{I}}(\boldsymbol{v}) = \boldsymbol{f}(\boldsymbol{v}, \mathbf{0}) - \boldsymbol{f}(\mathbf{0}, \mathbf{0}) + \frac{1}{\mathcal{Y}} \left\| \int\_{0}^{1} \frac{\partial \boldsymbol{\zeta}(\boldsymbol{v}, \boldsymbol{\sigma})}{\partial \mathbf{t}}(\boldsymbol{\omega}, \mathbf{0}) d\boldsymbol{\sigma} \right\|\_{H^{1}\_{0}(\Omega)}^{2} + \frac{1}{\mathcal{Y}} \left\| \int\_{0}^{1} \boldsymbol{\zeta}(\boldsymbol{v}, \boldsymbol{\sigma})(\boldsymbol{\omega}, \mathbf{0}) d\boldsymbol{\sigma} \right\|\_{L^{2}(\Omega)}^{2} \dots \tag{18}$$

#### **4. Characterizations**

In the recent section, we aim to find a full characterization for the averaged no-regret control and averaged low-regret control via optimality systems.

Theorem 1.1 There exists a unique averaged low-regret control *u<sup>γ</sup>* solution to (17), (18).

*Averaged No-Regret Control for an Electromagnetic Wave Equation Depending upon… DOI: http://dx.doi.org/10.5772/intechopen.95447*

**Proof.** We have for every *v*∈*L*<sup>2</sup> ð Þ Σ<sup>0</sup> : *J<sup>γ</sup>* ð Þ*v* ≥ � *J*ð Þ 0, 0 , this means that (17), (18) has a solution.

Let *v<sup>γ</sup> n* � �∈ *L*<sup>2</sup> ð Þ Σ<sup>0</sup> be a minimizing sequence such that

$$\lim\_{n \to \infty} J\_\gamma(v\_n^r) = J\_\gamma(u\_\gamma) = d\_\gamma. \tag{19}$$

We know that

$$J\_{\boldsymbol{\gamma}}(\boldsymbol{v}\_{n}^{\boldsymbol{r}}) = \boldsymbol{f}(\boldsymbol{v}\_{n}^{\boldsymbol{r}},\boldsymbol{0}) - \boldsymbol{f}(\boldsymbol{0},\boldsymbol{0}) + \frac{\mathbf{1}}{\boldsymbol{\gamma}} \left\| \int\_{0}^{1} \frac{\partial \boldsymbol{\zeta}(\boldsymbol{v}\_{n}^{\boldsymbol{r}},\boldsymbol{\sigma})}{\partial \boldsymbol{t}}(\boldsymbol{x},\boldsymbol{0}) d\boldsymbol{\sigma} \right\|\_{H\_{0}^{1}(\mathfrak{U})}^{2} + \frac{\mathbf{1}}{\boldsymbol{\gamma}} \left\| \int\_{0}^{1} \boldsymbol{\zeta}(\boldsymbol{v}\_{n}^{\boldsymbol{r}},\boldsymbol{\sigma})(\boldsymbol{x},\boldsymbol{0}) d\boldsymbol{\sigma} \right\|\_{L^{2}(\mathfrak{U})}^{2} \leq d\_{\boldsymbol{\gamma}} + \mathbf{1}. \tag{20}$$

This implies the following bounds

$$\left\| \left| \boldsymbol{v}\_{n}^{\prime} \right| \right\|\_{L^{2}(\Sigma\_{0})} \leq \mathsf{C}\_{\mathcal{T}},\tag{21}$$

$$\left\| \left| \int\_{0}^{1} \mathcal{Y}(v\_{n}^{\mathcal{I}}, \mathbf{0}, \sigma) d\sigma \right| \right\|\_{L^{2}(Q)} \leq \mathcal{C}\_{\mathcal{I}},\tag{22}$$

$$\left\| \frac{\partial^2 \zeta\_n}{\partial t^2} - \Delta \zeta\_n + p(\mathbf{x}, \sigma) \zeta\_n \right\|\_{L^2(Q)} \le \mathcal{C}\_{\mathcal{V}},\tag{23}$$

where *C<sup>γ</sup>* is a positive constant independent of *n*. Moreover, by continuity w.r.t. data and (21) we get

$$\left\|\left|\mathcal{y}\left(\nu\_{n}^{\mathrm{r}},\mathbf{0},\sigma\right)\right\|\right\|\_{L^{2}\left(Q\right)}\leq\mathsf{C}\_{\mathcal{T}}.\tag{24}$$

By similar way an by using (22) we obtain

$$\left\| \left| \zeta \left( v\_n^{\mathcal{I}}, \sigma \right) \right| \right\|\_{L^2 \left( 0, T; H\_0^1(\Omega) \right)} \le \mathcal{C}\_{\mathcal{I}}.\tag{25}$$

Then, from (21) we deduce that there exists a subsequence still denoted *v<sup>γ</sup> n* � � such that *v<sup>γ</sup> <sup>n</sup> \* uγ*weakly in *L*<sup>2</sup> ð Þ Σ<sup>0</sup> , and from (22) we get

$$
\gamma(v\_n^r, 0, \sigma) \rightharpoonup \mathcal{Y}\_\gamma \text{ weakly in } L^2(Q). \tag{26}
$$

Also, because of continuity w.r.t. data we have *y v<sup>γ</sup> <sup>n</sup>*, 0, *<sup>σ</sup>* � � *\* y u<sup>γ</sup>* , 0, *<sup>σ</sup>* � � weakly in *L*<sup>2</sup> ð Þ *<sup>Q</sup>* , by limit uniqueness *<sup>y</sup><sup>γ</sup>* <sup>¼</sup> *y u<sup>γ</sup>* , 0, *<sup>σ</sup>* � � solution to

$$\begin{cases} \frac{\partial^2 \mathcal{y}\_\gamma}{\partial t^2} - \Delta \mathcal{y}\_\gamma + p(\mathbf{x}, \boldsymbol{\sigma}) \mathcal{y}\_\gamma = \mathbf{0} & \text{in } Q, \\ \mathcal{y}\_\gamma = \begin{cases} u\_\gamma & \text{on } \Sigma\_0, \\ 0 & \text{on } \Sigma | \Sigma\_0, \end{cases} \\ \mathcal{y}\_\gamma(\mathbf{x}, \mathbf{0}) = \mathcal{y}\_0(\mathbf{x}); \frac{\partial \mathcal{y}\_\gamma}{\partial t}(\mathbf{x}, \mathbf{0}) = \mathcal{y}\_1(\mathbf{x}) \end{cases} \tag{27}$$

In other hand, use (24) and (22) to apply the convergence dominated theorem and, we have

$$\int\_{0}^{1} \mathcal{y}(v\_n^r, 0, \sigma) d\sigma \rightharpoonup \int\_{0}^{1} \mathcal{y}(u\_\gamma, 0, \sigma) d\sigma \text{ weakly in } L^2(Q). \tag{28}$$

From (25) we deduce the existence of a subsequence still be denoted by *ζ v<sup>γ</sup> <sup>n</sup>*, *<sup>σ</sup>* � �ð Þ *<sup>x</sup>*, 0 such that

$$
\zeta(v\_n^\gamma, \sigma)(\mathbf{x}, \mathbf{0}) \rightharpoonup \zeta(u\_\gamma, \sigma)(\mathbf{x}, \mathbf{0}) \text{ weakly in } H\_0^1(\Omega), \tag{29}
$$

then

$$\frac{\partial^2 \zeta\_n}{\partial t^2} - \Delta \zeta\_n + p(\varkappa, \sigma) \zeta\_n \to \frac{\partial^2 \zeta\_\gamma}{\partial t^2} - \Delta \zeta\_\gamma + p(\varkappa, \sigma) \zeta\_\gamma \text{ in} D'(Q), \tag{30}$$

where *D Q*ð Þ¼ *<sup>C</sup>*<sup>∞</sup> <sup>0</sup> ð Þ *Q* , and (23) leads to

$$\frac{\partial^2 \zeta\_n}{\partial t^2} - \Delta \zeta\_n + p(\varkappa, \sigma) \zeta\_n \rightharpoonup f \text{ weakly in} L^2(Q). \tag{31}$$

Again, by limit uniqueness *<sup>∂</sup>*<sup>2</sup>*ζγ <sup>∂</sup>t*<sup>2</sup> � <sup>Δ</sup>*ζγ* <sup>þ</sup> *p x*ð Þ , *<sup>σ</sup> ζγ* <sup>¼</sup> <sup>Ð</sup> <sup>1</sup> <sup>0</sup>*y u<sup>γ</sup>* , 0, *<sup>σ</sup>* � �*d<sup>σ</sup>* in *<sup>L</sup>*<sup>2</sup> ð Þ *Q* . Finally, *ζγ* is a solution to

$$\begin{cases} \frac{\partial^2 \zeta\_\gamma}{\partial t^2} - \Delta \zeta\_\gamma + p(\mathbf{x}, \boldsymbol{\sigma}) \zeta\_\gamma = \int\_0^1 p(u\_r, \mathbf{0}, \boldsymbol{\sigma}) d\boldsymbol{\sigma} & \text{in } Q\_\gamma \\ \zeta\_\gamma = \mathbf{0} & \text{on } \Sigma, \\ \zeta\_\gamma(\mathbf{x}, T) = \mathbf{0}, \frac{\partial \zeta\_\gamma}{\partial t}(\mathbf{x}, T) = \mathbf{0} & \text{in } \Omega. \end{cases} \tag{32}$$

The uniqueness of *u<sup>γ</sup>* follows from strict convexity and weak lower semicontinuity of the functional *<sup>J</sup><sup>γ</sup>* ð Þ*<sup>v</sup>* .■

After proving existence and uniqueness, we aim in the next theorem to give a full description to the average low-regret control for the electromagnetic wave equation.

Theorem 1.2 For all *γ* >0, the average low-regret control *u<sup>γ</sup>* is characterized by the following optimality system

*∂*2 *yγ <sup>∂</sup>t*<sup>2</sup> � <sup>Δ</sup>*y<sup>γ</sup>* <sup>þ</sup> *p x*ð Þ , *<sup>σ</sup> <sup>y</sup><sup>γ</sup>* <sup>¼</sup> 0, *∂*2 *ζγ <sup>∂</sup>t*<sup>2</sup> � <sup>Δ</sup>*ζγ* <sup>þ</sup> *p x*ð Þ , *<sup>σ</sup> ζγ* <sup>¼</sup> ð1 0 *y u<sup>γ</sup>* , 0, *σ* � �*dσ*, *∂*2 *ργ <sup>∂</sup>t*<sup>2</sup> � <sup>Δ</sup>*ργ* <sup>þ</sup> *p x*ð Þ , *<sup>σ</sup> ργ* <sup>¼</sup> 0, *∂*2 *qγ <sup>∂</sup>t*<sup>2</sup> � <sup>Δ</sup>*q<sup>γ</sup>* <sup>þ</sup> *p x*ð Þ , *<sup>σ</sup> <sup>q</sup><sup>γ</sup>* <sup>¼</sup> ð1 0 *ργ* <sup>þ</sup> *y u<sup>γ</sup>* , 0, *<sup>σ</sup>* � � � � *<sup>d</sup><sup>σ</sup>* � *yd* in *<sup>Q</sup>*, *<sup>y</sup><sup>γ</sup>* <sup>¼</sup> *<sup>u</sup><sup>γ</sup>* 0 ( on Σ<sup>0</sup> onΣnΣ<sup>0</sup> , *ζγ* ¼ 0, *ργ* ¼ 0, *q<sup>γ</sup>* ¼ 0 on Σ, *y<sup>γ</sup>* ð Þ¼ 0, *x* 0, *∂yγ <sup>∂</sup><sup>t</sup>* ð Þ¼ 0, *<sup>x</sup>* 0, *ζγ* ð Þ¼ *x*, *T* 0, *∂ζγ <sup>∂</sup><sup>t</sup>* ð Þ¼ *<sup>x</sup>*, *<sup>T</sup>* 0, *ργ* ð Þ¼� *<sup>x</sup>*, 0 <sup>1</sup> *γ* ð1 0 *∂ζ u<sup>γ</sup>* � � *<sup>∂</sup><sup>t</sup>* ð Þ *<sup>x</sup>*, 0 *<sup>d</sup>σ*, *∂ργ <sup>∂</sup><sup>t</sup>* ð Þ¼ *<sup>x</sup>*, 0 <sup>1</sup> *γ* ð1 0 *ζ u<sup>γ</sup>* � �ð Þ *<sup>x</sup>*, 0 *<sup>d</sup>σ*, *q<sup>γ</sup>* ð Þ¼ *T*, *x* 0, *∂qγ <sup>∂</sup><sup>t</sup>* ð Þ¼ *<sup>T</sup>*, *<sup>x</sup>* <sup>0</sup> in <sup>Ω</sup>, 8 >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>< >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>: (33)

*Averaged No-Regret Control for an Electromagnetic Wave Equation Depending upon… DOI: http://dx.doi.org/10.5772/intechopen.95447*

with

$$\mu\_{\gamma} = \frac{1}{N} \int\_{0}^{1} \frac{\partial p\_{\gamma}}{\partial \eta} d\sigma \,\mathrm{in}L^{2}(\Sigma\_{0}). \tag{34}$$

**Proof.** From the first order necessary optimality conditions, we have

$$J^{\prime}(\boldsymbol{u}\_{\boldsymbol{\tau}})(\boldsymbol{w}) = \left(\int\_{0}^{1} \boldsymbol{\chi}(\boldsymbol{u}\_{\boldsymbol{\tau}}, \mathbf{0}) d\sigma - \boldsymbol{\chi}\_{d}, \int\_{0}^{1} \boldsymbol{\chi}(\boldsymbol{w}, \mathbf{0}) d\sigma\right)\_{L^{2}(\Omega)} + N(\boldsymbol{u}\_{\boldsymbol{\tau}}, \boldsymbol{w})\_{L^{2}(\Sigma\_{0})} + \tag{1}$$

$$\left(\frac{1}{\mathcal{I}}\int\_{0}^{1} \boldsymbol{\zeta}^{\prime}(\boldsymbol{u}\_{\boldsymbol{\tau}})(\mathbf{0}) d\sigma, \int\_{0}^{1} \boldsymbol{\zeta}^{\prime}(\boldsymbol{w})(\mathbf{0}) d\sigma\right)\_{L^{2}(\Omega)} + \left(\frac{1}{\mathcal{I}}\int\_{0}^{1} \boldsymbol{\zeta}(\boldsymbol{u}\_{\boldsymbol{\tau}})(\mathbf{0}) d\sigma, \int\_{0}^{1} \boldsymbol{\zeta}(\boldsymbol{w})(\mathbf{0}) d\sigma\right)\_{L^{2}(\Omega)} = \mathbf{0} \tag{35}$$

for all *w* ∈ *L*<sup>2</sup> ð Þ Σ<sup>0</sup> . Now, let us introduce *ργ* <sup>¼</sup> *<sup>ρ</sup> <sup>u</sup><sup>γ</sup>* , 0 � � unique solution to

$$\begin{cases} \frac{\partial^2 \rho\_\gamma}{\partial t^2} - \Delta \rho\_\gamma + p(\mathbf{x}, \boldsymbol{\sigma}) \rho\_\gamma = \mathbf{0} & \text{in } Q, \\\rho\_\gamma = \mathbf{0} & \text{on } \Sigma, \\\rho\_\gamma(\mathbf{x}, \mathbf{0}) = -\frac{1}{\mathcal{Y}} \int\_0^1 \frac{\partial \zeta(\boldsymbol{u}\_\gamma)}{\partial t}(\mathbf{x}, \mathbf{0}) d\boldsymbol{\sigma}; \frac{\partial \rho\_\gamma}{\partial t}(\mathbf{x}, \mathbf{0}) = \frac{1}{\mathcal{Y}} \int\_0^1 \zeta(\boldsymbol{u}\_\gamma)(\mathbf{x}, \mathbf{0}) d\boldsymbol{\sigma} & \text{in } \Omega. \end{cases}$$

So that for every *w* ∈ *L*<sup>2</sup> ð Þ Σ<sup>0</sup> , we obtain

$$\begin{split} \mathcal{J}^{\prime}(\boldsymbol{u}\_{\mathcal{I}})(\boldsymbol{w}) &= \left( \int\_{0}^{1} \boldsymbol{y}(\boldsymbol{u}\_{\mathcal{I}}, \mathbf{0}) + \rho\_{\boldsymbol{r}} d\sigma - \boldsymbol{y}\_{\mathcal{I}}, \int\_{0}^{1} \boldsymbol{y}(\boldsymbol{w}, \mathbf{0}) d\sigma - \boldsymbol{y}(\mathbf{0}, \mathbf{0}) \right)\_{\boldsymbol{L}^{2}(\mathbb{Q})} + N(\boldsymbol{u}\_{\mathcal{I}}, \boldsymbol{w})\_{\boldsymbol{L}^{2}(\mathbb{Q})} \\ &= \mathbf{0} \end{split} \tag{37}$$

We finally define another adjoint state *q<sup>γ</sup>* ¼ *q u<sup>γ</sup>* � � as the unique solution of

$$\begin{cases} \frac{\partial^2 q\_\gamma}{\partial t^2} - \Delta q\_\gamma + p(\mathbf{x}, \boldsymbol{\sigma}) q\_\gamma = \int\_0^1 (\rho\_\gamma + \chi(\boldsymbol{u}\_\gamma, \mathbf{0}, \boldsymbol{\sigma})) d\boldsymbol{\sigma} - \chi\_d & \text{in } Q, \\ q\_\gamma = \mathbf{0} & \text{on } \Sigma, \\ q\_\gamma(\mathbf{x}, T) = \mathbf{0}, \frac{\partial q\_\gamma}{\partial t}(\mathbf{x}, T) = \mathbf{0} & \text{in } \Omega. \end{cases} \tag{38}$$

Then (35) becomes

$$
\mu\_{\gamma} = \frac{1}{N} \int\_0^1 \frac{\partial q\_{\gamma}}{\partial \eta} d\sigma \,\mathrm{in} L^2(\Sigma\_0). \tag{39}
$$

■ The previous Theorem gives a low-regret control characterization. For the noregret control, we need to prove the convergence of the sequence of averaged lowregret control to the averaged no-regret control. Then, we announce the following Proposition.

For some constant *C* independent of *γ*, we have

$$\left\|\left|u\_{\mathcal{I}}\right|\right\|\_{L^{2}(\Sigma\_{0})} \leq \mathsf{C},\tag{40}$$

$$\left\| \int\_{0}^{1} \mathcal{Y}(u\_{\gamma}, \mathbf{0}, \sigma) d\sigma \right\|\_{L^{2}(Q)} \leq \mathcal{C},\tag{41}$$

$$\left\|\left[\boldsymbol{\mathcal{y}}\left(\boldsymbol{u}\_{\boldsymbol{\gamma}},\mathbf{0},\sigma\right)d\sigma\right\|\_{L^{2}\left(Q\right)}\leq\mathsf{C},\tag{42}$$

$$\left\| \int\_{0}^{1} \frac{\partial \zeta(u\_{\gamma}, \sigma)}{\partial t}(\varkappa, 0) d\sigma \right\|\_{H^{-1}(\Omega)} \leq C\sqrt{\gamma},\tag{43}$$

$$\left\| \int\_{0}^{1} \zeta(u\_{\gamma}, \sigma)(\infty, 0) d\sigma \right\|\_{L^{2}(\Omega)} \leq C\sqrt{\gamma},\tag{44}$$

$$\left\|\rho\_{\mathcal{I}}\right\|\_{L^{\infty}\left(0,T;H^{1}\_{0}(\mathfrak{Q})\right)} \leq \mathsf{C},\tag{45}$$

$$\left\|\boldsymbol{q}\_{\boldsymbol{\mathcal{V}}}\right\|\_{L^{\infty}\left(0,T;H^{1}\_{0}(\Omega)\right)} \leq \mathbf{C}.\tag{46}$$

$$J\_{\mathcal{I}}(u\_{\mathcal{I}}) \le J\_{\mathcal{I}}(\mathbf{0}),\tag{47}$$

$$J(u\_{\mathcal{I}},\mathbf{0}) + \frac{1}{\mathcal{I}} \left\| \int\_{0}^{1} \frac{\partial \zeta(u\_{\mathcal{I}},\sigma)}{\partial \mathbf{t}}(\mathbf{x},\mathbf{0}) d\sigma \right\|\_{L^{2}(\Omega)}^{2} + \frac{1}{\mathcal{I}} \left\| \int\_{0}^{1} \zeta(u\_{\mathcal{I}},\sigma)(\mathbf{x},\mathbf{0}) d\sigma \right\|\_{L^{2}(\Omega)}^{2} \leq J(\mathbf{0},\mathbf{0}), \tag{48}$$

$$E\_{\rho\_{\gamma}}(t) = \frac{1}{2\alpha} \left[ \left| \frac{\partial \rho\_{\gamma}}{\partial t} \right|^2 + \left| \nabla \rho\_{\gamma} \right|^2 + q(\mathbf{x}, \sigma) \left| \rho\_{\gamma} \right|^2 \right] d\mathbf{x} = E\_{\rho\_{\gamma}}(\mathbf{0}) \le \mathbf{C},\tag{49}$$

$$u\_{\gamma} \rightharpoonup u \text{ weakly in} L^2(\Sigma\_0),\tag{50}$$

$$\left\| f(\mathsf{u}\_{\mathsf{T}}, \mathsf{g}) - f(\mathsf{0}, \mathsf{g}) - \gamma \right\| \left\| \mathsf{y}\_{0} \right\|\_{H^{1}\_{0}(\Omega)}^{2} - \gamma \left\| \mathsf{y}\_{1} \right\|\_{L^{2}(\Omega)}^{2} \leq \sup\_{\mathbf{g} \in H^{1}\_{0}(\Omega) \times L^{2}(\Omega)} (f(\mathsf{v}, \mathbf{g}) - f(\mathsf{0}, \mathbf{g})), \quad \text{(51)}$$

$$J(u, \mathbf{g}) - J(\mathbf{0}, \mathbf{g}) \le \sup\_{\mathbf{g} \in H\_0^1(\Omega) \times L^2(\Omega)} (f(v, \mathbf{g}) - f(\mathbf{0}, \mathbf{g})),\tag{52}$$

*Averaged No-Regret Control for an Electromagnetic Wave Equation Depending upon… DOI: http://dx.doi.org/10.5772/intechopen.95447*

Finally, we can present the following theorem giving a full characterization the average no-regret control.

Theorem 1.3 The average no-regret control *u* is characterized by the following optimality system

$$\begin{cases} \frac{\partial^2 \rho}{\partial t^2} - \Delta \mathbf{y} + p(\mathbf{x}, \boldsymbol{\sigma}) \mathbf{y} = \mathbf{0}, \\\\ \frac{\partial^2 \zeta}{\partial t^2} - \Delta \zeta + p(\mathbf{x}, \boldsymbol{\sigma}) \zeta = \int\_0^1 p(u, \mathbf{0}, \boldsymbol{\sigma}) d\boldsymbol{\sigma}, \\\\ \frac{\partial^2 \rho}{\partial t^2} - \Delta \boldsymbol{\rho} + p(\mathbf{x}, \boldsymbol{\sigma}) \boldsymbol{\sigma} = \mathbf{0}, \\\\ \frac{\partial^2 q}{\partial t^2} - \Delta q + p(\mathbf{x}, \boldsymbol{\sigma}) q = \int\_0^1 (\rho + \mathbf{y}(u, \mathbf{0}, \boldsymbol{\sigma})) d\boldsymbol{\sigma} - \boldsymbol{\chi}\_d & \text{in } Q, \\\\ \mathbf{y} = \begin{cases} \boldsymbol{u} & \text{on } \Sigma\_0 \\ \mathbf{0} & \text{on } \Sigma\_1 \end{cases}, \zeta = \mathbf{0} & \text{on } \Sigma, \\\\ \rho = \mathbf{0}; q = \mathbf{0} & \text{on } \Sigma, \\\\ \mathbf{y}(0, \mathbf{x}) = \mathbf{0}, \ \frac{\partial \mathbf{y}}{\partial t}(\mathbf{0}, \mathbf{x}) = \mathbf{0}, \\\\ \zeta(\mathbf{x}, T) = \mathbf{0}, \ \frac{\partial \zeta}{\partial t}(\mathbf{x}, T) = \mathbf{0}, \\\\ \rho(\mathbf{x}, 0) = \boldsymbol{\lambda}\_1(\mathbf{x}), \ \frac{\partial \sigma}{\partial t}(\mathbf{x}, \mathbf{0}) = \boldsymbol{\lambda}\_2(\mathbf{x}), \\\\ q(T, \mathbf{x}) = \mathbf{0}, \ \frac{\partial \mathbf{y}}{\partial t}(T, \mathbf{x}) = \mathbf{0} & \text{in } \Omega, \end{cases}$$

with

$$
\mu = \frac{1}{N} \int\_0^1 \frac{\partial p}{\partial \eta} d\sigma \,\ln L^2(\Sigma\_0),
\tag{54}
$$

and

$$\begin{aligned} \lambda\_1(\mathfrak{x}) &= \lim\_{r \to 0} -\frac{1}{r} \int\_0^1 \frac{\vartheta\_\gamma^\epsilon(\mathfrak{u}\_r)}{\mathfrak{d}\mathfrak{r}}(\mathfrak{x}, \mathbf{0}) d\sigma \text{ weakly in } H\_0^1(\mathfrak{Q}), \\\lambda\_2(\mathfrak{x}) &= \lim\_{r \to 0} \frac{1}{r} \int\_0^1 \mathfrak{f}\_\mathbf{0}^\epsilon(\mathfrak{u}\_r)(\mathfrak{x}, \mathbf{0}) d\sigma \text{ weakly in } L^2(\mathfrak{Q}). \end{aligned}$$

**Proof.** From (42) continuity w.r.t data, we can deduce that

$$\mathcal{Y}(u\_{\mathcal{I}},0,\sigma) \rightharpoonup \mathcal{Y}(u,0,\sigma) \text{ weakly in } L^2(\Omega),\tag{55}$$

solution to

$$\begin{cases} \frac{\partial^2 y}{\partial t^2} - \Delta y + p(\mathbf{x}, \boldsymbol{\sigma}) y = \mathbf{0} & \text{in } Q, \\ y = \begin{cases} u & \text{on } \Sigma\_0, \\ 0 & \text{on } \Sigma \backslash \Sigma\_0, \end{cases} \\ y(\mathbf{x}, \mathbf{0}) = y\_0(\mathbf{x}); \frac{\partial y}{\partial t}(\mathbf{x}, \mathbf{0}) = y\_1(\mathbf{x}) \end{cases} \tag{56}$$

Again, by (41) and dominated convergence theorem

$$\int\_{0}^{1} \mathcal{y}(u\_{\tau}, 0, \sigma) d\sigma \rightharpoonup \int\_{0}^{1} \mathcal{y}(u, 0, \sigma) d\sigma \text{weakly in} L^{2}(\Sigma\_{0}). \tag{57}$$

The rest of equations in (53) leads by a similar way, except the convergences of initial data *<sup>ρ</sup>*ð Þ *<sup>x</sup>*, 0 , *<sup>∂</sup><sup>ρ</sup> <sup>∂</sup><sup>t</sup>* ð Þ *x*, 0 which will be as follows.

From (43) and (44) we deduce the convergences of

$$-\frac{1}{\gamma}\frac{\partial \zeta(u\_{\gamma}, \sigma)}{\partial t}(\varkappa, \mathbf{0}) \rightharpoonup \lambda\_1(\varkappa) \text{ weakly in } H\_0^1(\Omega), \tag{58}$$

and

$$\frac{1}{\mathcal{Y}}\zeta(u\_{\mathcal{Y}},\sigma)(\mathfrak{x},\mathfrak{O}) \rightharpoonup \mathbb{A}\_2(\mathfrak{x}) \text{ weakly in } L^2(\mathfrak{O}).\tag{59}$$

#### **5. Conclusion**

As we have seen, the averaged no-regret control method allows us to find a control that will optimize the situation of the electromagnetic waves with missing initial conditions and depending upon a parameter. The method presented in the paper is quite general and covers a wide class of systems, hence, we could generalize the situation to more control positions (regional, punctual, … ) and different kinds of missing data (source term, boundary conditions, … ).

The results presented above can also be generalized to the case of other systems which has many biomedical applications. This problem is still under consideration and the results will appear in upcoming works.

#### **Acknowledgements**

This work was supported by the Directorate-General for Scientific Research and Technological Development (DGRSDT).

#### **Author details**

Abdelhak Hafdallah\* and Mouna Abdelli Laboratory of Mathematics Informatics and Systems (LAMIS), University of Larbi Tebessi, 12002 Tebessa, Algeria

\*Address all correspondence to: abdelhak.hafdallah@univ-tebessa.dz

© 2021 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

*Averaged No-Regret Control for an Electromagnetic Wave Equation Depending upon… DOI: http://dx.doi.org/10.5772/intechopen.95447*

#### **References**

[1] A. Hafdallah and A. Ayadi. Optimal control of electromagnetic wave displacement with an unknown velocity of propagation, International Journal of Control, DOI: 10.1080/00207179. 2018.14581, 2018.

[2] Savage LJ. The Foundations of Statistics. 2nd ed. New York: Dover; 1972

[3] Lions JL. Contrôle à moindres regrets des systèmes distribués. C. R. Acad. Sci. Paris Ser. I Math. 1992;**315**:1253-1257

[4] Lions JL. No-regret and low-regret control. Economics and Their Mathematical Models, Masson, Paris: Environment; 1994

[5] Hafdallah A. On the optimal control of linear systems depending upon a parameter and with missing data. Nonlinear Studies. 2020;**27**(2):457-469

[6] Jacob B, Omrane A. Optimal control for age-structured population dynamics of incomplete data. J. Math. Anal. Appl. 2010;**370**:42 48

[7] Baleanu D, Joseph C, Mophou G. Low-regret control for a fractional wave equation with incomplete data. Advances in Difference Equations. 2016. DOI: 10.1186/s13662-016-0970-8

[8] Mophou G, Foko Tiomela RG, Seibou A. Optimal control of averaged state of a parabolic equation with missing boundary condition. International Journal of Control. 2018. DOI: 10.1080/00207179.2018.1556810

[9] Mophou G. Optimal for fractional diffusion equations with incomplete data. J. Optim.Theory Appl. 2015. DOI: 10.1007/s10957-015-0817-6

[10] Lions JL. Duality Arguments for Multi Agents Least-Regret Control. Paris: College de France; 1999

[11] Nakoulima O, Omrane A, Velin J. On the Pareto control and no-regret control for distributed systems with incomplete data. SIAM J. CONTROL OPTIM. 2003;**42**(4):1167 1184

[12] Nakoulima O, Omrane A, Velin J. No-regret control for nonlinear distributed systems with incomplete data. Journal de mathématiques pures et appliquées. 2002;**81**(11):1161-1189

[13] Choudhury PK, El-Nasr MA. Electromagnetics for biomedical and medici- nal applications. Journal of Electromagnetic Waves and Applications. 2015;**29**(17):2275-2277. DOI: 10.1080/09205071.2015.1103984

[14] Zuazua E. Averaged control. Automatica. 2014;**50**(12):3077 3087

[15] Kian Y. Stability of the determination of a coefficient for wave equations in an infinite waveguide. Inverse Probl. Imaging. 2014;**8**(3): 713-732

[16] Nakoulima O, Omrane A, Velin J. Perturbations à moindres regrets dans les systèmes dis- tribués à données manquantes. C. R. Acad. Sci. Ser. I Math. (Paris). 2000;**330**:801 806

[17] A. Hafdallah and A. Ayadi. Optimal Control of a thermoelastic body with missing initial conditions; International Journal of Control, DOI: 10.1080/ 00207179.2018.1519258,2018.

[18] Aubin JP. Lanalyse non linéaire et ses motivations économiques. Paris: Masson; 1984

[19] Attouch H, Wets RJB. Isometries for the Legendre-Fenchel transform. Transactions of the American Mathematical Society. 1986;**296**(1): 33-60

Section 2
