Table Of ContentIEEE T R A N S A C T I 0 N S O N
MICROWAVE THEORY
AND TECHNIQUES
A PUBLICATION OF THE IEEE MICROWAVE THEORY AND TECHNIQUES SOCIETY
SEPTEMBER 1996 VOLUME 44 NUMBER 9 IETMAB (ISSN 0018-9480)
[email protected]
PAPERS
Low frequency leaky regime in covered multilayered striplines - F. Mesa ; R. Marques 1521 - 1525
Simultaneous time and frequency domain solutions of EM problems using finite element and CFH techniques -
M.A. Kolbehdari ; M. Srinivasan ; M.S. Nakhla ; Qi-Jun Zhang ; R. Achar 1526 - 1534
Application of a simple and efficient source excitation technique to the FDTD analysis of waveguide and microstrip circuits -
An Ping Zhao ; A.V. Raisanen 1535 - 1539
Analysis and design of feeding structures for microstrip leaky wave antenna - Yu-De Lin ; Jyh-Wen Sheen ; C.-K.C. Tzuang 1540 - 1547
Design and characterization of a 250-350-GHz fixed-tuned superconductor-insulator-superconductor receiver -
C.-Y.E. Tong ; R. Blundell ; S. Paine ; D.C. Papa ; J. Kawamura ; Xiaolei Zhang ; J.A. Stern ; H.G. LeDuc 1548 - 1556
A new model for microwave characterization of composite materials in guided-wave medium –
S. LeFrancois ; D. Pasquet ; G. Maze-Merceur 1557 - 1562
Complex media microstrip ridge structures: formulation and basic characteristics of ferrite structures - G.W. Hanson 1563 - 1568
Optimum mesh grading for finite-difference method - W. Heinrich ; K. Beilenhoff ; P. Mezzanotte ; L. Roselli 1569 - 1574
The eigenfunction expansion of dyadic Green's functions for chirowaveguides - Hon-Tat Hui ; A.K.N. Yung 1575 - 1583
Analytical nonlinear HEMT model for large signal circuit simulation - T. Tanimoto 1584 - 1586
Quasi-TEM analysis of coplanar waveguides with an inhomogeneous semiconductor substrate - Jean-Fu Kiang 1586 - 1589
Millimeter-wave dual-band microstrip patch antennas using multilayer GaAs technology –
D. Sanchez-Hernandez ; Q.H. Wang ; A.A. Rezazadeh ; I.D. Robertson 1590 - 1593
Planar millimeter-wave antennas using SiN/sub x/-membranes on GaAs - M. Stotz ; G. Gottwald ; H. Haspeklo ; J. Wenger 1593 - 1595
Application of the spatial finite-difference and temporal differential (SFDTD) formulation to cylindrical structure problems –
A.M.K. Chan ; Zhizhang Chen 1595 - 1600
Analysis of a double step microstrip discontinuity in the substrate using the 3D-FDTD method - Joong Chang Chun ; Wee Sang Park 1600 - 1602
New results using membrane-supported circuits: a Ka-band power amplifier and survivability testing –
T.M. Weller ; L.P.B. Katehi ; M.I. Herman ; P.D. Wamhof ; K. Lee ; E.A. Kolawa ; B.H. Tai 1603 - 1606
Axisymmetric modes of cylindrical resonators with cascaded inhomogeneous dielectrics - J.-F. Kiang 1606 - 1610
Precision broadband wavemeter for millimeter and submillimeter range –
Y.A. Dryagin ; V.V. Parshin ; A.F. Krupnov ; N. Gopalsami ; A.C. Raptis 1610 - 1613
A CAD-suitable approach for the analysis of nonuniform MMIC and MHMIC transmission lines –
A.H. Hamade ; A.B. Kouki ; F.M. Ghannouchi 1614 - 1617
Comments on "A new edge element analysis of dispersive waveguiding structures" - F.A. Fernandez ; Y. Lu ; G. Pan ; J. Tan 1618 - 1619
Author's reply - G. Pan ; J. Tan
(end)
IEEE TRANSACTIONS ON MICROWAVE THEORY AND TECHNIQUES, VOL. 44, NO. 9, SEPTEMBER 1996 152 1
Low Frequency Leaky Regime in
Covered Multilayered St riplines
Francisco Mesa, Member, IEEE and Ricardo MarquCs, Member, IEEE
Abstruct- This work studies the low frequency leakage of
power in covered multilayered striplines including the possibility
of several conductors on different interfaces. The appearance and
features of the lateral radiation is discussed, assuming that only
one parallel-plate waveguide mode is above cutoff. Comparison t - f
with quasi-TEM results shows that some leaky waves behave like
quasi-TEM modes with radiation losses. Lateral radiation and
ohmic conductor losses are also compared.
Zd
I. INTRODUCTION Fig. 1. Cross section of a covered noncoplanar multilayered stripline
0
WING to its versatility, noncoplanar and/or multilay-
ered striplines have been often proposed in the CAD the quasi-TEM solution. An additional comparison between
of microwave integrated circuit (MIC) and monolithic mi- typical values of leaky and conductor losses will show that
crowave/millimeter wave integrated circuits (MMIC) [ l]-[3]. the leaky losses can become as important as the ohmic losses
At low frequencies, these lines are usually analyzed using (or even more).
a quasi-TEM approach, which can account for losses due
to the use of lossy substrates [4] and/or finite conductivity 11. FULLW AVE ANALYSIS
metallizations [5]. Noncoplanar and/or multilayered stmctures
The structure under study is a covered multilayered
may also exhibit radiation losses due to the excitation of
and/or multiconductor stripline as shown in Fig. 1. An elec
surface and volume waves in the background waveguide ke.,
tromagnetic wave with the following dependence: E(r,t ) =
the remaining waveguide when all the strips are removed)
E(z,y)exp[-j(k,z - wf)],i s assumed to propagate in the
161-[8]. From a technological point of view, the presence of
line; where w is the angular frequency and IC, the complex.
these radiation losses can present serious problems for the propagation constant: IC, = ,Bz + jaz (with ,Bo > 0 and
proper performance of the systems. On the other hand, the
a, 5 0 being real numbers for a field propagating along:
leaky phenomenon could be eventually used for designing
the +z direction). The method of analysis employed in this
some microwave components (antennas, directional couplers,
work is an extension of the Galerkin method in the spectral1
etc.). In any event, a proper study of the propagation charac-
domain modified to take into account the changes in the
teristics of these structures at low frequency should consider
inversion-contour definition of the spectral Green’s dyad [7]
the possible existence of leakage either to prevent it or to take
For low frequency leaky waves, when only the dominanl
advantage of it. To our knowledge, little attention has been
waveguide mode is above cutoff, this inversion contour should1
still paid to its low frequency leaky regime, despite of the
be deformed to surround the Green’s dyad pole associated withi
possible practical significance of the lateral radiation at low
this dominant background waveguide mode [6], [7], and [9]
frequencies in these types of structures.
In this case, the leakage of power due to the excitation of the
In this paper we restrict ourselves to treat the low frequency
first waveguide mode is radiated at an angle 6’ from the strip
leakage of power by the dominant waveguide mode in covered
as a nonuniform mode with a complex propagation constant
noncoplanar and/or multilayered striplines. In this case, we
+
will not be concerned with the aspects related to the presence k = ,B,u ja,v (11’
of multiple poles near the real axis andor branch cuts of the
where ,Bg and as are real and positive numbers; and u andl
spectral Green’s function. The numerical results of the full-
v real and unit vectors in the (IC, z) plane. In lossless media.
wave analysis will be compared with those derived from a
this mode propagates along direction u and increases at the
quasi-TEM analysis. This comparison suggests that, under cer-
perpendicular one, v, namely: u . v = 0 [lo]. Therefore v
tain conditions, the leaky mode can be considered as a laterally
and u can be used to form the following right-handed se1
radiative mode whose near-field behavior is well described by
of unit and orthogonal vectors: {v,ay,u),s hown in Fig. 2.
Manuscript received February 13, 1995; revised May 24, 1996. This work The relationships between the line propagation constant andl
was supported by the DGICYT, Spain TIC95-0447.
the wavenumber of the excited waveguide mode can be
The authors are with the Grupo de Microondas, Departamento de
Electrhca y Electromagnetismo, Facultad de Fisica, Universidad de Sevilla, summarized as [ 101
Avenida Reina Mercedes s/n, 41012 Sevilla, Spain
Publisher Item Identifier S 001 8-9480(96)06381 -8
0018-9480/96$05,00 0 1996 IEEE
1522 IEEE TRANSACTIONS ON MICROWAVE THEORY AND TECHNIQUES, VOL. 44, NO. 9, SEPTEMBER 1996
V .3w 1.52 0.08Yt3
8
E
s 1.48 0
0.06.8
Qg) 1.44 .9*-.
E c
0.04g
4
W
Fig. 2. View of the orthogonal vectors: {v,a y,u } ..r,.,.l.(.. .. */...-* 0.02Rq
~~~ I.-*-
,,....*<. .”...“. ..... *.., E
:* I”.
n Z
aZ = -ag sin Q (3) 0 1 2 3 4 5 6 7 8 9 1
Freq(GHz)
Pi - a; = 7,2 (4) Fig. 3. Normalized propagation constants of a stripline with ans airgap
versus frequency. 6 = 3.5 mm, h = 4.45 mm, w = 6.35 mm, er = 2.6.
where 7,”is the squared wavenumber of the first waveguide Solid lines: Normalized phase constants of proper and improper stripline
modes. Dashed line: Normalized attenuation constant. Dotted line: Normalized
mode, which is a real and positive number for lossless media
wavenumber of the dominant parallel-plate mode.
u11.
The real part of the Poynting vector in the transverse
(z,x) plane, St, of the excited waveguide mode can be IV. NUMERICARL ESULTS
computed taking into account the LSM nature of this mode. Our numerical computations have been checked by com-
After a straightforward calculation (see Appendix in [7]), the paring with pertinent previously published results. We found
following expression is obtained a good agreement with the data reported in Fig. 2 of [SI. This
result was shown in [7] and now this structure is analyzed
versus frequency in Figs. 3 and 4. At low frequencies and
for the chosen dimensions, Fig. 3 shows that the normalized
(to ko = w a )ph ase constant of the quasi-TEM solution
where H, is the magnetic field component in the v direction; is higher than the normalized eigenvalue, y,/ko, of the first
that is, in the direction perpendicular to propagation. In a (and only above-cutoff) waveguide mode and very close to the
similar way, the imaginary part of St is found to be phase constant of a bound real proper mode (RPM). A complex
improper mode (CIM) (namely, a leaky mode) is also present
for frequencies above 1.6 GHz after the conjunction of two real
improper modes (RIM). When frequency increases, the quasi-
TEM data appear approximately at the same distance from the
These two equations show that a) there is no flux of reactive phase constants of the CIM and the RPM. This points out that
power along the direction of propagation, u, of the excited the quasi-TEM solution does not account properly for any of
waveguide mode; and b) there always exists a flux of active the full-wave solutions in the present case (as was expected
power in this direction, provided that 7,”> 0. in view of the h/X ratio of this structure). It is interesting
to note in Fig. 3 that the onset of the leaky mode occurs
with a phase constant greater than the wavenumber of the
111. QUASI-TEMA NALYSIS
first waveguide mode yg,m eaning that the solution lies in the
A quasi-TEM analysis has also been used in order to region that is usually called the spectra2 gap [12]. Fig. 4 shows
compare with the full-wave results. This analysis assumes the the normalized phase constant, Pg/k0, and the radiation angle,
usual quasi-TEM approach, that is, the longitudinal component 8, of the excited waveguide mode in addition to the normalized
of both the electric and magnetic fields, E, and B,, have phase constant, (P,/ko), of the stripline leaky mode. It can
been neglected, and the transverse fields have been expressed be seen that the leaky wave is fast with respect to that of
as the transverse gradients of some potential functions (the the excited waveguide mode for all frequencies (regarding the
propagation constant can be then expressed as a function of phase velocities: PZ < Pg), even in the spectral gap region
line capacitance and inductance per unit length). As is well (PZ > rg), and that the radiation angle takes nonzero values
known, this analysis can not evaluate the radiation losses from the onset of the leaky mode (where the leakage begins).
due to the leakage, but after some modifications (see for The leaky wave does not show any particularly peculiar
instance [4] and [5]), the ohmic losses in both the substrate behavior from its onset, even inside the spectral gap region,
and conductors can be evaluated with high accuracy in the and therefore the question of the excitability of such a leaky
low frequency range (i.e., when the line wavelength is much wave by a finite source (which relates to the physical meaning
greater than the transverse dimension of the line). The details of the wave) can not be answered in the frame of a purely (2D)
of the quasi-TEM analysis can be found in [5], and will not propagative analysis-as performed in this paper. However,
be presented here. considering the analysis of a simpler structure [13], it can be
MESA AND MARQUES: LOW FREQUENCY LEAKY REGIME IN COVERED MULTILAYERED STRIPLINES 1523
3
cs
8
2.9
8
8
2.8
i
.^
2.7
3 2.6
€
8
t l l 1 ‘1, 2.5
I I I I I > I I I I I ,
2 3 4 5 6 7 8 9
Freq(GHz)
Fig. 4. Dispersion curves and radiation angle for the leaky mode of Fig. 3. Fig. 5. Normalized propagation constants of an inhomogeneous stripline
versus the upper ground plane height at 2 GHz. h = 200 pm, w = 80 pm.
Solid lines: normalized phase constants of proper and improper stripline
modes. Dashed line: Normalized attenuation constant. Dotted line: Normalized
postulated that the degree of excitation of the leaky wave on
wavenumber of the dominant parallel-plate mode (yg/ ko). Dashed-dotted
the line by a finite source decreases as the solution enters into line: Quasi-TEM data.
the spectral gap region [ 131.
Next, we are interested in analyzing what happens when the
quasi-TEM solution for the phase constant becomes smaller small values suggest that the first waveguide mode is excited
than the wavenumber of the first waveguide mode yg. For this with a small amplitude. In other words, the contribution of this
purpose we have chosen the dimensions of the inhomogeneous waveguide mode to the near field is small and we can assume:
stripline and the frequency as shown in Fig. 5. A small ratio that the near field of the leaky mode is mainly a superpositioni
between the wavelength and the transverse dimensions of the of all the remaining nonuniform waveguide modes, which are
+
structure, (d h)/Xo 0.01, has been also imposed to assure well below cutoff at the assumed low frequencies. Therefore
N
both the validity of the quasi-TEM approach and that only one the eigenvalues of these latter modes are purely imaginary
waveguide mode is above cutoff. Also note that this structure numbers with absolute values much higher than the leaky
+
can be used to analyze one of the two fundamental modes (the mode propagation constant, (ys,nl >> lkzl, w@ (n
even mode) of a noncoplanar symmetrical broadside coupled 0). Under these circumstances, the differential operator of‘
line. In the inhomogeneous stripline, a leaky mode can be the Helmholtz equation satisfied by these waveguide modes,
+
induced varying the height of the upper ground plane. This d2/dy2- ki - k: w2pt, approaches d2/dy2 - k;, since
leaky mode turns into two RIM’S for d > 101 pm and the w2,w << k; and kz << k:. The Helmholtz operator can be then
bound mode (RPM) becomes a RIM after its phase constant approximated by the Laplace’s operator, d2/dy2- k;, for all
reaches the wavenumber of the dominant waveguide mode the meaningful near-field components of the leaky wave. On
(d 88 pm). A comparison with the quasi-TEM data clearly the other hand, when the quasi-TEM solution accounts for the
N
shows that the phase constant of the leaky mode coincides bound mode of Fig. 5, the contribution of the first waveguide
with the quasi-TEM results at the lower values of the upper mode to the expansion of the near line field does not seem
ground plane height. In this range cy, << ,&, and therefore to be negligible, making the line field less “TEM’ than in the
from (2) to (4): cyg << pg (these facts have been also tested leaky case.
numerically), implying that the first waveguide mode is excited A result showing the appearance of leaky modes in
as a quasiuniform mode. In this case (5) and (6) show that broadside-coupled structures is presented in Fig. 6. This figure
the flux of reactive power-at any direction-is negligible in shows the behavior of the normalized propagation constants
comparison with the flux of active power, and consequently the versus the value of the relative permittivity of the lowest
nature of this low frequency leaky mode should be considered substrate. Although this structure has two quasi-TEM modes,
radiative rather than reactive. The closeness of the different we have not plotted that mode corresponding to the highest
full-wave solutions to the quasi-TEM solution in the different value of the phase constant (even-like mode) because here this
d ranges suggest that, for d > 100 pm the real proper mode mode shows the common behavior of a bound mode. However,
has a field pattern very similar to the quasi-TEM solution, and the other quasi-TEM solution (QTS) appears very close to
for d < 80 pm it is the leaky wave field that resembles the either a leaky mode or a real proper mode. The phase constant
quasi-TEM field pattern. of the RPM runs very close to the dispersion curve for the first
An interesting feature observed in Fig. 5 is that the quasi- waveguide mode for higher values of E ~a,nd approaches the
TEM results appear closer to the phase constant of the leaky QTS at the low values of the relative permittivity. In the first
mode than to the phase constant of the bound mode in its region, the RPM seems to be a perturbation of the waveguide
respective ranges. This fact can be explained taking into mode caused by the presence of the strips. We can deduce from
account the small values of the normalized attenuation constant the curves of Fig. 6 that the appearance of the leaky mode is
of the leaky mode in Fig. 5 for all values of d, since these again related to the closeness of the quasi-TEM solution and
1524 IEEE TRANSACTIONS ON MICROWAVE THEORY AND TECHNIQUES, VOL. 44, NO. 9, SEFTEMBER 1996
20 50
Fig. 6. Normalized propagation constants of an broadside-coupled line as a ~~~~~~~~~
function of the relative permittivity of the lower substrate. h = w = 200 Fig. 7. Dispersion curves of the fundamental modes of a broadside-coupled
pm, Freq = 2 GHz. Solid lines: Normalized phase constants of proper microstrip line, h = 200 pm, w = 80 pm, d = 50 pm. Solid lines:
and improper modes. Dashed line: Normalized attenuation constant; Dotted normalized phase constants of the three fundamental full-wave modes. Dashed
line: normalized wavenumber of the dominant parallel-plate mode (y,/ko). lines: normalized attenuation constant. Dotted line: normalized wavenumber
Dashed-dotted line: Quasi-TEM data. of the dominant parallel-plate mode (y,/lc~).
waveguide mode (this G-mode seems to be a perturbation
the wavenumber of the dominant waveguide mode, although
of the waveguide mode). Fig. 7 also provides a comparative
the onset of the leaky regime cannot be easily predicted. At
analysis between the leakage losses of the E-mode and the
the highest permittivity values (with the QTS being lower
conductor losses induced by the presence of real copper strips
than the waveguide wavenumber), we find that the full-wave
solution that is closer to the QTS is the leaky mode. Since (0 = 4. lO'(nm)-l) with finite thickness (t = 5 pm). The
results for conductor losses have been computed using the
the attenuation constant of this mode is much lower than the
method developed in [5]. This comparison shows how ohmic
phase constant, the leaky mode can be considered again as a
losses, which are dominant at very low frequencies, can be
laterally radiative mode whose near- field behavior resembles
lower than lateral radiation losses, even at frequencies where
a quasi-TEM mode.
the phase constant of the leaky wave is still very close to the
Fig. 7 shows the propagation constants of a potentially
quasi-TEM solution.
leaky symmetrical broadside-coupled line as a function of
frequency. The suitability of a quasi-TEM approach could
V. CONCLUSION
be again expected in view of the reduced dimensions of this
structure; this is also confirmed by the nondispersive nature Low frequency leakage of power in covered noncoplanar
of the full-wave solutions for the phase constant plots of and/or multilayered stripline transmission lines has been ana-
both the odd mode (E-mode) and even mode (H-mode). lyzed. It has been verified that the low frequency leaky regime
In fact, a quasi-TEM analysis of this structure (assuming is expected to appear when the quasi-TEM values for modal
perfect conductors) yields almost identical results for the phase phase constants appears near (or below) the wavenumber
constants. Nevertheless, for small distances between strips, it of the dominant waveguide mode. In this circumstance, the
has been found that the E-mode of this structure becomes computed results for the phase constant of the leaky mode
leaky starting at very low frequencies, as can be seen in have been compared with the quasi-TEM results, showing a
Fig. 7. The attenuation constant of this leaky mode shows good agreement over a wide frequency range. This behavior
an almost quadratic dependence with the frequency (linear has been systematically found for small dimensions of the
dependence for the normalized attenuation constant, cr,/k~). structure when the quasi-TEM phase constant is smaller than
This dependence can be understood taking into account that the the wavenumber of the first waveguide mode. These results
attenuation constant must be an even ffunction of the frequency show that the low frequency leaky modes behave in these
[14]. This latter fact together with the linear dependence of the structures like a quasi-TEM mode with radiation losses.
phase constant with frequency suggests the validity of a first- Computed values for the power associated with the low
order approach in the series expansion of the fields in powers frequency leaky mode have also revealed that its nature is
of frequency. Therefore, for practical purposes, only one full- mainly radiative instead of reactive. A comparative analysis
wave computation would be required to obtain both the phase between radiation and conductor losses has also shown that,
and attenuation constants of a leaky mode for all frequencies in many cases, radiation losses can be higher than conductor
where the quasi-TEM approach yields very accurate results losses. Therefore, the characterization of the losses of the
for the phase constant. This accuracy seems to be assured for structure should include a full-wave analysis accounting for
small dimensions of the line when the quasi-TEM solution is the possible leakage of power. Nevertheless, when the phase
well below than the wavenumber of the first waveguide mode. constant is well approximated by the quasi-TEM results,
An additional bound mode, referred as G-mode in the figure, we have found a quadratic dependence of the attenuation
also appears very close to the dispersion curve of the first constant with frequency (this fact is somewhat expected after
MESA AND MARQUCS: LOW FREQUENCY LEAKY REGIME IN COVERED MULTILAYERED STRIPLINES 1525
a series expansion of the field in powers of frequency), which [lo] R. Marquts, F. Mesa, and N. K. Das, “Comments to the criterion of
implies that only one full-wave computation would provide leakage from printed circuit transmission lines,” ZEEE Trans. Microwave
Theory Tech., vol. 43, pp. 242, Jan. 1995.
the propagation constants over a wide frequency range.
[ 111 R. E. Collin, Field Theory of Guided Waves, 2nd ed. New York: 1EE.E
Press, 1991.
[121 H. Shigesawa, M. Tsuji, and A. A. Oliner, “The nature of the spectral-
ACKNOWLEDGMENT _ga_p b etween bound and leaky solutions when dielectric loss is mesent
The authors like to thank for providing in printed-circuit lines,” Raio Sci., 1993.
J’ [13] H. Ostner, J. Detlefsen, and D. R. Jackson, “Radiation from on(:-
the quasi-TEM data (including conductor losses) and also to dimensional dielectric leakv-wave antennas.” IEEE Trans. Anrennus
the referees for their helpful suggestions. Propagat., vol. 43, pp. 331-339, Apr. 1995.
[14] I. V. Lindell, “On the quasi-TEM modes in inhomogeneous multicon-
ductor transmission lines,” IEEE Trans. Microwave Theory Tech., vol.
REFERENCES
29, pp. 812-817, Aug. 1981.
M. Homo and F. Medina, “Multilayer planar structures for high-
directivity directional coupler design,” ZEEE Trans. Microwave Theory
Tech., vol. 34, pp. 1442-1449, Dec. 1986.
M. Tran and C. Nguyen, “Modified broadside-coupled microstrip lines
Francisco Mesa (M’95) was bom in CBdiz, Spain,
suitable for MIC and MMIC applications and a new class of broadside
on April 15, 1965. He received the Licenciadlo
coupled hand pass filters,”. IEEE Trans. Microwave Theory Tech., vol.
degree in June 1989 and the Doctor degree in De-
41, pp. 1336-1342, Aug. 1993.
cember 1991, both in physics, from the Universily
F. Masot, F. Medina, and M. Homo, “Analysis and experimental
of Seville, Spain.
validation of a type of three-Microstrip directional coupler,” to be
He is currently Associate Professor in the De-
published in IEEE Trans. Microwave Theory Tech., Sept. 1994.
partment of Applied Physics at the University of
M. Homo, F. Mesa, F. Medina, and R. MarquBs, “Quasi-TEM analysis
Sevilla, Spain. His research interests focus on elec-
of multilayered, multiconductor, coplanar structures with dielectric and
tromagnetic propagatiodradiation in planar linm
magnetic anisotropy including substrate losses,” ZEEE Trans. Microwave
with general anisotropic materials.
Theory Tech., vol. 43, pp. 1059-1068, Aug. 1990.
J. Aguilera, R. MarquCs, and M. Homo, “Quasi-TEM surface impedance
approaches for the analysis of MIC and MMIC transmission lines,
including both conductor and substrate losses,” ZEEE Trans. Microwave
Theory Tech., vol. 43, pp. 1553-1558, July 1995.
L. Carin and N. K. Das, “Leaky waves on broadside-coupled mi-
crostrip,” ZEEE Trans. Microwave Theory Tech., vol. 40, pp. 58-66, Ricardo MarquCs (M’94) was born in San Fer-
Jan. 1992. nando, CBdiz, Spain. He received the degree of
F. Mesa and R. Marques, “Integral representation of spatial Green’s Licenciado in physics in June 1983, and the degree
function and spectral domain analysis of leaky covered strip-like lines,” of Doctor in physics in July 1987, both from the
IEEE Trans. Microwave Theory Tech., vol. 43, pp. 828-837, Apr. 1995. University of Sevilla, Spain
D. N. Nghiem, J. T. Williams, D. R. Jackson, and A. A. Oliner, “Leakage Since January 1984, he has been with the De
of the stripline dominant mode produced by a small air gap,” in MTT-S partment of Electronics and Electromagnetism at
Dig., 1992, pp. 491494. the University of Sevilla, where he is currently
N. K. Das and D. M. Pozar, “Full-wave spectral-domain computation Associate Professor in Electricity and Magnetism.
of material, radiation and guided-waves losses in infinite multilayered His main fields of interest includes MIC deviccs
printed transmission lines,” IEEE Trans. Microwave Theory Tech., vol. design, wave propagation in anisotropic media, and
39, pp. 5443, Jan. 1991. electromagnetic theory.
1526 IEEE TRANSACTIONS ON MICROWAVE THEORY AND TECHNIQUES, VOL. 44, NO. 9, SEPTEMBER 1996
Simultaneous Time an ncv ai
Solutions of EM roblems
Finite Element and GF
M. A. Kolbehdari, Meera Srinivasan, Student Member, IEEE, Michel S. Nakhla, Senior Member, IEEE,
Qi-Jun Zhang, Senior Member, IEEE, and Ramachandra Achar, Student Member, IEEE
Abstruct- This paper describes an efficient method for both In this case, fast fourier transform (FFT) and its inverse
time- and frequency-domain solutions of electromagnetic (EM) can be used to move from one solution space to the other.
field problems. In this method EM field problems are formulated This can increase the computational time since in order to
using Laplace-domain finite element approach and are solved
achieve satisfactory resolution, FFT has to span longer time
using complex frequency hopping (CFH) technique. CFH is a
and frequency intervals with smaller increments [16]. In an
moment-matching technique which has been used successfully in
the circuit simulation area for solution of large set of ordinary attempt to improve the efficiency of FEM solution techniques,
differential equations. Problems consisting of Dirichlet, Neumann a Laplace domain FEM has been proposed in [17]. This
and combined boundary conditions can be solved using the method is an extension of the technique used for solution
proposed algorithm to obtain both time and frequency responses.
of heat transfer problems [18]-[19]. However, it is based
Several electromagnetic field problems have been studied using
on congruence transformation which involves computationally
the new technique and the speed-up advantage (one to three
orders of magnitude) compared to conventional finite element expensive process of determining all the eigenvalues and
technique is demonstrated. A good agreement between numerical eigenvectors of a large matrix.
results obtained using the proposed method and the previously From the conceptional point of view, the new technique pro-
published results has been found.
posed in this paper falls in the category of [ 171-[ 191. However,
the main difference is that it requires the computation of the
I. INTRODUCTION dominant natural modes only and thus eliminating the major
M ODELING and simulation based on electromagnetic computing cost.
field formulations are important for accurate analysis The new solution technique is based on complex fre-
and design of high-speed circuits and systems [l], [2]. One of quency hopping (CFH) [20]-[23], an expansion of asymptotic
the most common approaches used for the solution of electro- waveform evaluation (AWE) [24] recently developed in the
magnetic field problems is the finite element method (FEM) circuit simulation area, which yields a speed-up factor of
[3]-[7]. FEM formulations are either space/frequency-domain 10-1000 over conventional circuit simulators. It has been
or spacehime-domain. Space/frequency formulation leads to a extended to solution of static fields in VLSI interconnects
set of algebraic equations which have to be solved repeatedly [25], in ground/power planes [26] and thermal equations [27].
at many frequency points, while spacekime formulation leads CFH uses the concept of moment matching to obtain both
to a set of ordinary differential equations which have to be frequency- and time-domain responses of large linear networks
solved in the time-domain. The size of the equations to be through a lower order multipoint Pad6 approximation. It
solved is usually large and the conventional solution algo- extracts a relatively small set of dominant poles to represent
rithms are restricted by computing time and stability condition. a large network that may contain hundreds to thousands of
For example, frequently used time-stepping schemes [8]-[ 121 actual poles. CFH is particularly suitable for solving large
require the satisfaction of the stability condition, that is, the set of ordinary differential equations which make it a logical
ratio of spatial to temporal subdivision is to be greater than candidate for solving time-harmonic EM equations. The main
or equal to the speed of propagation [13]. This means the steps involved can be summarized as follows: first, the given
finer the mesh, the smaller the time step that must be chosen problem which is in the form of damped wave equation is
[ 141. Implicit variable time step integration algorithms [15] formulated using FEM and the resulting ordinary differential
can be used. However, they require solution of large set of equation is transformed to the Laplace domain; second, the
algebraic equations at each time point. Frequently, both time- Laplace domain output is expanded using Pad6 approximation
domain and frequency-domain results are of interest [ lo]-[ 11 1. at selected frequency points; third, information from each
expansion point is used to generate the output frequency
Manuscript received June 29, 1995; revised May 24, 1996. This work
response or alternatively a unified set of dominant polehesidue
was supported in part by the Natural Sciences and Engineering Research
Council of Canada, Micronet, a Canadian Network of Centres of Excellence on pairs; finally, the results are transformed to the time-domain
Microelectronics Devices, Circuits and Systems, and Bell-Northern Research. in either analytical or numerical form.
The authors are with the Department of Electronics, Carleton University,
The main advantages of the proposed method can be sum-
Ottawa, Ont, K1S 5B6 Canada.
Publisher Item Identifier S 001 8-9480(96)06397- 1. marized as follows:
0018-9480/96$05.00 0 1996 IEEE
KOLBEHDARI et al.: SIMULTANEOUS TIME AND FREQUENCY DOMAIN SOLUTIONS 1527
1) 10-1000 times faster than the conventional FEM solu- where !P, A, B, C, and q are defined as follows
tion techniques;
2) solution algorithm does not suffer from instability prob-
lems associated with the time-stepping methods;
3) produces simultaneously both the time- and frequency-
domain results. ep=l
The remaining part of this paper is organized as follows: Ne
in Section 11, the damped wave equations are derived from B = C = [TIep (7)
Maxwell’s equations and formulated using variational tech- e,=]
niques. Section I11 describes the AWE and CFH techniques
in the context of solving the FEM equations. To illustrate
the accuracy and efficiency of the proposed method, Section e,=l
IV presents numerical results for several electromagnetic field
where Ne denotes the total number of elements. A, B, and
problems. The efficiency of the proposed method is illustrated
C are N x N symmetric, positive definite matrices assembled
and discussed in Section V. Finally a brief conclusion of the
from [SIep and [TIepa nd N is the total number of nodes. q is
paper is presented in Section VI.
a vector of dimension N, assembled from [GI,,, containing
the forced terms attributed to the time or space excitation.
11. FORMULATION
[SIepa nd [TIepa re real symmetric square element matrices
For a homogeneous, isotropic, and a linear medium, starting and [GI,, is a column vector given by
from Maxwell’s equations the following scalar equation can
be derived as
where @(zy, , z,t ) represents either the electric or magnetic
field components, and E, p, cr are the permittivity, permeabil- (1111
ity, conductivity of the medium respectively. f(x,y , z,t ) is
related to the external excitation which could be function of
subject to the initial conditions
time, space or both. This expression can be applied to either a
diffusion problem where the second term in the left-hand-side
(12)
of (1) is dropped or to a wave equation where the third term
(13)
is dropped. Applying the finite element method [13] to (1)
where the interpolating functions are selected in exactly the where Rep denotes the domain of the finite elements. In (9) to
same fashion as compared to the time independent problems
(1 1) p refers to either 1, 2, or 3 corresponding to one, two, or
except now the nodal values are taken to be functions of time three-dimensional (3-D) finite element solution scheme while
rather than constants. Expanding the unknown function @ in ap is the corresponding shape function and ep represents the
the triangular finite element domain Re as
finite elements. Taking the Laplace transform of (4) results in
Ne, (Cs2+ Bs + A)X(s) = R(s) (1411
@(x, Y,2 , t)= @%(t)%Y,(2 ~) , (2)
n= 1 or
where Ne, is the total number of nodes in element Q,,, @i(t)
denotes the local set of unknown time dependent expansion Y(s)X(s)= R(s) (15)
coefficients and a,(z,y , z) are nodally based interpolation
where B = pcrB, C = e&, X(s) = C[9(t)],a nd R(s)
functions such as those regularly used in FEM triangular
is given by
formulations [7]. It is of interest to note that in (2), the spatial
+ + +
variables x, y and z are discretized whereas the temporal R(s) = sC(qo q;) BQo Q(s) (16)
variable t is not. An appropriate functional to be minimized
for (1) is where Q(s) = C[q(t)].
111. MOMENTM ATCHINGT ECHNIQUES
Moment matching techniques such as asymptotic waveform
where R represents the finite element region. By minimizing
evaluation (AWE) [24] and complex frequency hopping (CFH)
the functional of the problem and applying associated Dirichlet
[21] have been topics of intense research in the circuit sim-
and Neumann boundary conditions, a system of ordinary
ulation area in the recent years. They have been successfully
differential equations is obtained as
and efficiently applied for obtaining the solution of large set
+ +
Aq crpBq’ epCQ’’ = q (4) of ordinary differential equations.
1528 IEEE TRANSACTIONS ON MICROWAVE THEORY AND TECHNIQUES, VOL. 44, NO. 9, SEPTEMBER 1996
In general, moment matching technique approximates the less than full FFT analysis requiring hundreds or thousands of
frequency response of a Taylor series expansion in the com- frequency point analyzes.
plex s-plane. The cost of an expansion is approximately one
frequency point analysis. The moments (coefficients of the B. Pad6 Approximation
expansion) are matched to a lower order transfer function
Consider the response vector X(s) represented in (15).
using a rational Pad6 approximation. This transfer function
Taylor's series expansion of the output X(s) about a complex
can be used to obtain the output response. The moments are
frequency point s = so, is given by
generally taken from an s = 0 expansion (Maclaurin series).
Single Pad6 approximations are accurate near the point of 00
expansion in the complex s-plane and decrease in accuracy X(s) = M,(s - so), (17)
with increased distance from the point of expansion. Complex n=O
frequency hopping overcomes this problem and is summarized
where M, is the nth moment vector of the Taylor's expansion
in the following section.
about so and is given by
A. Complex Frequency Hopping
Complex frequency hopping is a method by which the
A recursive equation for the evaluation of the moments can
frequency response of a system is expanded in multiple Taylor
series expansions in the Laplace s-domain. The expansion be obtained in the form
points are chosen on or near the imaginary axis because n
poles that dominate the transient and frequency response of [Y(so)]Mn= - c[Y1rM,-T/l! (19)
a system are found there. The moments of the expansion are r=l
then matched to a rational Pad6 approximation. These Pad6
where
approximations have several useful properties, one of which is
the convergence of the poles of the approximation to the actual d '
[Y]'=- YI .
poles of a system. There are two approaches for generating 8s' s=s,
the response of a system.
In the first approach [20], [21], several expansion points are For each expansion point, the moments m, = [M,](,);
generated and the converged poles from each expansion are n = 0,1,2,. . . (2q - 1) of an output i are matched to a lower
compared. If two expansions have the same poles, then they order frequency-domain function in the form
are considered accurate within the radius of accuracy defined
by those poles. If no poles are found in common, then an H(s)= X[,](S)= -PL(S)
intermediate hop (expansion) is chosen. All poles within the QM(s)
radius of accuracy of each hop are then collected giving the
actual system poles near the imaginary axis. The frequency and
transient responses are then a closed-form function of these
poles and their residues. For given L and M, the coefficients of the numerator and
In the second approach [22], in order to obtain the frequency denominator of the transfer function are related to the moments
response, a set of rational transfer functions are generated at
by
a minimized set of expansion points. It is the value of these
transfer functions that is compared at points intermediate to rmL-M+l m ~ - ~ + 2...
the expansion points rather than a search for same poles. If
two transfer functions are found to give the same frequency
response at an intermediate frequency point between the two
generating hops then these transfer functions are considered
accurate. If this is not the case, then another hop is cho-
sen between the two expansion points and an expansion is
performed there. It was found empirically that this approach
1
generally requires lesser number of hops compared to the first L " M
approach. This means that the CPU time can be reduced. r
Further reduction was also noted due to the fact that no pole a, = mr-j b,
convergence was required at each expansion point. j=O
Summarizing, CFH ensures accuracy of an approximation
where r = 0,1,. . . ,L and m j = 0 if j < 0.
for a complete frequency range, using multiple expansion
points and corresponding Pad6 approximations at the fre-
C. Binary Search Algorithm
quencies of interest. Additional hops are generated at an
incremental CPU cost above the cost of the first. However, Pad6 approximation is very accurate near the point of
the number of hops typically needed ranges from 2-15, far expansion i.e., s = SO. However, the accuracy of Pad6
KOLBEHDARI et al.: SIMULTANEOUS TIME AND FREQUENCY DOMAIN SOLUTIONS 1.529
closed form function of the generated poles and residues or
tj”
alternatively by using inverse fast fourier transform (IFFT).
D. Moments Generation
To derive an expression for moments M,, we can rewrite
(14) using (16) and (17) at any arbitrary complex frequency
point s = so, and expanding the right hand side of (14) using
Taylor’s series, we get
[(s - s,)~C+ (S - s,)B +A + (S - s,)~s,C
M
+ +
Bs, Cs:] M,(s - so),
n=O
complex ‘s’ plane
Fig. 1. Generation of transfer function by CFH.
(231
approximation decreases as we move further away from the + +
point of expansion similar to the case of a Taylor’s series where R(s,) = B90 &(so), R’(s0) = C(9o @) ; i-
expansion. Q’(sO),R ”(s0) = Q”(so). .. and R”(s0) = Qn(so).
In order to check the accuracy of an approximation, two Equating the coefficients of the powers of (s- so) on both
sides, we get
expansion points are necessary. The accuracies of these two
expansion points can be verified by matching the poles gen- + +
[Cs: Bs, A] M, = R(s,) (24)
erated at these two expansion points [20], [21]. Alternatively,
+ + + +
the two expansion points can be verified for their accuracy by [CS: Bs, AIM1 = -[B 2soC]MO R’(so). (25)
finding the value of the transfer function at a point intermediate
to these two expansion points [22]. Number of expansions Generalizing we have
required to obtain a fairly accurate set of transfer functions
over a specified frequency range is controlled by a binary
search algorithm.
The steps involved in the binary search algorithm for the
second approach is summarized below.
for n 2 2.
Step 1) Set f~ = 0 and f~ = fmax (Fig. 1)
Step 2) Expand the system’s response at frequency f = f ~ .
Determine the coefficients of the corresponding IV. NUMERICARLE SULTS
transfer function HL(s)u sing (21) and (22). Example 1: For purposes of comparison the first problem
Step 3) Expand the system’s response at frequency f = chosen is an example reported in [17]. Consider a one-
fH. Determine the coefficients of the corresponding dimensional (1-D) diffusion problem shown in Fig. 2. The
transfer function HH( s). problem consists of a magnetic slot with bottom and side
i ( f+~
Step 4) Set f = f ~ )C.a lculate H~(j27r.f)a nd walls made of magnetic material with infinite permeability.
H~(j2x.f)I.f , IH~(j2x.f-) H~(j2rf)I< E, The time dependent magnetic field is illuminated on the top
where t represents pre-specified threshold for rel- surface of the magnetic slot producing a time dependent
ative error, GO TO Step 5. Otherwise expand at current flowing through a metal-filled slot of infinite length.
+
fmid = ( f ~f ~an)d d etermine Hmid (3). The excitation of the magnetic field f(t) in (1) on the top
Step 5) If no middle frequency fm;d is generated between surface of the magnetic slot is f(t) = H(t) = eat - ePt,
any two other frequency expansions, STOP, ELSE, which is an EMP type excitation with a = -4.0 x lo6 and
repeat Steps 2-4 between every two consecutive p = -4.76 x lo8. The other parameters are d = 1.0 x
frequency points (e.g., between f~ & finid and (m), ,LA = 47r x (Wm) and (T = 5.76 x lo5 (S/m).
fmid & fH). The boundary condition of the problem is H(t) = 0 at
A similar search algorithm can be used for the first approach z = 0. Fig. 3 shows the normalized time-domain magnetic
L211. field obtained by the proposed method at three different points
At the completion of the binary search algorithm, a set (x = -f , x = 2d , and 5 = a3d) .T he results are compared
of transfer functions are generated. The frequency response with the response obtained by solving the ordinary differential
at a particular frequency point is computed by choosing a equation (4) using conventional implicit numerical integration
transfer function valid for that region. This is repeated for all method [28]. A good agreement is observed between these:
other frequency points and the system’s frequency response two methods. These results are also in agreement with the:
is computed. The time-domain response is obtained as a analytical solution reported in [ 171.
