**3. The Helmholtz time harmonic model**

In 2-D, the geometry and boundary conditions do not vary along an axis (say the z-axis). Hence fields can be represented as a superposition of fields of two orthogonal polarizations using the linearity property. Any field can be decomposed into transverse magnetic (TM) and transverse electric (TE) parts for the z-variable. In the TM case (horizontal polarization), only the z-component of the electric field ð Þ **E** ¼ *a*^*zEz*ð Þ *x*, *y* and the x-y-components magnetic field ð**H** ¼ *a*^*xHx*ð Þþ *x*, *y a*^*yHy*ð Þ *x*, *y*Þ exist. However, in the TE case (vertical polarization), which is a dual of TM, only the z-component of the magnetic field ð Þ **H** ¼ *a*^*zHz*ð Þ *x*, *y* and the x-ycomponents electric field **<sup>E</sup>** <sup>¼</sup> *<sup>a</sup>*^*xEx*ð Þþ *<sup>x</sup>*, *<sup>y</sup> <sup>a</sup>*^*yEy*ð Þ *<sup>x</sup>*, *<sup>y</sup>*<sup>Þ</sup> � exist.

The material tensors in this case are defined as

$$\mathbf{e}\_{rc} = \begin{bmatrix} \boldsymbol{\varepsilon}\_{rc}^{\infty} & \boldsymbol{\varepsilon}\_{rc}^{\text{xy}} & \mathbf{0} \\ \boldsymbol{\varepsilon}\_{rc}^{\text{yx}} & \boldsymbol{\varepsilon}\_{rc}^{\text{y}} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \boldsymbol{\varepsilon}\_{rc}^{\text{xz}} \end{bmatrix} \quad \text{with} \quad (\mathbf{e}\_{rc})\_{sub} = \begin{pmatrix} \boldsymbol{\varepsilon}\_{rc}^{\text{xx}} & \boldsymbol{\varepsilon}\_{rc}^{\text{xy}} \\ \boldsymbol{\varepsilon}\_{rc}^{\text{yx}} & \boldsymbol{\varepsilon}\_{rc}^{\text{yy}} \end{pmatrix},\tag{9}$$

$$
\boldsymbol{\mu}\_r = \begin{bmatrix}
\boldsymbol{\mu}\_r^{\rm xx} & \boldsymbol{\mu}\_r^{\rm xy} & \mathbf{0} \\
\boldsymbol{\mu}\_r^{\rm xx} & \boldsymbol{\mu}\_r^{\rm yy} & \mathbf{0} \\
\mathbf{0} & \mathbf{0} & \boldsymbol{\mu}\_r^{\rm xx}
\end{bmatrix} \quad \text{with} \quad (\boldsymbol{\mu}\_r)\_{\rm sub} = \begin{pmatrix}
\boldsymbol{\mu}\_r^{\rm xx} & \boldsymbol{\mu}\_r^{\rm xy} \\
\boldsymbol{\mu}\_r^{\rm xx} & \boldsymbol{\mu}\_r^{\rm y}
\end{pmatrix}, \tag{10}
$$

The generalized homogeneous Helmholtz equation in TM and TE case can be written, respectively, in the following forms [15].

$$\nabla \cdot \left(\Lambda\_{\mu} \cdot \nabla E\_{\underline{x}}\right) + k\_0^2 \varepsilon\_{rc}^{xx} E\_{\underline{x}} = \mathbf{0},\tag{11}$$

$$\nabla \cdot (\mathbf{A}\_{\varepsilon} \cdot \nabla H\_x) + k\_0^2 \mu\_r^{xx} H\_x = \mathbf{0},\tag{12}$$

Where

$$\mathbf{A}\_{\mu} = \frac{1}{\left| (\mu\_r)\_{sub}^T \right|} (\mu\_r)\_{sub}^T,\tag{13}$$

$$\mathbf{A}\_{\varepsilon} = \frac{1}{\left| \left( \mathbf{e}\_{rc} \right)\_{sub}^{T} \right|} \left( \mathbf{e}\_{rc} \right)\_{sub}^{T}, \tag{14}$$

Where the superscript T indicates the transpose and j j*:* refers to the determinant of the corresponding submatrix; *k0* is the free-space wavenumber; *εrc* ¼ *ε<sup>r</sup>* � *jσ=ωε<sup>0</sup>* is the complex relative permittivity *εr*, and *μ<sup>r</sup>* are the relative permittivity and permeability, respectively; and *σ* is the conductivity. For isotropic mediums, Eq. (8) becomes

$$\nabla \cdot \left( \mu\_r^{-1} \nabla E\_x \right) + k\_0^2 \varepsilon\_n E\_x = \mathbf{0},\tag{15}$$

Or

$$
\frac{
\partial
}{
\partial\mathbf{x}
}
\left(\frac{\mathbf{1}}{\mu\_r}\frac{
\partial E\_x
}{
\partial\mathbf{x}
}\right) + \frac{
\partial
}{
\partial\mathbf{x}
}
\left(\frac{\mathbf{1}}{\mu\_r}\frac{
\partial E\_x
}{
\partial\mathbf{y}
}\right) + k\_0^2 \varepsilon\_n E\_x = \mathbf{0},
\tag{16}
$$

The same way Eq. (11) becomes

$$\nabla \cdot \left( \varepsilon\_{rc}^{-1} \nabla H\_x \right) + k\_0^2 \mu\_r H\_x = \mathbf{0},\tag{17}$$

Or

$$
\frac{
\partial
}{
\partial\mathbf{x}
} \left(\frac{1}{\varepsilon\_{\mathbf{r}}} \frac{
\partial H\_{\mathbf{z}}
}{
\partial\mathbf{x}
}\right) + \frac{
\partial
}{
\partial\mathbf{x}
} \left(\frac{1}{\varepsilon\_{\mathbf{r}}} \frac{
\partial H\_{\mathbf{z}}
}{
\partial\mathbf{y}
}\right) + k\_0^2 \mu\_r H\_{\mathbf{z}} = \mathbf{0},
\tag{18}
$$

#### **4. Scattering from a dielectric medium (soil)**

In [16], a new wave number model is proposed with the combination of the Peplinski principle and multiple scattering from particles in the soil medium. The new wave number is used in the computation of the path loss. In another recent

*Signal Propagation in Soil Medium: A Two Dimensional Finite Element Procedure DOI: http://dx.doi.org/10.5772/intechopen.99333*

work, sensitivity analysis of the Ku-band scattering coefficient to soil moisture was performed under single-polarized, dual-polarized, and dual-angular combinations [17]. Similarly, a model of parabolic equations for reflection and refraction in an environment with an obstacle where the area is decomposed into two different domains. The discrete mixed Fourier transformation is used to compute the field strength in the upper subdomain, and the finite difference method is used to calculate the field strength in the lower subdomain [18].

In our case, an infinitely large dielectric object of an arbitrary cross-section is considered and illuminated by an incident wave that is not a function of z. An illustration of the scattering problem for TM mode is shown in **Figure 1** where a general FEM is depicted. The problem is defined in TE mode by replacing the electric field with a magnetic field.

#### **4.1 The perfectly matched layer (PML)**

Since the domain extends to infinity, an artificial boundary or layer is used to truncate the computational domain. That is an absorbing boundary condition (ABC) or perfectly matched layer (PML) [15].

Locally-Conformal PML (LC-PML) is a powerful PML method whose implementation is straightforward. It is implemented by replacing the real coordinates ð Þ **ρ** ¼ ð Þ *x*, *y* inside the PML region with their corresponding complex coordinates ð Þ **ρ**~ ¼ ð Þ *x*~, ~*y* [15]. The computational domain is the union of the PML region, the free-space region, and the scatterer region.

Mesh generation is the process of representing the domain of interest as a collection of elements. The two most commonly used elements in 2-D problems are triangular and quadrilateral. **Figure 2** shows the mesh that is formed by triangular elements for a rectangular domain of the boundary value problem (BVP). The figure was generated in MATLAB.

For a computational domain with such curved boundaries, discretization errors will occur due to the inability to capture the exact geometry of the domain. Triangular elements are preferred due to their simplicity and the possibility of developing algorithms for automatic triangulation for a computational domain such as Delaunay triangulation [19, 20]. Although fewer elements are needed when quadrilateral elements are used, Triangular elements are well-suited for complex geometries and cause fewer numerical dispersion errors. Furthermore, the calculation of element matrices is easier in triangular elements. Discretization error might inevitably occur unless sufficient mesh density is used. This happens regardless of whether triangular or quadrilateral elements are used for meshing. One way to overcome this problem is to refine the mesh. Mesh refinement might be needed

**Figure 1.** *Electromagnetic scattering in soil: FEM modeling with PML.*

especially if the geometry has a curved boundary or some corners, sharp edges, small features, or discontinuities. This might be important especially if linear interpolation functions are used. Another approach that increases accuracy is to use high-order elements at the expense of increased computational load. To achieve this, we use extra nodes within an element and use high-order interpolation functions. **Figure 3** shows a simple mesh with six linear triangular elements having eight nodes. Local node numbers must follow the anticlockwise orientation in all elements to guarantee that the area of each element is obtained as a positive quantity. During the mesh generation, certain data arrays must be created [19].

An element connectivity matrix of size Mx3, where M is the number of elements vectors of node coordinates (*x* and *y*), each of size Nx1, where N is the number of nodes representation of the given functions (*px*, *py*, *q*, and *f* ) in each element (each of which is of size Mx1 Arrays containing special nodes and/or element. The connectivity matrix belonging to **Figure 3** is given in **Table 1**.

**Figure 2.** *Mesh formed by triangular elements for the rectangular domain of the BVP.*

**Figure 3.** *Mesh of a 2-D domain using linear triangular elements.*


**Table 1.** *Structure of element connectivity for triangular mesh.* *Signal Propagation in Soil Medium: A Two Dimensional Finite Element Procedure DOI: http://dx.doi.org/10.5772/intechopen.99333*

**Figure 4.** *Triangular mesh generated with Δh* ¼ *λ0=15.*

The Delaunay function is used to create the connectivity matrix and automatically enumerates the nodes. The triangular mesh is generated for a rectangular domain. The element size (els) = *Δh* ¼ *λ0=5*, *Δh* ¼ *λ0=10*, *Δh* ¼ *λ0=20* and the like. The higher the denominator, the finer the element. **Figure 4** shows the mesh generated in MATLAB for els = *Δh* ¼ *λ0=15*.

### **4.2 The scattered field formulation**

Fields in the presence of the scatterer can be decomposed into two parts:


We first assume the TM polarization case, where the incident field **<sup>E</sup>**inc <sup>¼</sup> *<sup>a</sup>*^*zE*inc *z* , and the scattered field **<sup>E</sup>**scat <sup>¼</sup> *<sup>a</sup>*^*zE*scat *z* , are expressed in terms of the z-component, the z-component of the total field can be expressed as *Ez* ¼ *E*scat *<sup>z</sup>* <sup>þ</sup> *<sup>E</sup>*inc *<sup>z</sup>* . For an isotropic case, the total field satisfies the differential Eq. [9].

$$-\frac{\partial}{\partial \mathbf{x}} \left( \frac{\mathbf{1}}{\mu\_r} \frac{\partial E\_x^{\text{scat}}}{\partial \mathbf{x}} \right) - \frac{\partial}{\partial \mathbf{x}} \left( \frac{\mathbf{1}}{\mu\_r} \frac{\partial E\_x^{\text{scat}}}{\partial \mathbf{y}} \right) - k\_0^2 \varepsilon\_r E\_x^{\text{scat}} = f(\mathbf{x}, \mathbf{y}), \tag{19}$$

Where the source term is given by

$$f(\mathbf{x}, \mathbf{y}) = \frac{\partial}{\partial \mathbf{x}} \left( \frac{\mathbf{1}}{\mu\_r} \frac{\partial E\_x^{\text{inc}}}{\partial \mathbf{x}} \right) + \frac{\partial}{\partial \mathbf{x}} \left( \frac{\mathbf{1}}{\mu\_r} \frac{\partial E\_x^{\text{inc}}}{\partial \mathbf{y}} \right) + k\_0^2 \varepsilon\_{\text{re}} E\_x^{\text{inc}},\tag{20}$$

This differential equation is a special form of Eq. (4a), where *<sup>u</sup>* <sup>¼</sup> *<sup>E</sup>*scat *<sup>z</sup>* ð Þ *x*, *y* , *px* <sup>¼</sup> *py* <sup>¼</sup> *<sup>p</sup>* <sup>¼</sup> <sup>1</sup> *<sup>μ</sup>* and *<sup>q</sup>* ¼ �*k<sup>2</sup> <sup>0</sup>εrc*.

The incident field can be arbitrary and usually chosen as a uniform plane wave since the incident field sources are sufficiently far away from the object. The incident field with the unit magnitude is given by

$$\mathbf{E}^{\text{inc}} = \hat{a}\_{\mathbf{z}} E\_{\mathbf{z}}^{\text{inc}} = \hat{a}\_{\mathbf{z}} \exp\left[jk\left(\mathbf{x}\cos\boldsymbol{\rho}^{\text{inc}} + \boldsymbol{\mathcal{y}}\sin\boldsymbol{\rho}^{\text{inc}}\right)\right],\tag{21}$$

where *φ*inc is the angle of incidence for the *x*-axis in cartesian coordinates.

Similar calculations can be performed in TE polarization mode by replacing the electric field with the magnetic field and permittivity with permeability. Hence, the differential equation in terms of the scattered magnetic field is given by

$$-\frac{\partial}{\partial \mathbf{x}} \left( \frac{\mathbf{1}}{\varepsilon\_{\text{re}}} \frac{\partial H\_{\text{x}}^{\text{scat}}}{\partial \mathbf{x}} \right) - \frac{\partial}{\partial \mathbf{x}} \left( \frac{\mathbf{1}}{\varepsilon\_{\text{re}}} \frac{\partial H\_{\text{x}}^{\text{scat}}}{\partial \mathbf{y}} \right) - k\_{0}^{2} \mu\_{r} H\_{\text{x}}^{\text{scat}} = f(\mathbf{x}, \mathbf{y}), \tag{22}$$

Where the source term is given by

$$f(\mathbf{x}, \boldsymbol{y}) = \frac{\partial}{\partial \mathbf{x}} \left( \frac{\mathbf{1}}{\varepsilon\_n} \frac{\partial H\_x^{\text{inc}}}{\partial \mathbf{x}} \right) + \frac{\partial}{\partial \mathbf{x}} \left( \frac{\mathbf{1}}{\varepsilon\_n} \frac{\partial H\_x^{\text{inc}}}{\partial \mathbf{y}} \right) + k\_0^2 \mu\_r H\_x^{\text{inc}},\tag{23}$$

This differential equation is a special form of Eq. (4a), where *<sup>u</sup>* <sup>¼</sup> *<sup>H</sup>*scat *<sup>z</sup>* ð Þ *x*, *y* , *px* <sup>¼</sup> *py* <sup>¼</sup> *<sup>p</sup>* <sup>¼</sup> <sup>1</sup> *<sup>ε</sup>rc* and *<sup>q</sup>* ¼ �*k<sup>2</sup> 0μ*.

For dielectric objects, the right-hand side of the matrix equation can be obtained by using the source term *f*(*x*, *y*) in each element. An alternative simpler approach is that the source term in Eq. (20) is just the differential operator applied to the known incident field within the object. The left-hand side of Eq. (19) is the same differential operator being applied to the unknown scattered field, which yields the left-hand side of the matrix equation ð Þ **Au** , where **u** refers to the vector of nodal values of the scattered field. Therefore, **Au** ¼ �**Au**inc must be satisfied within the scatterer region. After forming the global matrix, as usual, the right-hand side vector can be modified by just multiplying the incident field vector by the global matrix, i.e. **<sup>b</sup>** ¼ �**Au**inc � �, only for entries **<sup>b</sup>** corresponding to the nodes lying inside the object.

#### **4.3 Radar cross section**

The radar cross-section (RCS) of the scatterer is perhaps the most critical parameter that must be evaluated in the post-processing phase of FEM for the electromagnetic scattering problem [9]. RCS is the reflection of the scattering electromagnetic ware at an incident angle in a particular direction. In other words, it is the area capturing that amount of scattered power produced at the receiver in an isotropic medium. This is a density that is equal to that scattered by the actual target. It is a function of several parameters, such as operation frequency, polarization, illumination angle, observation angle, geometry, and material properties of the object. In 2-D, it is mathematically defined as

$$\sigma\_{2D} = \lim\_{r \to \infty} \, \_{\text{os}} 2\pi\rho \, \frac{\left| \mathbf{u}\_{\text{far}}^{\text{scat}} \right|^2}{\left| \mathbf{u}^{\text{inc}} \right|^2} \tag{24}$$

Where uscat far <sup>¼</sup> <sup>1</sup>ffiffi *<sup>ρ</sup>* <sup>p</sup> *f*ð Þ *φ* is the scattered electric or magnetic field at far-zone observed in a given direction, when *ρ* satisfies the inequality *ρ* ≫ *2D2=λ*, where *D* is the largest dimension of the scatterer and *λ* is the wavelength.

If the incident and observation directions are the same, the RCS is called monostatic or backscatter RCS; otherwise, it is referred to as bistatic RCS. In 2-D, RCS is usually normalized for *λ* (wavelength) or m (meter). The unit for *σ2D=λ* is dBw, and for *σ2D=*m is dBm.

*Signal Propagation in Soil Medium: A Two Dimensional Finite Element Procedure DOI: http://dx.doi.org/10.5772/intechopen.99333*
