Franchois, Ann

STIENS, J

Van den Bulcke, Sara

ZHANG, L

English

Ghent University Academic Bibliography

A 2.5D BiCGS-FFT Forward Solver to ModelScattering in a mm-Wave Imaging SystemS. Van den Bulcke1, A. Franchois1, L. Zhang2 and J. Stiens21 Department of Information Technology (INTEC-IMEC), Ghent University, Belgium,email: sara.vandenbulcke@ugent.be, ann.franchois@ugent.be2 ETRO-TW, Vrije Universiteit Brussel (VUB), Belgium, email: lzhang@etro.vub.ac.be, jstiens@etro.vub.ac.beAbstractWe present an exact forward solver to compute the three-dimensional electric field scattered by a two-dimensional inhomoge-neous dielectric object, illuminated with a given three-dimensional time-harmonic field. A volume integral equation approach isused and the resulting contrast source integral equation is discretized with the Method of Moments and solved iteratively with astabilized biconjugate gradient Fast Fourier Transform method. Simulation results are compared to analytic solutions for TE andTM polarizations and to measurement results, obtained with a free-space active millimeter wave imaging system.I. INTRODUCTIONA free-space active mm-wave imaging system is being developed at the VUB, in order to study the performance of variousimaging modes for applications in indoor security and non-destructive testing. The system consists of a mm-wave vectornetwork analyzer operating in the 75 to 300 GHz range, which measures the S-parameters in amplitude and phase with adynamic range of more than 80 dB. At the transmitting side a horn antenna emits an incident Gaussian beam, which isfocused by two parabolic reflectors at the object location. At the receiving side the scattered field is focused by two parabolicreflectors for image formation at a receiving horn antenna. This optical approach is a qualitative technique for image formation.Another imaging mode is a quantitative technique, which is derived from microwave imaging principles and is largely basedon the exact solution of Maxwell’s equations. In both cases, an exact model of the wave-field propagation and scattering isindispensable to carefully study system performance and imaging capabilities. Therefore we developed an exact forward solverfor two-dimensional (2D) inhomogeneous dielectric objects embedded in a homogeneous background medium. The objects areilluminated with a given three-dimensional (3D) time-harmonic incident field and the 3D scattered field is computed. This 2.5Dconfiguration is adopted since it reduces the computational burden while maintaining the capability of accurately studying thesystem performance. A volume integral equation approach is used, whereby a contrast source integral equation is discretizedwith the Method of Moments (MoM) and solved iteratively with a stabilized biconjugate gradient Fast Fourier Transform(BiCGS-FFT) method. The formulation of the algorithm is given in section II and comparisons of simulation results to analyticsolutions for TE and TM polarizations and to measurements are presented in section III.II. FORMULATIONThe 2.5D solver presented in this paper is based on the 2D-TE CG-FFT forward solver developed by [1]. We consider a 2Dinhomogeneous dielectric object embedded in free space with complex permittivity ²(r) = ²r(r)²0 = ²′(r)+ j²′′(r), with ²′(r)and ²′′(r) representing the real and imaginary part of ²(r) and r = (x, y) the position vector. The imaginary part of the relativecomplex permittivity ²r(r) is given by ²′′r (r) = σω²0 , with ω the angular frequency, ²0 the permittivity of vacuum and σ theelectric conductivity. The object is illuminated with a 3D time-harmonic incident field Ei(r, z) = [Ei1(r, z), Ei2(r, z), Ei3(r, z)]and the resulting scattered field is defined as Es(r, z) = E(r, z) − Ei(r, z) with E(r, z) the total field. The time dependenceexp(−jωt) is omitted in the following. The problem is formulated in terms of the unknown electric flux density D(r, z) =[D1(r, z), D2(r, z), D3(r, z)]. When Ei(r, z) is expanded as Ei(r, β)ejβz , the integral equation over the object domain S takesthe following formEi1(r, β)Ei2(r, β)Ei3(r, β) = 1²(r)D1(r, β)D2(r, β)D3(r, β)−k20 +∂2∂x2∂2∂x∂y jβ∂∂x∂2∂x∂y k20 +∂2∂y2 jβ∂∂yjβ ∂∂x jβ∂∂y k20 − β2A1(r, β)A2(r, β)A3(r, β) , r ∈ S, (1)where k0 = ω√²0µ0. The vector potential A(r, β) can be calculated asA(r, β) = F−1[F[G(r, β)] F[χ(r)D(r, β)²(r)]] (2)with F the forward 2D spatial Fourier transform and F−1 the inverse 2D spatial Fourier transform. The scalar 2D Green’sfunction is given by G(r, β) = j4H(1)0 (√k20 − β2|r|) and the normalized contrast function χ(r) is defined as χ(r) = ²(r)−²0²(r) .277The β-dependence is assumed in the following.The spatial coordinate r is discretized by means of a rectangular 2D grid with a cell size of ∆ in both directions and with N cellsin the x-direction and M cells in the y-direction. We introduce the following staggered grids for n = 0 . . . N and m = 0 . . .M :r(0)n,m = (n∆,m∆), r(1)n,m = ((n− 12 )∆, m∆), r(2)n,m = (n∆, (m− 12 )∆) and r(3)n,m = ((n− 12 )∆, (m− 12 )∆). In each cell withcenter r(0)n,m the complex permittivity ²(r) is assumed constant. Vectorial testing functions Ψ(p)(r(p)n,m−r) = ψ(p)(r(p)n,m−r)u(p)are introduced for p = 1, 2, 3. u(1), u(2) and u(3) represent the unit vectors in the x, y and z-direction. Multiplication of (1)with these vectorial testing functions and integration over S yields:∫Sψ(p)(r(p)n,m− r)Eip(r)dr =∫Sψ(p)(r(p)n,m− r)Dp(r)²(r)dr− k20∫Sψ(p)(r(p)n,m− r)Ap(r)dr+∫S∂pψ(p)(r(p)n,m− r)divA(r)dr,r(p)n,m ∈ S, p = 1, 2 (3a)∫Sψ(3)(r(3)n,m− r)Ei3(r)dr =∫Sψ(3)(r(3)n,m− r)D3(r)²(r)dr− k20∫Sψ(3)(r(3)n,m− r)A3(r)dr− jβ∫Sψ(3)(r(3)n,m− r)divA(r)dr,r(3)n,m ∈ S. (3b)In the derivation of this equations it was assumed that ψ(p)(r(p)n,m−r) = 0 for r ∈ ∂S. ∂1 and ∂2 represent the partial derivativeswith respect to x and y. The Gauss theorem is used for every subdomain where ψ(p)(r(p)n,m−r)divA is continuous differentiableand the continuity of the normal component of ψ(p)(r(p)n,m − r)divA through the interfaces between the subdomains is used.The electric flux density D(r), the incident field Ei(r) and the vector potential A(r) are expanded as:Dp(r) = ²0N∑k=1M∑l=1dp;k,l ψ(p)(r− r(p)k,l ), (4a)Eip(r) =N∑k=1M∑l=1eip;k,l ψ(p)(r− r(p)k,l ), (4b)Ap(r) =N∑k=1M∑l=1ap;k,l ψ(p)(r− r(p)k,l ). (4c)When inserting (4) in (3) and carrying out the divergence operator as well as interchanging the order of integration andsummation, the weak form of the integral equation is obtained for every r(p)n,m ∈ S (p = 1, 2, 3). The functions ψ(p)(r), p =1, 2, 3 are defined as:ψ(1)(r) = Λ(x; 2∆)Π(y; ∆), (5a)ψ(2)(r) = Π(x; ∆)Λ(y; 2∆), (5b)ψ(3)(r) = Λ(x; 2∆)Λ(y; 2∆), (5c)in which Λ(x; 2∆) represents a one-dimensional triangle function in the x-direction of support 2∆ and Π(y; ∆) represents aone-dimensional pulse function in the y-direction of support ∆.In (2) the field values of D(r) are replaced by the discrete values dp;n,m and 2D spatial discrete Fourier transforms are used.Therefore the Green’s function is integrated over a circle with center r(0)n,m and radius ∆2 . The resulting value is then dividedby the surface of that circle:gn,m ={( j∆k0 )J1(k0∆2 )H(1)0 (k0∆√n2 +m2)( j∆k0 )[H(1)1 (k0∆2 ) +4jpi∆k0]n = m = 0.(6)This procedure gives no problem in handling the singularity of the Green’s function at r = 0 [2]. This leads toap;n,m = ∆2DFT−1[DFT [gn,m]DFT [χ(p)n,mdp;n,m]], p = 1, 2, 3. (7)in whichχ(1)n,m =12(χ((n− 1)∆,m∆) + χ(n∆,m∆)), χ(2)n,m =12(χ(n∆, (m− 1)∆) + χ(n∆,m∆)), (8a)χ(3)n,m =14(χ((n− 1)∆, (m− 1)∆) + χ((n− 1)∆,m∆) + χ(n∆x, (m− 1)∆x) + χ(n∆x,m∆x)). (8b)The system of equations (3) and (7) is iteratively solved for the unknowns dp;n,m using a BiCGS-scheme [3].278III. NUMERICAL AND EXPERIMENTAL RESULTSA. Comparison to analytic solutions for TM and TE polarizationsIn the simulations with our 2.5D BiCGS-FFT forward solver the initial estimate for dp ( p = 1, 2, 3) is chosen zero andthe BiCGS iterations are stopped when the relative error drops below 10−8. The object is a dielectric cylinder with relativepermittivity ²r = 2 and with a radius equal to one wavelength (λ0 = 1 mm, f = 300 GHz). We first consider TM polarization.The incident field Ei(r) = Ei3(r) is a line source at a distance of 10λ0 from the cylinder. The scattered fields Es3(r) arecalculated for 256 points on a circle with a radius of 10λ0. A good agreement with the analytical solutions is observed in Fig.1. Fig. 2 shows the normalized root mean square error (NRMSE) for five different discretizations: 8, 10, 20, 100 and 500 cellsper λ0. The error decreases as O(∆). Next we consider TE polarization. The incident field Ei(r) = Ei1(r) is a plane wavepropagating in the y-direction, which is polarized along the x-direction. The scattered fields Es1(r) are calculated in the samepoints as before. In Fig. 3, we observe again a good agreement with the analytical solutions. Fig. 4 shows the NRMSE for thesame discretizations as in the TM case.0 50 100 150 200 250 30000.0050.010.0150.020.0250.03index detector pointAmplitude Es3  exact20 cells per wavelengthFig. 1. Amplitude of scattered field for TM simulation100 101 102 10310−310−210−1100number of cells per wavelengtherrorFig. 2. Error in TM simulation0 50 100 150 200 250 30000.20.40.60.81index detector pointAmplitude Es1  exactsolver: 20 cells per wavelengthFig. 3. Amplitude of scattered field for TE simulation100 101 102 10310−410−310−210−1100number of cells per wavelengtherror Es1  Fig. 4. Error in TE simulation279B. Comparison to measurements for TE-polarizationFig. 5 represents the measurement setup used to validate our model. A W-band BWO (Backward Wave Oscillator) emits a94GHz wave which is modulated by a lock-in amplifier with a 20 kHz signal and which is focused into a gaussian beam. Theobject is positioned along the beam axis to scatter the field. At the detecting side, an open WR10 waveguide probe is usedas the receiving antenna. To detect the received signal power, a W-band planar detector with a sensitivity of 550 mV/mW isapplied. Its output voltage is fed back into the lock-in amplifier. A 30dB Faraday isolator is placed between the probe and thedetector to eliminate the reflection from the planar detector side. The probe, isolator and detector are mounted together on acomputer-controlled sled, which is able to move both parallel and perpendicular to the beam axis in a horizontal plane.The object is a metal cylinder with a diameter of 37.74 mm (which is about 11.5λ0, f = 94 GHz). In our model we choose²′r = 1 and ²′′r = 100, which corresponds to σ = 523. No higher values for ²′′r are chosen because the convergence of theBiCGS algorithm significantly slows down. The object is discretized using 10 cells per λ0, yielding in a grid size of 256x256cells. The incident field is a gaussian beam which illuminates the object completely. The distance between the source and theobject is 45 cm. The total field amplitude is measured on a line of length 8 cm at a distance of 22.5 cm away from the centerof the cylinder in 81 points, spaced 1 mm. Fig. 6 shows a comparison of a simulation with a measurement. There is a goodagreement in the center of the line. At the edges the difference is propably due to the low value chosen for ²′′r = 100 in ourmodel.Lock-inAmplifierComputerW-bandBWOPlanarDetectorFaradayIsolatorLensHornAntennaObjectMotorStage OpticalTable20 kHz Modulation SignalOpenWaveguideProbeAbsorberChamberFig. 5. Schematic diagram of scattering measurement setup Fig. 6. Amplitude of total field: gaussian beamIV. CONCLUSIONSA weak formulation of the BiCGS-FFT method for dielectric scatterers is presented. The simulated field curves are in verygood agreement with the analytical solutions for TM and TE polarizations. The convolution structure in the vector potentialallows the application of FFT’s, yielding a fast and efficient method for solving scattering problems. A first comparison tomeasurements indicates that PEC materials are demanding test objects. In the future, we will also do measurements on variousdielectric objects. Moreover, measurements will be done at the VUB with a millimeter wave vector network analyzer so thatamplitude and phase can be measured and compared to simulations.REFERENCES[1] P. Zwamborn and P. M. van den Berg, “A Weak Form of the Conjugate Gradient FFT Method for Two-Dimensional TE Scattering Problems,” IEEETrans. on Microwave Theory and Techniques, vol. 39(6), pp. 953-960, June 1991.[2] J. H. Richmond, “Scattering by a dielectric cylinder of arbitrary cross section,” IEEE Trans. on Antennas and Propagation, vol. AP-13, pp. 334-341,May 1965.[3] X. M. Xu, Q. H. Liu and Z. Q. Zhang, “The stabilized biconjugate gradient fast Fourier transform method for electromagnetic scattering,” J. Appl.Computat Electromag. Soc., vol. 17(1), pp. 97-103, March 2002.280                                                   MEDITERRANEAN     MICROWAVE     SYMPOSIUM 2006       September 19-21, 2006     Starhotel President     Genova, Italy                   Conference Proceedings                                       www.elettromagnetismo.it/conf/MMS2006                 

A 2.5D BiCGS-FFT forward solver to model scattering in a mm-wave imaging system

https://biblio.ugent.be/publication/358089/file/572691

A 2.5D BiCGS-FFT forward solver to model scattering in a mm-wave imaging system

Abstract

Similar works

Full text

Available Versions

Ghent University Academic Bibliography