Computational electromagnetic analysis of partially-ﬁlled rectangular waveguide at X-band frequencies

. In this paper, the full-wave analysis of the direct scattering problem in the rectangular waveguide is derived. The numerical calculation results of the scattering characteristics are presented. To show the advantage of our proposed model, we compared direct problem results with three-dimensional simulation software Ansys HFSS calculations. An excellent agreement is observed when compared these two approaches.


Introduction
The response of materials to electromagnetic (EM) fields is a key physical parameter to understanding the material's behavior and implementing the dielectric materials in many applications. The two parameters, which determine electromagnetic field propagation in any material or any space, are the electrical permittivity and magnetic permeability of the material. Typical measurements of these material parameters involve measurements of the complex relative permittivity (ε * = ε − ıε ) and/or complex relative permeability (µ * = µ − ıµ ), which consist of real and imaginary parts [1]. There is a number of dielectric characterization techniques developed over half a century [24]. Nevertheless, the demand for accurate dielectric material measurements at radio frequencies (RF) and microwave frequencies have increased over the years. This is mainly due to the increased development of communication technology, microwave circuit design, and materials sciences. Additionally, scientists are interested not only in technological implementation of telecommunications, but the dielectric characterization of biological materials also could provide valuable information about the impact of electromagnetic fields on biological tissues [5,36]. Each technique has its limitations associated with the frequency bands, the materials to be used, and its geometry.
Many methods exist for the assessments of complex dielectric permittivity and complex magnetic permeability of the materials. The majority of these methods are only derivative variations of several approaches to the scientific-engineering problem. Based on the implemented measurement structure, these methods can be classified into freespace methods [6,20], transmission/reflection (TR) line methods [11,15,19,21,30,32], resonant cavities techniques [12,13,18], and open-ended coaxial methods [7,10,31]. The free-space method is commonly used to characterize materials in order to estimate their influence on telecommunications despite the undesirable reflection and diffraction causing lower accuracy of the measurements. The accuracy of dielectric measurements is much higher for cavity resonator techniques, but they are applicable only over a narrow frequency band. Besides, it requires a new design of the resonator and a new size of sample for each frequency range. The open-ended coaxial methods require no matching of the sample thus lightening the preparation of the samples. Unfortunately, it is critical to have the samples polished smoothly in order to avoid air gaps for solid specimens that affect the accuracy of the methods. The requirements of the surface of the sample are the main disadvantages of the open-ended coaxial methods.
The transmission/reflection line is a popular broadband measurement method of a specialized structure designed to conduct electromagnetic waves in a contained manner. Transmission line techniques, including reflection techniques, utilize transmission line concepts, where a piece of dielectric material is placed inside a transmission line, and an electromagnetic wave is directed at the sample. There are different types of transmission lines, such as two-pair wires lines, coaxial lines, strip lines, microstrip lines, waveguides, etc. The most prevalent types of transmission lines are waveguides and the coaxial transmission line. In this method, only the fundamental modes (TEM mode in coaxial line and TE mode in waveguides) are assumed to propagate.
Many electromagnetic (EM) characterization techniques based on rectangular waveguides have been developed. Typically, the Nicholson-Ross-Weir (NRW) technique [21,34] has been used in which a sample is loaded to fill the cross-section of the rectangular waveguide, and the transmission and reflection coefficients are measured. Nicolson-Ross-Weir is a methodology used to obtain the permittivity and permeability from the S-parameters. This method combines the values of S 11 and S 21 to develop an equation system that allows the estimation of the reflection and transmission coefficients. The advantage of the NRW technique is the ability to directly solve for the complex ε and µ in closed form. The limitation of this method is the special requirements associated with the test sample: the test samples must be linear, homogenous, and isotropic. Additionally, the rectangular waveguide method requires sample preparation so that the samples fit tightly into the structure, typically cross-section of the waveguide. Nevertheless, the sample requirement to completely fill the guide's cross-section is extremely difficult to implement for certain cases, and, in practice, air gaps between the sample and conducting boundaries may exist [8], leading to the excitation of higher-order modes. Several techniques have also been developed to reduce the size requirement of the test sample. Alternative techniques to these approaches reduce the width of the sample so that the guide is only partially filled, avoiding, in most cases, the air gaps effect [4,23]. A method for measuring samples using a partially-filled waveguide has been presented in [3]. A mode-matching (MM) technique has been developed to accommodate the high-order modes excited at the discontinuity of the waveguide, and Newton's method is employed to extract the electromagnetic properties of samples with low-loss and high-loss parameters. The basis for the mode-matching technique has been founded by Wexler [35]. He stated in full-wave analysis the amplitudes of normal modes were chosen to satisfy boundary conditions thus solving waveguides discontinuities problems. This method is also useful for multimode propagations. Since it conforms closely to physical reality, it has the widest applications. This is the main advantage of the MM technique over other methods. There are also numerous other techniques for solving waveguide discontinuity problems. Direct and inverse scattering problems for electromagnetic waves occur in partially-filled waveguides. The direct problem R, T (ε * , µ * ) is stated as that of finding an appropriate functional form for the reflection factor R and/or transmission factor T (ε * is the complex relative permittivity, µ * is the complex permeability). To solve such a problem for a rectangular waveguide, Maxwell's equations are usually solved in the rectangular coordinate system. There are several methods proposed to solve this problem. Propagation of modes in some rectangular waveguides using the finite-difference time-domain (FDTD) method has been proposed in [9]. An advantageous finite element method (FEM) for the rectangular waveguide problem has been developed [33]. Various commercial software allowing threedimensional simulation of complex geometries in the waveguide or other transmission line (CST, HFSS, FEKO, COMSOL, etc.). However, analytical solutions of the direct problem R, T (ε * , µ * ) require a large amount of CPU resources as well as time. It limits its uses in the solution of the direct problem. Compared with the numerical simulation with commercial software, a proper theoretical model can accurately analyze the transversally inhomogeneous waveguides. The characterization model that we present can be divided into two main parts. The first part is the full-wave analysis of the loaded waveguide regions. In the second part, the solutions corresponding to loaded regions are "stitched" at the proper boundaries. Here the mode-matching technique is employed. It involves dividing an irregular structure into several regular parts for which the electrodynamics equations can be solved analytically in terms of propagation modes. The mode amplitudes are then obtained by matching electric and magnetic fields on the interfaces between the parts. The general idea is that only a small number of lowest-order modes are often required to obtain sufficient accuracy, thus reducing the number of calculations. This technique leads to the theoretical value of the reflection factor R and/or transmission factor T .
The inverse problem ε * , µ * (R, T ) consists of the determination of the complex relative permittivity (ε * ) and/or complex relative permeability (µ * ). The method should involve the measurement of the reflected (S 11 and S 22 ) and transmitted (S 21 and S 12 ) signals. The relevant scattering parameters relate closely to the complex permittivity and permeability of the material. There are various approaches for obtaining permittivity and/or permeability from S-parameters. Each of conversion has different advantages and limitations. Therefore, we are interested to have the most efficient technique. Some commercial simulation software has an optimization function that uses various optimization https://www.journals.vu.lt/nonlinear-analysis algorithms in order to obtain predefined S-parameters. If the predefined S-parameters are the same as measured by VNA, thus we can solve the inverse problem -to assess the electrical properties of the materials. However, such a solution requires considerable time.
As mentioned earlier, in order to accurately determine dielectric constant and loss tangent, or another characteristic, it is necessary to have a mathematical model that relates the measured values of the S-parameters of the network analyzer to the characteristic of the material under test (MUT). However, there have been significant recent improvements in these techniques due to advances in computational electromagnetic (CEM) commercial simulation tools and numerical nonlinear solvers [2,13,14,16,17,22,[27][28][29]. These commercial computational tools have relaxed the geometrical restrictions on the specimen, allowing for the characterization of specimens with irregular shapes as well as complex anisotropy. Despite these advances, commercial software demands accurate results as well as the huge amount of CPU time and memory resources. CEM can be used in solving EM problems by numerical methods using models that are adopted for solving a specific problem. Hence, the main goal is to present our approach for inferring for inferring the complex permittivity from S-parameters, measured in a partially-filled waveguide, and to demonstrate the accuracy and numerical efficiency.
The remainder of the paper is organized as follows. In Section 2, we introduce the solution of the direct problem. In Section 3, the numerical calculations are carried out using the theoretical model, and the results are compared with HFSS simulation results. Section 4 concludes the paper.

Direct problem
To obtain the relative permittivity ε * = ε − ıε and permeability µ * = µ − ıµ of the dielectric material in the frequency range from ∼10 GHz to ∼100 GHz, we measure the complex reflection coefficient R and the complex transmission coefficient T of the dominant TE 10 electromagnetic wave from the sample of the material in the shape of a rectangular rod placed in the center of the rectangular waveguide ( Fig. 1). It is assumed that the width a and the height b of the waveguide are such that only the main (TE 10 ) mode can propagate at the measuring frequencies.
To obtain accurate relations between ε * , µ * and R, T , we consider the measurement setup as consisting of three regions: a -empty waveguide to the left of the sample (towards the generator); b -a waveguide of the length d, partially filled in the center with the measured sample of width w, filling the entire height of the waveguide b so that the system is symmetric along the x-axis and uniform along the y-axis; c -empty waveguide to the right of the sample (towards the matched load).
The reflection factor R and transmission factor T are determined by the scattering from the sample inside the waveguide.
There are several methods to solve this electrodynamics problem. A relatively simple one would be the frequency-domain finite-difference (FD) method to solve the Maxwell equations on rectangular coordinate system with proper boundary conditions. While easy to implement, the FD method is quite resource (both memory and processor) intensive (depending on the required precision), and the direct problem R, T (ε * , µ * ) takes considerable time to complete. Since we are interested in the reverse problem ε * , µ * (R, T ), which in nontrivial cases can only be solved by iteration, the efficiency of the method is very important. Therefore, we are interested to have a most efficient method albeit more complex analytically.
Here the mode-matching approach comes to mind. It involves dividing an irregular structure into several regular parts (we already have done that) for which the electrodynamic equations can be solved analytically in terms of propagation modes. The mode amplitudes are then obtained by matching electric and magnetic fields on the interfaces between the parts. The general idea is that only a small number of lowest-order modes is often required to obtain sufficient accuracy, thus reducing the amount of calculations.

Regions a and c: regular (empty) waveguide
The only propagating mode is TE 10 . Additional higher-order modes can appear close to the sample as a result of scattering. However due to the particular symmetry of our system only TE n0 modes can ever appear, where n is an odd number. The TE n0 mode has only three nonzero electric field (E y ) and magnetic field (H x and H z ) components: where h = k 2 − g 2 , k = 2π/λ = ω √ ε 0 µ 0 , g = nπ/a. The propagation constant h is real only for the TE 10 mode and imaginary for the higher-order modes. C is the mode amplitude.
Thus each mode in the region a is characterized by the mode amplitude C a i and propagation constant h i = k 2 − g 2 i , where g i = (1 + 2ı)π/a. We assign indices i = 0, 1, 2 . . . in the order of increasing |g i |. For the region c, we relabel the mode amplitudes accordingly: C c i .

Region b: partially-filled rectangular waveguide
For the empty part (0 x (a − w)/2 and (a + w)/2 x a) the field equations are similar to those of the empty waveguide (regions a and c): However, this time g 1 = ω 2 ε 0 µ 0 − h 2 is not equal to nπ/a, as in the case of an empty waveguide, but is to be determined from special boundary conditions, described below. For the filled part ((a − w)/2 r (a + w)/2), we can write (taking symmetry into account): The fields E y and H z must be continuous on the boundary between empty and filled parts, so the following boundary conditions must be satisfied: Substituting the above expressions for the fields, and taking into acount, that these equations must hold for arbitrary z, thus the propagation constant h of a particular mode must be the same for both layers, we obtain: A nonzero solution for the mode amplitudes exists only if Since g 1 and g 2 depend on the same propagation constant h via the relations g 2 1 = k 2 − h 2 and g 2 2 = k 2 ε * µ * − h 2 , the solution of this equation gives a set of discrete values of h i , g 1i , and g 2i each i corresponding to a particular TE mode. It can be obtained by iteration. The solution also defines the relation between the amplitudes C 1 i and C 2 i : . so that the ith mode is fully defined by h i and C 1 i , which we relabel as h bi and C b i to identify their relevance to the region b.

Mode matching
The reflection coefficient R we are trying to calculate is a ratio of amplitudes of the incident and reflected TE 10 modes in the waveguid. Similarly, the transmition coefficient T is a ratio of amplitudes of the incident and transmitted (passed-through) TE 10 modes. These amplitudes can be obtained by solving a general problem of obtaining the amplitudes of all modes involved in the scattering process. Once it is possible to define the electromagnetic field components in terms of the waveguide modes, the mode amplitudes can be calculated from the field continuity requirement: where z = 0 is the boundary between regions a and b, which we select as the phase reference plane, and z = d is the boundary between the regions b and c. Assuming the amplitude of the incident TE 10 wave to be equal to 1, for region a, we expand the fields in terms of the modes: Here f a i are shorthands for the known functions of x described in the previous sections: In the region b, each mode travels both ways, thus Here for 0 x (a − w)/2 and (a + w)/2 x a, and for (a − w)/2 r (a + w)/2, In the region c, all modes travel toward the matched load (positive z): The field continuity requirement at z = 0 becomes and at z = d, Equalities (1)-(4) are between functions. To transform them to the system of linear equations, we use the orthogonality property of the waveguide modes:

In our case, it is equivalent to
If we multiply both sides of (1) and (3) by f b j , j = 0, . . . , m and also those of (2) and (4) by f a j , j = 0, . . . , n, and take an integral from 0 to a, we get a system of 2(n+1)+ 2(m + 1) linear equations with the same number of unknown mode amplitudes: The known coefficients are the overlap integrals . . , n, j = 0, . . . , m, which all have analytical expressions i = 0, . . . , n, j = 0, . . . , m.
https://www.journals.vu.lt/nonlinear-analysis Solving this system of linear equations gives us R and T along with all the other mode amplitudes. The numbers n an m, needed to obtain the required accuracy, have to be estimated empirically. To match the accuracy of the measurements, only a few (3 to 5) modes are normally required.
3 Results of numerical simulations and discussions

Verification of solution of the direct problem
A rectangular waveguide filled with the rectangular rod-shaped specimen method is a variant of the partially-filled waveguide method. A cylindrical rod model was previously used for dielectric permittivity calculations [8]. In order to verify our solution of the direct problem, the nonmagnetic (its permeability is equal to unity) dielectric material was chosen for the simulations. It will be studied the configuration where the rectangular rod-shaped dielectric specimen is partially filling the cross-section of the rectangular waveguide. Theoretical transmittance/reflectance characteristics are calculated with the three-dimensional simulation software HFSS, which is performed by using the finiteelement model (FEM) for EM analysis. To study the rectangular specimen partiallyfilled rectangular waveguide more specifically, characteristics of the same parameters are numerically calculated by using the proposed theoretical model and compared to HFSS simulations results.
The proposed model can be applied for the frequency band from ∼10 GHz to ∼100 GHz. But for simplicity, it will be focused on TE n0 modes in subsequent analyses, and the X-band frequencies regime (8−12 GHz) of the waveguide is chosen to verify the validity of the verified model and illustrate physical phenomena in the theoretical analysis. The structural parameters are set as follows: a = 23 mm, w = 0.6 mm, d = 1.35 mm (see Fig. 1). For the evaluation of the calculation accuracy of our model, it was calculated the dependence of complex reflectance R * on the number of modes taken into account. The calculations were performed when the rectangular waveguide is filled with a baseline material with small magnetic permeability, which ε = 100, ε = 10, and complex permeability µ = 1, µ = 0 at a frequency 10 GHz. The results are shown in Fig. 2.
As a matter of fact, the change of complex reflection coefficient R * is 0.001 when 20 or 30 modes are taken into account. The obtained result shows that using our model with 30 modes taken into account is more than enough to obtain accurate results. For this reason, further calculations are performed with 20 modes taken into account. The accuracy of HFSS simulation depends on the mesh cell dimensions. Iterations are performed by employing adaptive meshing. The mesh is automatically tuned to give the most accurate and efficient mesh possible subdividing the mesh in the regions where the largest gradients in the E-fields occur. After each adaptive pass, HFSS compares the S-parameters from the current mesh to the results of the previous mesh (the criteria DeltaS ). If the answers have not changed by the user-defined value or DeltaS , then the solution has converged, and the current or previous mesh can be used to perform a frequency sweep. If the solution has converged, then technically, the previous mesh is as good as the current mesh. In this case, Ansoft HFSS will use the previous mesh.  According to the parameter DeltaS , we can judge the accuracy of the calculations. In our calculations, it was defined DeltaS < 0.001. Performing HFSS simulations with different values of parameters at multiple frequencies proves to be very time-consuming. The duration of calculation of R * with 20 modes taken into account takes considerably less time with the same accuracy in comparison to the calculation with HFSS by using the finite element method. The simulations with HFSS parameters, optimization goals, and useful tips are discussed more detailed in the thesis by Rudys [25] and in [26].
For material with complex permittivity ε = 100, ε = 10 and complex permeability µ = 1, µ = 0, the modulus of reflection and transmission coefficients obtained using the proposed theoretical model and the Ansys HFSS program are plotted in Fig. 3. We observe an excellent agreement between two different computational methodologies.
The proposed partially-filled rectangular waveguide model accurately incorporates only a small number of lowest-order and obtains sufficient accuracy, thus reducing the amount of calculations of the direct problem. It was separately evaluated relative errors of the reflection and transmission coefficients, and it can seen that errors are a bit higher for the transmission coefficients but remain below 1.1% in the worst case. It was established that errors in the reflection coefficients remain low and do not exceed even 0.05. This result validate our full-wave analysis in the case of dielectric material (complex permittivity ε = 100, ε = 10 and complex permeability µ = 1, µ = 0). Typically, realistic experimental verification of the mathematical model can be provided, but the previous considerations make you think that the differences between the presented model and HFSS simulations are considerably lower than the possible measurement uncertainty. This problem arises due to the many error factors during experimental measurements, i.e. vectorial network analyzer errors, parasitic parameters of the real system, errors in sample geometry, etc. To estimate the accuracy of the realistic measurements, it will be determined the possible measurement uncertainty due to the network analyzer and calibration kit as it is provided by the manufacturer in the specifications of any VNA. It is possible to download a free uncertainty calculator from manufacturer's web page to generate accuracy curves for any VNA. We have chosen a vectorial network analyzer that shows high accuracy at 10 GHz frequency -keysight technologies PNA N5242B. The technical datasheet of the selected PNA states that the S-parameters magnitude dynamic accuracy is 0.01 dB, and phase dynamic accuracy is 0.8 • at 10 GHz. Differences in magnitude between all suggested model values and HFSS simulation did not exceed 0.003 dB in magnitude and 0.1 • in phase. The verification by using commercial software HFSS based on the finite element method showed consistent results with considerably higher accuracy. Thus, the experimental measurement results are not presented as they cannot provide information on the accuracy of the presented model.

Sensitivity analysis
The performed direct EM analysis of material under test (MUT) plays an important role in the inverse problem developed to find the complex relative permittivity and permeability from the S-parameters measured with the vectorial network analyzer (VNA). Therefore, the applicability of the presented partially-filled waveguide technique depends on the sensitivity of S-parameters to small changes in the material properties. A comprehensive sensitivity analysis should be performed to determine how the predefined parameters of the sample influence the sensitivity of the technique and how it impacts overall measurement uncertainty. Generally, the goal of sensitivity analysis is to describe how much model output values are affected by changes in model input values. As well as, sensitivity analysis can provide information about possible error sources. The sensitivity analysis of the partially-filled rectangular waveguide at the X-band for extracting the permittivity and permeability was performed by a series of simulations.
Reflectance at 8, 10, and 12 GHz. In Fig. 4 presented baseline material (ε = 100, ε = 10, µ = 1, µ = 0) was selected to analyze the impacts of slight variations in the real and imaginary parts of the permittivity. In order to investigate the sensitivity of the (a) (b) Figure 4. Sensitivity analysis of the (a) amplitude and (b) phase to the slight variation of the imaginary and real parts of the permittivity, respectively. The real part of the permittivity ε was varied by ±1% and ±5%. The imaginary part of the permittivity ε was varied by a factor of 1/2, 2, 1/4, and 3/4. The baseline material was selected with the following parameters: ε = 100, ε = 10, µ = 1, µ = 0. measurement technique to changes in the loss factor ε , the imaginary component of the permittivity was varied by a factor of 1/2, 2, 1/4, and 3/4, i.e. ε = 5, 20, 2.5, 7.5, while ε = 100, respectively. Similarly, the real part of the permittivity ε was varied by ±1% and ±5%, while maintaining ε as a constant value, such that ε = 101, 99, 105, 95, while ε equals 10, respectively. Simulations were performed with constant permeability value (µ = 1, µ = 0) at these complex permittivity values. Contrary to fully-filled waveguide techniques, machining imperfections do not limit the accuracy of this technique, therefore, the dimensions of the specimen are not considered as huge error factors, and we do not incorporate structural parameters into the performed sensitivity analysis. Therefore, the specimen geometrical dimensions of 0.6 mm × 1.35 mm (w × d, schematically illustrated in Fig. 1) were held constant for these simulations. The shift of the reflection coefficients and phases after sensitivity analysis with respect to the slight variation of the imaginary and real parts of the permittivity are shown in Figs. 4(a), 4(b), respectively. The changes in the amplitude and phase directly correlate to changes in ε and ε , respectively. The sensitivity to changes in the real ε and imaginary ε parts of the permittivity should be directly compared to the theoretical measurement repeatability of the network analyzer. As a conservative threshold, we consider the measurement repeatability of the network analyzer to be 0.02 dB and 0.1 • in amplitude and phase, respectively. For this simulated material, a change of the imaginary part of dielectric permittivity of specimen results in a change of amplitude of reflection coefficient in the range 0.16−0.75 dB as shown in Fig. 4(a), which is greater than the measurement uncertainty of the network analyzer, 0.02 dB. Therefore, the presented theoretical model overcomes the conservative network analyzer tolerances used for this analysis even at small enough changes in ε . For both cases (±1% and ±5%), the measurement sensitivity to the real part of the permittivity is greater than the measurement repeatability of the network analyzer in phase. For the change in 1% of the real part of the permittivity resulted in a ∼0.40 • change in phase, as shown in Fig. 4(b), which is over 4 times the measurement repeatability of the network analyzer. It should be noted that sensitivity studies can provide a general assessment of model precision; however, it does not include additional parasitic uncertainties due to the grooves in the waveguides or placement of the rod in the center of the waveguide, etc. In general, performed sensitivity analysis showed that the method presented in this work is robust in terms of variances in material properties to ensure accurate characterization.

Considerations about the magnetic field
Distribution of the electric field within the specimen is characterized by the dielectric constant. One may observe that the electric field inside the sample is discontinuous (Fig. 5). The discontinuity of the electric and the magnetic fields is observed even within the samples with much lower dielectric permittivity (for instance, ε = 100, ε = 10). The magnetic field in the sample is stronger than in the excitation wave (Fig. 5). Thus, we can expect the system scattering parameters to be sensitive to the magnitude of the magnetic permeability of the sample.  It will be calculated the reflection and transmission coefficients when µ = 1, µ = 0 and when the magnetic permeability is slightly higher -µ = 1.5, µ = 0. The magnitude of the S-parameters is shown in Fig. 6.
It is clearly seen that in these cases the values of the S-parameter modules may differ more than the measurement error, so if the sample has even a small magnetic permeability (compared to the dielectric permittivity), distorted dielectric permeability measurement results can be obtained. Thus, when using the partially-filled waveguide method for measurements, we must take into account the possible, even small, magnetic permeability of the sample.

Conclusions
In this work, an X-band partially-filled rectangular waveguide measurement method that permits the characterization of the complex permittivity for dielectric material with magnetic properties specimen was presented. The main objective of this paper was to develop the techniques that enable us to solve problems with dielectric materials in the crosssection of the rectangular waveguide. The proposed mathematical model is based on mode-matching of a partially-filled rectangular waveguides. It shows excellent agreement with calculations performed with commercial HFSS software. It was demonstrated that the sensitivity of the technique is considerably higher than the repeatability of the vectorial network analyzer. The suggested model accounts for the value of complex magnetic permeability, and it can be evaluated simultaneously (this is an advantage of our model).