Structure of element connectivity for triangular mesh.
A two-Dimensional Finite Element Method of electromagnetic (EM) wave propagation through the soil is presented in this chapter. The chapter employs a boundary value problem (BVP) to solve the Helmholtz time-harmonic electromagnetic model. An infinitely large dielectric object of an arbitrary cross-section is considered for scattering from a dielectric medium and illuminated by an incident wave. Since the domain extends to infinity, an artificial boundary, a perfectly matched layer (PML) is used to truncate the computational domain. The incident field, the scattered field, and the total field in terms of the z-component are expressed for the transverse magnetic (TM) and transverse electric (TE) modes. The radar cross-section (RCS), as a function of several other parameters, such as operating frequency, polarization, illumination angle, observation angle, geometry, and material properties of the medium, is computed to describe how a scatterer reflects an electromagnetic wave in a given direction. Simulation results obtained from MATLAB for the scattered field, the total field, and the radar cross-section are presented for three soil types – sand, loam, and clay.
- electromagnetic signal in soil
- finite element method
- radar cross-section
- underground wireless communication
Signal propagation underground and in soil medium constitute the backbone of the Internet of underground things (IoUT), which power many applications such as precision agriculture, border monitoring for intrusion detection, pipeline monitoring, etc. [1, 2, 3]. In the recent past, the study of signal propagation in soil medium for underground wireless communication has focused mainly on empirical techniques [4, 5, 6, 7, 8]. The most commonly used modeling techniques for implementing IoUT include electromagnetic waves, magnetic induction, and acoustic waves [9, 10]. Electromagnetic field analysis is performed in this chapter using numerical modeling with the finite element method to examine the signal strength in the soil medium.
The process of numerical modeling of how electromagnetic fields propagate and interact with physical objects and the environment is usually referred to as computational electromagnetics (CEM), numerical electromagnetics. The primary motivation of this process is to develop efficient approximations to Maxwell’s equations through numerical schemes for cases where closed-form analytical solutions of Maxwell’s equation cannot be obtained due to the complexity of geometries, material parameters, and boundary conditions. Therefore, several real-life problems that are not analytically computable, such as electromagnetic scattering, antenna radiation, electromagnetic wave propagation, electromagnetic compatibility, etc., can effectively be solved by numerical techniques. The mathematical model of the electromagnetic problem is usually obtained in terms of partial differential equations, integral equations, or integro-differential equations derived from Maxwell’s equations and a set of a priori constraints of the problem such as boundary and initial conditions material parameters and geometry. The problem is ideally defined on an infinite-dimensional function space. Numerical methods apply a discretization to the continuum to reduce infinite degrees of freedom to a finite degree of freedom. In other words, the solution of an infinite dimension-dimensional problem is projected into a finite-dimensional space. Hence, the solution to the problem becomes amenable on a digital computer. The main philosophy in most of the numerical methods is to apply the divide-and-conquer strategy. The idea is to divide an intractable continuous problem into smaller pieces (divide), express the solution over each small piece (conquer), and then combine the piecewise solutions to obtain a global solution. In this chapter, FEM will be applied to the two-dimensional boundary value problem in EM wave propagation through the soil to evaluate the signal strength of the wave propagation in soil. This evaluation will be based on the incidence angle of the transmitted wave. The radar cross-section of the scatterer will be used to evaluate the direction of the wave.
2. Time-harmonic EM model
When deriving the wave equation for electromagnetic wave propagation through a vacuum or a dielectric medium such as soil or pure water, the free charge density and the free current density . However, in the case of conductors like seawater or metals, we cannot control the flow of charges and, in general, is certainly not equal to zero. With this, Maxwell’s equations for linear media assume the form .
Apply the curl to (iii) and (iv), and we obtain modified wave equations for the electric field
where is the angular frequency,
where is the unknown function to be determined, and are given functions. , and are given functions in boundary conditions (BCs); and refers to boundaries where Dirichlet, Neumann, and mixed BCs are imposed, respectively; is the unit vector normal to the boundary in the outward direction.
where is the weight function.
The weak form is applied in each element domain, and element matrices are formed by expressing the unknown function as a weighted sum of nodal shape functions. The sum of line integrals of two neighboring elements cancels out while combining the element matrices. Therefore, the line integral can be omitted for interior elements and should be considered only for elements adjacent to the boundary. Due to its special form, the line integral makes easier the imposition of mixed types of BCs. The mesh generation will be discussed first, and the shape functions will be given, and then the finite element solution of the time-harmonic problem of electromagnetic scattering in a dielectric medium (soil) will be presented.
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 and the x-y-components magnetic field exist. However, in the TE case (vertical polarization), which is a dual of TM, only the z-component of the magnetic field and the x-y-components electric field exist.
The material tensors in this case are defined as
The generalized homogeneous Helmholtz equation in TM and TE case can be written, respectively, in the following forms .
Where the superscript T indicates the transpose and refers to the determinant of the corresponding submatrix; is the free-space wavenumber; is the complex relative permittivity , and are the relative permittivity and permeability, respectively; and is the conductivity. For isotropic mediums, Eq. (8) becomes
The same way Eq. (11) becomes
4. Scattering from a dielectric medium (soil)
In , 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 work, sensitivity analysis of the Ku-band scattering coefficient to soil moisture was performed under single-polarized, dual-polarized, and dual-angular combinations . 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 .
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) .
Locally-Conformal PML (LC-PML) is a powerful PML method whose implementation is straightforward. It is implemented by replacing the real coordinates inside the PML region with their corresponding complex coordinates . 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 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 .
An element connectivity matrix of size Mx3, where M is the number of elements vectors of node coordinates (
|Element (||Node 1||Node 2||Node 3|
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) = , , and the like. The higher the denominator, the finer the element. Figure 4 shows the mesh generated in MATLAB for els = .
4.2 The scattered field formulation
Fields in the presence of the scatterer can be decomposed into two parts:
The first one is the incident field which is produced without the scatterer and
The second one is the scattered field produced by an equal amount of current induced on the scatterer or the surface enclosing the scatterer.
We first assume the TM polarization case, where the incident field , and the scattered field , are expressed in terms of the z-component, the z-component of the total field can be expressed as . For an isotropic case, the total field satisfies the differential Eq. .
Where the source term is given by
This differential equation is a special form of Eq. (4a), where , and .
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
where is the angle of incidence for the
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
Where the source term is given by
This differential equation is a special form of Eq. (4a), where , and .
For dielectric objects, the right-hand side of the matrix equation can be obtained by using the source term
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 . 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
Where is the scattered electric or magnetic field at far-zone observed in a given direction, when satisfies the inequality , where
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 is dBw, and for is dBm.
5. Scattering inside the dielectric medium
For each such element, the far-zone field is computed and then superposed over all elements, as follows:
Since FEM is formulated in terms of the scattered electric field, computation of the derivatives of the scattered field requires additional effort. This can simply be achieved by using the weighted sum of derivatives of shape functions in terms of nodal scattered fields.
5.1 Scattering inside the dielectric medium, TM case
For the dielectric object, the TM case, the interior part of the object should be included. Since the magnetic current density becomes nonzero, the integral containing the magnetic current density should be evaluated. The magnetic current density has
Finally, the scattered electric field can be obtained as follows:
Here, can be determined as the average of nodal field values connected to the boundary.
5.2 Scattering inside the dielectric medium, TE case
For dielectric object, TE case: The magnetic current density is nonzero and has only
6. Simulation results and discussion
The operating frequency used in the simulation is 300 MHz, and the free-space wavelength is . The mesh size is approximatively . The dielectric medium is a rectangular computational domain, the relative permittivity for sandy soil , the relative permittivity for loamy soil , and the relative permittivity for clay soil . The angle of incidence is , that is, the incident plane wave travels in the direction. The wavelength inside the dielectric medium (soil) decreases in terms of the wavelength in free-space and relative permittivity, i.e. , and hence, mesh size resolution corresponds to resolution inside the dielectric. Therefore, smaller elements are used in dielectric media to preserve the level of accuracy. All simulation results are obtained using MATLAB.
Figure 5 shows scattering from a dielectric medium for TE mode, with sandy soil characteristics. In (a) and (b), the scattered and total fields are shown respectively. We observe a minimal electric field scattering and an intense total electric field inside the dielectric medium of sandy soil. This is because this medium is porous and allows for better signal propagation.
Figure 6 shows scattering from a dielectric medium for TM mode, with sandy soil characteristics. In this case, we observe a minimal magnetic field scattering in (a) and an intense total magnetic field in (b) inside the dielectric medium of sandy soil.
Figure 7 shows scattering from a dielectric medium for TE mode, with loamy soil characteristics. In (a) and (b), the scattered and total fields are shown respectively. We observe an intense electric field scattering and a minimal total electric field inside the dielectric medium of loamy soil. This is because this medium is less porous and presents some challenges in signal propagation.
Similarly, Figure 8 shows scattering from a dielectric medium for TM mode, with loamy soil characteristics. In this case, we observe a high magnetic field scattering in (a) and a low total magnetic field in (b) inside the dielectric medium of loamy soil.
Figure 9 shows scattering from a dielectric medium for TE mode, with clay soil characteristics. In (a) and (b), the scattered and total fields are shown respectively. We observe very high electric field scattering and a low total electric field inside the dielectric medium of clay soil. This is because this medium is non-porous and presents a very poor signal propagation.
Similarly, Figure 10 shows scattering from a dielectric medium for TM mode, with clay soil characteristics. In this case, we observe a very high magnetic field scattering in (a) and a low total magnetic field in (b) inside the dielectric medium of clay soil.
Figures 11–13 show the bistatic RCS profiles to describe how scatterers reflect the incident electromagnetic wave in a given direction. This is the area intercepting that amount of power which, when scattered in the soil medium, produces at the receiver a density which is equal to that scattered by the actual target. The RCS is a function of several parameters, such as operation frequency, polarization, illumination angle, observation angle, geometry, and properties of the soil medium. It is shown in the TE modes (a) and TM modes (b) for sandy soil, loamy soil, and clay soil, respectively.
In a two-Dimensional Finite Element Analysis of EM wave Propagation through the soil, a boundary value problem (BVP) used to solve the time-harmonic electromagnetic problem in 2-D, has been expressed in its generic form. In TM and TE cases, the Helmholtz model has considered an infinitely large dielectric object of an arbitrary cross-section for scattering from a dielectric medium and illuminated by an incident wave. Since the domain extends to infinity, an artificial boundary, an absorbing boundary condition (ABC), or a perfectly matched layer (PML), has been used to truncate the computational domain. The incident field, the scattered field, and the total field in terms of the z-component are expressed for the TM and TE modes. The radar cross-section (RCS), a function of several parameters, such as operation frequency, polarization, illumination angle, observation angle, geometry, and material properties of the medium, has been computed to describe how a scatterer reflects an incident electromagnetic wave in a given direction. Simulation results for the scattered field, the total field, have been presented for soil types, and the radar cross-section for different element refinements have also been presented.
Salam A, Subsurface MIMO. A Beamforming Design in the Internet of underground things for digital agriculture applications. MDPI Journal of Sensors and Actuators. 2019. DOI: 10.3390/jsan8030041
Vuran MC, Salam A, Wong R, Irmak S. Internet of underground things in precision agriculture: Architecture and technology aspects. Ad Hoc Netw. 2018
Akyildiz IF, Sun Z, Vuran MC. Signal propagation techniques for wireless underground communication networks. Phys. Commun. J. 2009
Jabbar S, Asif Habib M, Minhas AA, Ahmad M, Ashraf R, Khalid S, et al. Analysis of factors affecting energy-aware routing in wireless sensor network. Wirel. Commun. Mob. Comput. 2018
Konda A, Rau A, Stoller MA, Taylor JM, Salam A, Pribil GA, et al. Soft microreactors for the deposition of conductive metallic traces on planar, embossed, and curved surfaces. Adv. Funct. Mater. 2018
Vuran MC, Akyildiz IF. Channel model and analysis for wireless underground sensor networks in soil medium. Phys. Commun. 2010
Sun Z, Akyildiz I. Magnetic induction Communications for Wireless Underground Sensor Networks. IEEE Trans. Antennas Propag. 2010
Engmann F, Adu-Manu K S, Abdulai J-D, Katsriku F. Application of Prediction Approaches in Wireless Sensor Networks. In: Wireless Sensor Networks – Design, Deployment, and Application. Intech Open. 202. DOI: 10.5772/intechopen.94500
Banaseka FK, Katsriku FA, Abdulai J-D, Adu-Manu KS, Engmann F. Signal propagation models in soil medium for the study of wireless underground sensor networks: A review of current trends. In: Hindawi Journal of Wireless Communication and Mobile Computing. 2021
Saeed N, Alouini M-S, Al-Naffouri TY. Toward the internet of underground things: A systematic survey. IEEE Communication Survey and Tutorials. 2019; 21(4):3443-3466
David JG. Introduction of Electrodynamics. 2nd Ed. ISBN No 0–13–481367-7. Prentice – Hall: Inc; 1989
Arshad K, Katsriku F, Lasebae A. Radiowave VHF propagation modeling in forest using finite elements. In proceedings – 2006 international conference on information and communication technologies, pp. 2146–2149, Damascus. Syria. 2005
Arshad K, Katsriku F. Lasebae a. Modeling obstructions in straight and curved rectangular tunnels by finite element approach. Journal of Electrical Engineering–Bratislava. 2007; 59(1):09-13
Arshad K, Katsriku F, Lasebae A. Effects of different parameters on attenuation rates in circular and arch tunnels. Piers Online. 2006; 3(5):607-611
Ozgun O, Kuzuoglu M. Matlab-based finite element programming in electromagnetic Modeling. CRC Press is an imprint of Taylor and Francis Group Informa. 2019
Banaseka FK, Franklin H, Katsriku F, Abdulai J-D, Ekpezu A, Wiafe I. Soil medium electromagnetic scattering model for the study of wireless underground sensor networks. Hindawi Journal of Wireless Communication and Mobile Computing. 2021
Li T, Hajnsek I, Chen K-S. Sensitivity analysis of Bistatic scattering for soil moisture retrieval. MDPI Remote Sens. 2021; 13:188. DOI: 10.3390/rs13020188
Li A-Q, Yin C-Y, Wu W. Accurate Prediction of Radio Wave Propagation in an Environment of Dielectric Ground and Obstacles Based on the Principle of Domain Decomposition. Antennas and Propagation: IET Microwave; 2021. DOI: 10.1049/mia2.12130
Delaunay B. Sur la sphere vide. Bulletin de l’Academie des Sciences de l’URSS. Classe des sciences mathematiques et naturelles. Vol. 1934; 6:789-800
Shewchuk JR. Lecture notes on Delaunay mesh generation. University California at Berkeley. 2012