FREE-SPACE CHARACTERIZATION OF CONDUCTOR-BACKED ABSORBING MATERIALS USING AN APERTURE SCREEN By Korede Oladimeji A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements for the degree of Electrical and Computer Engineering – Doctor of Philosophy 2022 ABSTRACT Absorbing materials are often adhered to conducting surfaces for the purpose of controlling the electromagnetic field scattered by objects. Many of these materials have both electric and magnetic properties, and these properties may degrade with time, thereby decreasing the effectiveness of the coatings. Because of this, it is important to accurately assess the health of the coatings by interrogating them with an electromagnetic wave and analyzing the inter- action of the wave with the coating materials. Ideally, the permeability and permittivity of the coatings would be measured and compared to baseline values. Since these measurements must be done in the field, the materials cannot be removed from the underlying conducting surfaces. One convenient way to measure the permittivity and permeability of coatings is to il- luminate a coated surface using an antenna placed at a certain standoff distance from the coated object. Standoff techniques do not involve physical contact with the coatings and thus reduce the possibility of damaging the coating during the measurement process. Because both complex permeability and complex permittivity are desired, two sufficiently different measurements of the complex reflected field are required. Previous studies have shown that varying the polarization or incidence angle of the interrogating field does not provide enough variation on the reflected field for robust measurements of the material parameters. These studies have also shown that applying a material layer in front of the coating does not alter the information about the coating available from measurements of the reflected field, and is thus ineffective. In the technique proposed here, one measurement is made by illuminating the coated surface with a plane wave and a second measurement is made by illuminating the coated surface with a conducting screen containing an aperture placed immediately on top. This approach has proven effective with a waveguide contact probe, and the purpose here is to assess its viability as a free-space technique. The specific case of a narrow rectangular aperture is emphasized due to its simplicity of analysis compared to other aperture shapes. The constitutive parameters are extracted by comparing the measured reflected field in the presence of the aperture to the reflected field obtained from a numerical model that has been developed. The model is based on plane wave excitation of an infinite layered medium with the reflected field found by numerically solving a magnetic field integral equation. The numerical solution is validated by using a radiation problem with a line source placed in the aperture. Error analysis is used to compare the efficacy of the proposed aperture method to that of the two-thickness method (which, although effective, cannot be applied in the field). Calibration of the approach is also considered and measured results are described. For Anna, Ella and the little one that’s on the way. iv ACKNOWLEDGMENTS I would like to express my most sincere gratitude to Dr Rothwell for his support. I am thankful that he chose to continue on as my advisor through periods that were challenging in getting my work complete. I would also like to thank my committee members: Dr Frasch, Dr Chahal and Dr Shanker. They provided invaluable advice and showed great patience that I am truly grateful for. Also, thanks to Dr Pierre for planting the seed in my heart as an undergraduate student that led me to pursue a graduate degree. My family provided constant encouragement and believed even in moments when I had doubt. Thank you Mom for your constant checking in. Thank you Temi, Bowo and Iyin for all the help proof-reading a document that was only slightly meaningful to you. Thanks for the uplifting phone calls as well. Last but most certainly not least, thank you Anna for your perseverance through this process. I’m grateful for your understanding over the span of time when we had to focus on various aspects of this work. I couldn’t not have done this without you. And thanks for always making us smile Ella. v TABLE OF CONTENTS CHAPTER 1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1. Motivation and Background . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2. Free-space Aperture Technique . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.3. Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 BIBLIOGRAPHY . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 CHAPTER 2 Theoretical Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.1. Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.2. Plane-wave reflection from a conductor-backed slab . . . . . . . . . . . . . . 10 2.3. MFIE analysis of a slotted conductor above a conductor backed slab . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.4. Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 BIBLIOGRAPHY . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 CHAPTER 3 Numerical Solution to Integral Equations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 3.1. Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 3.2. General MFIE Formulation (Version 1 MFIE) . . . . . . . . . . . . . . . . . 28 3.3. Hallen-Type MFIE Formulation (Version 2 MFIE) . . . . . . . . . . . . . . . 62 3.4. Considerations for Numerical Accuracy . . . . . . . . . . . . . . . . . . . . . 85 3.5. Solution Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91 3.6. Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 BIBLIOGRAPHY . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 CHAPTER 4 Extraction of Material Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 4.1. Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 4.2. Considerations for Calibration . . . . . . . . . . . . . . . . . . . . . . . . . . 99 4.3. Effect of  and µ Variation on the Scattered Field . . . . . . . . . . . . . . . 101 4.4. Constitutive Parameter Extraction . . . . . . . . . . . . . . . . . . . . . . . 105 4.5. Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109 BIBLIOGRAPHY . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110 CHAPTER 5 Aperture Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 5.1. Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 5.2. Inversion Problem Solution . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112 5.3. Error Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113 5.4. Aperture Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 129 5.5. Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137 vi BIBLIOGRAPHY . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 138 CHAPTER 6 Experiments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139 6.1. Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139 6.2. MSU Measurement Setup and Implementation . . . . . . . . . . . . . . . . . 140 6.3. External Measurements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 155 6.4. Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172 BIBLIOGRAPHY . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 173 CHAPTER 7 Conclusions and Further Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 174 7.1. Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 174 APPENDIX A: GREEN’S FUNCTION DERIVATION . . . . . . . . . . . . . . . . . 176 APPENDIX B: STANDARD ERROR PROPAGATION METHOD . . . . . . . . . . 187 vii CHAPTER 1 Introduction 1.1 Motivation and Background 1.1.1 Material Characterization The accurate determination of the constitutive parameters ( and µ) of materials is of par- ticularly great importance in various fields ranging from bioengineering to agriculture and medicine [1] [2] [3]. Material characterization involves the employment of appropriate mea- surement and extraction techniques to obtain the permittivity () and permeability (µ) of materials of interest. Knowledge of the values of these parameters and changes in their properties provide valuable insight that is useful in informing design choices and monitoring process quality. This plays a major role in applications ranging from communication to microelectronics, manufacturing, and food processing, just to name a few. For instance, the degradation in performance of an antenna radome over time may be monitored based on changes in the permittivity of the radome material. The problem considered here is that of the characterization of material coatings on the surface of aircraft. The very nature of this problem makes it necessary that whatever mea- surement approach is taken be a reflection-only approach. This is because transmission measurements are precluded since the materials in question are conductor-backed. Also, techniques that require access to different material thicknesses or displacements of the ma- terial relative to the conductor backing are impractical since the coatings are affixed to the aircraft surface and cannot be excised. What is desired is an in situ conductor-backed material characterization approach. In this work, a free-space aperture characterization method is proposed. The proposed method is non-destructive, non-contact and does not require the availability of transmission measurements. Before delving into the details, some general background and current state of the art in material characterization techniques are briefly discussed. 1 1.1.2 Material Characterization Techniques Material characterization is founded on the measurement of the response of a material in the presence of an applied electromagnetic field. The constitutive parameters of a sample material can be extracted upon considering two implications of the interaction between an electromagnetic field and the sample: (1) a robust and repeatable means of exposure may be realizable, and (2) an accurate and tractable mathematical model describing this mechanism can be developed where possible. Furthermore, there must be access to two distinct complex measurements since two parameters ( and µ) are to be extracted. This is analogous to the solution of an equation with two variables: two independent equations are required in order to solve for them. All materials characterization techniques employ a measurement compared against a mathematical model with the constitutive parameters obtained by obtaining agreement be- tween measurements taken in the lab and the theoretical model. Important considerations are the complexity of the model and the feasibility of the material exposure method. Broadly speaking, the simpler the model and the more flexible the exposure method the better. As such, what is desired is a technique that combines a well-understood theoretical model with a suitable method for exposing the material sample to an interrogating electromagnetic field. Additionally, a strong response is desired as it may difficult to successfully extract the pa- rameters otherwise. Material characterization techniques can be broadly classified into the following ap- proaches: (1) Guided-wave technqiues, (2) Resonant Cavity techniques, and (3) Free-space techniques. Some factors to be cognizant of when choosing which technique is appropriate for a particular application include: frequency range of interest, temperature, desired level of accuracy, MUT (Material-Under-Test) size restrictions and MUT form (solid, liquid or gas). Also of importance is whether a destructive or non-destructive measurement is required and whether contact can be made with the MUT. The choice of which technique to employ usually involves trade-offs between these. 2 1.1.2.1 Guided-Wave Techniques Guided-wave techniques can be further classified into open waveguide, filled waveguide and waveguide probe techniques. In the case of open waveguide and filled waveguide techniques, sample excision and preparation is necessary to ensure that the MUT is ready to be “held.” As such, neither of these methods are well suited for the problem being considered since they are not non-destructive. Waveguide probe techniques, on the other hand, are accommodating to in-situ interro- gation of the surface of materials [4] [5]. There is a need in the case of these techniques, however, for the MUT to be planar as the waveguide aperture must make proper contact with the surface of the material being interrogated. Hence, they are also not applicable to the problem here. Filled guide techniques include coaxial, stripline and rectangular waveguide techniques [6] [7]. Coaxial waveguide techniques are applicable in scenarios where the materials can be fashioned into a coaxial shape. The constitutive parameters can be obtained by measuring S-parameter data from a network analyzer connected to the coaxial waveguide. The material properties can then be extracted from processed S-parameter data. Rectangular waveguide techniques are best suited to applications where the MUT is solid and can be readily machined/fashioned into a rectangular shape to fit into a waveguide fixture. S-parameter measurements through the waveguide are obtained and the constitutive parameters are extracted by processing the S-parameter data with known algorithms [8] [9]. Rectangular waveguide techniques are not suitable for the problem being considered. This is because they require the excision of the MUT and an imposition of sample size restrictions. Additionally, reflection measurements and transmission measurements (which are unavailable for a conductor-backed material) are needed for implementing the widely used Nicolson-Ross-Weir algorithm. Work has been done to extend the waveguide technique to the characterization of conductor-backed materials using a two-iris method [10]. But the need to fashion the MUT sample into sizes that can be fitted within a waveguide appropriately 3 still remains. Other guided-wave implementations include dual-ridged waveguides[11], stepped flange waveguides [12], and the triaxial applicator system [13]. These different implementations have advantages over the other approaches that arise such as their applicability to magnetic materials. However, in certain cases the corresponding mathematical models may become intractably complex. 1.1.2.2 Resonant Cavity Techniques Resonant cavity methods usually boast great accuracy and support for high temperature measurements [14] [15]. However, they are very narrow-banded and require special sample preparation that may be difficult and costly [16]. They are usually applied in determining the permittivity and loss tangent of low-loss materials. The dielectric properties are determined by first taking a measurement of the resonant frequency and Q of an empty cavity. The same measurements are carried out with the material present and the constitutive parameters can be computed from the resulting shift in frequency, volume and Q-factor. Resonant cavity techniques are by definition not non-destructive, are only applicable at discrete frequencies and also require special sample preparation. Therefore, they are not suitable for the conductor-backed measurement desired. 1.1.2.3 Free-space Techniques Free-space techniques are applicable over a wide frequency range with practical MUT sample size constraints imposed at lower frequencies. They are also applicable to a wide range of materials and useful in cases where measurements at high temperatures are desirable [17]. A figure illustrating a free-space measurement system is shown in Figure 1.1. Reflection and/or transmission measurements are obtained from a measurement system with the MUT illuminated by a transmit antenna and received by another antenna. A focused beam system can be utilized to reduce the effects of diffraction and time gating can be employed to minimize the impact of multiple reflections. In contrast to guided-wave and resonant cavity methods previously discussed, free-space 4 Figure 1.1: Free-Space Setup illustration. techniques by definition can be implemented to satisfy the non-contact and non-destructive conditions desired for this application. Furthermore, the presence of a conductor backing narrows down the available options to reflection-only free-space methods since there is a lack of access to transmission measurements. As will be discussed in the next section, there are constraints and challenges that make the available reflection-only free-space methods inapplicable or troublesome to implement. Hence the need for the proposed technique. 1.2 Free-space Aperture Technique This technique has the the advantage of being non-contact, robust and flexible. Before delving further into the details we first consider free-space techniques where only reflection measurements are available and/or obtainable. These restrictions necessitate the development of a novel free-space approach which can be implemented in-situ in a non-destructive, non-contact manner. Additionally, the ap- proach must have a mathematical model that is computationally tractable and realizable. A survey of the literature revealed the following reflection-only methods:the Two-Thickness 5 Method [18], the Air/Conductor Backed Method [19], the Layer-Shift Method [20], the Two- Polarization Method, and the Frequency Varying Method [21]. The two-thickness methods, as its name suggests, involves the carrying out/execution of two separate measurements with a conductor-backed MUT. The two measurements are taken with two samples of the same MUT having different thicknesses [19]. The air/conductor backed method similarly involves the performance of two measure- ments: the first with the MUT by itself without a conductor backing and the other with an identical sample of the MUT with a condcutor backing. The layer-shift method involves one measurement of the MUT in a conductor-backed configuration with the conductor placed right behind and in contact with the sample, and a second measurement with conductor backing displaced a certain distance away from the surface of the sample. These reflection-only techniques outlined are not well suited for the problem being con- sidered for various reasons including this need for contact, inapplicability due to the curved nature of the aircraft surface, and lack of access to implement any displacement of the MUT. Also, the two-polarization method is highly sensitive to errors in the knowledge of aspect angle. To recap, the technique proposed in this dissertation is a free-space technique that can be implemented in situations where transmission measurements are precluded given that the sample being considered is conductor backed and adhered to the conductor. The permittivity and permeability are determined by interrogating the target under test in two different distinct scenarios. One scenario involves the illumination of the material and measurement in the specular direction. The second scenario involves the illumination of the target with a metallic aperture adjacent to it and measurement of the reflection in a non-specular direction. After two measurements are taken under the aforementioned scenarios, a root search method is employed to minimize the difference between the measurements and theoretical results so that the permittivity and permeability can be extracted. 6 1.3 Summary The proposed material characterization technique is a reflection-only, in situ, nondestructive and non-contact approach. The upcoming chapter provides the theoretical framework for the free-space aperture technique. The choice of this approach is motivated by the unique circumstance under which mea- surements are being taken for the problem being considered here. The approach satisfies the desire for a tractable mathematical model in conjunction with an appropriate physical method for holding the MUT and capturing the interaction with an electromagnetic field. The upcoming chapters detail the proposed free-space aperture technique. Chapter 2 delves into the foundations of the formulation of the mathematical model. In Chapter 3, the numerical solution to the mathematical model is framed in terms of integral equations with the details of the solution provided. The extraction approach is detailed in Chapter 4 along with details on the implementation for the extraction of  and µ. Chapter 5 describes the implementation of the optimization that is employed for determining the slot dimensions that lead to the optimal extraction parameters. The measurements carried out to validate the approach are discussed in Chapter 6. In Chapter 7, concluding remarks and future work are highlighted. 7 BIBLIOGRAPHY [1] A. R. Von Hippel. Dielectric Materials and Applications: Papers by Twenty-two Con- tributors. Technology Press books in science and engineering. M.I.T. Press, 1954. [2] H. E. Bussey. Measurement of rf properties of materials a survey. Proceedings of the IEEE, 55(6):1046–1053, 1967. [3] L. F. Chen, C. K. Ong, C. P. Neo, V. V. Varadan, and V. K. Varadan. Microwave Electronics: Measurement and Materials Characterization. Wiley, 2004. [4] M Hyde, J. W. Stewart, M. Havrilla, W. P. Baker, E. J. Rothwell, and D. Nyquist. Nondestructive electromagnetic material characterization using a dual waveguide probe: A full wave solution. Radio Science, 44, 2009. [5] G. D. Dester, E. J. Rothwell, and M. J. Havrilla. An extrapolation method for improving waveguide probe material characterization accuracy. IEEE Microwave and Wireless Components Letters, 20(5):298–300, 2010. [6] H. Yue, K. L. Virga, and J. L. Prince. Dielectric constant and loss tangent measurement using a stripline fixture. IEEE transactions on components, packaging, and manufac- turing technology. Part B, Advanced packaging, 21(4):441–446, November 1998. [7] W. Barry. A broad-band, automated, stripline technique for the simultaneous mea- surement of complex permittivity and permeability. IEEE Transactions on Microwave Theory and Techniques, 34(1):80–84, 1986. [8] A. M. Nicolson and G. F. Ross. Measurement of the intrinsic properties of materials by time-domain techniques. IEEE Transactions on Instrumentation and Measurement, 19(4):377–382, 1970. [9] W. B. Weir. Automatic measurement of complex dielectric constant and permeability at microwave frequencies. Proceedings of the IEEE, 62(1):33–36, 1974. [10] G. D. Dester. Electromagnetic Material Characterization of a Conductor-backed Ma- terial Using the Two Layer, Two Thickness, and Two Iris Waveguide Probe Methods: Error Analysis, Simulation, and Experimental Results. Michigan State University. De- partment of Electrical Engineering, 2008. [11] M. W. Hyde and M. J. Havrilla. Simple, broadband material characterization using dual-ridged waveguide to rectangular waveguide transitions. IEEE Transactions on Electromagnetic Compatibility, 56(1):239–242, 2014. [12] J. P. Massman, M. J. Havrilla, K. W. Whites, and M. W. Hyde. A stepped flange waveguide technique for determining tapered r-card sheet impedance. In 2010 Asia- Pacific Microwave Conference, pages 1769–1772, 2010. [13] S. Karuppuswami, E. J. Rothwell, P. Chahal, and M. J. Havrilla. A triaxial applicator for the measurement of the electromagnetic properties of materials. Sensors, 18(2), 8 2018. [14] J. Sheen. Microwave measurements of dielectric properties using a closed cylindrical cavity dielectric resonator. IEEE Transactions on Dielectrics and Electrical Insulation, 14(5):1139–1144, 2007. [15] Keysight Application Note. Split Post Dielectric Resonators for Dielectric Measure- ments. Keysight Technologies, 2017. [16] M. S. Venkatesh and V. Raghavan. An overview of dielectric properties measuring techniques. Canadian Biosystems Engineering / Le Genie des biosystems au Canada, 47:15–30, 01 2005. [17] B. Clarke. Measurement of the dielectric properties of materials at RF and microwave frequencies, pages 409–458. 01 2007. [18] J. Baker-Jarvis, E. J. Vanzura, and W. A. Kissick. Improved technique for determining complex permittivity with the transmission/reflection method. IEEE Transactions on Microwave Theory and Techniques, 38(8):1096–1103, 1990. [19] R. Fenner, E. J. Rothwell, and L. Frasch. A comprehensive analysis of free-space and guided-wave techniques for extracting the permeability and permittivity of materials using reflection-only measurements. Radio Science, 47, 01 2012. [20] A. A. Kalachev, I. V. Kukolev, S. M. Matitsin, L. N. Novogrudskiy, K. N. Rozanov, A. Sarychev, and A. V. Seleznev. The methods of investigation of complex dielectric permittivity of layer polymers containing conductive inclusions. MRS Proceedings, 214, 01 2011. [21] S. Wang, M. Niu, and D. Xu. A frequency-varying method for simultaneous measure- ment of complex permittivity and permeability with an open-ended coaxial probe. IEEE Transactions on Microwave Theory and Techniques, 46(12):2145–2147, 1998. 9 CHAPTER 2 Theoretical Model 2.1 Introduction This chapter details the theoretical work that was carried out to describe the proposed aperture technique. As mentioned in the previous chapter, a mathematical model is needed to describe the interaction between the interrogating electromagnetic field and the MUT. As a starting point, reflection from a coated conductor surface is considered. As will be shown, this is a problem for which closed form expressions relating the wave impedance to the permittivity and permeability can be derived. Following this, a theoretical model is developed for analysis of a slotted conductor above a conductor-backed slab. 2.2 Plane-wave reflection from a conductor-backed slab Before the theory for the aperture problem is formulated, the field reflected by a plane wave incident on a conductor-backed material layer without an aperture screen is considered. Figure 2.1 shows an illustration of the general problem. A plane was is incident at an angle θ0 on a conductor-backed material with thickness t having constitutive parameters  and µ. The field can be polarized either in a parallel or perpendicular sense to the plane of incidence. Given this scenario, global reflection coefficients can be expressed as follows depending on the polarization [1] The Fresnel (interfacial) reflection coefficients above are Z − Z0⊥ Γ⊥ = ⊥ (2.1) Z⊥ + Z0⊥ and Zk − Z0k Γk = , (2.2) Zk + Z0k 10 Figure 2.1: Conductor-backed material layer with slotted conductor directly on top. 11 where the various impedances are given by the expressions kη Z⊥ = , (2.3) kz η0 Z0⊥ = , (2.4) cos θ0 kz η Zk = , (2.5) k Z0k = η0 cos θ0 (2.6) and kz2 = k 2 − k02 sin2 θ0 , k 2 = ω 2 µ, k02 = ω 2 µ0 0 , η 2 = µ/, η02 = µ0 /0 . The term P is the propagation factor for the wave traversing the material layer, and is given by P = e−jkz t . It should be noted that the product of the permittivity and permeability appears in the wavenumber k, while the ratio of these quantities appears in the intrinsic impedance η. This behavior will be exploited later in the inversion algorithm. 2.3 MFIE analysis of a slotted conductor above a conductor backed slab A theoretical model for the field scattered by a conductor-backed material with a slotted conductor placed above it is considered next. As mentioned earlier, in this case a closed form expression cannot be obtained. An integral equation is formulated for the field in the aperture and then the scattered field can be computed with the aperture field used as an equivalent source. Consider a slot in a ground plane over a conductor-backed material region as shown in Figure 2.1. The material region (Region 2) has a permittivity 2 , permeability µ2 and thickness t while the half space above it (Region 1) has parameters 1 and µ1 . 12 The slot has a length L and a width w and is excited by a line current located at the point x = xo in the aperture or by an incident plane wave. The plane wave incidence case is of interest for the purpose of material characterization. The line current case is useful for validation since solutions obtained for the radiated field can be compared with results in the literature, as will be shown in the next section. A magnetic field integral equation can be obtained by applying the boundary conditions on the electric and magnetic field across the slotted conductor. As illustrated in the figure, the slot is assumed to be “thin” with the assumption that ko w  1 and w  L. Given this assumption, the electric field within the slot may be approximated as E ~ s = ŷEs . Let E A y be the aperture field immediately above the slot (i.e. in the free space region) and EyP P be the aperture field immediately below the slot (i.e. in the region bounded by the parallel conducting plates). In general, a magnetic current can be expressed as ~ m = −n̂ × E. K ~ (2.7) From (2.7) magnetic currents can be expressed for the two regions as   K~ A = −ẑ × ŷE A = x̂E A , (2.8) m y y   K P P = ẑ × ŷE P P = −x̂E P P , ~m (2.9) y y with K ~ A being the magnetic current in the half space immediately above the slotted con- m ductor and K ~ P P that immediately below the slotted conductor in the material region. m Using the Hertz potential representation of the fields, the electric and magnetic fields can be written as E ~ = −jωµ∇× Π ~ m and H~ = k2Π ~ m +∇(∇· Π ~ m ), respectively [2]. Furthermore, the corresponding wave equation for the Hertz potential is, as derived in the next section,   ~ 2 ∇ +k Π 2 ~ m = − Jm . (2.10) jωµ 13 The solution for Π ~ m produced by the magnetic current K ~ m on the aperture surface is Z ↔ ~ m (~r 0 ) K ~ m (~r) = Π G (~r|~r 0 ) · ds0 , (2.11) jωµ SA ↔ where G (~r|~r 0 ) is the dyadic Green’s function for either the free-space or parallel plate regions. The details of the Hertz potential formulation will be briefly considered in the following section. 2.3.1 Vector Potentials If it is assumed that no electric charges are present (ρe = 0, where ρe is the electric volume charge density) then the divergence of the electric field is zero. This implies that the electric field is the curl of some vector quantity since the divergence of a curl is zero. Therefore, the electric displacement can be expressed as D~ = ∇ × F~ (2.12) where F~ is the electric vector potential. Substituting D ~ from (2.12) into Ampere’s law, ∇×H ~ = jω D ~ (2.13) leads to ∇×H ~ = jω∇ × F~ , (2.14) which after rearranging becomes   ∇× H ~ − jω F~ = 0. (2.15) 14 Since the curl of the gradient of some quantity is zero, the quantity in parenthesis in (2.15) can be rewritten as jω F~ − H ~ = −∇φm . (2.16) Next (2.12) can be substituted back into Faraday’s law (assuming homogeneity), ∇×D ~ = −jωB~ − J~m , (2.17) to give ∇ × ∇ × F~ = −jωB ~ − J~m . (2.18) Furthermore, (2.16) can be rearranged so that ~ = jω F~ + ∇φm . H (2.19) Recalling that B~ = µH ~ gives   ∇ × ∇ × F = −jωµ jω F + ∇φm − J~m . ~ ~ (2.20) Utilizing (2.19) in the Magnetic Gauss’ law, ~ = ρm , µ∇ · H (2.21) it can be seen that   µ∇ · ∇φm + jω F~ = ρm , (2.22) 15 which can be rearranged to give 1 ∇2 φm + jω∇ · F~ = ρm . (2.23) µ Equation (2.20) can be expanded and simplified using the vector identity ∇ × ∇ × F~ = ∇∇ · F~ − ∇2 F~ yielding   ∇∇ · F~ − ∇2 F~ = −J~m + −jωµ jω F~ + ∇φm , (2.24) which upon expansion and rearrangement gives   −∇2 F~ − k 2 F~ + ∇ ∇ · F~ + jωµφm = −J~m (2.25) where k, the wavenumber, is given by the relation k 2 = ω 2 µ. Employing the Lorentz gauge condition, the quantity in parenthesis in (2.25) is set to zero so that ∇ · F~ = −jωµφm . (2.26) Therefore, ∇ · F~ φm = − , (2.27) jωµ and (2.25) becomes ∇2 F~ + k 2 F~ = J~m (2.28) which is a wave equation for the potential F~ . Similarly, the wave equation for the scalar 16 potential can be obtained by substituting (2.26) into (2.23) resulting in the expression ρm ∇2 φm + jω (−jωµφm ) = (2.29) µ which can be rewritten as ρm ∇2 φm + k 2 φm = . (2.30) µ From (2.12), the electric field can be expressed in terms of F~ as ~ = 1 ∇ × F~ . E (2.31)  The magnetic field can be similarly found by substituting (2.27) in (2.19) so that ! ~ = jω F~ + ∇ −∇ · F~ H , (2.32) jωµ or ! ~ ~ = jω F~ + ∇∇ · F H . (2.33) k2 Alternatively, F~ can be expressed in terms of the Hertzian vector potential Π~ m as F~ = jωµΠ~ m. (2.34) Substituting (2.34) in (2.28) , the inhomogeneous Helmholtz equation for the Hertzian vector potential is found to be J~m ∇2 Π ~ m + k2Π ~m = . (2.35) jωµ 17 A solution to (2.35) can be obtained by convolving the Green’s function with the source term on the right side of the equation. Going back to considering the scattering problem, in the half-space region above the slotted conductor, the space-domain dyadic Green’s function is given by [3] ↔ ↔ G = I G. (2.36) In this region, the Hertz potential is therefore Kmx (~r 0 ) 0 Z ~ m (~r) = x̂ Π G(~r|~r 0 ) ds (2.37) jωµ1 SA where G is the free-space Green’s function [3], Z∞ Z∞ −p|z−z 0 | G(~r, ~r0 ) = 1 e ~ ρ−~ ej k·(~ ρ 0 ) dk dk . (2.38) x y (2π)2 2p −∞ −∞ In the region bounded by the parallel conducting plates, ↔ ↔P P G =G = x̂GP P PP PP xx x̂ + ŷGyy ŷ + ẑGzz ẑ, (2.39) ↔P P where G represents the parallel plate Green’s function [4]. Geometrical considerations lead to the conclusion that GP P PP xx = Gyy since the underlying material is infinite in extent in the x and y directions. Because the slot is assumed to be narrow, K ~ m is assumed to be only x-directed, then Kmx (~r 0 ) 0 Z ~ m (~r) = x̂ Π P P Gxx (~r|~r ) 0 ds (2.40) jωµ2 SA 18 and the magnetic field is given as ~ = k2Π H ~ m + ∇(∇ · Π ~ m ). (2.41) The slot magnetic field is primarily x-directed, and ! ∂ 2 Hx = k2 + Πmx , (2.42) ∂x2 with the Hertzian potentials being K A (~r 0 ) 0 Z ΠAmx (~r) = G(~r |~r 0 ) mx ds , z > 0 (2.43) jωµ1 SA and Z PP 0 ΠP P GP P 0 Kmx (~r ) ds0 , z > 0 mx (~r) = xx (~r |~r ) jωµ (2.44) 2 SA in the free-space and material regions respectively. Referring back to the expressions for the currents in equations (2.8) and (2.9), from (2.42) the magnetic fields can be expressed as 2EyA (~r 0 ) 0 !Z ∂2 0 HxA = k12 + G(~r |~r ) ds , z > 0 (2.45) ∂x2 jωµ1 SA and −EyP P (~r ) 0 !Z ∂2 P P 0 HxP P = k22 + Gxx (~r |~r ) ds , z < 0. (2.46) ∂x2 jωµ2 SA Invoking the boundary conditions at the interface, the tangential magnetic field is discon- 19 tinuous by a jump in the current so ẑ × (H~ −H ~ ) = K. ~ (2.47) 1 2 Substituting the values of the magnetic fields at the interface and the value of the surface current, ẑ × (x̂HxA + x̂Hzi − x̂HxP P ) = −ŷIg(x − xo ) (2.48) and HxA − HzP P = −Hxi − Ig(x − xo ), (2.49) where g(x) is a current distribution chosen subject to the condition Z L 2 g(x − x0 )dx = 1. (2.50) −L 2 Also, given that the tangential electric field is continuous EyA = EyP P = Ey . (2.51) Substituting (2.49) and (2.51) in (2.45) and (2.46) gives 2Ey (~r 0 ) 0 Ey (~r 0 ) 0 !Z !Z ∂2 0 ∂2 0 k12 + G(~r |~r ) ds + k22 + P P Gxx (~r |~r ) ds ∂x2 jωµo ∂x2 jωµ SA SA (2.52) = −Hxi − Ig(x − xo ), z > 0, ∀(x, y) ∈ SA . This is an MFIE for the aperture field Ey . 20 2.3.2 Field Expansion in Terms of Slot Voltages An alternate expression of (2.52) can be written in terms of the voltage across the slot. Remembering that the slot is narrow, the voltage across it exhibits only x-dependence and the electric field and voltage are related by Ey (x, y) = V (x)f (y) (2.53) where −w Z2 V (x) = − Ey dy. (2.54) w 2 Consequently, (2.52) can be rewritten as L w Z2 Z2 G(x, y|x0 , y 0 ) ! ∂2 k12 + 2V (x0 )f (y 0 )dy 0 dx0 + ∂x2 jωµ1 x0 =− L y 0 =− w2 2 L w (2.55) ! Z2 Z2 0 0 ∂2 GP P xx (x, y|x , y ) V (x0 )f (y 0 )dy 0 dx0 k22 + ∂x2 jωµ2 L 0 w x0 =− 2 y =− 2 = −Hxi (x, y) − Ig(x − xo ). The electric field, Ey, can be substituted from (2.53) and (2.54) becomes w Z2 V (x) = V (x) f (y)dy (2.56) −w 2 21 which implies that the above integral of f (y) is unity i.e. w Z2 f (y)dy = 1. (2.57) −w 2 In order to solve (2.55), a description of the magnetic current distribution is required. The function A f (y) = r (2.58) y 2   1− w/2 is chosen so as to satisfy the quasi-static edge singularities. Substuting (2.58) into (2.57), w Z2 dy 2A r  = 1. (2.59) y 2  0 1− w/2 y Setting u = , (2.59) becomes w/2 Z1 w du 2A p =1 (2.60) 2 1 − u2 0 and upon integration of (2.60) and some algebraic manipulation, it is found that 2 A= . (2.61) πw The value obtained for A in (2.61) be substituted into (2.59) and finally 2 f (y) = r  . (2.62) y 2  πw 1− w/2 22 The equation in (2.55) is multiplied by a weighting function W (y) and integrated over y, L w w Z2 Z2 Z2 G(x, y|x0 , y 0 ) ! ∂2 k12 + 2V (x0 )f (y 0 )W (y)dy 0 dx0 dy+ ∂x2 jωµ1 y 0 =− w w x0 =− L2 2 y=− 2 L w w (2.63) ! Z2 Z2 Z2 0 0 ∂2 GP P xx (x, y|x , y ) V (x0 )f (y 0 )W (y)dy 0 dx0 dy k22 + ∂x2 jωµ2 0 L 0 w y=− w x =− 2 y =− 2 2 = −Hxi (x, y) − Ig(x − xo ). After matching the field at the center of the slot (by choosing ω(y) = δ(y)), the MFIE can be written as L L ! Z2 ! Z2 ∂2 ∂2 k12 + G1 (x, x0 )V (x0 )dx0 + k22 + G2 (x, x0 )V (x0 )dx0 ∂x2 ∂x2 (2.64) −L2 −L 2 L L = −H(x) − Iδ(x − x0 ), − ≤x≤ 2 2 where w Z2 2G(x, y = 0|x0 , y 0 ) G1 (x, x0 ) = f (y 0 )dy 0 , (2.65) jωµ1 −w 2 w Z2 0 0 GP P xx (x, y = 0|x , y ) f (y 0 )dy 0 , G2 (x, x0 ) = (2.66) jωµ2 −w2 23 and H(x) = Hxi (x, 0). (2.67) 2.3.3 Alternative Form for Low-loss Materials An alternate approach which is found to be more amenable to application with low-loss media is also employed to solve the problem. L 2 Adding and subtracting the term k22 G1 (x, x0 )V (x0 )dx0 , (2.64) can be rewritten as R −L 2 L ! Z2 ∂2 G1 x, x0 + G2 x, x0 V (x0 )dx0 h    i k22 + ∂x2 −L2 L   Z2 (2.68) = k22 − k12 G1 (x, x0 )V (x0 )dx0 − H(x) − Iδ(x − x0 ), −L2 L L − ≤x≤ 2 2 which takes the form ! ∂2 k2 + F (x) = G(x), (2.69) ∂x2 which has a solution [1] Zx 1 F (x) = G(u) sin k(x − u)du + C1 sin kx + C2 cos kx (2.70) k x1 where C1 and C2 are arbitrary constants and x1 is arbitrary. Choosing x1 = −L 2 leads to 24 a Hallen-type MFIE as follows L Z2 h G1 (x, x ) + G2 (x, x ) V (x0 )dx0 0 0 i −L 2   L ! Zx  Z2 k22 − k12  −  0 0 0  G1 (u, x )V (x )dx  sin k2 (x − u)du  k22 (2.71)    −L2 − L 2 Z x 1 = H(u) sin k2 (x − u)du − C1 sin k2 x − C2 cos k2 x, k2 −L 2 L L − ≤x≤ . 2 2 Equations (2.64) and (2.71) fully describe the problem and numerical solutions will be pre- sented in the next chapter. 25 2.4 Summary The MFIE for the slot voltage has been derived in this chapter. In the next chapter the details of the solution of the MFIE using the Method of Moments are outlined both in the original case and for the Hallen alternative. 26 BIBLIOGRAPHY [1] E. J. Rothwell and M. J. Cloud. Electromagnetics. CRC Press, 2018. [2] E. A. Essex. Hertz vector potentials of electromagnetic theory. American Journal of Physics, 45(11):1099–1101, 1977. [3] R. E. Collin. Green’s Functions, pages 55–172. 1991. [4] C. Tai and P. Rozenfeld. Different representations of dyadic green’s functions for a rectangular cavity. IEEE Transactions on Microwave Theory and Techniques, 24(9):597– 601, 1976. 27 CHAPTER 3 Numerical Solution to Integral Equations 3.1 Introduction Now that the MFIE has been derived, what follows is: (a) the derivation of the numerical solution to the problem, and (b) validation of the solution obtained. The main concerns in the implementation of a problem solution are convenience in expressing the results and the management of computational expense. Additionally, close attention has to paid to numerical accuracy in the implementation. The approach utilized in the solution of the MFIE developed in the previous chapter is discussed in the sections that follow. Two solution approaches are considered. The first approach is applicable to the problem in general but difficulties are encountered when it is employed with low-loss materials. This approach is therefore utilized when the MUT is lossy. A second approach, which is better suited to dealing with low-loss materials, is also developed and is employed when the slot is in free space or over a low-loss MUT above a conducting screen. Finally, validation of solutions obtained are carried out to ensure accuracy before pro- ceeding to the optimization of the slot geometry. Numerical accuracy is discussed briefly as this is important in obtaining an optimal slot design. 3.2 General MFIE Formulation (Version 1 MFIE) An initial formulation is developed for the problem and presented in this section. The MFIE obtained in the previous chapter is considered and the two parts of it are examined separately. 28 Recall the MFIE L L ! Z2 ! Z2 ∂2 ∂2 + k12 G1 (x, x0 )V (x0 )dx0 + + k22 G2 (x, x0 )V (x0 )dx0 ∂x2 ∂x2 −L2 −L 2 (3.1) = −Hxi (x, 0) − Iδ(x − xo ), L L − ≤x≤ . 2 2 The free-space Green’s function is utilized in the region above the slot (free-space) whereas in the material region a parallel plate Green’s function is employed, i.e. w   Z2 2G x, y = 0|x0 , y 0 G1 (x, x0 ) = f (y 0 )dy 0 , (3.2) jωµ1 −w 2 w   Z2 GP P x, y = 0|x0 , y 0 xx G2 (x, x0 ) = f (y 0 )dy 0 . (3.3) jωµ2 −w 2 The free-space Green’s function is given by Z Z∞ G(x, y|x0 , y 0 ) = 1 1 j~k·(~ e ρ−~ρ 0 ) dk dk (3.4) x y (2π)2 2p1 −∞ where ~k = x̂kx + ŷky , (3.5) p21 = kx 2 + k2 − k2, y 1 (3.6) 29 and ρ~ = x̂x + ŷy = x̂x. (3.7) The Green’s function in the parallel plate region is Z Z∞ e−p2 t " # GP P 0 0 1 1 ~ ρ−~ ej k·(~ ρ 0 ) dk dk xx (x, y|x , y ) = 1+ x y (3.8) (2π)2 p2 sinh(p2 t) −∞ where p22 = kx 2 + k2 − k2. y 2 (3.9) Complete derivations for these Green’s functions are shown in Appendix 7.1. Using the appropriate Green’s functions from (3.4) and (3.8), G1 and G2 in (3.1) are expressed as w   Z2 Z Z∞ G1 (x, x0 ) = 2  1 1 j~k·(~ e ρ 0 ) dk dk  f (y 0 )dy 0 ρ−~ (3.10) x y jωµ1 (2π)2 2p1  −w2 −∞ and w   Z2 Z Z∞ e−p2 t " # G2 (x, x0 ) = 1  1 1 1+ ~ ρ−~ ej k·(~ ρ 0 ) dk dk  × x y jωµ2 (2π)2 p2 sinh(p2 t)  (3.11) −w2 −∞ f (y 0 )dy 0 respectively. The portion from (3.11) in square brackets is rewritten by taking advantage of 30 the definition of the sinh function so that e−p2 t e−p2 t 1+ =1+2 (3.12) sinh(p2 t) ep2 t − e−p2 t leading to e−p2 t 2 1+ =1+ . (3.13) sinh(p2 t) e2p2 t − 1 This alternative version will be used in expanding (3.11) later. The expressions in (3.10) and (3.11) are rearranged so that they can be expressed as convenient integrals over dkx and dky . Substituting ~k and ρ~ from (3.5) and (3.7) respectively into (3.10) gives w   Z Z∞  Z2 0 1 2 1 jk x −jk x 0 −jky y 0 0 0  G1 (x, x ) =  e x e x e f (y )dy  dkx dky . (3.14)  2 (2π) 2jωµ1 p1    −∞ −w2 Defining the following integral w Z2 Q(ky ) = e−jky u f (u)du, (3.15) −w 2 allows (3.14) to be written as Z Z∞ jkx (x−x0 ) 1 1 e G1 (x, x0 ) =  Q ky dkx dky . (3.16) 2 (2π) jωµ1 p1 −∞ Recalling that f is a distribution representing a magnetic current and is chosen to satisfy 31 the quasi-static edge singularities as done in (2.62), f can be substituted in Q resulting in w Z2 4 cos(ky u) Q(ky ) = r 2 du, (3.17) πw  u 0 1− w/2 which after employing the substitution x = u/(w/2) becomes Z1 cos ky w x   4 w 2 Q(ky ) = p dx. (3.18) πw 2 1 − x2 0 The integral in (3.18) can be evaluated using (3.753.2) of [1] Z1 cos(ax) π p dx = J0 (a), (3.19) 1 − x2 2 0 so that  w Q(ky ) = J0 ky (3.20) 2 where J0 is the Bessel function of the first kind and order zero. Therefore, Z Z∞ jkx (x−x0 )  1 1 e w G1 (x, x0 ) = J0 ky dkx dky . (3.21) (2π)2 jωµ1 p1 2 −∞ Substituting from (3.13) and similarly taking advantage of (3.20), (3.11) can be rewritten as Z Z∞ jkx (x−x0 )    1 1 e 2 w G2 (x, x0 ) = 1+ J0 ky dkx dky . (3.22) 2 (2π) jωµ2 p2 e2p2 t − 1 2 −∞ 32 This makes it possible to express the functions G1 and G2 concisely as Z Z∞ 0 G1,2 (x, x0 ) = W1,2 (kx , ky )ejkx (x−x ) dkx dky (3.23) −∞ where 1 1  w W1 = J0 ky (3.24) (2π)2 jωµ1 p1 2 and    1 1 2 w W2 = 1+ J0 ky . (3.25) (2π)2 jωµ2 p2 e2p2 t − 1 2 After carrying out the derivatives with respect to x, the MFIE in (3.1) can be rewritten in the form Z Z∞ L/2 0 Z W1 (kx , ky )(k12 − kx2) V (x0 )ejkx (x−x ) dx0 dkx dky + −∞ −L/2 Z Z∞ L/2 0 Z W2 (kx , ky )(k22 − kx2) V (x0 )ejkx (x−x ) dx0 dkx dky (3.26) −∞ −L/2 = −H(x) − Iδ(x − x0 ), − L/2 ≤ x ≤ L/2, 3.2.1 Method of Moments (MoM) Solution The Galerkin method [2] is employed for solving the integral equation in (3.26) by choosing the same expansion and testing functions. The functions are chosen so that the integrals obtained can eventually be computed in closed form. The voltage is expanded in terms of piecewise sinusoidal functions. The voltage across the slot is expanded as a set of of piecewise 33 sinusoidal functions XN V (x) = an fn (x) (3.27) n=−N with the basis function given by fn (x) = f0 (x − n∆) (3.28) where L ∆= (3.29) 2N + 2 and   sin k (∆ + x) −∆≤x≤0  0 f0 (x) = sin k0 (∆ − |x|) = . (3.30)  sin k0 (∆ − x)  0≤x≤∆ The voltage expression from (3.27) is substituted into (3.26) to give N Z Z∞ L/2 0 Z W1 (kx, ky)(k12 − kx2) fn (x0 )ejkx (x−x ) dx0 dkx dky + X an n=−N −∞ −L/2 N Z Z∞ L/2 0 Z W2 (kx, ky)(k22 − kx2) fn (x0 )ejkx (x−x ) dx0 dkx dky X an (3.31) n=−N −∞ −L/2 = −H(x) + Iδ(x − x0 ), − L/2 ≤ x ≤ L/2. 34 The expression in (3.31) is multiplied by fm (x) and integrated over dx which gives N Z Z∞ h    i 2 2 2 2 X an k1 − kx W1 (kx , ky ) + k2 − kx W2 (kx , ky ) × n=−N −∞ L L Z2 Z2 (3.32) 0 fn (x0 )e−jkx x dx0 fm (x)ejkx x dxdkx dky −L 2 −L 2 = bm , m = −N, −N + 1, ..., N. where L L Z2 Z2 bm = − fm (x)H(x)dx − I fm (x)δ(x − x0 )dx. (3.33) −L 2 −L 2 It will be now be shown that the result in (3.32) can be conveniently expressed as a system of linear equations. First, a new quantity Tn defined as L Z2 Tn = fn (x)ejkx x dx (3.34) −L2 is considered. Recalling the definition of fn (x) in (3.30), (3.34) can be expanded as n∆ Z Tn = sin k0 (x − [n − 1]∆)ejkx x dx+ (n−1)∆ (3.35) (n+1)∆ Z sin k0 (−x + [n + 1]∆) ejkx x dx. n∆ Substituting u = (x − [n − 1]∆) and u = (−x + [n + 1]∆) in the first and second portions 35 of the right side of (3.35), respectively, gives Z∆ Z∆ Tn = ejkx u ejkx (n−1)∆ sin k0 udu + e−jkx u ejkx (n+1)∆ sin k0 udu. (3.36) 0 0 Evaluating the integrals in (3.36) yields jkx u ∆ jk (n−1)∆ e Tn = e x [jkx sin k0 u − k0 cos k0 u] + k02 − kx2 0 (3.37) −jkx u ∆ jk u(n+1)∆ e e x [−jkx sin k0 u − k0 cos k0 u] . k02 − kx2 0 Upon evaluation, (3.37) becomes 1 n Tn = ejkx n∆ [jkx sin k0 ∆ − k0 cos k0 ∆] + ejkx (n−1)∆ k0 2 k0 − kx2 (3.38) o +ejkx n∆ [−jkx sin k0 ∆ − k0 cos k0 ∆] + ejkx (n+1)∆ k0 . Collecting terms and simplifying further, 1 jk x n∆   Tn = e − 2k0 cos k0 ∆ + 2k0 cos kx ∆ . (3.39) k02 − kx 2 Therefore, L Z2 k0 Tn = fn (x)ejkx x dx = 2 ejkx n∆ (cos kx ∆ − cos k0 ∆) . (3.40) 2 k0 − kx 2 −L2 36 By inspection, it is seen that L Z2 0 k0 Tn∗ = fn (x0 )e−jkx x dx0 = 2 e−jkx n∆ (cos kx ∆ − cos k0 ∆). (3.41) 2 k0 − kx 2 −L2 Substituting (3.40) and (3.41) into (3.32) gives N Z Z∞ h    i k12 − kx 2 W (k , k ) + k 2 − k 2 W (k , k ) × X an 1 x y 2 x 2 x y n=−N −∞ 4k02 e−jkx (m−n)∆  2 (cos kx ∆ − cos k0 ∆)2 dkx dky k02 − kx2 (3.42) L L Z2 Z2 =− fm (x)H(x)dx − I fm (x)δ(x − x0 )dx, −L 2 −L2 m = −N, −N + 1, ..., N. Utilizing Euler’s relation and taking advantage of the fact that the integrands are even, (3.42) can be written as N ZZ∞ X an W1 (kx , ky )Smn (k1 , kx )dkx dky + n=−N 0 N ∞ (3.43) X ZZ an W2 (kx , ky )Smn (k2 , kx )dkx dky n=−N 0 = bm , m = −N, −N + 1, ..., N 37 where k02 Smn (ki , kx ) = 16(ki2 − kx 2) (cos kx ∆ − cos k0 ∆)2 cos kx (m − n)∆, 2 2 (k0 − kx ) 2 (3.44) i = 1, 2. Finally, considering (3.43), the MFIE can be expressed as a matrix equation in the form XN [Amn + Bmn ] an = bm . (3.45) n=−N where ZZ∞ Amn = W1 (kx , ky )Smn (k1 , kx )dkx dky , (3.46) 0 ZZ∞ Bmn = W2 (kx , ky )Smn (k2 , kx )dkx dky , (3.47) 0 and L L Z2 Z2 bm = − fm (x)H(x)dx − I fm (x)δ(x − x0 )dx. (3.48) −L2 −L2 3.2.2 Evaluation of Matrix Elements This section details the calculation of the entries in the matrix equation in (3.45). 38 3.2.2.1 Evaluation of Amn Examining Amn alone in more detail, Z∞ 4 1 Amn = F (k1 , kx )Smn (k1 , kx )dkx (3.49) π 2 jωµ1 0 where Z∞ J0 (ky w2 ) dk , β 2 = k 2 − k 2 . F (k, kx ) = q y x (3.50) ky2 + β 2 0 Using the relation (6.552.1) in [1] Z∞     dx βy βy J0 (xy) q = I0 K0 , (3.51) 2 2 2 2 0 x +β (3.50) becomes  q   q  w 2 2 w 2 2 F (k, kx ) = I0 kx − k K0 kx − k (3.52) 4 4 where I0 and K0 are the modified Bessel functions of the first and second kind respectively, and order zero. The integral to infinity over dkx in (3.49) can be split into two pieces about a value K as follows: ZK 4 1 Amn = F (k1 , kx )Smn (k1 , kx )dkx + π 2 jωµ1 0 Z∞ (3.53) 4 1 F (k1 , kx )Smn (k1 , kx )dkx π 2 jωµ1 K = A1mn + A2mn 39 This is useful because a value of K can be chosen such that A1mn can be computed directly in closed form. K is chosen such that K > k0 . The singularity at kx = k0 which occurs in A1mn is removable. In order to compute A2mn , Smn is considered in more detail. For the purpose of further simplification, the trigonometric portion of the expression in (3.44) is considered by itself as ρ = (cos kx ∆ − cos k0 ∆)2 cos kx (m − n)∆. (3.54) Expanding (3.54) leads to ρ = cos2 kx ∆ cos kx (m − n)∆ − cos k0 ∆ cos kx ∆ cos kx (m − n)∆ (3.55) + cos2 k0 ∆ cos kx (m − n)∆ Using the cosine square identity, (3.55) becomes   1 1 ρ= + cos 2kx ∆ cos kx (m − n)∆ − 2 cos k0 ∆ cos kx ∆ cos(m − n)∆ 2 2 (3.56) + cos2 k0 ∆ cos kx (m − n)∆. Rearranging (3.56) gives   1 2 1 ρ= + cos k0 ∆ cos kx (m − n)∆ + cos 2kx ∆ cos kx (m − n)∆− 2 2 (3.57) 2 cos k0 ∆ cos kx ∆ cos kx (m − n)∆. Finally, using the trigonometric identity 1 cos A cos B = [cos (A − B) + cos (A + B)] , (3.58) 2 40 the expression in (3.57) can be rewritten as   1 2 ρ= + cos k0 ∆ cos kx (m − n)∆ 2   1 + cos kx (m − n − 2)∆ + cos kx (m − n + 2)∆ (3.59) 4   − cos k0 ∆ cos kx (m − n − 1)∆ + cos kx (m − m + 1)∆ . An integral Z∞ 2 2 2) k0 (k − kx Ξl (k) =  2 F (k, kx ) cos kx l∆dkx (3.60) 2 2 K k0 − kx can be defined so that upon inserting ρ, and consequently the trigonometric portion of Smn into (3.53), A2mn becomes   2 4 1 1 2 Amn = + cos k0 ∆ Ξm−n + π 2 jωµ1 2  (3.61) 1   Ξ − Ξm−n+2 − cos k0 ∆ Ξm−n−1 + Ξm−n+1 . 4 m−n−2 It should be pointed out again that although there is a singularity at kx = k0 , it is not encountered since K > k0 For l > 0, the QDAWF routine is employed for evaluating the integral given the form Z∞ f (x) cos(ωx)dx. (3.62) a QDAWF uses an adaptive scheme to integrate such Fourier integrals over a semi-infinite interval [3]. It carries out the integration over a number of subintervals and invokes an extrapolation procedure to estimate the integral [4]. However, for the case where l = 0, Ξl (k) converges slowly since the cosine product is no longer present. An alternative expression of the integral that converges faster can be obtained by some mathematical manipulation. 41 Substituting in l = 0, (3.60) becomes Z∞ 2 2 2) k0 (k − kx Ξ0 (k) =  2 F (k, kx )dkx . (3.63) 2 2 K k0 − kx Furthermore, considering the case with a free-space overlay, k = k0 , and (3.63) can be rewritten as Z∞ k02 Ξ0 (k0 ) = F (k0 , kx )dkx . (3.64) k02 − kx 2 K The entire integrand above can be expressed as a function f (kx ), so that (3.64) becomes Z∞ Ξ0 (k0 ) = f (kx )dkx . (3.65) K As kx approaches infinity, the integrand above goes to zero since this is necessary for con- vergence. The approximation (9.7.5) in [5] 1 1 I0 (z)K0 (z) ∼ {1 + + ...} (3.66) 2z 8z 2 can be employed so that as kx approaches infinity, k02 k02  q   q  w w 1 f (kx ) = I0 2 2 kx − k0 K0 kx − k0 ∼ f A (kx ) = 2 2 . k02 − kx2 k02 − kx2 w q 4 4 2 2 4 kx − k0 2 (3.67) 42 The large-argument asymptotic form of the integrand in (3.64) can thus be expressed as k02 1 2 1 A f (kx ) = − = − k02  3/2 . (3.68) 2 − k2 w q kx 2 k 2 − k2 w 2 2 0 4 x 0 kx − k0 This implies that for large kx , f (kx ) ∼ O(kx −3 ). The expression for Ξ0 in (3.65) can be reconsidered and written in the form Z∞h i Z∞ Ξ0 (k0 ) = f (kx ) − f A (kx ) dkx + f A (kx )dkx . (3.69) K K This is advantageous because the integral of the difference in square brackets converges faster than the original integral of f (kx ) since the integrand decays more rapidly with kx . Additionally, the integral of f A (kx ) can be computed in closed form as follows Z∞ Z∞ A 2 2 dkx f (kx )dkx = − k0 3/2 (3.70) w  2 2 K K kx − k0 This integral is equal to Z∞ ∞   " # 2 kx 2 K f A (kx )dkx = − k02 − q  = 1− . (3.71) w 2 2 k0 kx − k0 K 2 w K 2 − k02 K Finally, substituting f (kx ) and f A (kx ) from (3.67) and (3.71) respectively into (3.69) gives Z∞   −k0 2  q   q  w 2 2 w 2 2 I(k) =    I0  kx − k0 K0 kx − k0 +  kx2 − k2 4 4 K 0  " # (3.72)  2 2 1  2 K k w 0 2 3/2  dkx + w 1 − 2 2 . K − k kx − k02  0 43 3.2.2.2 Evaluation of Bmn After substituting W2 from (3.25) into (3.47), Bmn is given by ∞ J0 (ky w2) 1+ ZZ   1 1 2 Bmn = Smn (k2 , kx )dkx dky , (3.73) 2 (2π) jωµ2 p2 e2p2 t − 1 0 which can be rewritten as Z∞ 1 1 Bmn = F (k2 , kx )Smn (k2 , kx )dkx 2 4π jωµ2 0 (3.74) Z∞ 1 1 + F̄ (k2 , kx )Smn (k2 , kx )dkx . 4π 2 jωµ2 0 Here F (k, kx ) is given by (3.52) and Z∞ J0 (ky w 2) e−2p2 t F̄ (k2 , kx ) = dky (3.75) p2 1 − e−2p2 t 0 with q p2 = kx 2 + k2 − k2, <{p2 } > 0. (3.76) y 2 The expression for Bmn can thus be considered in two parts, Bmn = Bmn 1 + B2 , (3.77) mn as done previously for Amn in (3.53) with ∞   J k w 1 = 1 1 ZZ 0 y 2 Bmn Smn (k2 , kx )dkx dky (3.78) 4π 2 jωµ2 p2 0 44 and ∞   J k w " −2p2 t # 2 = 1 1 ZZ 0 y 2 2e Bmn Smn (k2 , kx )dkx dky . (3.79) 4π 2 jωµ2 p2 1 − e−2p2 t 0 The portion of (3.79) in square brackets can be renamed Υ so that e−2p2 t Υ=2 . (3.80) 1 − e−2p2 t Rearranging,   1 Υ = 2e−2p2 t . (3.81) 1 − e−2p2 t The portion in brackets above can be written as a geometric series [6] ∞ Υ = 2e−2p2 t e−2lp2 t X (3.82) l=0 which can be rewritten as ∞ e−2lp2 t . X Υ=2 (3.83) l=1 Substituting Υ in (3.79) leads to   ∞ ∞  Z∞ J0 (ky w   2 ) e−2lp2 t dk Z   2 = 1 1 X  Bmn 2 y Smn (k2 , kx )dkx . (3.84) 4π 2 jωµ2  p2  l=1kx =0 ky =0    45 Using Z∞ q −α x2 +β 2  q  dx β e Jν (γx) q = I1 α2 − γ 2 − α 2 2 2 ν 2 0 β +x (3.85)  q  β 2 2 × K1 α −γ −α 2 ν 2 [Re α > 0, Re β > 0, γ > 0, Re ν > −1] from (6.637.1) in [1], (3.84) can be rewritten as Z∞ q  ∞ kx2 − k2 2 = 1 1 X 2 Bmn 2 I0  [δl − 2lt] 2 4π jωµ2 2 l=1 kx =0 q  (3.86) kx 2 − k2 2 × K0  [δl − 2lt] Smn (k2 , kx )dkx 2 where s w2 δl = 4l2 t2 + . (3.87) 4 Finally, a more concise way of expressing (3.86) is given by ∞ Z∞ 2 1 1 X Bmn = 2 Fl (k2 , kx )Smn (k2 , kx )dkx . (3.88) 4π 2 jωµ2 l=1 kx =0 where q  q  2 kx − k 2 2 − k2 kx Fl (k, kx ) = I0  [δl − 2lt] K0  [δl − 2lt] . (3.89) 2 2 46 Thus, Bmn can be expressed as ∞ Z∞ 1 1 X Bmn = l Fl (k2 , kx )Smn (k2 , kx )dkx , (3.90) 4π 2 jωµ2 l=0 kx =0 where l = 1 for l = 0 and l = 2 for l > 0. It should be noted that in the case where k2 is complex Smn has no singularity. On the other hand, when k2 is real, Bmn is computed in the same manner as Amn . 3.2.3 Evaluation of right-hand side elements The computation of the right-hand side in equation (3.45) is examined in more detail in the sections that follow. 3.2.3.1 Plane-wave Excitation Considering the case of plane-wave illumination, bm from (3.45) can be computed in closed form. It is assumed that the plane wave incident on the slot is y-polarized with a propagation vector in the x − z plane. This assumption is made since the incident field couples well into the slot direction in this case. For a y-polarized plane wave incident on the surface of the slot, as shown in Figure 3.1, the incident wave vector and incident electric field can be expressed as ~k i = −k x̂ sin θ − k ẑ cos θ (3.91) 0 0 0 0 and ~ i = E ŷejk0 x sin θ0 ejk0 z cos θ0 E (3.92) 0 47 Figure 3.1: Figure showing y-polarized plane-wave excitation. respectively where θ0 is the incident angle. Thus the incident magnetic field which is given by ~i ~ i H~i = k ×E (3.93) η0 is H~ i = −ẑ E0 ejk0 x sin θ0 ejk0 z cos θ0 sin θ + η0 0 (3.94) E x̂ 0 ejk0 x sin θ0 ejk0 z cos θ0 cos θ0 . η0 48 Similarly, given that ~k r is ~k r = −k x̂ sin θ + k ẑ cos θ , (3.95) 0 0 0 0 the reflected electric field is E~ r = −E ŷejk0 x sin θ0 e−jk0 z cos θ0 (3.96) 0 and the reflected magnetic field is ~r ~ r ~r = k ×E . H (3.97) η0 Substituting (3.96) into (3.97) leads to H~ r = ẑ E0 ejk0 x sin θ0 e−jk0 z cos θ0 sin θ η0 0 (3.98) E +x̂ 0 ejk0 x sin θ0 e−jk0 z cos θ0 cos θ0 . η0 The total excitation field which is a superposition of the incident and reflected fields is ~i+H H ~ r = −ẑ2j E0 ejk0 x sin θ0 sin θ sin (k z cos θ ) η0 0 0 0 (3.99) E + x̂2 0 ejk0 x sin θ0 cos θ0 cos (k0 z cos θ0 ) . η0 The total field in the plane of the aperture (z = 0) becomes H~i+H ~ r = 2x̂ E0 ejk0 x sin θ0 cos θ . (3.100) η0 0 In the case of normal incidence (θ0 = 0), H~i+H ~ r = 2x̂ E0 . (3.101) η0 49 Since L = ∆(2m + 2) from (3.29), (3.48) can be rewritten so that (m+1)∆ Z bm = − H(x)f0 (x − m∆)dx. (3.102) (m−1)∆ Here, H is given by (3.100) and f0 by (3.30). Substituting gives m∆ Z E bm = − 2 0 cos θ0 ejk0 x sin θ0 sin k0 (∆ + x − m∆)dx η0 (m−1)∆ (3.103) (m+1)∆ Z E − 2 0 cos θ0 ejk0 x sin θ0 sin k0 (∆ − x + m∆)dx. η0 m∆ If u is chosen as u = x − (m − 1)∆ and u = (m + 1)∆ − x for the first and second integrals in (3.103), respectively, then Z∆ E bm = −2 0 cos θ0 ejk0 u sin θ0 ejk0 (m−1)∆ sin θ0 sin k0 udu η0 0 . (3.104) Z∆ E −2 0 cos θ0 e−jk0 u sin θ0 ejk0 (m+1)∆ sin θ0 sin k0 udu η0 0 An alternative expression E bm = −2 0 cos θ0 ejk0 m∆ sin θ0 × η0 Z∆ h (3.105) ejk0 u sin θ0 e−jk0 ∆ sin θ0 + e−jk0 u sin θ0 ejk0 ∆ sin θ0 sin k0 u i 0 is obtained by pulling out the terms that have no dependence on the integration variable. 50 This can then be rewritten as E bm = −2 0 cos θ0 ejk0 m∆ sin θ0 × η0 Z∆ h (3.106) ejk0 (u−∆) sin θ0 + e−jk0 (u−∆) sin θ0 sin k0 udu i 0 which, after invoking Euler’s relation, gives   Z∆ E bm = −2 0 cos θ0 ejk0 m∆ sin θ0 2 cos [k0 (u − ∆) sin θ0 ] sin k0 udu . (3.107)   η0 0 After integration, (3.107) becomes E bm = −4 0 cos θ0 ejk0 m∆ sin θ0 × η0  ∆ (3.108) cos [k0 (u − ∆) sin θ0 − k0 u] cos [k0 (u − ∆) sin θ0 + k0 u] − 2 [k0 sin θ0 − k0 ] 2 [k0 sin θ0 + k0 ] 0 which reduces to E0 bm = −2 cos θ0 ejk0 m∆ sin θ0 × η0 k0   (3.109) cos[k0 ∆] cos[k0 ∆] cos[k0 ∆ sin θ0 ] cos[k0 ∆ sin θ0 ] − − + sin θ0 − 1 sin θ0 + 1 sin θ0 − 1 sin θ0 + 1 after evaluation. Pulling like terms together (3.109) can be rewritten as E0 bm = −2 cos θ0 ejk0 m∆ sin θ0 × η0 k0 " # (3.110) 2 2 cos[k0 ∆] − cos[k0 ∆ sin θ0 ] sin2 θ0 − 1 sin2 θ0 − 1 51 and finally E0 1 bm = −4 ejk0 m∆ sin θ0 [cos(k0 ∆ sin θ0 ) − cos(k0 ∆)]. (3.111) η0 k0 cos θ0 For the case of normal incidence (θ0 = 0) E0 bm = −4 [1 − cos(k0 ∆)] (3.112) η0 k0 3.2.3.2 Line-current Excitation It is necessary to examine the case of excitation with a line current as the resulting antenna problem is important for validation of the numerical solution. This will be shown in more detail later in the chapter. Considering the case of a line excitation within the slot, bm can be expressed as x Zm 1 bm = − H(u) sin k2 (xm − u)du (3.113) k2 −L 2 with H(x) defined as H(x) = H i (x, 0) + Ii(x, 0). (3.114) Since H i (x, 0) = 0 and i(x, 0) = δ(x − x0 ), x Zm I bm = − δ(u − x0 ) sin kx (xm − u)du. (3.115) k0 −L 2 52 Evaluating (3.115) which leads to   0  xm < x0 bm = . (3.116)  − I sin(k2 xm − x0 )  xm > x0 k2 For a center-fed slot, x0 = 0 and therefore   0  xm < 0 bm = . (3.117)  − I sin(k2 xm )  xm > 0 k2 3.2.4 Scattered Field Calculation 3.2.4.1 Using Sommerfeld Integral Green’s function The scattered field can be calculated using the Sommerfeld Integral Green’s function. The electric field can be described as [7] E~ = −jωµ∇ × Π ~m (3.118) where Z ~ m (~r 0 ) K ~ m (~r) = Π G(~r|~r 0 ) ds0 (3.119) jωµ SA is the magnetic Hertzian potential. Here K ~ m , the magnetic current across the slot, is K~ m (~r) = x̂2Ey (x, y) = x̂2V (x)f (y) (3.120) and Z Z∞ −pz G(~r, ~r 0 ) = 1 e ~ ρ−~ ej k·(~ ρ 0 ) dk dk (3.121) x y (2π)2 2p −∞ 53 where ~k = x̂kx + ŷky (3.122) and p2 = kx 2 + k2 − k2. y 0 (3.123) The y and z components of E ~ are ∂Πmy Ey = −jωµ (3.124) ∂z and ∂Πmy Ez = jωµ (3.125) ∂y respectively. Therefore, Z Z∞ −pz ρ 0 ) dk dk dx0 dy 0 Z 1 e ~ ρ−~ Ey = 2V (x0 )f (y 0 ) ej k·(~ x y (3.126) (2π)2 2 SA −∞ which can be rewritten as w   Z Z∞  Z2 1 0 −jk y 0 0 Ey =  f (y )e y dy   (2π)2    −∞ y 0 =− w 2   (3.127) L Z2 0 0  jk x ky y −pz   ×  0 V (x )e −jk x x dx  e x e e dkx dky .     x0 =− L 2 54 Considering the integral over x0 in (3.127) separately, a new quantity L Z2 0 I= V (x0 )e−jkx x dx0 (3.128) −L 2 can be defined. As done earlier in (3.27) – (3.30), V can be represented by a linear combi- nation of piecewise sinusoids so that (3.128) becomes N Z∆ an e−jkx (n−1)∆ sin k0 ue−jkx u du+ X I= n=−N 0 (3.129) N Z∆ an e−jkx (n+1)∆ sin k0 ue+jkx u du. X n=−N 0 After evaluating the integrals in (3.129), N X e−jkx n∆ I= an × k02 − kx2 n=−N h (3.130) − jkx sin k0 ∆ − k0 cos k0 ∆ + k0 ejkx ∆ + jkx sin k0 ∆ − k0 cos k0 ∆] + k0 e−jkx ∆ . i The integral over y 0 in (3.127), which has been encountered earlier in (3.15) is simply w Z2 0  w f (y 0 )e−jky y dky = J0 ky . (3.131) 2 −w2 55 Substituting (3.130) and (3.131) into (3.127) yields ∞ Z Z∞ −jkx n∆ 2k0 X e Ey = an [cos kx ∆ − cos k0 ∆] (2π)2 2 k02 − kx n=−N −∞ (3.132)  w × J0 ky ejkx x ejky y e−pz dkx dky 2 which can be alternatively expressed in the form N Z∞ 2k Ey = 0 X an F (kx )Sn (kx )dkx (3.133) π2 n=−N kx =0 where Z∞ q  w z 2 +k 2 −k 2 kx F (kx ) = J0 ky cos ky ye y 0 dk (3.134) y 2 ky =0 and cos kx ∆ − cos k0 ∆ Sn (kx ) = cos kx (x − n∆). (3.135) k02 − kx 2 Considering at first a case where kx < k0 with γ 2 = k02 − kx 2 > 0, Zγ q  w −jz γ 2 −ky2 F (kx ) = J0 ky cos ky ye dky 2 0 (3.136) Z∞  q w  −z ky2 −γ 2 + J0 ky cos ky ye dky . 2 γ The split in the integral is implemented to ensure the square root of a positive number is considered over each integration region. On the other hand, if kx > k0 , again letting 56 2 − k2 γ 2 = kx 0 Z∞  q w −z ky2 +γ 2 F (kx ) = J0 ky cos ky ye dky (3.137) 2 0 Invoking Euler’s relation, (3.136) can be expanded to give Zγ   w  q F (kx ) = J0 ky 2 cos ky y cos −jz γ − ky dky2 2 0 Zγ   w  q −j J0 ky cos ky y sin −jz γ 2 − ky2 dky (3.138) 2 0 Z∞  q w −z ky2 −γ 2 + J0 ky cos ky ye dky . 2 γ An expression for Ey can finally be obtained by substituting (3.137) and (3.138) into (3.133) so that 57 N 2k Ey = 0 X an × π2 n=−N q  2 2 ( Zk0  kZ0 −kx  q    w 2 − k 2 − k 2 dky  Sn (kx )dkx     J 0 ky cos k y y cos z k 0 x y   2   kx =0 ky =0 q  k −k 2 Zk0  Z0 x   q    w  2 2 2  −j  J0 ky cos ky y sin z k0 − kx − ky dky   Sn (kx )dkx   2  kx =0 ky =0    q  Zk0  Z∞ 2 2 2  − z ky +kx −k0    w  + J0 ky cos ky y e dky  Sn (kx )dkx   2    kx =0  q ky = k02 −kx 2    q   ∞ ∞ 2 2 2  − z ky +kx −k0 Z Z )   w  +  J0 ky cos ky y e dky   Sn (kx )dkx  2 kx =k0 ky =0 (3.139) An alternative approach to defining the scattered field is covered next. 3.2.4.2 Using Free Space Green’s function Using the free space Green’s function, the scattered field can be defined as Z  ~ r) = 1 n̂ × E × ∇0 ψdS 0 0  E(~ ~ (3.140) 2π SA where 1 + jkR −jkR ∇ψ = e R̂, (3.141) R2 58 with R being the distance between a field point being considered and the center of the slot. As is done conventionally, the source points are denoted by prime coordinates while the field points are denoted without the prime. The y-directed E ~ field can be expressed as Ey (x, y) = V (x)f (y) (3.142) where −w Z2 V (x) = − Ey (x, y)dy. (3.143) w 2 Substituting (3.142) in (3.140) leads to L w Z2 Z2 ~ r) = 1 ẑ × ŷV (x0 )f (y 0 ) × R̂ 1 + jkR (jkR) 0 0 h i E(~ e dx dy (3.144) 2π R2 x0 =− L y 0 =− w2 2 The distance between the field point ~r = xx̂ + y ŷ + z ẑ (3.145) and the source point r~0 = x0 x̂ + y 0 ŷ (3.146) is given by ~ = (x − x0 )x̂ + (y − y 0 )ŷ + z ẑ. R (3.147) 59 Assuming that R >> W (the width of the slot), ~ ≈ (x − x0 )x̂ + y ŷ + z ẑ R (3.148) and the distance q R≈ (x − x0 )2 + y 2 + z 2 . (3.149) Substituting R~ from (3.148) in (3.144) L W Z2 Z2 ~ r) = 1 E(~ [−x̂ × (y ŷ + z ẑ)] V (x0 )f (y 0 ) 1 + jkR (jkR) 0 0 e dx dy (3.150) 2π R3 x0 =− L 0 2 y =− 2 W which leads to L Z2 ~ r) = 1 (z ŷ − yẑ) E(~ V (x0 ) 1 + jkR −jkR 0 e dx (3.151) 2π R3 L 2 Again, as done previously in (3.27) – (3.30), V can be represented by a linear combination of piecewise sinusoids so that (3.151) becomes L N Z2 z 1 + jkR −jkR 0 f0 (x0 − n∆) X Ey = an e dx . (3.152) 2π R3 n=−N −L2 60 q Using R = (x − x0 )2 + y 2 + z 2 and inserting f0 from (3.30) N n∆ Z z  1 + jkR sin k0 x0 − (n − 1)∆ e(−jk0 R) dx0 X  Ey = an 2π R 3 n=−N (n−1)∆ (3.153) N (n+1)∆ Z z X  1 + jkR sin k0 (n − 1)∆ − x0 e(−jk0 R) dx0  + an 2π R 3 n=−N n∆ Letting u = −x0 + (n + 1)∆ (which implies that x0 = −u + (n + 1)∆) and substituting x0 in (3.153) leads to N Z∆ z X 1 + jk0 Rn (−jkRn ) Ey = an sin k0 u e du 2π Rn3 n=−N 0 (3.154) N Z∆ z 1 + jk0 R¯n −jk R¯n X  + an sin k0 e 0 du 2π ¯3 n=−N Rn 0 q q where R¯n = (x + u − [n + 1]∆)2 + y 2 + z 2 and Rn = (x − u − [n − 1]∆)2 + y 2 + z 2 . Rearranging (3.154) gives N Z∆ " # z X 1 + jk0 Rn −jk Rn 1 + jk0 R¯n −jk Rn ¯ Ey = an sin k0 u e 0 + e 0 du. (3.155) 2π n=−N Rn3 R¯n 3 0 It should be noted that an alternative way of expressing (3.155) is possible by taking advan- tage of the fact that (1 + jk0 Rn ) e−jk0 Rn = (1 + jk0 Rn ) (cos k0 Rn − j sin k0 Rn ) (3.156) 61 leading to (1 + jk0 Rn ) e−jk0 Rn = [cos k0 Rn + k0 Rn sin k0 Rn ] (3.157) + j [k0 Rn cos k0 Rn − sin k0 Rn ] . 3.3 Hallen-Type MFIE Formulation (Version 2 MFIE) Recalling that the Hallen-type MFIE obtained in the previous chapter in Section 2.3.3 from equations (2.68) through (2.71) for V (x) is L Z2 h G1 (x, x0 ) + G2 (x, x0 ) V (x0 )dx0 i −L 2   L ! Zx  Z2 k22 − k12  −  G (u, x0 )V (x 0 )dx 0   sin k (x − u)du 1 2  k22 (3.158)     −L2 −2 L Zx 1 =− H(u) sin k2 (x − u)du − C1 sin k2 x − C2 cos k2 x, k2 −L 2 L L − ≤x≤ , 2 2 the details of expanding G1 and G2 for the second approach employed are outlined. The expression for G1 from (3.21) is rewritten using the definition of p1 from (3.6) so that   Z∞ Z∞ 1 2 0  J0 (ky w2) G1 (x, x0 ) = ejkx (x−x )   dky   dkx (3.159) (2π)2 jωµ1 q ky2 + (kx2 − k2)  kx =−∞ ky =0 1 62 as this provides the advantage that the integral contained in square brackets can be expressed in closed form. The aforementioned integral can be defined separately as Z∞ J0 (ky w2) Ω(kx ) = q dky . (3.160) 2 2 ky + (kx − k1 ) 2 ky =0 For real values of k1 , the approach used to compute the integral depends on whether k1 < kx or k1 > kx . For the case where k1 < kx , using (6.552.1) from [1] Z∞ J (xy)  ay   ay  p0 dx = I0 K0 , a>0 (3.161) x 2 + a2 2 2 0 leads to  q   q  w 2 2 w 2 2 Ω(kx ) = I0 kx − k1 K0 kx − k1 . (3.162) 4 4 In order to evaluate the integral when k1 > kx , Ω(kx ) is rewritten in the form Zβ Z∞ J0 (ky w 2 ) J0 (ky w 2 ) dk Ω(kx ) = q dky + q y (3.163) j β 2 − ky2 ky2 − β 2 0 β where β 2 = k12 − kx 2 > 0. (3.164) ky Letting u = , (3.163) can be expressed as β Z1 Z∞ J0 (β w 2 u) J0 (β w u) Ω(kx ) = p βdu + p 2 βdu. (3.165) jβ 1 − u2 β u2 − 1 0 1 63 The closed form integrals, Z1 J (xy) π h  y i2 p0 dx = J (3.166) 1 − x2 2 0 2 0 and Z∞ J (x, y) π y  y  p0 dx = − J0 Y , (3.167) x2 − 1 2 2 0 2 1 from (6.552.4) and (6.552.6) in [1] are used to evaluate the integral in (3.165) so that Ω(kx ) becomes βw 2 π         π βw βw I(kx ) = −j J − J Y0 (3.168) 2 0 4 2 0 4 4       π βw βw βw = −j J0 J0 − jY0 (3.169) 2 4 4 4  q   q  π w 2 2 (2) w 2 2 = −j J0 k1 − kx H0 k1 − kx , (3.170) 2 4 4 (2) since H0 (z) = J0 (z)−jY0 (z). Using (9.6.3) from [1] and remembering that Bessel functions of even order are even functions of their argument, I0 (jz) = J0 (z) (3.171) and π (2) K0 (jz) = −j H0 (jz), (3.172) 2 64 Therefore, (3.162) is valid for all kx . Taking (3.162) and (3.170) and substituting them in for Ω(kx ), (3.159) can be expressed as Z∞ 1 4 G1 (x − x0 ) = cos[kx (x − x0 )]Ω(kx )dkx (3.173) 2 (2π) jωµ1 kx =0 where   q   q  π w 2 2 (2) w 2 2  −j 2 J0 4 k1 − kx H0 4 k1 − kx , kx < k1   Ω(kx ) =  q   q  . (3.174) w 2 2  I0 4 kx − k1 K0 4 kx − k1 ,   w 2 2 kx > k1 Similarly, G2 can be considered in more detail so that it can be expressed in a form more convenient for computation. Examining (3.22), G2 can be divided into two portions as G2 (x, x0 ) = GA 0 B 2 (x, x ) + G2 (x, x ) 0 (3.175) where Z Z∞ jkx (x−x0 )  0 1 1 e w GA 2 (x, x ) = (2π)2 jωµ J0 ky dkx dky (3.176) 2 p2 2 −∞ and Z Z∞ jkx (x−x0 ) " −p t #  0 1 1 e e 2 w GB2 (x, x ) = (2π)2 jωµ J0 ky dkx dky . (3.177) 2 p2 sinh(p2 t) 2 −∞ Following a procedure similar to that just outlined for G1 leads to Z∞ GA (x − x 0) = 1 4 cos[kx (x − x0 )]Ω̄(kx )dkx (3.178) 2 (2π)2 jωµ2 kx =0 65 where  q   q  w 2 2 w 2 2 Ω̄(kx ) = I0 kx − k2 K0 kx − k2 . (3.179) 4 4 Considering GB 2 further, Z∞ 0 1 1 0 GB2 (x − x ) = (2π)2 jωµ ejkx (x−x ) × 2 kx =−∞ (3.180) Z∞ J0 (ky w −p2 t " # 2 ) e dky dkx , p2 sinh(p2 t) ky =−∞ where as before p22 = kx2 + k2 − k2. y 2 (3.181) As done in (3.13), the ratio within square brackets can be rewritten using the definition of the sinh function so that e−p2 t e−2p2 t =2 (3.182) sinh(p2 t) 1 − 2e−2p2 t which considering the sum to infinity of a geometric series with common ratio r for −1 < r < 1, ∞ 1 rk = X , (3.183) 1−r k=0 can be represented in the form of the geometric series ∞ ∞ 2e−2p2 t e−2qp2 t = 2 e−2qp2 t . X X (3.184) q=0 q=1 66 Substituting from (3.184) into (3.180) and taking advantage of the even nature of the inte- grand over ky leads to ∞ Z∞ 0 1 2 0 GB ejkx (x−x ) × X 2 (x − x ) = (2π)2 jωµ 2 q=1 kx =−∞   (3.185) ∞   Z J ky w  0 2 −2qp2 t   √ e dky  dkx  p2  ky =0 Consider the integral q  Z∞ ky2 +β 2 t J0 (ky w2) e −2q Iq (kx ) = q dky (3.186) ky2 + β 2 0 where β 2 = kx 2 − k2. 2 (3.187) Using (6.637.1) from [1] Z∞ q −α x2 +β 2  q   q  dx β β e J0 (γx) q = I0 α 2 + γ 2 − γ K0 α2 + γ 2 − γ , 2 2 2 2 0 β +x Re{α} > 0, Re{β} > 0, γ > 0, (3.188) gives q  q  2 2 kx − k2  2 2 kx − k2    Iq (kx ) = I0  δq − 2qt  K0  δq + 2qt  (3.189) 2 2 67 where s w2 δq = 4q 2 t2 + . (3.190) 4 Going back to (3.185), GB 2 can finally be written as ∞ Z∞ B 0 1 8 X h 0 i G2 (x − x ) = cos kx (x − x ) Iq (kx )dkx . (3.191) (2π)2 jωµ2 q=1k =0 x 3.3.1 MoM Solution For the MoM solution, pulse functions are employed for the expansion and point matching is implemented. The slot voltage can be expressed as a linear combination of pulse functions so that the voltage is expanded in the form XN V (x) = an Pn (x) (3.192) n=1 where   1 xn − ∆ 2 ≤ x ≤ xn + 2 ∆  Pn (x) = (3.193)  0  elsewhere and   L 1 xn = − + n − ∆. (3.194) 2 2 68 After substituting the expanded voltage in (3.158) and point-matching at x = xm , the following system of linear equations is obtained xn + ∆ N Z 2h G1 (xm , x ) + G2 (xm , x ) dx0 0 0 X i an n=1 xn − ∆2   x + ∆ N ! xZm  nZ 2 k22 − k12  0 0 X   − an G1 (u, x )dx  sin k2 (xm − u)du  k 2  (3.195) n=1 2   −L2 x − n 2 ∆ Zx 1 =− H(u) sin k2 (x − u)du − C1 sin k2 xm − C2 cos k2 xm , k2 −L2 m = 1, 2, ..., N. The expression in (3.195) can be written as a system of N equations in N + 2 unknowns thus XN an (Amn + Bmn ) + C1 sin k2 xm + C2 cos k2 xm = bm (3.196) n=1 where xn + ∆ Z 2h G1 (xm , x0 ) + G2 (xm , x0 ) dx0 , i Amn = (3.197) xn − ∆ 2   x + ∆ ! xZm  nZ 2 k22 − k12  Bmn = −  0 0  G1 (u, x )dx   sin k2 (xm − u)du (3.198)  k22    −L2 x − n 2 ∆ 69 and Zx 1 bm = − H(u) sin k2 (x − u)du. (3.199) k2 −L2 The two additional equations required to provide a system of N + 2 equations in N + 2   unknowns are obtained by implementing the boundary condition on V (x), i.e. V − L 2 =   V L 2 = 0. It is important to note that if k1 = k2 , Bmn vanishes and dealing with (3.195) is greatly simplified. On the other hand, if k1 6= k2 Bmn is difficult to compute. However, here it is used specifically for the free space case where k1 and k2 are equal. 3.3.2 Calculation of matrix elements 3.3.2.1 Calculation of Amn The computation of Amn is considered when the slot is immersed in free space i.e. the backing material is air. For the free-space case, k1 = k2 = k0 and Amn from (3.197) can be written in the form Amn = A1mn + A2mn + A3mn (3.200) Splitting G2 into two portions as done in (3.175), Amn is written as xn + ∆ Z 2h G1 (xm , x ) + G2 (xm , x ) + G2 (xm , x ) dx0 0 0 0 i Amn = A B (3.201) xn − ∆2 = A1mn + A2mn + A3mn . (3.202) 70 Upon substituting G1 , GA B 2 and G2 from (3.173), (3.176) and (3.177) respectively, the terms in Amn are xn + ∆ Z 2 Z∞ 1 4 cos kx (xm − x ) Ω(kx )dkx dx0 , 0 h i A1mn = (3.203) 2 (2π) jωµ1 xn − ∆2 kx =0 xn + ∆ Z 2 Z∞ 1 4 cos kx (xm − x0 ) Ω̄(kx )dkx dx0 , h i A2mn = (3.204) 2 (2π) jωµ2 xn − ∆2 kx =0 and xn + ∆ Z 2 ∞ Z∞ 1 8 X cos kx (xm − x ) Iq (kx )dkx dx0 0 h i A3mn = (3.205) 2 (2π) jωµ2 ∆ q=0k =0 xn − 2 x where   q   q   π w 2 2 (2) w 2 2 −j 2 J0 4 k1 − kx H0 4 k1 − kx kx < k1   Ω(kx ) =  q   q  (3.206)  w 2 2 w 2 2 4 kx − k1 K0 4 kx − k1 kx > k1 , I0    q   q  w 2 − k2 K w 2 2 Ω̄(kx ) = I0 4 kx 2 0 4 kx − k2 (3.207) and q  q  kx2 − k2 k 2 − k2 2  x 2  Iq (kx ) = I0  δq − 2qt  K0  δq + 2qt  , (3.208) 2 2 71 with s w2 δq = 4q 2 t2 + . (3.209) 4 Inspecting (3.203), (3.204) and (3.205), it can be seen that if k1 = k2 then A1mn = A2mn . For the slot in free-space, Amn = 2A1mn . Considering A1mn alone, if the integral in (3.203) is assumed to be uniformly convergent then the integrals can be rearranged with the expression becoming   ∆ Z∞  xnZ+ 2  1 1 4  0 0  Amn = Ω(kx )  cos kx (x − xm )dx   dkx . (3.210) (2π)2 jωµ1    kx =0 xn − 2 ∆ After computing the integral in the square brackets in (3.210), Z∞ 1 4 sin kx2∆ A1mn = Ω(kx ) cos [kx (m − n)∆] dkx . (3.211) 2π 2 jωµ1 kx kx =0 Let l = m − n. Substituting Ω(kx ) from (3.206) into (3.211), k 1 4  π  Z 1 sin kx ∆  q w  1 Amn = −j 2 cos kx l∆J0 2 2 2 k1 − kx dkx + 2π 2 jωµ1 2 kx 4 0  Zk1 sin kx2∆   q  1 4 −π w 2 2 cos kx l∆J0 k1 − kx × 2π 2 jωµ1 2 kx 4 (3.212) 0  q  w Y0 k12 − kx 2 dkx + 4 Z∞ sin kx2∆  q   q  1 4 w 2 2 w 2 2 cos kx l∆I0 kx − k1 K0 kx − k1 dkx . 2π 2 jωµ1 kx 4 4 k1 72 As done in equations (3.69) through (3.72), the third integral in (3.212) can be expressed in ∞ R h i ∞ f (kx ) − f A (kx ) dkx + f A (kx )dkx where f A (kx )dkx is the large-argument R the form K K asymptotic form of the integrand. Again, this is done because the difference in the integral of the difference in square brackets converges more rapidly than the original integral. The third integrand from (3.212) can be rewritten so that sin kx2∆  q   q  w 2 2 w 2 2 cos kx l∆I0 kx − k1 K0 kx − k1 kx 4 4       1 1 1 (3.213) = sin kx l + ∆ − sin kx l − ∆ × 2kx 2 2  q   q  w 2 − k2 K w 2 2 I0 4 kx 1 0 4 kx − k1 which can be approximated as sin kx2∆  q   q  w 2 2 w 2 2 cos kx l∆I0 kx − k1 K0 kx − k1 kx 4 4       (3.214) 1 1 1 2 1 ≈ sin kx l + ∆ − sin kx l − ∆ q 2kz 2 2 w k2 − k2 x 1 as kx approaches infinity. Simplifying (3.214) further leads to sin kx2∆  q   q  w 2 2 w 2 2 cos kx l∆I0 kx − k1 K0 kx − k1 kx 4 4    1 ∆ − sin k l − 1 ∆  (3.215) 1 sin kx l + 2 x 2 ≈ . w kx2 Thus, the third integral in (3.211) can be computed by taking advantage of (3.215) and using the relation Z∞ sin(ax) dx = sin(a) − a Ci(a) (3.216) x2 1 73 where Z∞ cos t Ci(x) = − dt (3.217) t x is the cosine integral. The expression in (3.216) is obtained since it is known from (3.721.1) in [1] that Z∞ cos(ax) dx = − Ci(a) (3.218) x 1 and by employing integration by parts. Therefore, A1mn becomes k 1 4  π  Z 1 sin kx ∆  q w  1 Amn = −j 2 cos kx l∆J0 2 2 2 k1 − kx dkx 2π 2 jωµ1 2 kx 4 0  Zk1 sin kx2∆   q  1 4 −π w 2 2 + cos kx l∆J0 k1 − kx × 2π 2 jωµ1 2 kx 4  0q  w 2 2 Y0 k1 − kx dkx 4 Z∞" sin kx2∆  q  1 4 w 2 2 + cos kx l∆ I0 kx − k1 × 2π 2 jωµ1 kx 4 (3.219) k1  q  ! # w 2 − k2 − 2 K0 kx 1 dkx 4 wkx   1 4 1 ∆ + 2 cos [k1 l∆] sin k1 2π 2 jωµ1 k1 w 2      ! 1 1 − k1 l + ∆ Ci k1 l + ∆ 2 2        1 4 1 1 1 + k1 l − ∆ Ci k1 l − ∆ . 2π 2 jωµ1 k1 w 2 2 A2mn is computed in a similar fashion as A1mn . Moving on to A3mn alone, following the 74 same procedure as in (3.210) A3mn can be rewritten as ∞ Z∞ sin kx2∆ 3 1 8 X Amn = Iq (kx ) cos kx (m − n)∆dkx . (3.220) 2π 2 jωµ2 kx q=1k =0 x Again, let l = m − n. Substituting from (3.208) into (3.220) leads to k ∞ Z 2 sin kx ∆ 1 8  π  A3mn = 2 cos k l∆× X −j x 2 2π jωµ2 2 kx q=1k =0 x q  q  k22 − kx2   k 2 − k2 2 x  J0  δq − 2qt  J0  δq + 2qt  dkx 2 2 k  X ∞ Z 2 sin kx ∆ 1 8 −π 2 cos k l∆× + x 2π 2 jωµ2 2 kx q=1k =0 (3.221) x q  q  k22 − kx2   k 2 − k2 2 x  J0  δq − 2qt  Y0  δq + 2qt  dkx 2 2 ∞ Z∞ sin kx ∆ 1 8 X 2 cos k l∆× + x 2π 2 jωµ2 kx q=1k =k x 2 q  q  2 − k2 kx k 2 − k2 2  x 2  I0  δq − 2qt  K0  δq + 2qt  dkx . 2 2 3.3.2.2 Calculation of Bmn Recall from (3.195) that   x + ∆ Zm  nZ 2 x k22 − k12  0 0     Bmn = G1 u, x dx   sin k2 (xm − u) du. (3.222)  k22    L − 2 xn − 2 ∆ 75 Let a new function I be set as xn + ∆ Z 2 G1 u, x0 dx0 .   I(u) = (3.223) xn − ∆2 Inserting the expression for G1 from (3.173) xn + ∆ Z 2 Z∞ −1 8 0 I(kx )dkx dx0 h  i I(u) = cos kx u − x (3.224) 2 (2π) 2jωµ1 xn − ∆ 2 kx =0 Rearranging,   ∆ Z∞  xnZ+ 2  1 2 h  0 i 0 I(u) = −   I(kx )  cos kx x − u dx   dkx (3.225) 2π jωµ1    kx =0 ∆ xn − 2 and letting y = x0 − u   ∆ −u Z∞  xn +Z2  1 2 I(u) = −   I(kx )   cos kx ydy   dkx (3.226) 2π 2 jωµ1   kx =0 xn − ∆2 −u gives Z∞   ∆   1 2 sin kx 2 I(u) = − I(kx ) 2 cos kx (u − xn ) dkx (3.227) 2 2π jωµ1 kx kx =0 after evaluating the integral in brackets. 76 Rearranging (3.227) leads to Z∞   1 4 sin kx ∆ 2 I(u) = − I(kx ) cos (kx (u − xn )) dkx . (3.228) 2π 2 jωµ1 kx kx =0 Inserting (3.228) into (3.222) gives x Zm " Z∞ k22 − k12 sin kx ∆ # 1 4 2 cos k (u − x ) Bmn = − − I(kx ) x n k22 2π 2 jωµ1 kx −L kx =0 (3.229) 2 × sin k2 (xm − u)du. Let V be defined as x Zm V(kx ) = cos (kx (xn − u)) sin (k2 (xm − u)) du (3.230) −L 2 which becomes xm cos [kx xn + k2 xm − u(kx + k2 )] V(kx ) = 2(kx + k2 ) −L 2 xm (3.231) cos [kx xn − u − k2 (xm − u)] − 2(kx − k2 ) −L 2   after carrying out the integration. Since xn = − L 2 + n − 1∆ , 2  1 V(kx ) = cos (kx xn + k2 xm − kx xm − k2 xm ) 2(kx + k2 )       L L − cos kx xn + k2 xm + kx + k2 (3.232) 2 2      1 L L + cos (kx xm − xm )) − cos kx xm + − k2 xm + . 2(kx − k2 ) 2 2 77 Rearranging results in 1 V(kx ) = 2(kx + k2 )      L L × cos kx xn + k2 xm − kx xm − k2 xm − cos kx xn + k2 xm + kx + k2 (3.233) 2 2       1 L L − cos[kx (xn − xm )] − cos kx xn + − k 2 xm + , 2(kx − k2 ) 2 2 and simplifying leads to 1 V(kx ) = 2(kx + k2 )        L L × cos kx xn + k2 xm − kx xm − k2 xm − cos kx xn + − k2 xm + (3.234) 2 2       1 L L − cos[kx (xn − xm )] − cos kx xn + + k 2 xm + . 2(kx − k2 ) 2 2 Expanding (3.234) leads to   1 1 1 V(kx ) = cos (kx (m − n)∆) − 2 kx + k2 kx − k2          1 1 1 − cos kx n − ∆ cos kx m − ∆ 2 2 2         1 1 1 − sin kx n − ∆ sin kx m − ∆ (3.235) 2 2 kx + k2          1 1 1 + cos kx n − ∆ cos kx m − ∆ 2 2 2         1 1 1 + sin kx n − ∆ sin kx m − ∆ 2 2 kx − k2 Invoking the trigonometric identity, 1   1 1  V(kx ) == cos kx (m − n)∆ − 2 kx + k2 kx − k2   1  1   L  1 1 − cos kx (n − )∆ cos k2 (m − )∆ − 2 2 2 kx + k2 kx − k2   1  1   L  1 1 + sin kx (n − )∆ sin k2 (m − )∆ + (3.236) 2 2 2 kx + k2 kx − k2 78 3.3.3 Calculation of Right Hand-Side Elements 3.3.3.1 Plane Wave Incidence (Scattering Problem) Considering plane wave excitation as depicted in Figure 3.1, the total excitation field is given by ~ r = −ẑ2j E0 e(jkx sin θ0 ) sin θ sin k z cos θ   H~i+H η0 0 0 0 E   +x̂2 0 e(jk0 x sin θ0 ) cos θ0 cos k0 z cos θ0 (3.237) η0 In the plane of the aperture where z = 0, ~i+H H ~ r = 2x̂ E0 cos θ e(jk0 x sin θ0 ) (3.238) η0 0 which in the case of normal incidence, is ~i+H H ~ r = 2x̂ E0 . (3.239) η0 Recalling that x Zm 1 bm = H(u) sin k2 (k2 (xm − u)) du (3.240) k2 −L2 gives for normal incidence x Zm 1 E bm = 2 0 cos θ0 e(jk0 u sin θ0 ) sin k2 (xm − u)du (3.241) k2 η0 −L2 79 after substituting H(u). Setting a variable w = xm − u, leads to x Zm E0 bm = −2 cos θ0 eu[jk0 sin θ0 ] sin k2 (xm − u)du, (3.242) η0 k2 −L2 or   m− 12 ∆ Z E cos θ0 bm = −2 0 e(xm jk0 sin θ0 ) e(−wjk0 sin θ0 ) sin k2 wdw, (3.243) η0 k2 0   L 1 given that xn = − 2 + n − 2 ∆. Evaluating the integral in (3.243) leads to E cos θ0 jk xm sin θ bm = −2 0 e 0 0 η0 k2 (3.244)   ( −jk u sin θ m− 1 ∆) e 0 0 2 × [−jk0 sin θ0 sin k2 u − k2 cos k2 u] 2 2 k2 − k0 sin θ02 0 which, after inserting the limits, gives E cos θ0 jk xm sin θ bm = −2 0 e 0 0 η0 k2 (3.245) e−jk0 u sin θ0 ( ) × [−jk0 sin θ0 sin k2 (m − 1)∆ − k2 cos k2 (m − 1)∆] . k22 − k02 sin2 θ0 3.3.3.2 Line-Current Excitation (Radiation Problem) As mentioned earlier, the radiation problem is also considered because it is employed later for validation purposes. In doing this, a line source is placed at the center of the slot and the resulting field is calculated. Figure 3.2 illustrates the problem. 80 Figure 3.2: Figure showing line-current excitation in the plane of the slot. 81 Recall that x Zm 1 bm = − H(u) sin k2 (xm − u)du (3.246) k2 −L 2 Since H(x) = H i (x, 0) + Ii(x, 0) (3.247) and i(x, 0) = δ(x − x0 ) (3.248) (3.246) becomes x Zm I bm = − δ(u − x0 ) sin k2 (xm − u)du (3.249) k2 −L 2 which implies that   0,  xm < x0 bm = . (3.250)  − I sin k2 (xm − x0 ),  xm > x0 k2 3.3.4 Calculation of Scattered Field for Plane Wave Incidence 3.3.4.1 Using Free Space Green’s Function Consider Figure 3.2. The electric field can be expressed as [7] Z  ~ r) = 1 n̂0 × E ~ × ∇0 ψds0  E(~ (3.251) 2π SA 82 where 1 + jkR −jkR ∇0 ψ = e R̂ (3.252) R2 with R being the distance between the source and field points and ψ some field quantity. Recalling from (3.142) that Ey (x, y) = V (x)f (y), (3.253) results in L w Z2 Z2 ~ r) = 1 0 0 1 + jkR −jkR dkx0 dy 0 . h i E(~ ẑ × ŷV (x )f (y ) × R̂ e (3.254) 2π R 2 x0 =− L y 0 =− w2 2 If it is assumed that the radial distance from the field point to the center of the slot is considerably larger than the width of the slot (R  w), then ~r = x̂x + ŷy + ẑz, (3.255) and ~r0 ≈ x̂x0 . (3.256) ~ being approximated as This results in R R~ = ~r − ~r 0 ≈ (x − x0 )x̂ + y ŷ + z ẑ (3.257) 83 so that q R≈ (x − x0 )2 + y 2 + z 2 . (3.258) Substituting R into (3.254) gives L w Z2 Z2 ~ r) = 1 1 + jkR −jkR 0 0 E(~ [−x̂ × (y ŷ + z ẑ)] V (x0 )f (y 0 ) e dx dy (3.259) 2π R3 x0 =− L y 0 =− w2 2 which becomes L Z2 ~ r) = 1 = (z ŷ − yẑ) 1 + jkR −jkR 0 E(~ V (x0 ) e dx (3.260) 2π R3 −L 2 after considering that the integral over y 0 is unity. Choosing the voltage V (x0 ) as a linear combination of pulse functions as done earlier, the y-component of the E-field is xn + ∆ N Z 2 z X 1 + jkR −jkR 0 Ey = an e dx . (3.261) 2π R3 n=1 xn − ∆ 2 After replacing the exponential function using Euler’s relation, Ey becomes xn + ∆ N Z 2 z X 1 Ey = an [cos kR − j sin kR + jkR cos kR + kR sin kR] dx0 (3.262) 2π R3 n=1 xn − ∆ 2 84 and separating real an imaginary parts, N z X Ey = an 2π  n=1  ∆ ∆  xnZ+ 2 xn + 2  (3.263) − Z  kR sin kR + cos kR 0 kR cos kR sin kR 0  × dx + j dx .   R3 R3  ∆ xn − 2 ∆ xn − 2 3.4 Considerations for Numerical Accuracy It is important to understand the accuracy of the solution for the scattered field. This makes it possible for error analysis to be carried out to determine the optimal slot design. There are several variables to be considered so far as the numerical accuracy of the results are concerned. First, the numerical integration involved in computing the MoM entries must be done with prescribed accuracy. Second, the number of partitions in the MoM solution must be set appropriately. Lastly, the number of terms used in the expansion of the Green’s function must also be carefully chosen. Investigations are carried out to assess the effects of the number of partitions and the number of terms in the Green’s function expansion on the accuracy of the MoM solution. A slot with dimensions 20 mm by 1 mm cut into a metallic screen above a layer of the commercial MagRAM FGM-125 (r = 7.32−j0.0464 and µr = 0.576−j0.484) in the normal incidence case is considered at 2.88 GHz. In both cases, the parameter numbers required for stabilization are examined. As used here, “stabilization” denotes when the difference between successive terms is less than the accuracy of the network analyzer used for measurements (an Agilent 8510C VNA) which is 0.04 dB and 2◦ for the magnitude and phase, respectively. Figures 3.3 through 3.6 show the results. From the figures, it can be seen that the scattered field satisfies the required criteria when the number of Green’s functions terms and number of slot partitions used in computation are 25 and 41, respectively. The accuracy attained at this point for the scattered field is comparable to that of the measurement system used in 85 the laboratory. Both numbers are used in subsequent investigations into the error analysis described later. Similar results have been obtained with different slot dimensions. It is important to note that it is possible to change the parameters to suit computations with different slot dimensions. 86 Figure 3.3: Magnitude of the scattered field computed for a slot above a conductor-backed layer of FGM-125 as a function of the number of terms used in the Greens function expansion. Number of partitions is N = 41. 87 Figure 3.4: Phase of the scattered field computed for a slot above a conductor-backed layer of FGM-125 as a function of the number of terms used in the Greens function expansion. Number of partitions is N = 41. 88 Figure 3.5: Magnitude of the scattered field computed for a slot above a conductor-backed layer of FGM-125 as a function of the Number of partitions. Number of terms used in the Greens function expansion is 25. 89 Figure 3.6: Phase of the scattered field computed for a slot above a conductor-backed layer of FGM-125 as a function of the number of partitions. Number of terms used in the Greens function expansion is 25. 90 3.5 Solution Validation In order to validate the solutions outlined in the preceding sections, the results must be compared against known numerical results. No results were found in the literature for comparison for the case of the field scattered by a slot in an aperture screen above a material layer. However, the literature search yielded references to compare against the radiation of a slot antenna. Assuming a line source excitation, the input impedance of the slot antenna can be compared to published results. The results in [8] allow for comparison against slot antenna radiation above a dielectric half-space and above a conductor-backed dielectric. The slot antenna radiation characteristics are considered in free-space, above a half-space of dielectric material and above a conductor-backed free-space region. Given a line source with amplitude I placed at a position x = x0 , as shown in Figure 3.2, the input impedance can be found once the voltage across the slot has been computed as XN V (x) = an f (x − n∆). (3.264) n=−N The input impedance, which can be expressed as Zin = R + jX, (3.265) is given by N P an f (x − n∆) n=−N Zin = . (3.266) I 91 3.5.1 Slot Antenna above a Half-Space Material Layer and above a Conductor- backed Material Layer The input impedance of a center-fed slot antenna above a half-space of dielectric material is computed. The material considered has a permittivity of r = 2.55. The results obtained, shown in Figure 3.7, are compared with the work done by Kominami et al [9]. The slot impedance plot in Figure 4 of [9] agrees closely with that in Figure 3.7. The discrepancy in the results is likely attributable to the use of EDB (entire-domain basis) functions by the authors of the paper for the representation of the slot voltage, whereas PWS (piecewise sinusoidal) functions are used in this work. The case of a slot above a conductor-backed layer of the same dielectric material (r = 2.55) is also shown in Figure 3.7, although this case is not considered by the authors. The thickness of the material layer considered is t = 0.25λ0 . It is instructive to note that the input resistance is considerably higher in the vicinity of the resonance and exhibits a narrower Q-curve for the conductor-backed case. This suggests that the inclusion of the conductor backing accentuates the resonance. Finally, a slot antenna placed above a layer of MagRAM FGM-125 is examined and the input impedance is computed as a function of slot normalized length. First, a case of a half-space of FGM-125 is considered and then a case with a 1/8 inch thick layer of FGM-125 is also considered. The results obtained are shown in Figure 3.8. The Q of the conductor backed antenna is lower than that of the half-space case. This is due to dissipation of energy in the material. The bandwidth of the antenna is quite significant when compared against the lossless case in Figure 3.7 and this bodes well for wideband characterization of absorbing materials. 92 Figure 3.7: Input impedance of a slot antenna located in a conducting screen above a half space with  = 2.550 and µ = µ0 and above a conductor-backed layer with the same properties. The slot length is L and the slot width is w = 0.02L. The layer thickness for conductor-backed case is t = 0.25λ0 . N = 41 piecewise sinusoids were used to represent the slot voltage. 93 Figure 3.8: Input impedance of a slot antenna located in a conducting screen above a half space of FGM-125 with  = (7.32−j0.0464)0 and µ = (0.576−j0.484)µ0 above a conductor- backed layer with the same properties. Slot length is L, slot width is w = 0.002L and layer thickness for conductor-backed case is t = 0.125 inches. N = 41 piecewise sinusoids were used to represent the slot voltage. 94 3.5.2 Slot Antenna above Free Space Layer with Varying Thicknesses The effect of the thickness of the material layer on the slot impedance is investigated by varying the thickness, t, of a material layer underneath a conducting screen with the mate- rial set as free space. This is compared to work published by Ostner [8]. The results, shown in Figure 3.9, are consistent with expectations. The antenna becomes shorted out as the ground plane is moved close to the slot. On the other hand, as the distance between the ground plane and the slot is increased, the radiation resistance rises until it even exceeds that of the unbacked case at about 0.3λ. This effect may be accounted for by ducting of energy away from the slot in the form of TEM guided wave modes [8]. This is possibly advantageous for material characterization since it is useful in interrogating the underlying material. There is good agreement between the results generated from the solution and Os- tner’s results although there are some slight discrepancies. It is possible that the differences can be accounted for by differences in the feed models and the computational techniques employed. 95 Figure 3.9: Input impedance of a slot antenna located in a conducting screen with a con- ducting plate located a distance t behind. The region between the slot and the plate is free space. The slot length is L = 50 mm and the slot width isw = 1 mm. N = 41 piecewise sinusoids were used to represent the voltage. 96 3.6 Summary In this chapter, numerical results are given for solution of the MFIE presented in Chapter 2. The results have been shown to be consistent with theoretical expectations and have been validated against results published in the literature for slot antenna radiation. The material parameter extraction procedure is examined in detail in the next chapter. 97 BIBLIOGRAPHY [1] I. S. Gradshteyn and I. M. Ryzhik. Table of integrals, series, and products. Else- vier/Academic Press, Amsterdam, seventh edition, 2007. Translated from the Russian, Translation edited and with a preface by Alan Jeffrey and Daniel Zwillinger. [2] R. Harrington. Origin and development of the method of moments for field computation. IEEE Antennas and Propagation Magazine, 32(3):31–35, 1990. [3] IMSL, Inc. The IMSL libraries: FORTRAN subroutines for mathematics and statistics. Math/Library for mathematical applications, Stat/Library for statistical analysis, Sfun/ Library. IMSL, Houston, TX, USA, 1989. [4] Jacob Coetzee. Efficient spectraldomain analysis of open asymmetrical multiconductor transmission lines. Microwave and Optical Technology Letters, 8:206 – 209, 03 1995. [5] Milton Abramowitz and Irene A. Stegun. Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables. Dover, New York, ninth dover printing, tenth gpo printing edition, 1964. [6] George Arfken. Mathematical Methods for Physicists. Academic Press, Inc., San Diego, third edition, 1985. [7] Constantine A Balanis. Antenna theory: analysis and design. Wiley-Interscience, 2005. [8] H. Ostner and E.M. Biebl. Planar slot antennas backed by a ground plane. In Proceedings of IEEE Antennas and Propagation Society International Symposium, pages 612–615 vol.2, 1993. [9] M. Kominami, D. Pozar, and D. Schaubert. Dipole and slot elements and arrays on semi-infinite substrates. IEEE Transactions on Antennas and Propagation, 33(6):600– 607, 1985. 98 CHAPTER 4 Extraction of Material Parameters 4.1 Introduction The feasibility of the extraction procedure for the proposed technique hinges on the assump- tion that the slot provides a perturbation that produces an observable difference between the conductor-backed measurement and the slotted conductor measurement. A measurement technique that provides “enough information” to extract  and µ is desirable. In this chap- ter, investigations into whether the scattered field has sufficient dependence on the material parameters are presented. The specifics of the extraction procedure are also detailed in this chapter. First, there is a short discussion on peculiarities in calibration for this technique and the need for particular care in the calibration approach. 4.2 Considerations for Calibration Consider a plane wave incident at an angle θ0 on a slotted conductor above a conductor- backed material sample as shown in Figure 4.1. The following parameters are chosen: an incident angle θ0 = 45 degrees, a slot length L = 20 mm, a slot width w = 1 mm, a layer thickness t = 0.125 inches, a frequency of 2.88 GHz and an observation point 50 cm from the center of the slot. N = 100 piecewise sinusoids are used to represent the slot voltage in the Galerkin’s method solution for the scattered field. An incident field strength of E0 = 1 V/m is assumed. The sample is assumed to be MagRAM FGM-125 which has representative parameters r = 7.319669 − j0.046408 and µr = 0.575582 − j0.484231. Figures 4.2 and 4.3 show the computed scattered field in amplitude and phase respectively vs. scattering angle θ measured from the normal to the slot for both FGM-125 and free-space [1]. For the dimensions considered here, the amplitude of the scattered field in both cases are quite close. The small size of the slot accounts for this. The phase values, on the other hand, are markedly different. Calibration is required before the implementation of the proposed technique and per- formance of measurements on the conductor-backed material. In free-space measurement 99 Figure 4.1: Figure showing plane wave incidence on slotted conductor surface. techniques, it is typical to carry out calibration using a conducting plate [2]. However, such an approach is not applicable here because the field incident on a shorting plate is only re- flected in the specular direction. This specular reflection makes it unfeasible to measure the response in the backscatter direction. An alternative approach involves calibration against a slotted conductor in free-space. This is achieved by computing the ratio between the scat- tered field amplitude for a slot above a conductor backed layer of material and the scattered field for a slot in a conductor in free space. Figures 4.4 and 4.5 show the computed scattered ratios in amplitude and phase for the same slot dimensions and parameters above for FGM- 125 and free-space. The figures show that the ratios remains fairly flat despite changes in the scattering angle. This implies that there is some flexibility in the choice of the incident 100 Figure 4.2: Scattered field amplitude for a slot above a conductor-backed layer of FGM-125 and for a slot in a screen immersed in free space. Incident field strength is 1 V/m. Incident field angle is θ0 = 45◦ . Slot length is L = 20 mm , slot width is w = 1 mm, and layer thickness for conductor-backed case is t = 0.125 in. N = 100 piecewise sinusoids were used to represent the slot voltage. Observation distance is R = 50 cm and frequency is f = 2.88 GHz. angle. However, it may be desirable to avoid choosing the specular angle because of the dominance of the specular reflection at that angle. Angles close to grazing should also be avoided because of the weakened signal strength and diffraction effects. An incident angle close to θ0 = 45◦ is chosen and an observation angle at θ = 45◦ is similarly chosen. 4.3 Effect of  and µ Variation on the Scattered Field The viability of the proposed technique can be assessed by varying notional values of  and µ and observing what effect the variations have on the scattered field. Figures 4.4 and 4.5 show results obtained when both  and µ are varied by 10%. As can be observed, there is a reduction in the scattered ratio of about 5 dB and an increase in the phase of about 101 Figure 4.3: Scattered field phase for a slot above a conductor-backed layer of FGM-125 and for a slot in a screen immersed in free space. Incident field strength is 1 V/m. Incident field angle is θ0 = 45◦ . Slot length is L = 20 mm , slot width is w = 1 mm, and layer thickness for conductor-backed case is t = 0.125 in. N = 100 piecewise sinusoids were used to represent the slot voltage. Observation distance is R = 50 cm and frequency is f = 2.88 GHz. 102 Figure 4.4: Ratio of scattered field amplitudes for a slot above a conductor-backed layer of FGM-125 to the scattered field for a slot in a screen immersed in free space. Incident field strength is 1 V/m. Incident field angle is θ0 = 45◦ . Slot length is L = 20 mm, slot width is w = 1 mm, and layer thickness for conductor-backed case is t = 0.125 in. N = 100 piecewise sinusoids were used to represent the slot voltage. Observation distance is r = 50 cm and frequency is f = 2.88 GHz. 4.5◦ . These are measurable changes that suggest that the proposed technique is viable for extracting the constitutive parameters of a material sample. 103 Figure 4.5: Ratio of scattered field phases for a slot above a conductor-backed layer of FGM- 125 to the scattered field for a slot in a screen immersed in free space. Incident field strength is 1 V/m. Incident field angle is θ0 = 45◦ . Slot length is L = 20 mm, slot width is w = 1 mm, and layer thickness for conductor-backed case is t = 0.125 in. N = 100 piecewise sinusoids were used to represent the slot voltage. Observation distance is r = 50 cm and frequency is f = 2.88 GHz. 104 Figure 4.6: Plane wave incident on a conductor-backed material sample. 4.4 Constitutive Parameter Extraction 4.4.1 Extracting  and µ from S-Parameter Measurements Consider Figure 4.6 where a plane wave is incident on a conductor-backed layer of material. The propagation factor P is given by P = e−jkz ∆ , (4.1) 105 where kz2 = k02 − k 2 sin2 θ with k 2 = ω 2 µ. The global reflection coefficient can be expressed for perpendicular and parallel polarization as Γ − P2 R⊥ = ⊥ , (4.2) 1 − Γ⊥ P 2 and Γk − P 2 Rk = , (4.3) 1 − Γk P 2 respectively. Here, the perpendicular interfacial reflection coefficient Γ⊥ is given by Z − Z0⊥ Γ⊥ = ⊥ (4.4) Z⊥ + Z0⊥ where kη Z⊥ = (4.5) kz is the wave impedance in the material region while η0 Z0⊥ = (4.6) cos θ √ is the wave impedance in free space. Here, η = µ is the intrinsic impedance of the medium √ and η0 = µ0 0 is the intrinsic impedance of free space. Similarly, the parallel interfacial reflection coefficient is given by Zk − Z0k Γk = (4.7) Zk + Z0k 106 where kz η Zk = (4.8) k is the wave impedance for parallel polarization within the material and Z0k = η0 cos θ (4.9) is the free-space wave impedance. Substituting (4.5) and (4.6) into (4.4) gives kη cos θ − kz η0 Γ⊥ = (4.10) kη cos θ + kz η0 while substituting (4.8) and (4.9) into (4.7) gives kz η − kη0 cos θ Γk = . (4.11) kz η + kη0 cos θ In order to extract  and µ, two distinct measurements are required: 1) a reflection measurement of a conductor-backed material without a slot, and 2) a reflection measurement of a conductor-baked material with a slot. The two different reflection measurements are denoted by S11o and S s respectively. For the first measurement scenario, the global reflection 11 coefficient can be expressed as o = Γo − P 2 S11 , (4.12) 1 − Γo P 2 from (4.12), Γo can be determined as o + P2 S11 o Γ = 1 + S11o P2. (4.13) This demonstrates that the interfacial reflection coefficient can be computed from the mea- 107 sured global reflection coefficient. Consider the case of perpendicular polarization. The expression in (4.10) can be rear- ranged to give η 1 + Γ⊥ kz = . (4.14) η0 1 − Γ⊥ k cos θ Similarly for the case of parallel polarization, (4.11) can be rearranged so that η 1 + Γk k cos θ = . (4.15) η0 1 − Γk kz p Since η/η0 = µr /r , the ratio of the permeability and permittivity can be computed if o is measured. Recalling from k 2 = ω 2 µ and k 2 = k 2 − k 2 sin θ that only products of  S11 z 0 r and µr appear in k and kz , the ratio µ/r can be written in terms of measured S-parameters and the product µr r . The theoretical expressions for the scattered field in the second measurement scenario (with the slot present) are also functions of ratios and products of r and µr . The computed ratios can be substituted so that only the product remains unknown. Defining a relative wave number kr = k/k0 , S11 o can be used to eliminate η = η/η and r 0 µr = kr ηr . This leaves the complex permittivity, r , as the only remaining unknown. It can be determined by solving the equation T (k ) − S s = 0 S11 (4.16) r 11 where S11s is the measured reflection coefficient with the slot in place as mentioned earlier T (k ) is the theoretical global reflection coefficient with the slot in place. Once k and S11 r r is found the complex permittivity and permeability can be computed from  = kr /ηr and µ = kr ηr respectively. As will be discussed in the next chapter, the equation in (4.16) is solved by utilizing a root solver. 108 4.5 Summary In this chapter, special calibration considerations that arise due to the nature of the proposed slotted conductor technique have been highlighted. Based on results presented, the technique appears to be able to capture changes in the complex constitutive parameters. Additionally, the extraction approach following S-parameter measurements has been described in detail. In the next chapter, sensitivity of the extraction procedure is considered and verified. 109 BIBLIOGRAPHY [1] K.J. Vinoy and R.M. Jha. Radar Absorbing Materials: From Theory to Design and Characterization. Springer US, 1996. [2] Philip G. Bartley and Shelley B. Begley. A new free-space calibration technique for materials measurement. In 2012 IEEE International Instrumentation and Measurement Technology Conference Proceedings, pages 47–51, 2012. 110 CHAPTER 5 Aperture Optimization 5.1 Introduction Before an experimental implementation of the proposed technique was undertaken, investi- gations were carried out to understand the robustness of the technique. The sections that follow detail error analysis performed to assess the sensitivity of the technique to uncer- tainties in the values of the measured S-parameters. Following that, assessments of the performance of aperture screen technique are presented. Finally, optimizations on the aper- ture screen parameters are carried out. This is important as the aperture screen design must be carefully chosen to facilitate the extraction process. Ultimately, there is a need to identify optimal slot parameter combinations that ensure that the constitutive parameters are accurately extracted when the aperture screen technique is implemented in the field. 111 5.2 Inversion Problem Solution As covered in the previous chapter, the material constitutive parameters ( and µ) can be extracted by solving the transcendental equation T (k ) − S s = 0. S11 (5.1) r 11 s is the measured reflection with the aperture screen in place and S T (k ) is the Here, S11 11 r T is obtained by taking theoretical global reflection coefficient with the screen in place. S11 the ratio between the theoretical scattered field with the aperture screen in place, Ey , and the incident field amplitude, E0 . The equation in (5.1) can be solved by utilizing a root search approach. However, a solution may be difficult to obtain in the presence of noise or measurement error when using the root search approach. A solution can alternatively be found by minimizing the following error E = S11 T (k ) − S s 2 . = 0. (5.2) r 11 An inversion code was written in MATLAB for implementing the minimization. A solution is found given a set of design parameters that is, for a given aperture length and width the complex permittivity and permeability of a material are extracted. The forward problem is implemented by making calls to an executable generated after compiling Fortran code for the numerical solutions. The hybrid MATLAB-Fortran approach makes it possible for the forward problem to be computed quickly while allowing for access to the various solvers available in MATLAB for the minimization. Newton’s method is used to perform the mini- mization. The code is validated by “going around in a circle.” S parameters are generated for assumed values of  and µ using the forward problem and then inserted into the inverse problem to get back the original values of the constitutive parameters. 112 Error Analysis Parameters Frequency, f (GHz) 2.88 Slot length, l (mm) 20 Slot width, w (mm) 1 Material thickness, t (mm) 3.175 Incidence Angle (◦ ) 90 Observation Angle ( ) ◦ 90 Observation Distance (mm) 500 Table 5.1: Parameters utilized for Error Analysis. 5.3 Error Analysis In order to assess the robustness of the proposed technique, the material extraction can be tested by investigating the sensitivity of the resulting extracted constitutive parameters to uncertainties in the S-parameter measurement. Two approaches were explored for this sensitivity analysis, namely the Monte Carlo method and the standard error propagation method. The aim of this analysis is to provide some insight into how noise propagates through the inversion process. For a non-robust technique, the error would be greatly amplified. Unless otherwise stated, the parameters shown in Figure 5.1 are considered for the error analysis in the later sections. The Monte Carlo method and the standard error propagation method are discussed briefly before delving into the error analysis results. 5.3.1 Monte Carlo Method The Monte Carlo method is the primary approach considered for the error analysis because of its simplicity and the insight it provides into the statistical behavior of processes. Monte- Carlo analysis is used for assessing the viability of processes with the with evaluation of a large number of tests. Example of scenarios where Monte Carlo type analyses have been employed with great success include design of robots, tuning of antennas, design of mechanical fixtures and optimization of chemical plants [1]. In general, Monte-Carlo analyses find application in situations where the execution of multiple experiments is not physically feasible or may be exorbitantly expensive. The ability to carry out multiple simulations in relatively quick fashion may be more practical. This 113 makes it possible to assess the behavior of processes in ways that would have been impossi- ble before widespread improvements and access to computational resources allowed for the running of a great number of numerical tests [1]. For the purposes here, Monte Carlo sim- ulation provides insight into the robustness of the proposed technique. The simulations are carried out by perturbing S-parameter measurements with Gaussian noise and observing the propagation of error through the extraction process. 5.3.2 Standard Error Propagation Method The standard error propagation method is implemented as an alternative to the Monte Carlo Method. The error propagation method involves the computation of derivatives and provides the advantage of quicker computation [2]. The mathematical details of the approach are covered in Appendix 7.1. In implementing this approach, the Cauchy-Reimann equations can also be leveraged to reduce the number of required computations. The S-parameters can be expressed in polar form as S = Aejφ where A is the amplitude and φ is the phase, thus ∂00 0 r = A ∂r (5.3) ∂φ ∂A and ∂0r ∂00 = −A r . (5.4) ∂φ ∂A Since the complex permittivity and permittivity are extracted from the inversion proce- dure, only one derivative needs to be computed with respect to the amplitude for complex  and another for complex µ. The derivatives with respect to the phase can be obtained from the Cauchy-Reimann relations in (5.3) and (5.4). 5.3.3 Monte Carlo Simulation Results A 125 mil thick conductor-backed layer of FGM-125 with constitutive parameters r = 7.319099 − j0.046408 and µr = 0.575582 − j0.484231 is considered at 2.88 GHz. The forward problem for a layer without the aperture screen is used to obtain the reflection S-parameter 114 for the case of normal incidence. The forward problem for a layer with the aperture screen (dimensions: 20 mm by 1 mm) in place is also used to obtain reflection S-parameters for the normal incidence case. As covered in the previous chapter, in practice it is disadvantageous to chose the specular angle as the observation angle due to the dominance of the specular reflection. However, from the numerical calculations the scattered field can be computed as a standalone quantity so this problem is avoided here. Random Gaussian noise is added to the reflection S-parameters. The noise parameters are chosen to correspond to worst-case uncertainty values for an Agilent 8510C network analyzer which are 0.02 dB in amplitude and 1◦ in phase [3]. The noisy S-parameters generated are taken as “measured” data from which  and µ are extracted. The extraction process is performed N = 100 times and the resulting statistics are generated. The standard deviation of the extracted parameters corresponds to the propagated error. It should be noted that a larger value for N leads to greater confidence in the statistical distribution obtained. The standard deviation obtained for  is 0.136 + j0.357 while that for µ is 0.0155 + j0.0361. These values correspond to 2%, 3% and 7% errors in 0 , µ0 and µ00 , respectively. There is a very large value in the error of 00 but this is likely attributable to the fact that the extracted value is very small and cannot be properly captured by the proposed technique. In Figures 5.1 through 5.4, histograms showing the real and imaginary parts of the extracted parameters are presented. Although a larger choice for the value N will lead to even more statistically meaningful results, the results show that the extracted results form a Gaussian- like distribution as expected. Examining the mean values of  and µ (7.30 − j0.253 and 0.561 − j0.500 respectively) reveal that they are close to their nominal values but suggests the need for a larger number of trials. The appropriate number of trials for the analysis is determined by plotting the extracted parameters against increasing values of N. The results which are presented in Figures 5.5 through 5.8 show the errors in the constitutive parameters normalized against their nominal values. It can be seen that the results become fairly stable beyond N = 500 trials. Since 115 these results are normalized, they also represent an estimate of the percent error in the extracted parameters due to VNA uncertainty. It can be seen that the real and imaginary part of µ exhibit about a 5% and 7% error respectively while the real part of  exhibits error of about 3%. The imaginary part of  experiences an extremely large error of 1300% which, as mentioned earlier, is due to the very small nominal value. 5.3.3.1 Comparison against Two-Thickness Techniques The results obtained for the propagated error using the aperture technique seem quite promis- ing. However, it is important to compare these results against another material extraction technique to properly assess the performance of the proposed technique. The two-thickness method is considered for this comparison. It is chosen since it is also a reflection-only technique involving the presence of a conductor backing as is the case with the proposed technique. The two-thickness technique is described in detail in [4]. A similar Monte Carlo error analysis is carried out using the two thickness method and the ratios of the standard deviations of the constitutive parameters generated from the aperture technique and the two-thickness are computed and plotted. The resulting ratios are shown in Figures 5.9 through 5.12 for 0 , 00 , µ0 and µ00 . A ratio less than unity indicates that the proposed technique performs better then the two-thickness method. This is the case for all four constitutive parameters. The results thus demonstrate that the proposed technique is indeed promising for accurate extraction of the parameters. The two thicknesses considered for the two-thickness implementation are 125 mils and 250 mils (that is, one layer of FGM-125 for the first thickness and two layers of FGM-125 for the second thickness). 5.3.4 Monte Carlo-Standard Error Propagation Comparison Table 5.2 shows results from a comparison of the Monte Carlo method and the error propa- gation method for the proposed technique. Again, a 3.175 mm thick layer of FGM 125 with the following parameters is used: slot length = 20 mm, slot width = 1 mm, frequency = 2.88 GHz, r = 7.319099 − j0.046408 and µr = 0.575582 − j0.484231. The uncertainty in the 116 Figure 5.1: Real part of  extracted from noisy data for N = 100 Monte Carlo trials. 117 Figure 5.2: Imaginary part of  extracted from noisy data for N = 100 Monte Carlo trials. 118 Figure 5.3: Real part of µ extracted from noisy data for N = 100 Monte Carlo trials. 119 Figure 5.4: Imaginary part of µ extracted from noisy data for N = 100 Monte Carlo trials. 120 Figure 5.5: Error in the real part of  normalized to the nominal value of the real part of , as a function of the number of Monte Carlo trials. 121 Figure 5.6: Error in the imaginary part of  normalized to the nominal value of the imaginary part of , as a function of the number of Monte Carlo trials. 122 Figure 5.7: Error in the real part of µ normalized to the nominal value of the real part of µ, as a function of the number of Monte Carlo trials. 123 Figure 5.8: Error in the imaginary part of µ normalized to the nominal value of the real part of µ, as a function of the number of Monte Carlo trials. 124 Figure 5.9: Ratio of the error in the real part of  extracted using the aperture method to the error in the real part of  extracted using the two thickness method. 125 Figure 5.10: Ratio of the error in the imaginary part of  extracted using the aperture method to the error in the imaginary part of  extracted using the two thickness method. 126 Figure 5.11: Ratio of the error in the real part of µ extracted using the aperture method to the error in the real part of µ extracted using the two thickness method. 127 Figure 5.12: Ratio of the error in the imaginary part of µ extracted using the aperture method to the error in the real part of µ extracted using the two thickness method. 128 Nominal Values µ Monte Carlo σ Monte Carlo σ Standard Error Propagation µ0r = 0.575582 0.571117589435909 0.029018113230405 0.044402395382609 µ00 r = 0.484231 0.481316041647672 0.030315888555527 0.013964580830703 0 r = 7.319669 7.310639156841075 0.221753603346064 0.326685227301673 00 r = 0.046408 0.051783551099407 0.563802891858627 0.283864841444781 Table 5.2: Comparison of Monte-Carlo and Standard Error Propagation Method. S-parameters is assumed to be 0.02 dB in amplitude and 1◦ in phase. N = 60 trials were used for the number of trials in the Monte Carlo method. The excitation field was assumed to be a plane wave normally incident and polarized across the slot, and the scattered field was computed 500 mm from the center of the slot along an angle 45◦ from the normal. For both approaches, there are huge error values in the extraction of 00 . This is, however, not surprising since as stated earlier this technique is not well suited to the extraction of the very small values of the imaginary part of the permittivity. 5.4 Aperture Optimization Several parameters including the dimension of the aperture, the direction and polarization of the incident wave, the distance and direction of observation must be carefully chosen in the design of the aperture screen. The observation angle is considered first. Monte Carlo simulations are performed to assess the effect of varying the observation angle with the other parameters set as follows: frequency = 2.88 GHz, slot length = 20 mm, slot width = 1 mm, material thickness = 3.175 mm and N = 500 runs. The observation angle is varied while normal incidence is assumed. The results obtained are shown in Figures 5.13 and 5.14. It should be immediately apparent that the results are not symmetric about zero. This is because for each Monte Carlo run a new set of random numbers are generated. As such, it not unexpected that the results are asymmetric. An observation angle of 45◦ is chosen since there does not appear to be a huge dependence on the observation angle. This angle is chosen and fixed for the other parameter investigations. The slot dimensions are considered next. A Monte Carlo simulation is carried out with 129 Figure 5.13: Effect of observation (scattering) angle on propagated error for the real part of , the real part of µ and the imaginary part of µ. 130 Figure 5.14: Effect of observation (scattering) angle on propagated error for the imaginary part of . 131 Figure 5.15: Effect of slot length on propagated error for the real part of . 132 Figure 5.16: Effect of slot length on propagated error for the imaginary part of . 133 Figure 5.17: Effect of slot length on propagated error for the real part of µ. 134 Figure 5.18: Effect of slot length on propagated error for the imaginary part of µ. 135 Optimal Slot Design f (GHz) 2.88 l (mm) 20 w (mm) 1 t (mm) 125 Table 5.3: Chosen optimal aperture screen design parameters. the slot width kept fixed at 1 mm and the slot length varied. For assessing the effect of changing the slot dimensions, the two-thickness method is considered for comparison again. The results are presented in Figures 5.15 through 5.18. The results show that the error in 0 decreases markedly as the slot length is increased from 10 mm to 24 mm. Around 22 mm, which corresponds to the resonant length in free-space, a minimum occurs. On the other hand, for 00 , the error does not seem to be affected to the same extent relative to the two-thickness method. For µ0 and µ00 , an increase in the slot length leads to increase in the error. 5.4.1 Aperture Design Choice These results obtained in Figures 5.15 through 5.18 show that a “good” choice in the slot length can be determined by striking a balance between the error propagated in  and µ. It should be noted that different length combinations may be chosen to obtain the best results for the various parameters. However, that approach has not been considered in detail in this work. A major determinant of the choice in the aperture design parameters is the practicality of implementation of the design in an experimental setup. An experimental implementation can be carried out by taping a layer of copper tape over a layer of the material-under-test with a conductor backing. How achievable this process is has a major impact on what parameters are chosen in the design and testing of the aperture technique. 136 5.5 Summary In this chapter, sensitivity analysis carried out to assess the robustness of the proposed aper- ture technique are presented. The results show that the technique demonstrates reasonable error values in the extraction of the real part of the permittivity (0 ), the real and imaginary parts of the permeability (µ0 and µ00 ) whereas large error values are obtained for the imagi- nary part of the permittivity (00 ). The results obtained for the constitutive parameters other than the imaginary part of the permittivity suggest that the proposed technique is viable for experimental implementation. In cases where an accurate extraction of  is desired, it is recommended that a different method be utilized. 137 BIBLIOGRAPHY [1] Samik Raychaudhuri. Introduction to monte carlo simulation. In 2008 Winter Simulation Conference, pages 91–100, 2008. [2] K.J. Vinoy and R.M. Jha. Radar Absorbing Materials: From Theory to Design and Characterization. Springer US, 1996. [3] Agilent Technologies. Agilent Technologies 8510C Network Analyzer System Operating and Programming Manual. Agilent Technologies Inc.., Santa Rosa, CA, U.S.A., 2001. [4] R. Fenner, E. J. Rothwell, and L. Frasch. A comprehensive analysis of free-space and guided-wave techniques for extracting the permeability and permittivity of materials using reflection-only measurements. Radio Science, 47, 01 2012. 138 CHAPTER 6 Experiments 6.1 Introduction Upon the establishment and theoretical validation of the proposed aperture technique, ex- perimental implementation followed. As covered in the previous chapter, the error analysis results predict a robust technique having performance that compares favorably against the two-thickness technique. A question arises: does this robust performance hold true in the field? Experimental validation should involve the measurement of the field reflected by a conductor-backed material layer and the field scattered by a conductor-baked material cov- ered by an aperture screen. To that end, experiments were set up in the range at Michigan State University (MSU) to assess the proposed technique and compare measured results to theoretical predictions. The experiments at MSU provided encouraging results for measurements of conductor- backed material layers (without the inclusion of the aperture screen). These results will be presented in this chapter. Challenges were encountered that precluded the implementation of the aperture screen measurements in the range at MSU. Beyond the initial conductor-backed material measurements, further experiments were carried out externally in collaboration with Dr Lydell Frasch at Boeing. As emphasized earlier, unlike many of the other characterization techniques mentioned in Chapter 1, the proposed technique provides the promise of a nondestructive and non-contact means of characterization. It should be noted that the technique does not strictly meet this criteria in its current experimental implementation. It is better described only as being non- destructive since the illuminating mechanism makes no direct contact with the material being interrogated. However, contact has to made since it is necessary to place an aperture adjacent to the sample being measured for the experiments. It is anticipated that with the future design of an appropriate applicator and modification of the numerical solution approach it can be truly non-destructive and non-contact. The different measurements carried out are 139 detailed in the sections that follow. The purpose of the experiments are to demonstrate the ability of the proposed aperture technique in successfully extracting the constitutive parameters of a given material in the field. 6.2 MSU Measurement Setup and Implementation The reflectivity arch range at MSU pictured in Figures 6.1 and 6.2 was used for the first series of measurements that were carried out. As shown in the figures, the range consists of a semi-circular metallic rail structure with a pedestal located at its center for holding samples to be measured. The range has a diameter of 6.096 m. Transmitting and receiving horn antennas with converging lenses are mounted on the range at a height of 1.219 m above the floor. The antennas (American Electronics Laboratory model H-1498 with an operating bandwidth between 2 and 18 GHz) can be slid along the rails to attain different angular positions for the plane-wave illumination of the target to be interrogated. This way, the desired incidence and reflection angles can be chosen. The lenses provide a focused beam centered on the sample with low field strength at the edges of the sample. Absorbers were strategically placed in and around the range to prevent unwanted reflections during the measurement process. A figure illustrating the setup used in the measurement procedure is shown in Figure 6.3. Ports 1 and 2 of an Agilent E5071C vector network analyzer (VNA) were connected to the transmitting and receiving antennas, respectively, through coaxial cables with the target to be measured placed at the center of the arch range as shown. Although both antennas can operate in either transmit or receive mode, only reflection measurements are required for the aperture technique. The reflection measurements were obtained in a bistatic manner with both antennas placed close to each other as shown in 6.4. A bistatic angle of 18◦ was used in the arrangement; the two antennas were placed along adjacent radials 9◦ from the normal to the surface of the material sample being measured. The reflections were measured by recording the transmission S-parameter, S21 from the transmitting to the receiving port of the VNA. Proper placement and alignment of the material samples was checked with a 140 laser level. 6.2.1 Measurement Procedure The procedure for one complete measurement set involves a background measurement, a measurement of a metallic plate, a measurement of the MUT backed by the metallic plate and finally another measurement of the metallic plate. The purpose of the background measurement is to provide a means for subtracting out clutter and isolating the response of the MUT. As such, it should be noted that any modification of the environment after the background measurement would require going back to the beginning of this work-flow as remeasuring the background would be necessitated. The metallic plate measurement is used as a calibrator with the calibration process detailed in [1] employed. For a series of measurements, two of such metallic plate measurements are carried out prior to the first measurement and after the final measurement. This is done in order to have a reference for keeping track of the positional stability of the MUTs and drift across the span of time over which the multiple measurements are performed. The stability and proper alignment of the MUT within the mounting frame used was confirmed using a laser level to check for vertical alignment of the sample through the course of the measurement procedure. Unless otherwise stated, the following analyzer settings were employed during the mea- surements: 16 averages with the averaging turned on, an IF bandwidth of 1 KHz, a start frequency of 2 GHz and an end frequency of 18 GHz and 1601 points. A standard 2-port cali- bration was carried out at the ends of the coaxial cables connected to the horn antennas with the ECal (electronic calibration) functionality on the VNA. These parameters were estab- lished after testing out various combinations and striking a balance that enabled reasonable measurement times while still ensuring sufficient dynamic range in the measurements. The measured results are obtained in the form of S-parameter measurements from the VNA in a complex (real and imaginary) format. The corresponding magnitude and phase results are computed in a MATLAB post-processing script. 141 6.2.2 Measurement Setup Validation The measurement setup was initially validated by taking a few measurements of known materials. Acrylic and Poly Vinyl Chloride (PVC) samples were measured in the conductor- backed scenario. The results obtained were compared with the theoretical values predicted by equations in Chapter 2 (Section 2.2) for a conductor-backed slab of material. The measured phase results obtained for a 2 ft by 2 ft Acrylic sample layer having a one inch (2.54 cm) thickness are shown in Figure 6.5. These results are compared with theoretical phase values geneated with a dielectric constant of 2.31 which is representative for Acrylic [2] and a loss tangent of zero with the material being assumed as lossless. There is good agreement between the measured and theoretical phase between 2 and 18 GHz but they seem to begin to deviate at higher frequencies. Similarly, Figure 6.6 shows measured phase results obtained for a 2 ft by 2 ft polyvinyl chloride (PVC) material layer having a one inch (2.54 cm) thickness. The phase values computed using the theoretical expressions from Section 2.2 are again used for comparison. These results are generated with the dielectric constant set to 3.3 [2] and the loss tangent set to zero (lossless assumption). Figure 6.6 shows that there is reasonable agreement between the measured and theoretical results. It should be noted that the discrepancies in results at higher frequencies in both cases are likely attributable to artifacts of the processing steps employed in the measurement calibration as discussed in [1]. 6.2.3 RAM Measured Results The amplitudes of the reflection measurements for Acrylic and PVC are not shown because there was no accurate knowledge of the loss tangent values of the samples. For compar- ing amplitude measurements, a 0.125-inch thick sample of MT-26 C-RAM was measured. C-RAM is a lossy, carbon filled foam that was chosen because it was readily available and convenient to work with [3]. The availability of C-RAM with adhesive backing allowed for convenience in attaching the sample to a conductor backing. Also, it served as a good start- ing point for measurements because of its non-magnetic properties. A Nicolson-Ross-Weir 142 (NRW) waveguide technique implementation [4], [5] was also employed in taking measure- ments to characterize available C-RAM samples and compare against data published by its manufacturer [3]. A waveguide setup similar to that shown in Figure 6.7 was used. An X-band waveguide was used with a C-RAM sample of thickness 0.125 inches cut to fit in the cross-section of a waveguide spacer of the same thickness. Reflection amplitude and phase results from free-space measurements of a 2 ft by 2 ft, 0.125-inch thick C-RAM sample are compared against the theoretical results generated using constitutive parameters measured by the NRW technique. A comparison of the measured reflection amplitude and theoretical amplitude is shown in Figure 6.8. As can be seen, there is good agreement between the measured and theoretical amplitude results with discrepancies less than 0.5 dB. On the other hand, the measured and theoretical reflection phase shown in Figure 6.9 do not show good agreement. In the case of the phase results, there are large discrepancies at higher frequencies. Furthermore, the behavior exhibited by the phase measurements is not physically possible as the measured results show an upward trend with an increase in frequency. Over several measurement trials the phase showed marked variation while the amplitude remained relatively unchanged. The inability to properly capture the phase precluded the measurement of the aperture screen with any confidence. Further experiments were carried out on the baseline conductor-backed material setup to investigate the possible sources of measurement error. 6.2.4 Alignment Investigation In order to determine the source of the measurement error, a 2 ft by 2 ft and 0.225-inch (0.5715 cm) thick sheet of Plexiglass was measured five separate times. The average and standard deviations of the reflection phase and amplitude were calculated. The system was calibrated between each measurement. The reflection amplitude results are shown in Figure 6.10 while Figure 6.11 shows the reflection phase results. In each figure the black line shows the average of the measurements, while the green and red lines show the 95% confidence interval (mean ±2 standard deviations). The blue line shows the theoretical 143 curves computed assuming a dielectric constant of 2.31 and a loss tangent of zero. The results show a significant variation between measurements. This variability may account for the poor performance of the C-RAM measurements. It is hypothesized that the variation in measurement is due to a slight misalignment of the sample with the transmitting and receiving horns. It is difficult to accurately re- position the sample in the pedestal mount after it has been removed, due to slight rotations and difficulty achieving a true vertical stature. Note that at 18 GHz, a shift in distance of about 0.023 mm produces a 2-way phase change of 1 degree, and thus very accurate positioning is crucial to producing reliable measurements. A much more secure measurement system, with the sample placed against a rigid vertical structure, is needed. Alternatively, a focused-beam reflection system, where the sample is placed on a rigid horizontal surface and the antennas are positioned above, would be appropriate. A system that meets the aforementioned requirements was available to Dr Lydell Frasch of Boeing and experiments carried out by him are presented in the sections that follow. It is instructive to note that simulations can also be carried out using a tool such as HFSS or CST to further investigate the alignment problems. However, that path was not consid- ered as the computational resources available were not amenable to a full-wave simulation approach. 144 Figure 6.1: Reflectivity arch range at MSU. Figure 6.2: Reflectivity arch range at MSU – Alternate view showing sample target. 145 Figure 6.3: Measurement Setup in MSU Arch Range. 146 Figure 6.4: Bistatic Geometry Setup. 147 Figure 6.5: Measured and theoretical reflection phase of a 1-inch thick layer of conductor- backed Acrylic. 148 Figure 6.6: Measured and theoretical reflection phase of a 1-inch thick layer of conductor- backed PVC. 149 Figure 6.7: Waveguide Material Measurement Setup. 150 Figure 6.8: Measured and theoretical reflection amplitude of a 0.125-inch thick layer of MT-26 C-RAM. 151 Figure 6.9: Measured and theoretical reflection phase of a 0.125-inch thick layer of MT-26 C-RAM. 152 Figure 6.10: Measured and theoretical reflection phase of a 0.225-inch thick layer of Plexi- glass. 153 Figure 6.11: Measured and theoretical reflection amplitude of a 0.225-inch thick layer of Plexiglass. 154 6.3 External Measurements As mentioned earlier, issues with alignment accuracy made it unfeasible to proceed with the aperture screen measurements in the range at MSU. In collaboration with Dr Lydell Frasch, additional experiments were performed at the Boeing facility in Seattle. The setup utilized for these measurements is a table-top system which ensures better stability in the position of the sample. This is in contrast to the initial measurement setup in the MSU arch range where the vertical orientation of the MUT made accurate alignment difficult to achieve. The configuration that was employed in these measurements is depicted in 6.12. The setup included transmit and receive horn antennas that were moved by positioners to obtain desired incident and scattering angles. The horns were oriented so that the direction of the electric field is perpendicular to the page in the figure shown, which is along the y- direction described by the coordinate system included. An absorber baffle was also present to reduce coupling between the transmit and receive antennas. For all the measurements, calibration was carried out against a metallic plate placed at the same horizontal position as the material samples. It should be noted that no focusing lenses were incorporated into the measurement setup. 6.3.1 Measurement Validation 6.3.1.1 Conductor-Backed Material Measurements For validation purposes, measurements of a conductor coated MagRAM sample that were taken in the absence of an aperture screen were compared against results generated from the conductor-backed theory shown in Section 2.2. The material properties of the MagRAM sample were determined previously at Boeing utilizing a hybrid approach that combined results obtained from waveguide and free-space implementations of the NRW technique. Figures 6.13 and 6.14 show the material properties of the MagRAM sample utilized in these measurements. The measured results that were obtained are compared against theoretical results and shown in Figures 6.15 through 6.18. Figures 6.15 and 6.16 show the reflection magnitude and 155 phase response of the conductor-backed MagRAM, respectively, measured at a bistatic angle of 20◦ . The results are seen to show agreement except for deviations that can be observed where resonance occurs in the magnitude response. Figures 6.17 and 6.18 show similar results when a bistatic angle of 60◦ was considered. Again, there is agreement between the theoretical and measured results for both the reflection magnitude and phase except for discrepancies in the vicinity of the resonance. There a couple of factors that may account for these discrepancies in the results. The thicknesses across the material sample surface were not uniform as MagRAM materials are non-rigid. Also, it is not unusual for differences to exist in the constitutive parameters of MagRAM samples across different material batches. Similar results are obtained for measurements of a foam material and the material AN- 74 which is a non-magnetic absorbing material. The permittivity properties of AN-74 are shown in Figure 6.19. The phase results for the 60◦ bistatic case of the AN-74 measurements are shown in Figure 6.20. As was observed in the case of the MagRAM measurements, the magnitude and phase results generally show similar behavior but begin to deviate as resonances in the magnitude responses are approached. 6.3.1.2 Aperture Screen Adjacent Conductor-Backed Material Measurements Measurements results for two aperture screens were considered. Firstly, a “small” slot hav- ing dimensions l = 20 mm and w = 1 mm was measured and, secondly, a “large” slot with dimensions l = 254 mm and w = 12.7 mm was measured. The aperture screen was implemented by placing a 5 mil thick layer of Aluminum foil with a slot cut out over the conductor backed material layer. Recalling earlier discussions regarding the scattering angle choice, there was a need to refrain from choosing specular scattering angles in order to avoid scenarios where the slot response is dominated by the much larger specular reflection. As such, the slot responses are considered in the case of normal incidence and a scattering angle of 45◦ . The measured reflection magnitude and phase results without a slot, with the inclusion of the small slot, and with the inclusion of the large slot are obtained for conductor-backed 156 MagRAM and conductor-backed AN-74 samples. These results are shown in Figures 6.21 through 6.24. It is notable that the results with the slots differ markedly from the conductor- backed results for the case of AN-74, particularly for the large slot. AN-74 is a more lossy material and is less susceptible to reflections of parallel-plate waveguide modes from the edge of the material sample, which are not accounted for in the theory where the material layer is assumed to be infinite in extent. The proposed calibration approach where a ratio is taken between the conductor-backed material response in the presence of the slotted conductor and the slotted conductor response in free space is not easily achievable with a table top system. An alternative approach was utilized for validation of the measurements with the aperture screen present. The approach establishes a comparison between the the measurements and the theoretical model by taking ratios of the material response for the conductor-backed MagRAM and the conductor-backed AN-74. The equality of the ratios M agRAMm M agRAMth = . (6.1) AN 74m AN 74th should hold if the theoretical and measured results are in complete agreement. Using the described approach, it is possible to compare measured results against the theoretical model. The results obtained are shown in Figure 6.25. 157 Figure 6.12: Depiction of measurement system employed for experiments at Boeing. 158 10 8 r 6 r r 4 2 0 -2 2 4 6 8 10 12 14 16 18 Frequency (GHz) Figure 6.13: Relative permittivity for MagRAM sample utilized in reflection measurements. 159 2.2 r 2 r 1.8 1.6 r 1.4 1.2 1 0.8 0.6 2 4 6 8 10 12 14 16 18 Frequency (GHz) Figure 6.14: Relative permeability for MagRAM sample utilized in reflection measurements. 160 0 -5 -10 |S21| (dB) -15 -20 S21 from and -25 Measured S21 -30 -35 4 6 8 10 12 14 Frequency (GHz) Figure 6.15: MagRAM reflection magnitude validation for 20◦ bistatic measurement (inci- dence angle = 10◦ and scattered angle = 10◦ ). 161 200 150 100 50 Phase (°) 0 S21 from and -50 Measured S21 -100 -150 -200 4 6 8 10 12 14 Frequency (GHz) Figure 6.16: MagRAM reflection phase validation for 20◦ bistatic measurement (incidence angle = 10◦ and scattered angle = 10◦ ). 162 0 S21 from and Measured S21 -5 -10 |S21| (dB) -15 -20 -25 4 6 8 10 12 14 Frequency (GHz) Figure 6.17: MagRAM reflection magnitude validation for 60◦ bistatic measurement (inci- dence Angle = 30◦ and scattered angle = 30◦ ). 163 200 150 100 Phase of S21 (deg) 50 0 -50 -100 S21 from and Measured S21 -150 -200 4 6 8 10 12 14 Frequency (GHz) Figure 6.18: MagRAM reflection phase validation for 60◦ bistatic measurement (incidence angle = 30◦ and scattered angle = 30◦ ). 164 3.5 r 3 r 2.5 r 2 1.5 1 0.5 2 4 6 8 10 12 14 16 18 Frequency (GHz) Figure 6.19: Relative permittivity for AN-74 sample utilized in reflection measurements. 165 200 150 100 Phase of S 21 (deg) 50 0 -50 -100 -150 2 4 6 8 10 12 14 16 18 Frequency (GHz) Figure 6.20: AN-74 reflection phase validation for 60◦ bistatic measurement (incidence angle = 30◦ and scattered angle = 30◦ ). 166 0 No Slot Small Slot -10 Large Slot -20 |S21| (dB) -30 -40 -50 -60 -70 2 4 6 8 10 12 14 16 18 Frequency (GHz) Figure 6.21: S21 magnitude response of conductor-backed MagRAM without an aperture screen present, with the small slot, and with the large slot. 167 200 No Slot 150 Small Slot Large Slot 100 Phase of S 21 (deg) 50 0 -50 -100 -150 -200 2 4 6 8 10 12 14 16 18 Frequency (GHz) Figure 6.22: S21 phase response of conductor-backed MagRAM without an aperture screen present, with the small slot, and with the large slot. 168 0 No Slot Small Slot -10 Large Slot -20 |S21| (dB) -30 -40 -50 -60 2 4 6 8 10 12 14 16 18 Frequency (GHz) Figure 6.23: S21 magnitude response of conductor-backed AN-74 without an aperture screen present, with the small slot, and with the large slot. 169 200 No Slot 150 Small Slot Large Slot 100 Phase of S 21 (deg) 50 0 -50 -100 -150 -200 2 4 6 8 10 12 14 16 18 Frequency (GHz) Figure 6.24: S21 phase response of conductor-backed AN-74 without an aperture screen present, with the small slot, and with the large slot. 170 3 Theoretical Measured Escattered-MagRAM/Escattered-AN-74 2.5 2 1.5 1 0.5 0 2 4 6 8 10 12 14 16 18 Frequency (GHz) Figure 6.25: Comparison between measured and theoretical results for MagRAM and AN-74. 171 6.4 Summary In this chapter, experiments carried out to implement the proposed aperture screen technique are detailed and the results obtained are presented. Experiments carried out in the MSU arch range revealed the sensitivity of the measurement setup to inaccuracies in alignment and the importance of mechanical stability in ensuring that the reflection phase responses are accurately measured. An alternative approach for executing the measurements was needed for the implementation of the aperture screen technique. Access to external measurements carried out by Dr Lydell Frasch allowed for further experimentation. The results obtained from these measurements compare favorably against theory for the conductor-backed material measurements in the absence of the aperture. However, there was a need for the developments of a different calibration approach before measurement results and theoretical results could be compared with the aperture screen. The comparisons reveal that differences in the theoretical and experimental implementations of the aperture technique have to be explored further. Discussions of these differences are covered in the next chapter. 172 BIBLIOGRAPHY [1] Bradley T. Perry, Edward J. Rothwell, and Leo C. Kempel. A comparison of the measured pulse response of layered materials using time- and frequency-domain systems [measure- ments corner]. IEEE Antennas and Propagation Magazine, 49(5):117–123, 2007. [2] Edward J. Rothwell. Extraction of the wideband dielectric properties of a material layer using measured natural frequencies. IEEE Transactions on Antennas and Propagation, 58(2):620–623, 2010. [3] Cuming Microwave Corporation. C-RAM MT-26 .125” x 24” x 24” PSA “Pres- sure Sensitive Adhesive.”. https://stores.cumingmicrowave-online-store.com/ 2-c-ram-mt-26-125-x-24-x-24-psa-pressure-sensitive-adhesive/. [4] W.B. Weir. Automatic measurement of complex dielectric constant and permeability at microwave frequencies. Proceedings of the IEEE, 62(1):33–36, 1974. [5] A. M. Nicolson and G. F. Ross. Measurement of the intrinsic properties of materials by time-domain techniques. IEEE Transactions on Instrumentation and Measurement, 19(4):377–382, 1970. 173 CHAPTER 7 Conclusions and Further Work 7.1 Summary The purpose of the work in this dissertation is to demonstrate the viability of a proposed aperture screen material characterization technique for the extraction of the constitutive parameters of conductor-backed materials. First, a theoretical formulation of the aperture screen problem is presented. Then, a theoretical demonstration is described with a numerical solution of the problem implemented and its robustness shown by considering the sensitivity of the technique to perturbations in the amplitude and phase values of an incident plane wave on an aperture screen adjacent to a conductor-backed material layer. The technique has been shown to exhibit performance that compares favorably against the two-thickness technique. An experimental implementation in the MSU reflectivity arch range revealed align- ment challenges that hindered immediate complete validation of the aperture screen theory through measurements. As such, only measurements of conductor-backed materials in the absence of the aperture screen are described. External measurements obtained from Dr Ly- dell Frasch made further experimentation possible. The results obtained and presented in this work reveal that the proposed technique shows great promise as an approach for the non-destructive characterization of the material properties of materials. However, further work is necessary to make the technique implementable in the field. There are a few distinctions between the proposed theoretical model and a realizable experimental implementation that need to be explored in further detail. These include the consideration of a non-zero aperture screen thickness in the theoretical formulation. Also, non-uniformities in the thickness of the material layer being measured, which are not captured in the theoretical implementation, should be considered in more detail. A drawback to the proposed method is that the aperture screen needs to be placed in direct contact with the MUT. This may not be desirable as it may lead to contamination of 174 the material. Also, it may be difficult to eliminate air gaps when contact is made with the MUT. In order to make the technique truly non-contact, an offset between the screen and material can be introduced, possibly in the form of cushioned stand-offs. This will require modification of the theoretical model to include the stand-off layer between the aperture screen and the MUT. Ultimately, this work has has introduced a new free-space technique which can pave the way for the realization of truly non-destructive and non-contact characterization of conductor-backed lossy materials. 175 APPENDIX A: GREEN’S FUNCTION DERIVATION The problem considered is shown Figure 7.1. The setup consists of an infinitesimally thick PEC layer with a thin aperture above a conductor-backed material region. The aperture has dimensions − L L w w 2 < x < 2 and − 2 < y < 2 with the material layer having dimensions for thickness −t < z < 0. The problem is solved by employing Hertzian potential Green’s functions of which the details of their derivation are shown here. Parallel Plate Configuration Green’s function The Hertzian potential wave equation given by J~m ∇2 Π~ m + k2Π ~m = − (A.1) jωµ ~ and H can be solved so that the dyadic Green’s function is identified. In turn, E ~ can be computed using E~ = −jωµ∇ × Π ~m (A.2) and   ~ 2 ~ H = k Πm + ∇ ∇ · Πm . ~ (A.3) Equation (A.1) can be solved by the method of superposition combining the homogeneous solution and a particular solution as ~m = Π Π ~ pm + Π~h (A.4) m where Π~ pm and Π ~ h satisfy m ~ ∇2 Π ~ pm = − Jm ~ pm + k 2 Π (A.5) jωµ 176 Figure A.1: Figure showing general problem geometry. 177 and ∇2 Π~h 2~ h = 0 (A.6) m+k Π m respectively. Dealing with (A.5) and (A.6) can be simplified by decomposing then into three separate scalar equations for each rectangular coordinate so that p p Jm ∇2 Πm + k 2 Πm = − i , (A.7) i i jωµ ∇2 Πh 2 h mi + k Πmi = 0 (A.8) where i = x, y or z. Solutions to these equations can be found by using the Fourier transform method. Since the configuration is invariant along the x and y directions, a 2-dimensional Fourier transform is carried out using the Fourier transform pair Z∞ Z∞ ~ f˜(~λ, z) = f (~r)e−j λ·~r dxdy (A.9) −∞ −∞ and Z∞ Z∞ 1 ~ f (~r) = f˜(~λ, z)ej λ·~r dkx dky (A.10) (2π)2 −∞ −∞ where ~r = x̂x + ŷy + ẑz, and ~λ = x̂kx + ŷky . Taking advantage of the Fourier transform differentiation property, (A.7) and (A.8) become   ∂2 J˜mi ~λ, z p p     Π̃m ~λ, z − p2 Π̃m ~λ, z = − (A.11) ∂z 2 i i jωµ 178 and ∂ 2 h ~    Π̃m λ, z − p2 Π̃h mi ~λ, z = 0 (A.12) ∂z 2 i where p2 = λ2 − k 2 and Z∞ Z∞ p p ~ Π̃m (~λ, z) = Πm (~r)e−j λ·~r dxdy, (A.13) i i −∞ −∞ Z∞ Z∞ Π̃h ~ Πh −j~λ·~r dxdy mi (λ, z) = mi (~r)e (A.14) −∞ −∞ and Z∞ Z∞ J˜mi (~λ, z) = Jmh (~r)e−j~λ·~r dxdy. (A.15) i −∞ −∞ Solutions to the second-order partial differential equations (A.11) and (A.12) are obtained in the next sections. In order to get a solution to (A.11) a Fourier transform is performed over z and using the Fourier tranform differentiation property as done earlier, Π̃ ˜ p is found mi to be ˜ p   J˜˜mi (~λ, kz ) ~ Π̃m λ, kz = . (A.16) i  jωµ kz2 + p2 Using Cauchy’s Integral Theorem, the inverse transform is taken and an expression for the Hertzian potential in the complex λ-plane is obtained as  Z  J˜˜m (~λ, z) p 0 dz 0 .   ~ Π̃m λ, z = G̃ λ; z|z p ~ i (A.17) i jωµ z0 179 The spectral domain Green’s function in (A.17) is given by  e−p|z−z 0 | 0  p ~ G̃ λ; z|z = (A.18) 2p where z and z 0 are the source and field points respectively and p is the spectral domain propagation factor. The solution to (A.12) can expressed in the standard form ~λ, z = W + ~λ e−pz + W − ~λ epz       Π̃h mi i i (A.19) where Wi+ and Wi− represent amplitude coefficients for upward traveling and downward traveling reflected waves respectively. Putting (A.17) and (A.19) together, the complete spectral-domain Hertzian potential is  Z e−p|z−z 0 | J˜˜m (~λ, z) dz 0 + Wi+ ~λ e−pz + Wi− ~λ epz .      ~ Π̃mi λ, z = i (A.20) 2p jωµ z0 A more concise expression can be obtained by taking advantage of the fact that the distance in the exponential depends on the field location being less than or greater than the source location as given by the relationship   z − z0,  z > z0 |z − z 0 |= (A.21)  z 0 − z,  z < z0. This makes it possible to rewrite (A.20) as   V + e−pz + W + e−pz + W − epz ,  z > z0 Πmi = i i i (A.22) −  V e +W e  pz + −pz − + Wi epz , z < z0 i i 180 where   ±pz 0 J˜ ~λ, z 0 Z e mi Vi± (~λ) = dz 0 (A.23) 2p jωµ z0 correspond to upward and downward traveling waves launched by the source. The spectral coefficients are computed next. Determination of the Spectral Coefficients In order to determine the spectral coefficients Wi+ and Wi− , the boundary conditions at z = 0 and z = −t are enforced on the spectral-domain Hertzian potential. Ex and Ey are expanded from (A.2) giving ∂Πmy   ∂Πmz Ex = −jωµ − (A.24) ∂y ∂z and   ∂Πmx ∂Πmz Ey = −jωµ − . (A.25) ∂z ∂x Since the tangential electric field vanishes on the surface of a PEC, the tangential boundary conditions on the spectral-domain Hertzian potential are found to be ∂Πmi = 0 ... i = x, y (A.26) ∂z and the normal boundary condition Πmz = 0. (A.27) As will be shown in the sections that follow, implementation of these boundary conditions can be used to compute the spectral coefficients. 181 Tangential Boundary Conditions Substituting (A.22) into (A.26), the tangential boundary condition at z = 0 can be applied leading to the expression p −Vi+ − Wi+ + Wi− = 0 ... i = x, y   (A.28) which can be rearranged to give Wi− = Vi+ + Wi+ ... i = x, y. (A.29) Similarly enforcing the boundary condition at z = −t by substituting (A.22) into (A.26), leads to the expression p Vi− e−pt − Wi+ ept + Wi− e−pt = 0 ... i = x, y   (A.30) which can be rewritten in terms of Wi+ as Wi+ = e−2pt Vi− + Wi−   ... i = x, y. (A.31) The expression for Wi− can be found by substituting (A.31) into (A.29) and solving so that − Vi− e−pt + Vi+ ept Wi = (A.32) ept − e−pt is obtained. Substituting (A.32) into (A.29) and solving for Wi+ gives + e−pt (Vi− + Vi+ ) Wi = . (A.33) ept − e−pt Now that Wi+ and Wi− have been found they can be substituted into (A.20) with the resulting expression manipulated so that the tangential spectral-domain Hertz potential can 182 be obtained. A final form of    Z  J˜m ~λ, z i Π̃mi ~λ, z = G̃t ~λ; z|z 0 dz 0   (A.34) jωµ z0 is desired for the Hertz potential where G̃t is the tangential spectral-domain Green’s function. This can be further expressed as a combination of a principal component and a scattered component thus G̃t ~λ; z|z 0 = G̃pt ~λ; z|z 0 + G̃st ~λ; z|z 0 .       (A.35) Substituting (A.32) and (A.33) into (A.20) gives   0 ˜ 0 ˜ 0 ˜ e−p|z−z | Jmi 0 e−pt · e−pz  e−pz Jmi 0 epz Jmi 0  Z Z Z Π̃mi = dz + dz + dz  2p jωµ epd − e−pd 2p jωµ 2p jωµ  z0  z 0 z0  0 ˜ 0 ˜ epz  e−pz Jmi 0 epz Jmi 0  Z Z + dz + dz  (A.36) ept − e−pt 2p jωµ 2p jωµ  z0 z0 from which G̃pt and G̃st are found to be e −p|z−z 0 | pt 0 G̃ (~λ; z|z ) = (A.37) 2p and e−p(z−z 0 +t) + e−p(z+z 0 +t) + e−p(−z+z 0 +t) + e−p(−z−z 0 −t) st 0 G̃ (~λ; z|z ) =   (A.38) pt 2p e − e −pt respectively. 183 Boundary Condition on the Normal Field The normal boundary condition is enforced at z = 0 by substituting (A.22) into (A.27) and setting z to zero so that Wz− = −Vz+ − Wz+ (A.39) Substituting (A.22) in (A.27) and enforcing the boundary condition at z = −t yeilds the expression Wz+ = −e−2pt Vz− + Wz− .   (A.40) Following a procedure similar to that outlined in the previous section, Wz+ and Wz− are found to be −   e −pt −Vz + Vz + Wz+ = (A.41) ept − e−pt and − Vz− e−pt − Vz+ ept Wz = (A.42) ept − e−pt The normal spectral-domain Hertz potential can be obtained by substituting (A.41) and (A.42) into (A.20) and comparing the resulting expression to the desired form   Z  J˜mz ~λ, z G̃n ~λ; z|z 0 dz 0 .  Π̃mz = (A.43) jωµ z0 Again, following the same steps as in the previous section G̃n is obtained as G̃n ~λ, z|z 0 = G̃pn ~λ, z|z 0 + G̃sn ~λ, z|z 0       (A.44) 184 where e−p|z−z 0 | pn 0 G̃ (~λ; z|z ) = (A.45) 2p and 0 0 0 0 sn ~ 0 e−p(z−z +t) − e−p(z+z +t) + e−p(−z+z +t) − e−p(−z−z −t) G̃ (λ; z|z ) =   . (A.46) 2p ept − e−pt Dyadic Green’s function The spectral-domain Green’s functions expressed in (A.35) and (A.44) can be manipulated to produce versions that are practical and that are employed in the MFIEs used in this work. Collecting terms in (A.37), (A.38), (A.45) and (A.46) and using the definitions of the hyperbolic functions the spectral domain Green’s functions can be written as      cosh p t − |z − z 0 | + cosh p t + z + z 0 G̃t ~λ, z|z 0 =  (A.47) 2p sinh(pt)      cosh p t − |z − z 0 | − cosh p t + z + z 0 G̃n ~λ, z|z 0 =  . (A.48) 2p sinh(pt) Considering only the tangential spectral-domain Green’s function, for z = z 0 = 0, (A.47) becomes t  ~ 0  2 cosh(pt) G̃ λ, z|z = (A.49) 2p sinh(pt) which can be rewritten as  sinh(pt) + e−pt G̃t ~λ, z|z 0 =  . (A.50) p sinh(pt) 185 Simplifying (A.50) yeilds −pt " #  1 e 0  t G̃ ~λ, z|z = 1+ (A.51) p sinh(pt) Finally, the space-domain Green’s function employed is obtained by taking an inverse Fourier transform as follows Z∞ Z∞ 1 0 G̃t ~r|~r 0 = G̃t ~λ; z|z 0 ejλ·(~r−~r ) dξdη     (A.52) (2π)2 −∞ −∞ 186 APPENDIX B: STANDARD ERROR PROPAGATION METHOD Let S1 = Aejφ1 and S2 = Aejφ2 represent measured quantities where S = S11 or S = S21 as is appropriate depending on the S-parameter being measured. It is desired that (S1 , S2 ) = 0r (S1 , S2 )+j00 0 00 r (S1 , S2 ) and µ(S1 , S2 ) = µr (S1 , S2 )+jµr (S1 , S2 ) be extracted. The propagated error can be computed using the standard formulas v u  0 2  0 2  0 2  0 2 u t 2 ∂r 2 ∂r 2 ∂r 2 ∂r σ 0 = σA +σ + σA +σ (B.1) r 1 ∂A1 φ 1 ∂φ1 2 ∂A2 φ 2 ∂φ2 and v u  00 2  00 2  00 2  00 2 u t 2 ∂r 2 ∂r 2 ∂r 2 ∂r σ 00 = σA +σ + σA +σ (B.2) r 1 ∂A1 φ1 ∂φ1 2 ∂A2 φ2 ∂φ2 for the real and imaginary parts of the constitutive parameters respectively. Similar expres- sions can be written out for µ. It should be noted that the following partial derivatives are ∂0r ∂0r ∂0r ∂0r called the ”amplification factors”: , , , . ∂A ∂A ∂A ∂A If r is differentiable then there are relationships between the amplification factors implied by the Cauchy-Riemann equations. These relationships are ∂0r 1 ∂00 r, = (B.3) ∂A A ∂φ 1 ∂0r ∂00 =− r, (B.4) A ∂φ ∂A ∂00 0 r = A ∂r , (B.5) ∂φ ∂A and 187 ∂00 r = − 1 ∂r 0 (B.6) ∂A A ∂φ Computing the derivative of  with respect to A gives ∂ ∂0 ∂00 = r +j r. (B.7) ∂A ∂A ∂A Considering (B.7), the real and imaginary parts can be split into ∂0 ( ) ∂ Re = r (B.8) ∂A ∂A and ∂00 ( ) ∂ Im = r, (B.9) ∂A ∂A respectively. The expressions in (B.4) and (B.5) can be substituted into (B.1) and (B.2) to give v  0 2 2  0 2 2 ∂00 ∂00 u   u t 2 ∂r 2 r 2 ∂r 2 r φ 0 = σA +σ A1 + σA +σ A2 (B.10) r 1 ∂A1 φ1 ∂A1 2 ∂A2 φ2 ∂A2 and v  00 2 2  00 2 2 ∂0r ∂0r u   u t 2 ∂r 2 2 ∂r 2 φ 00 = σA +σ A1 + σA +σ A2 (B.11) r 1 ∂A1 φ1 ∂A1 2 ∂A2 φ2 ∂A2 Similar expressions can be obtained for µ It should be noted that A1 and A2 from (B.10) and (B.11) are expressed in terms of linear amplitudes. If A1 and A2 are given in decibels then the chain rule can be applied so 188 that ∂0r ∂0 ∂A = r . (B.12) ∂AdB ∂A ∂AdB Here, AdB = 20 log10 A (B.13) which can be rearranged to give A ln 10 A = e dB 20 . (B.14) Taking the derivative of (B.14) with respect to AdB leads to ∂A ln 10 AdB ln 10 = e 20 (B.15) ∂AdB 20 which is simply ∂A A ln 10 = . (B.16) ∂AdB 20 Therefore, (B.12) can be rewritten as ∂0r ∂0r ln10   = A . (B.17) ∂AdB ∂A 20 And similarly for 00 , ∂00 00   r = ∂r ln10 A . (B.18) ∂AdB ∂A 20 The expressions in (B.17) and (B.18) can be substituted in (B.10) and (B.11) so that the 189 propagated error terms become v  !2  !2 0 2 00 u u u σ 2 ∂r ln 10 ∂r A A1 + σφ 2 (A1 )2 dB1 ∂A1 20 1 ∂A1 u u σ0 = u (B.19) u r !2  !2  u ∂r0 ln 10 2 ∂r00 t + σ2 2 (A2 )2  u AdB A 2 + σφ2 ∂A 2 ∂A 2 20 2 and v  !2  !2 ∂00 ∂0r u 2 u r ln 10 2 u σ 2 AdB A1 + σφ (A1 )2 ∂A1 20 1 ∂A1 u u 1 σ 00 = u (B.20) u r !2  !2  u ∂00 ln 10 2 ∂r0 u t+ σA2 r A 2 + σφ (A2 )2  . dB2 ∂A2 20 2 2 ∂A2 Expressing the error terms in the fashion shown in (B.19) and (B.20) provides the advantage ! ∂ , ∂ , ∂µ ∂µ that only four complex derivatives are needed for the computing the error , . ∂A1 ∂A2 ∂A1 ∂A2 Again, similar error expressions can be obtained for the permeability. 190