close

Вход

Забыли?

вход по аккаунту

?

Wire grid modeling for microwave heating and thermal runaway

код для вставкиСкачать
r a a
1i t B
National Library
or Canada
Canadian T h e s e s S ervice
Biblioih&que nationale
du Canada
S ervice d e s th e s e s can ad ien n es
Ottawa, Canada
K1A 0N4
NOTICE
AVIS
The quality of this microform is heavily dependent upon the
quality of the original th esis submitted for microfilming.
Every effort h as been m ad e to en su re the highest quality of
reproduction possible.
La quality de cette microforme depend grandement de la
quality de la t l i s e so u m ise au microfilmage. Nous avons
tout fait pour assurer une qualit6 sup6rieure de reproduc­
tion.
If pages are missing, contact th e university which granted
the degree.
S'il m anque des p a g es, veuillez communiquer avec
I'universit6 qui a conf6r6 le grade.
S om e p a g e s may h ave indistinct print especially if the
original p a g e s were typed with a poor typewriter ribbon or
if the university sent u s an inferior photocopy.
La quality d'impression de certaines p a g es peut laisser a
d^sirer, surtout si !es p a g e s originales ont 6t6 dactylogra­
p h i e s k I'aide d'un ruban u s6 ou si I’universit6 n ou s a fait
parvenirune photocopie d e quality in irieu re.
Reproduction in full or in part of this microform is governed
by the Canadian Copyright Act, R.S.C. 1970, c. C-30, and
subsequent am endm ents.
La reproduction, m em e partielle, de cette microforme est
soum ise k la Loi ca n a d ien n e sur le droit d'auteur, SRC
1970, c. C-30, et s e s am endem ents subsSquents.
NL-339
It.
88/04) c
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
Canada
W IR E GRID M ODELING FOR M ICRO W A V E HEATING
AND THERM AL RUNAW AY
by
Tien D. Pham
A thesis submitted in conformity with the requirements
for the Degree of Master of Applied Science in the
University of Toronto
Department of Electrical Engineering
University of Toronto
Toronto, Ontario
© Tien D. Pham
1990
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
eg
H
S . B
National Litxary
of Canada
Bibliotheque nationals
du Canada
Canadian Theses Service Service des theses canadienoes
O ttaw a. Canada
XIA ON 4
The author has granted an irrevocable non­
exclusive licence allowing the National Library
of Canada to reproduce, loan, distrfoute or sell
cop ies of his/her thesis by any means and in
any form or format, making this thesis available
to interested persons.'
L’auteur a accorde une licence irrevocable et
non exclusive permettant & la Bibfiothdque
nationale du Canada de reproduire, prSter.
distribuerou vendre des copies d e sa th ese
de quekjue manidre et sous quelque forme
que c e soit pour mettre d es exemptaires d e
cette these a la disposition des personnes
interessees.
The author retains ownership of the copyright
in his/her thesis. Neither the thesis nor
substantial extracts from it may be printed or
otherwise reproduced without his/her per­
mission.
L’auteur conserve la propriety du droit d’auteur
"qui protege sa these. Ni la these ni d es extraits
substan tiets de celle-ci ne doivent etre
imprimes ou autrement reproduits sans son
autorisation.
ISBf!
0 - 3 1 5 -6 5 4 4 8 -1
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
THE UNrVERSTTY OF TORONTO LIBRARY
MANUSCRIPT THESIS
AUTHORITY TO DISTRIBUTE
NOTE: T h e AUTHOR will sig n in o n e o f the tw o places indicated. It is the intention of th e University that
there b e NO RESTRICTION on the distribution of the publication of t h e s e s s a v e in exceptional
ca ses.
(a)
Im m ediate publication in microform by th e National Library is authorized.
Author’s signature ____________________ ^ C
Date--------- ------------— lA — (JLj-ii;--
-O R -
(b)
Publication by the National Library is to be postponed until_______________________
(normal maximum d ela y is tw o years).
Author's sig n a tu r e __________________________________________
19
D ate
This restriction is authorized for r e a s o n s which seem to m e, as Chairman of th e G raduate D epartm ent of
______________________________________________ , to be sufficient.
Signature of G raduate D epartm ent C h a irm a n __________________________________________________________
D a t e ________________________________________
B O R R O W E R S undertake to give proper credit for any u se m ade of the th e s is , and to obtain the co n se n t of
the author if it is proposed to m a k e e x te n siv e quotations, or to reproduce th e th e sis in w hole or in part.
A ddress
Signature of borrower
I
I
i
I
-
I
I
I
I
I
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
Date
|
ABSTRACT
The wire grid modeling method based on the thin wire program by J.H. Richmond
has been modified to model both conducting and dielectric bodies. The method is vali­
dated bv comparing the numerical and analytical solutions for the field and power distri­
butions in a dielectric-filled terminated waveguide model. The method is then used to
invest gate the microwave heating of some dielectric materials under various cir­
cumstances. Computation of cower dissipation in thin dielectric bodies with high loss
factors simulates some commonly observed phenomena in microwave heating like edge
heating and thermal runaway, which can be satisfactorily explained by simple surfacewave theory. One-dimensional microwave heating in a half-space has also been studied
for several dielectric materials to give some insights into the thermal runaway effect. A
numerical method is utilized to solve the coupled wave-thermal equation, which appears
in the problem as a result of the temperature dependence of the dielectric constant and
loss factor o f many dielectric materials.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
ACKNOWLEDGMENTS
I wish to express my sincere gratitude to my supervisor Professor Keith G. Balmain
for his expert guidance and enthusiastic support during the course of this work. I must
also thank Dr. Colin C. Bantin and Dr. Mark A. Tilston for many stimulating discussions
and suggestions on antennas and electromagnetic theory. The assistance of Dr. Claude
Lorenson and Christine Gallemeault of Alcan International Ltd. is also greatly appreci­
ated.
I would like to thank everyone in the lab who has helped me in many ways. They
have provided the atmosphere and knowledge that have made studies at this university an
enjoyable experience.
Finally, I would like to thank Professor Balmain and Alcan International Ltd. for
providing financial support and technical materials that have made this research project
possible.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
TABLE OF CONTENTS
INTRODUCTION................................................................................................................
1
1. FUNDAMENTAL PR IN C IPLES...................................................................................
3
1.1. Basic theory of microwave heatin g.............................................................
3
1.1.1. Dielectric lo s s ...................................................................................
3
1.1.2. Power dissipation.............................................................................
5
1.1.3. Dielectric properties..................................................
6
1.2. Thin-wire implementation of the moment method ...................................
9
1.3. Fields in rectangular waveguides and co n tain ers......................................
10
1.3.1. Fields in rectangular w aveguides....................................................
10
1.3.2. Fields in rectangular containers ......................................................
14
1.3.3. Excitation of TE w av es....................................................................
16
2. WIRE GRID MODELING OF CONDUCTING AND
DIELECTRIC BODIES .........................................................................................
19
2.1. Modeling of continuous media using wire g r id s ........................................
19
2.2. Dielectric heating in a small, deep rectangular co n tain er.........................
26
2.3. Dielectric heating in other rectangular containers ....................................
34
3. ONE-DIMENSIONAL MICROWAVE HEATING THERMAL R U N A W A Y .......................................................................................
38
3.1. Wave propagation in inhomogeneous media .............................................
38
3.2. Numerical method for coupled wave-thermal equation............................
41
3.2.1. Boundary condition considerations.................................................
42
3.2.2. Numerical approach .........................................................................
44
3.3. Thermal runaway in one-dimensional microwave h eatin g .......................
48
3.3.1. Example 1: Cooked beef ju ic e ........................................................
49
3.3.2. Example 2: Nylon 66 .......................................................................
51
3.3.3. Example 3: Thawing of frozen water ............................................
53
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
4. TRANSVERSE MICROWAVE HEATING - THERMAL RUNAWAY
59
4.1. Surface wave theory ......................................................................................
59
4.1.1. TE w aves...........................................................................................
61
4.1.2. TM w aves..........................................................................................
63
4.2. Edge heating in infinitely long dielectric strips..........................................
64
4.2.1. The TE case ................................................................................. —-
64
4.2.2. The TM case .....................................................................................
67
4.3. Edge heating in rectangular dielectric s la b s ................................................
4.4. Thermal runaway in transverse microwave heating...................................
68
73
CONCLUSIONS.............................................................................................................
78
.APPENDICES
A IMPLEMENTATION OF ELECTRIC AND MAGNETIC
SYMMETRIES IN THE THIN-WIRE PROGRAM .......................................
79
B COMBINATION OF SYM M ETRIES.....................................................................
83
REFERENCES ................................................................................................................
86
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
1
INTRODUCTION
A widely-used application of the electromagnetic moment method is the wire grid
modeling of two-dimensional or three-dimensional bodies in electromagnetic scattering
problems. By this method, a solid body is simulated by a grid of wires, from which many
electromagnetic characteristics of the body can be determined [1], [2]. The method has
been widely used to compute, for example, radar cross sections of aircraft and fields of
parabolic reflectors. In general, the main use of the wire grid method has been the model­
ing o f metallic-surface bodies, followed by computations of far field radiated from the
bodies.
The application of this method can be extended significantly if one ca^ model
bodies of other materials, rather than only perfect conductors. The motivation for this
idea has its origin in the electromagnetic problems of microwave heating, where one
wants to determine the power dissipation in food when exposed to microwave energy.
Basically, microwave heating problems are similar to scattering problems involving con­
ducting bodies. There are, however, several differences between these two classes of
problems. First, in microwave heating problems, one is interested in finding the near field
or the induced current that gives rise to power dissipation in a lossy dielectric body. Usu­
ally, near field computation is more demanding than far field computation with respect to
numerical accuracy [3]. Second, scattering problems involving conducting bodies are
essentially simpler than microwave heating problems, because the zero tangential electric
field boundary condition can be applied on the surface of a conductor in most practical
cases. The fact that no comparably simple boundary condition exists on the surface of a
dielectric body makes the field computation inside it a relatively difficult task. Neverthe­
less, it will be shown that a dielectric body can be successfully modeled by a wire grid
properly loaded with lumped impedances; furthermore, the wire grid modeling method
offers the unique flexibility of using wire-, surface- and volume elements together to sup­
port the modeling of most practical conducting and dielectric bodies. Therefore the
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
2
method is suitable for the analysis of microwave heating problems, including uneven
heating and thermal runaway in food materials, which are commonly observed in domes­
tic and commercial microwave ovens. Since microwave heating is a relatively new and
unexplored research area, some basic investigations on the runaway effect will be first
conducted for the one dimensional case with the hope that this will provide some insight
into the phenomenon.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission .
3
1. FUNDAMENTAL PRINCIPLES
This chapter contains a review of the theory of microwave heating along with a
short summary of the thin-wire implementation of the moment method. Also, the theory
of a waveguide filled with lossy dielectric will be introduced. These fundamental princi­
ples serve as a basis for the discussions taking place in the chapters that follow.
1.1. Basic theory of microwave heating
It has long been known that an insulating material can be heated by applying
energy to it in the form of electromagnetic waves. The origin of this heating lies in the
ability of the electric field to polarize the charges in the material and the inability o f this
polarization to follow the rapid reversals of the electric field. Coupled w; h the polariza­
tion effects, a dielectnc material can be heated through direct conduction effects due to,
for example, the redistribution of charged particles under the influence o f the externally
applied electric field. The theory of such polarization and conduction effects has been
extensively studied and summarized in the past [4]. In this chapter, however, only a short
and simple formulation will be presented to give the reader a basic understanding of the
theory of microwave heating.
1.1.1 Dielectric loss
At microwave frequencies (300 MHz < / < 30 GHz), dielectric loss in a food
material is due mainly to dipolar polarization in which the permanent dipoles contained
in the material reorient under the influence of a changing electric field. Another source of
polarization arises from charge build-up in interfaces between components in heterogene­
ous systems, termed Maxwell-W agner polarization. Other kinds o f polarization are elec­
tronic and atomic polarizations, which dominate at v e r v high frequency (optical) and are
not involved in microwave heating. For all kinds of polarization, dipole charges (bound
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
4
(bound charges) are induced and the field vector that relates to the bound charges is
termed the polarization P. The electric fiux (displacement) density vector D is given by
D = EoE + P
(1-1)
D = £r EoE
(1-2)
P = (er - l ) e o E
(1-3)
By definition,
which together with eq. (1) yields
In general, cr is a complex constant due to the phase difference between P and E. Furth­
ermore, the polarization vector P lags the electric field vector E, so the relative complex
permittivity must be expressed by
£. = er - j e p
(1-4)
where £r and £p are both positive and real. It should be noted that no conduction effect
has been taken into account so far. For a dielectric having the conductivity a , the total
current density is given by M axwell’s equation
VxH = J = o E + - ^ dr
where c E isthe conduction current density and
3D
(1-5)
isthe displacement current density.
For sinusoidal electric field variations with an em time variation, eq. (5) can be rewrit­
ten as
J = a E + y'o)£o£rE
= ;cd£o (e'r - j Z p ’ - } - £ - ) E
GjEO
(1-6)
M
. . .
In eq. (6), a/co£o represents conductivity loss and ep polarization loss.Since, with most
dielectric measuring techniques, it is difficult to separate the losses due toconduction and
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
5
those due to polarization, it is most convenient to group them together as the loss factor
which, in general, replaces ep in eq. (4).
The ratio of the loss factor to the dielectric constant is the loss tangent
( 1-8)
1.1.2 Power dissipation
Microwave heating involves the conversion of electromagnetic energy into heat.
Energy is transported through space or any medium, by means of electromagnetic waves.
The power flow through a closed surface S can be calculated by the complex Poynting
theorem
s
V
(1-9)
The surface integral on the left hand side represents the rate of flow of energy inward
through the surface enclosing the volume, the first term in the right hand side represents
the stored magnetic and electric energy in the volume and the last term the average dissi­
pated power. The dissipated power density is given by
p = ^coeoei-'E’ -E = y taeo e" | E | 2
(1 10)
or in terms of current density J given in eq. (6)
( 1- 11)
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
6
1.1.3 Dielectric properties
AT
The power dissipated in a dielectric body depends on the loss factor er as well as
the electric field E, which is indirectly dependent on the complex permittivity £r .
Knowledge of the dielectric properties of the material to be processed by microwave
heating is essendal for prediction o f the heating pattern inside the body when it is
exposed to microwave illumination. It has long been known that both the dielectric con­
stant er and the loss factor e r are frequency and temperature dependent A considerable
amount of data for inorganic, organic and biological materials has been tabulated on the
variation of £r with / and T. For dielectric materials with polar structures, the dependence
of £r on frequency / can be explained from the dielectric relaxation of the material.
Debye first advanced a theory which expresses £r and Er as functions o f / i n closed forms
'
**
due to a given kind of polarization [4], The variation of £r and £r w ith /in the entire fre­
quency range is then a combined effect of all kinds of polarization with each of them
being dominant in a particular frequency band. In the context o f this thesis, we are
mainly interested in the dependence of Er and £r on temperature because it has a direct
influence on the dynamic heating process and sometimes can lead to the thermal runaway
phenomenon. Generally, the chemical bonds in a dielectric material are dependent on
temperature and hence change the polarizability and the conductivity o f the material
which in turn change the dielectric constant and the loss factor. O f special interest are
food materials that usually contain a significant amount o f water. Liquid water is strongly
polar in its structure, and therefore able readily to absorb microwave energy and convert
it to heat. Therefore water absorbed in food materials has a marked effect upon its dielec­
tric properties. Bengtsson and Risman [5] showed that there are large differences in £r
and £r between ice and liquid water and the dielectric properties o f most foodstuffs fol­
low a similar pattern to those of water in bulk when plotted as a function o f temperature,
as shown Fig. 1-1. The sharp increase of £r around 0°C suggests a tendency to thermal
runaway when the food materials are defrosted. As soon as a few droplets of liquid water
are produced, the microwave energy is preferentially dissipated in the liquid phase
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
7
causing uneven thawing with damaging effect. The temperature dependence of dielectric
properties of some typical food materials above freezing point [6] is shown in Fig. 1-2.
Because the power dissipated in the dielectric materials also depends on the magnitude of
the electric field, thermal runaway can happen in many cases, where the loss factor
remains constant but the field distribution is shifted toward the regions of relatively high
temperature due to the inhomogeneity in the material. Although the dielectric properties
of a material are definitely o f great importance, other factors such as specific heat, ther­
mal conductivity, mass transfer and especially, the latent heat associated with the phase
changes in a thawing process may also affect the heating process significantly. Due to the
quick heating process in most cases where there is no phase change, it is practically
enough to consider the effects caused by the behaviour of the dielectric properties only.
D istilled
w a te r
80
iked
7Q
[ M ash ed
p o ta to
60
I
I
w
1
50
R aw b eet
A
40
(
30
Cooked
20
C ooked
ham \
^
10
, y ,.—
-20
______
25 -
20
T CC)
20 .
Cooked
v
%j
C ooked / i
ham
/ .
20
-10
0
20
T CC )
60
Fig. 1-1 Dielectric properties of some food materials around freezing point [5]
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
8
•
Cooked
bee!
a
Raw
♦
Cooked b e ef
3 0 0 MHz
9 1 5 MHz
2 4 5 0 MHz
15 0
80
juice
beef
60
100
u
r *—
■
50
.oA —4—-4*~
T
25
45
T ( °C )
65
•
5
C o o k e d t u r k e y j ui c e
Raw turkey
Cooked turkey
25
T
45
( 5C)
65
* > .' •
100
5
25
45
T (°C )
65
45
T CC )
65
Fig. 1-2 Dielectric properties of some food materials above freezing point [6]
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
9
1.2 Thin-wire implementation of the moment method
The numerical solutions to electromagnetic scattering problems can be formulated
in terms of currents induced in the scattering bodies. A solution for these currents is first
sought and then the scattered field is calculated from these induced currents. There are
many methods of solving for the induced currents and some of these methods have been
summarized by seven 1 authors [7], [8]. Generally, these methods are based on formulat­
ing an integral equv ;on • : r th'* induced currents. By expanding the current distribution as
the weighted sum of the . ter; hers o' some basic set, the integral equation can be used to
create a system jF lin
• .’ quatior.s from which the current distribution can be solved
algebraically. This method o f solution is generally known as the moment method.
A special implementation of this method using thin wires was developed and pro­
grammed by J.H. Richmond [2] to analyze antennas and structures of interconnected wire
segments. One of the motivations for the development o f the program was wire grid
modeling of conr.mous surfaces, from which the scattering properties of the surfaces can
be determined e.vily. Such questions as "How dense should the wire grid be to properly
model a given surface ?" or "what is the best wire radius to grid spacing ratio ?" arise
automatically in c:-'section with the modeling of wire grids and have been investigated
by several authors [9], [10], [11]. As a matter of fact, the program’s effectiveness has led
to the publication of many papers concerning the modeling aspects of scatterers. The use­
fulness of the program is due mainly to the straightforward approach it offers the user
when modeling scatterers with bodies of arbitrary shape. In effect, the physical
boundary-value problem relating to the electromagnetic scattering problem is reduced to
a network problem in which the network’s parameters must be found as the impedance
matrix of the scatterer in order to solve for the induced currents. It must be realized that
the range of problems solvable on any real computer extends to bodies not exceeding
several wavelengths in size. This limitation is due to both the cost of computing time and
the limited memory capacity of computers. In an effort to reduce the number o f unknown
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
10
currents in a wire grid model, electric and magnetic symmetries have been utilized in the
program. The concept of ehctxic and magnetic symmetries is best illustrated in Jordan
and Balmain [12]. Hovvr er, a short treatment of these symmetries together with some
practical aspects of t'tei use ivt the program can be found in Appendices A and B. As a
step in the evolution o t h e program and its use, volume modeling of three dimensional
dielectric bodies will be presented in Chapter 2. The use of this new approach in the rest
of the thesis will be concentrated on the computation of power dissipation or near field in
various wire grid models.
As a final note, the original version of the thin-wtire program by Richmond has
been reformulated by M. Tilston [13] to suppress the computation of nonphysical asym­
metric fields and to allow analysis of multiradius thin wire structures.
1.3 Fields in rectangular waveguides and containers
It is seen from eq. (10) that the power density is proportional to the loss factor and
the squared magnitude of the electric field. A knowledge o f the field distribution inside
the dielectric body is therefore essential in determining the microwave heating pattern.
An understanding of the basic principles o f microwave heating can be gained by studying
the electromagnetic field distribution inside a rectangular metallic container filled with
lossy dielectric. This configuration offers the unique advantage o f having a simple
geometry and boundary conditions. In fact, the well established theory of rectangular
waveguides can be applied directly in some cases to give the solution for the electric field
distribution inside the container.
1.3.1 Fields in rectangular waveguides.
The theory of rectangular waveguide can be found in any standard text about elec­
tromagnetics [12], so only a short summary of the theory will be presented here to facili­
tate the discussions that follow.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
11
In general, Maxwell’s curl equations and the wave equations are
VxH = ;coeE
V2E = y 2E
VxE = ~ M i H
(1-12)
V2H = t2H
(1 -1 3 )
e = £o (e'r-jE r)
(1-14)
where
y = V -r to2i j r
Consider the waveguide filled with lossy dielectric shown in
tions in the z direction may be expressed
Fig. 1-3.Assume that varia­
as e -Y' z, where ingeneral j 2= oq + j $ z is the
propagation constant in the z direction, yet to be determined. M axwell’s curl equations
written in rectangular co-ordinates for the dielectric region within the guide become
x
Fig 1-3 Geometry of perfectly conducting rectangular waveguide
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
and the wave equations for Ez and Hz are
+ 72z E z = - co2)i££7
ox2
(1-16)
dy2
(1-17)
Equation (15) can be rewritten to express Hx, Hy , Ex and Ey in terms of Ez and Hz alone.
There must exist either Ez or Hz so that other components do not all vanish. It is con­
venient to divide the possible field configurations within the guide into two sets,
transverse magnetic (TM) waves for which Hz = 0, and transverse electric (TE) waves for
which Ez = 0. Assuming perfect conductivity for the walls of the guide, one can set the
tangential electric fields at the walls equal to zero:
Ex = Ez = 0
at y = 0 and y = b
Ey = Ez - 0
at x = 0 and x - a
Using the separation of variables technique and applying the specified boundary condi­
tions to solve for Ez and Hz from eqs. (15) and (16), we get the complete solution for TM
and TE waves as follows:
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
13
TM waves.
—y2E2oa
_yz
Ex = ------ =— cosAx siaBy e
h
- y zEz0B
_jyZ
Eyy = -------=
h i --- sinAx cosBy e
Ez - E z o sinAx sinBy e Y'z
j(3iCEzQB
Hx = -----h2
sinAc cosBy e
_v
,
(1-18)
-jox.Ez
qA
.
-77
____
p in
//„y = --------=
----- cosAxv sinBy
e"
l2
h:
Hr= 0
n
•1
TE waves.
jw\iHzQB
.
_TjZ
= -------=----- cosAx sinBy e
h2
„
-J' w H zoA
. , n _v
£ yv = ---------=
h2 ----- sinAx cosBy e
Ez = 0
Tz//zoA .
-yz
/ / , = -----=— sinAx cosBy e “
h2
yzHznB
_YZ
/ / v = ----- ;— cosAx sinBy e
(1-19)
, ,
Hz = HzQ cosAx cosBy e
where
A =^ ~
a
£ = —
b
m,n = 0,1,2,... and /n+rz > 0
yz = ^ A 2 + B l + f
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
( 1-20 )
( 1- 2 1 )
14
(1-22)
h 1 =7^ + C02p£
For waveguides filled with lossless dielectric, yz is either imaginary or real, which
corresponds to propagating or cutoff (attenuating) modes, respectively. For waveguides
filled with lossy dielectric, yz is always a complex constant, which means every mode
propagates and attenuates simultaneously and there is no sharp division between pro­
pagating and cutoff modes.
1.3.2 Fields in rectangular containers
A rectangular container can be considered as a rectangular waveguide terminated at
z = l as shown in Fig. 1-4. Due to reflection at the bottom wall, a TM/TE wave travelling
in the positive z direction will result in a corresponding wave travelling backward in the
negative z direction. The field distribution within the container can be found as sum of
these two waves. Recall that in the derivation of TM/TE waves above we have made the
assumption that the z-variation of the field within the guide can be expressed as e " , so
z
x
a
Fig. 1-4 Geometry of rectangular container: side and bottom walls
are perfectly conducting, while the top (z = 0) is open.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
15
that backward travelling waves can be represented by assuming a variation o f the form
e lz~. It is easy to see that the transverse (x- and y-) variations must be identical for for­
ward and backward travelling waves, at least to satisfy the zero tangential electric field
condition at z = / for all values of x and y (0 <x <a, 0 <y< b). With this in mind, the zvariation of every E or H component can be expressed as
Z(z) = e~y' z + T e y’z
(1-23)
The reflection coefficient T can be found by applying the boundary condition Ex = 0 at
z = /. Then we have
Z ( l)^ -e ~ y‘l + T e y‘l = 0
which gives
T = -e ~ lyi1
(1-24)
Thus, thefield distribution within the container due to a particular mode is the same as
given in eq. (18) or (19) with the z-variation e~y,z beingreplaced by Z(z) given
in eqs.
(23) and (24). The total field is then a combination of all the modes existing within the
container. When the container is illuminated by some external electromagnetic source,
different modes will be excited depending on the geometry of the container, the material
filling the container and the illumination of the container. Determination o f the magni­
tude for each mode would be a complicated task. It is well known from the theory of
waveguides that only a few lowest order modes can propagate in a given waveguide due
to the cutoff state o f the higher order modes. In a container, which can be regarded as a
terminated and usually short waveguide, the evanescent field due to the cutoff modes can
be significant at the opening of the container and thus cannot be ignored. It turns out,
however, that for containers with small cross-sectional dimensions, a properly chosen
excitation source actually can generate only the lowest (dominant) mode.
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
16
1.3.3 Excitation o f T E waves
From eq. (19), it is seen that TEm0 wave has the sole electric component Ey given
by
Ey = Ey o sinAx e~^
(1-25)
where Ey o is a constant If the origin of the coordinate system is placed at the center of
waveguide’s opening, we have
Ey =EyQ cosAx e~y‘
(1-26)
For lossless dielectric, y2 = _/'pz, where Pz is a real constant, and eq. (26) becomes
Ey = —^-[exp(-yAx - j $ zz) + expfjAx - ypzz)]
(1-27)
Furthermore, it follows from eq. (21) that,
p 2= c o V
so eq. (27) can be
= > i2 + P z
C1'2 8 )
rewritten as
£
Ey = —^-[exp(-jP(nrsin0 + zcosG)) + exp(-jP(-xsin8 + zcosQ))]
(1-29)
where
0 = ±sin_1(-^-) = ±sin_1( - ^ - )
P
flP
(1-30)
Eq. (29) represents two uniform plane waves propagating in the dielectric medium with
propagation constant P and direction angles 9 given by eq. (30). Thus, any TE^o wave
can be generated by exciting two plane waves propagating in the given directions. These
plane waves can be excited by some incident field coming from outside the dielectric
medium as in practical microwave heating configurations. Consider a plane wave
incident upon the air-dielectric boundary at angle 9,- shown in Fig. 1-5. The plane wave
enters the dielectric medium at angle 9 determined by Snell’s law
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited without perm ission.
17
Fig. 1-5
sin0 ;
c
v1
sine =n = 7 7 = aj7p
( 1 ' 3 1 )
where v \ and v 2 are the velocity of plane wave in air and dielectric, respectively.
From eqs. (30) and (31), we get
-
. A
c mn
cm 7t
sinG; = — ——-— =
c o /(3
p a
,, ^
(1 -32)
co a
To generate a particular TEOTo mode, the incident angle 0,- is thus not determined by the
permittivity but only by the dimension of the container and the frequency. Strictly speak­
ing, Snell’s law for ray optics only applies for an infinitely extended boundary. For a
finite boundary as in the actual case of a rectangular container, the refracted field is not
only the plane wave with refraction angle 0 but rather a spectrum of plane waves entering
the dielectric medium at different angles [14].
For lossy dielectrics, a TE,„o wave can still be expressed as sum of two plane
waves. However, these plane waves are non-uniform, that is, the field magnitudes are not
uniform but exponentially decaying in one direction along the wave fronts. The refrac­
tion angle 0 is no longer given by eq. (30) but by a more complicated expression [15].
For most practical lossy dielectrics with relatively small loss tangent, however, eq. (30)
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
18
still gives a good approximation for 0 .
Despite these limitations, as we will see in Chapter 2, the excitation method dis­
cussed above is fairly effective in generating a given mode in a given container although
it is in general more difficult to generate a pure high order mode than to do so for a lower
one. Finally, it can be shown that any T E ^, mode can be expressed as sum of four plane
waves and hence can be generated in the same manner as above.
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
19
2. WIRE GRID MODELING OF CONDUCTING AND DIELECTRIC BODIES
The idea o f using a wire grid model to represent a dielectric body has its origin in
the established concept that a continuous perfectly conducting surface can be success­
fully modeled using a surface wire grid [1]. When illuminated by some electromagnetic
source, the scattered field produced by the wire grid approaches the true field produced
by the continuous surface if the wire grid spacing is chosen to be small enough. Since no
field or current can exist inside a solid perfectly conducting body, it is electromagnetically equivalent to an enclosed perfectly conducting surface and hence can be modeled
approximately as a surface wire grid. On the other hand, if the body is dielectric, fields
and currents vtithin the body are nonzero which in turn necessitates a volume wire grid
model for the body as will be discussed in this chapter.
2.1 Modeling o f continuous media using wire grids
Consider the two-dimensional case where a perfectly conducting plate, shown in
Fig. 2 -la, is excited by some external source. The induced surface current is in general a
vector quantity with x- and y-components, so a wire structure representing the plate must
have both x- and y-oriented wire segments, as shown in Fig. 2 -lb. A rectangular patch of
the plate is thus represented by a rectangular wire loop of the wire grid. As mentioned
before, the choice of wire grid spacing d and wire radius r has led to many discussions in
the past [9 ]-[ll], [30]. Intuitively, one can expect that the smaller the wire grid spacing
becomes, the better the wire grid should represent the continuous surface. However, in
practical modeling, the number of wire segments must be limited and hence there is a
maximum allowable wire grid spacing at which the wire grid will still represent the sur­
face satisfactorily. It has been generally found that d = \I 10, where X is the free space
wavelength (at the actual frequency), is adequate to give good results [11], [30]. In many
cases, however, a larger wire grid spacing can be tolerated, practically up to X/4 in the
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
20
b)
a
Fig. 2-1 P-rfecly conducting plate and wire grid model
cases where the current distribution is not undergoing quick changes throughout the sur­
face. On the other hand, no definite answer has been given regarding the choice of the
wire radius. Most authors, however, agree that a radius o f r ~ d H k would be the best
choice. Ludwig [11] suggests that the rule r —d t2 k must always oe followed and he
shows quite cleariy that too large a radius is just as bad as a too small a radius.
Consider the solid three dimensional dielectric body shown in Fig. 2-2a. The body
can be thought of as being made up o f elementary dielectric cubes, one of which is out­
lined in the lower left part of the body. When an electromagnetic field is incident upon
the body, the induced volume currents require a volumetric wire structure having x-, yand z-oriented segments: such a structure with equal-length segments is shown in Fig. 22b, in which part of the structure outlines tne elementary cube. At mis stage, however, it
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
21
is still unclear how the wire grid should be modified to model the body with a given per­
mittivity. To answer this question, let us look at the logical build-up o f the model as fol­
lows. Consider an elementary dielectric cube shown in Fig. 2-3a, which can be
represented by the elementary wire cube shown in Fig. 2-3b. Since the body is in general
a lossy dielectric, conduction and displacement currents actually flow within the body
and the dielectric cube cannot be correctly represented by the wire cube unless it is
loaded with lumped impedances to take into account the conduction and displacement
currents. The idea of using lumped impedances to represent the loss and phase shift asso­
ciated with these currents is similar to using lumped impedances to model a transmission
b)
a)
Fig. 2-2 Dielectric body and wire grid model
d
a)
b)
Fig. 2-3 Elementary dielectric cube and its wire cube model
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
22
line approximately.
If the material has the complex permittivity e = (er - j e r )eo, the effective conduc­
tivity is defined by
a e^ = E r, e 0 CL)
(2 -1 )
Furthermore, if the x-, y- and z-components of the volume current are uniform over the
respective cross sections, the resistance related to a current component, say Ix, is
R=
d
<5efjA
1
Geffd
1
£ r EgOki
(2-2 )
Similarly, the capacitance is given by
C = ( e ; - l ) e o — = (er - D M
a
(2-3)
Note that the elementary wire grid is physically located in free space with relative per­
mittivity equal to unity in comparison with er of the real medium so, the "additional"
capacitance must be given as above.
The resistance and capacitance must be distributed to the four segments oriented in
the x direction, which gives a total load for each wire segment to be 4R in parallel with
C /4. In the thin-wire program, it is only possible to load a segment at the ends and in
order to avoid asymmetry, the total lumped impedances must be distributed equally at the
two ends of each segment as shown in Fig. 2-4. It follows that
4R
HH
C/4
R.
HH
C,
Fig. 2-4 Total load and equivalent end Lads for an edge segment
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
Re
HH
C,
23
(2-4)
and
(2-5)
The above resistance and capacitance are connected in parallel because the conduction
and displacement currents add in Maxwell’s equations. Thus a lumped impedance is gen­
erated and this impedance is inserted at each end of each segment in the elementary wire
cube. The dielectric body shown in Fig. 2-2a consists of five elementary dielectric cubes
and it can be modeled by joining five elementary wire cubes together as shown in Fig. 25. Consequently, groups o f two, three or four segments o f adjacent elementary wire grids
are merged into single segments which we will call surface -, comer - and volume seg­
ments, respectively, corresponding to their positions in the resultant model. Those seg­
ments that are not merged with any other segments are termed edge segments. The
lumped loads for the first three kinds of segment are the parallel combinations of the
edge impedances given in eqs. (4) and (5). Of course, one can only tell the differences
among these kinds of segment by looking at their respective lumped loads. Without the
lumped loads, all segments will be physically identical. In terms of edge impedances, the
surface -, comer - and volume impedances are given by
( 2 -6)
a
dielectric body is thus completely represented by a wire grid model loaded with
lumped impedances, where the wire segments allow calculation of scattered field from
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
24
/
\
Surface segment
w ith/?,, C,
Edge segment
with/?,, C,
Com er segment
w ith Rc, Cc
Volume segment
w ith/?,, Cy
Fig. 2-5 Building up the wire grid model of the dielectric body shown in Fig. 2-2a
by joining elementary wire cubes.
the body and the lumped impedances allow calculation of loss and phase shift within the
body. We will concentrate on modeling bodies using cubic elementary wire grids in this
thesis. Non-cubic but still rectangular building blocks can be modeled easily following
the same principles. Further study of non-rectangular building blocks is required before
bodies of arbitrary shape can be modeled.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission .
25
(0.5,0.5, 1.0)
— source dipole
6
/
( 0 .01 , 0.01 , 0 .01 )
/
■
• / /
center
/
/
/
?. = 0.12 m at / = 2.45 GHz
(0.00, 0.00, 0.00)
Fig. 2 -6
-10
-15
-20
tu
S.
f=2.45 Ghz
-25
iu
d=10 mm
-30
-35
^0
-45
0.4
0.6
0.8
2.4
1.2
r[mm]
Fig. 2-7 Normalized field magnitude vs. wire radius
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
26
Before going over to computation with various models, let us try to confirm the
r -dl2.il rule as the best choice of wire radius in volume wire grid modeling. It should be
noted that an unloaded elementary wire cube represents a perfectly conducting block,
within which the field should be equal to zero. Consider the configuration shown in Fig.
2-6 where the wire cube is illuminated by the sinusoidal field from a source dipole placed
some distance away. The field magnitude at the center of the wire cube is computed and
plotted against the wire radius as shown in Fig. 2-7. Note that the field has been normal­
ized to the incident field, that is, the field that would exist without the wire cube, to give
some idea of the effective shielding of the grid. The normalized field takes a minimum at
r = 0.0012 m, which gives dir = 8.33 which is somewhat larger than the expected ratio
of 2k. Aversa [3] conducted a similar experiment with a larger hollow box and he found
the best spacing-to-radius ratio to be about 7. Since the radius step in his experiment was
big, the minimum could as well be somewhat larger or smaller. However, doubling or
halving the "best" radius definitely gives worse results. As a conclusion about the choice
of wire radius, the rule d =2kr should be followed but only approximately. In general, the
ratio d/r = 2 k will be used through the modeling work contained in this thesis.
2.2. Dielectric heating in a small, deep rectangular container
Having established the fundamentals for wire grid modeling of materials, we wish
to validate the approach by modeling and computing the field distribution within some
rectangular containers filled with lossy dielectrics. Besides having the practical value of
objects that could be used in microwave heating, rectangular containers seem to offer a
particularly useful way to validate the method through comparison of numerical results
with analytical solutions known from the theory of rectangular waveguides.
Consider the rectangular container filled with lossy dielectric shown in Fig. 2-8.
The container must be chosen so that only the lowest modes can propagate while higher
modes are suppressed. Ideally, we wish to have only one mode within the container so
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
air-dielectric boundary
b/2.
perfectly conducting container
z
Fig. 2-8 Geometry of container filled with dielectric
that the analytical solution is exactly known and can easily be compared with the result
deduced from the wire grid modeling approach. It is generally known from the theory of
waveguides that waveguides with sufficiently small cross section can support only the
lowest TEjo mode, while being cut off for the higher modes. This can be seen from eq.
(1-27), where the propagation constant yz is essentially a real constant for high order
modes, which means these modes attenuate rather than propagate along the z direction.
From this point o f view, it is desired to choose a and b as small as possible; a should be
large enough to support only the T E 10 mode and b can be chosen arbitrarily small.
Mathematically, we have
v -
I,
V a
s?— TTrnr?----T T -’— "7 T
b
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
28
U b < a, the cutoff condition for all modes except TEio requires
( - ) 2 - cd2|1£^.£o < 0 < (— )2 - c o V ; e o
a
a
(2-8)
or
k
2k
— — — <a< — — —
wV(i£r eo
“ VliCrEo
(2 -9 )
For butter with er = 3.5—j 0.9 at / = 2.45 GHz, a can thus be chosen in the interval
0.0327 m < a < 0.0654 m. We choose, however, a = 0.07 m in order to have more data
points for the evaluation of the numerical result. With this value o f a, the TE 30 mode will
be cut off but the TE 20 mode will not. By using symmetrical excitation, however, all odd
TEmo modes, that is TE2o, TE 40, ... will be suppressed and the dominance of the TE 10
mode is ensured. The depth I of the container must be large enough to show the field
variation along the z direction. For the TE 10 mode, the propagation constant is deter­
mined from eq. (7) with m — 1, n = 0 to be yz = 13.79 + j 86.05 m-1, which gives the
wavelength of the TE 10 wave Xz — 0.073 m. We choose / = 0.05 m. The wavelength o f a
plane wave in an infinite dielectric medium is given by
X = y ,
(5 = l m ( y } = I m { ^ J I ^ }
(2 - 10 )
which gives X = 0.065 m. The grid spacing can be chosen to be d = 0.005 m, which is
about X/13. Now, the last dimension b can be chosen as small as 0.02 m. In fact, a
smaller value of b could be chosen but this would not allow us to see the field variation
along the y direction, which is necessary to determine the existence of certain modes
within the container. With the chosen values of a, b and /, the size of the wire grid model
for the container becomes (14x4x1 OK As the origin of the new coordinate system is
placed at the center of the container’s opening, the container and the filling dielectric
material are fully symmetric about the x and y axes, so that electric and magnetic sym­
metries can be utilized to reduce the number of wire segments and dipole current modes.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
29
The fundamental theory of these symmetries is treated in Appendix A. Hence, only a
quarter o f the full wire grid model needs to be set up. Fig. 2-9 shows the quarter wire grid
models for the container and the filling dielectric separately. Note that the segments at
the edges, on the surfaces and inside the dielectric model are loaded with edge-, surfaceand volume impedances respectively and are determined from eqs. (4)-(6):
Re = 3260 Q.
Ce = 0.55 pF
Rs = 1630 Q
Cs = 1.10 pF
Rv = 815 f t
Cv = 2.20 pF
The effect of enclosing the dielectric with a perfectly conducting container is that the
lumped impedances on the segments shared by both models are shorted. This is a direct
consequence o f the fact that wire segments are themselves perfectly conducting and as
building elements for the container, they are unloaded. The complete dielectric-in­
container wire grid model is shown in Fig. 2-10. Also shown in this figure is the source
monopole oriented in the y direction. Using ME-symmetry in the XY-plane, the source
configuration becomes two dipoles placed symmetrically about the X = 0 plane at a dis­
tance R = lm from the center of the container’s opening. This distance is large enough so
that the radiation field from the dipoles approximates two y-polarized plane waves
incident at angle 0,- to the z axis. To give a proper excitation for the T E ^ mode, the
incident angle 0/ is chosen according to eq. (1-38):
„
. _i,cm it,
. _i,
3T08T
n_ , c , o
-------- ) = sin ------------ 3------- ) - 6 1
co a
2-2.45T0 -0.07
0 , = sin
The wire grid model of a dielectric in a container together with the excitation source
form a complete wire configuration as input for the thin- wire program. A fter the induced
currents on each segment have been computed, the near field at a particular location can
be determined as the sum of the radiating fields from all the segments. In general, the
field values are most reasonably sampled at the centers of the wire cubes since these are
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
30
//
t
4
V/
/ f/
A/
f,/
V/
V*
V/
Vt
■
dielectric
-
T /T
7
7
7
7
L a L 1 7 i i l J z. 7
i \ I 2 ? 2 i 2 2 I 7 2 7_7
i I 2 2 7 \ 2I I I 27 7
2I I 7 17 7
t I i. 2 I 2 7
7
I
I 7 77
AA 2 t t
1I I l ? 7
L AI 2 1 I
2I 2 l ?7
LA 2 I I
AC
A 1 1 1 J 2I J 2 7
7
r
t t
* / /
/ / /
/
/ /
I
I
II I I
II
I
I I I
I
I?
I I
5
.A 1 E EEEE
/.J
container
d=0.005 m
a=0.07 m
b=0.02 m
1=0.05 m
Fig. 2-9 Quarter wire grid model for container and dielectric
source monopole
R
6,
M
a/2/
f - -
0 ; = 61 °
IL ACr '
a <- AL
i l l t j l V- I I I I I I
Q l
in
H l V J
11 11
■H
R=1 m
IV
V
Fig. 2-10 Complete quarter wire grid model for dielectric in container
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
31
the building blocks of the model. Fig. 2-11 shows the variation of the \Ey | components
along the x axis for y = 0.0025 m and z = 0.0025m; the field values are normalized to the
value at y = 0 in order to be compared easily with the cosine distribution shown in the
solid curve. The numerical results follow the cosine distribution closely, indicating that
no higher order modes of m > 2 are present in the actual field. It follows, therefore, that
only the TE10, TEH and TM U modes are contributing to the total field. Fig. 2-12 shows
| Ey | and | Ez | as functions of x for different values of y and z = 0.0025 m. The Ex com­
ponent is negligible and therefore not plotted. Clearly, the Ey is the dominant component,
indicating that only a small portion of the TMu mode is present. Furthermore, the small
variation o f Ey along the y direction shows that the T E n mode is also insignificant.
Together, these modes represent no more than 1/10 of the total Ey component. Fig. 2-13
shows the variation of normalized Ey component along the z direction. It appears that the
numerical result displays a stronger attenuation than the theoretical standing wave pattern
for the TEio mode shown in the dashed curve, although the deviation is not large. A mix­
ing o f TEio and T E n modes by the ratio 10:1 and phase difference <{>= 7t produces a pat­
tern matching much better with the numerical result. This is no surprise since the T E n
mode is cut off and hence attenuates more than the TEio mode. A contribution of the
T E U mode to the total field is therefore likely to intensify the attenuation. The numerical
result obtained is in fact very good since the propagation and the attenuation presented in
the standing wave pattern is very sensitive to any change of the complex pennittivity e.
Tests have been done with other dielectric materials and the numerical results agree quite
closely with the analytical ones.
Power dissipation in each wire cube can be determined by simply summing the
ohmic losses on the lumped resistances belonging to the ware cube
P = Z T ( «i l / i I +*«2/?2 )
(2-11)
i=l 1
where —R[\ l}\ and —rt /2
are the ohmic losses at the two ends of i ’th segment of the
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
32
0.9
lEyl/lEyOl
0.7
0.6
0.5
: cos(x*pi/0.07)
0.4
0.3
0.2
(x): numerical result
y =0.0025m , z=0.0025m
0.005
0.01
0.02
0.015
0.025
0.03
0.035
X [m]
Fig. 2-11 Normalized field distribution Ey at the opening of
rectangular container filled with lossy dielectric
50
45
Field Magnitude (V/m)
40
35
30
25
20
15
10
(x): lEyl for y=O.025m. z=0.025m
(o): lEyl for y=0.075m, z=0.025m
(*): lEzl for y=0.025m. z=0.025m
(+): lEzl for y=0.075m. z=0.025m
5
0
0.005
0.01
0.02
0.015
0.025
0.03
0.035
X [ml
Fig. 2-12 Distributions o f some field components in container
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission .
33
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
- : Analytical solution : TEIO
: Analytical solution : TEIO + 0.1*TE11
(o): Numerical result for v=0.0025m , z=9.0025m
Z 1m
Fig. 2-13 Standing wave pattern in rectangular container filled with
lossy dielectric er = 3.5 —j 0.9
y ^ .0 0 2 5 m , r=0.0025ra
0.9
0.8
(x): numerical result
0.7
: (cos(x*pt/0.07))
0.6
o
0.5
0.4
03
0.2
0.1
0.005
0.01
0.02
0.015
0.025
0.03
0.035
X [ml
Fig. 2-14 Power distribution at the opening of the container
I
I
\
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
34
block, and c, = 1, 1/2, 1/3 and 1/4 for edge-, comer-, surface-, and volume segments
respectively. This factor is to take into account the ohmic losses shared by adjacent wire
cubes with common segments and lumped resistances.
Fig. 2-14 shows the distribution of normalized power dissipation along the x axis.
The numerical result agrees well with the theoretical cosine-squared distribution
although slight deviations are observed for small values of PIP o- The most probable
explanation is that the contribution of the Ez component to the total power dissipation
becomes apparent when the dominant component Ey is small. There is likewise good
agreement between numerical ar.d analytical results for the power distribution along the z
axis but the result is not shown.
2.3. Dielectric heating in o th e r rectangular containers
The container considered in the previous section is only a special case, in which the
dimensions have been chosen together with a special source configuration to support only
the lowest TE mode. In practice, containers can have various shapes and sizes, and are
usually exposed to a very intricate microwave field such as the one in a commercial
microwave oven, which can have many polarizations and magnitude distributions. The
power dissipation in the dielectric can therefore have a very complicated pattern due to
the complex field distribution within the container. Consider for example the container
and source dipoles shown in Fig. 2-15. The dimensions of the container have been
chosen to support several modes beyond the lowest one. The power density distributions
in the XY plane for different z values are shown in Fig. 2-16. It is noted that the plots are
of different scales to demonstrate the shape variation of the power absorption along the z
axis. The power distribution at z = 0.25 cm is seen to have two main lobes at x = 0 and
four smaller lobes at the comers. It appears the shape is formed by a mixing of T E ^ ,
TE/TMU , TE/TM12) TE/TM 31 and TE/TM 31 modes with TE/TM n being the dominant
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
35
2cm source dipole
10cm
E -s-
lm
lm
7cm
.5cm
Container
Dielectric
fU ^ rv J J
Fig. 2-15 Heating of dielectric (er=3.5-/0.9) in (20xl4x3)d container
with d = 0.5 cm being the grid spacing for the model
Fig. 2-16 Power density distributions at different depths inside the
dielectric in container shown in Fig. 2-15
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
36
components. As different modes have different propagation and attenuation constants,
the shape of the power distribution changes along the z axis. Symmetries are obtained
about the x and y axes because only one source configuration ME is used in this case.
Non-symmetric
power distributions may result by combining
different
source
configurations (see App. B).
One of the motivations behind microwave heating research is to obtain heating uni­
formity in food materials. If one could generate any mode with desired magnitude and
phase in a given container, it would be possible to synthesize a uniform heating pattern
by superimposing different modes with properly chosen magnitudes and phases. The
difficulty is, o f course, how to generate these modes. Even for a container large enough
to support a certain high order mode, it is difficult to control the mode excitation because
of the diffraction over the edges of the container. Consider the container shown in Fig.
2-17. The model and source configuration is similar to that shown in Fig. 2-10. The
source dipoles are at distance i? = l m from the center of the container’s opening and
several incident angles 0,- are chosen in the interval [0°,90°]. The heating patterns along
the x axis at y = 0.00375 m and z = 0.00375 m for different incidence angles are shown
in Fig. 2-18. It should be noted that 0/ = 17° and 0,- = 61° are incidence angles for TEio
and TE 30 modes, respectively, according to eq. (1-38). The heating patterns for these
angles exhibit approximately single mode field distributions as expected. At other angles,
the heating patterns show a mixing of T E ^ and TE 30 modes with different amplitude
ratios and relative phases. It is interesting to note that for 0X
- = 0°, the contribution of the
TE 30 mode is quite significant while it completely vanishes at 0,- = 17°. For 0; = 90°, the
power absorption in the dielectric arises totally from diffraction of the incident field over
the edges of the container. Further study of field transmission in containers made of con­
ducting and dielectric materials seems necessary to optimize the heat distribution.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
37
M
2cm source dipole
R=lm '
21cm
E
3cm
2.25cm
dielectric
Fig. 2-17 Heating of dielectric (er=3.5- j 0.9) in (20x4x3)d container
with d = 0.75 cm being the grid spacing for the model
0.02
0.018
0.016
0.014
m
0.012
E
•Si
a.
0.01
0.008
0.006
0.004
0.002
-0.15
-
0.1
0.05
-0.05
0.15
X [m|
Fig. 2-18 Power density distribution inside the container shown in Fig. 2-17
for different incidence angles.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
38
3. ONE-DIM ENSIONAL M ICROW AVE HEATING - TH ERM A L RUNAWAY
Non-uniform distribution of heat in dielectric materials when processed by
microwave power may be caused by the limited penetrability of the electromagnetic
wave as well as by standing wave or other effects. This non-uniformity is sometimes
intensified by the inhomogeneity of the material where the above effects are reinforced.
Due to the temperature dependence o f the relative permittivity, the heat generated within
the material can change the dielectric structure in such a manner as to concentrate the
power dissipation more and more into the hottest regions, which leads to thermal runa­
way with damaging effects for the material. Besides the problem of determining the field
distribution in an inhomogeneous material, the mutual dependence of the power dissipa­
tion and temperature through the dielectric properties obviously contributes to the com­
plexity of the microwave heating analysis. It seems therefore natural first to investigate
the problem of one-dimensional microwave heating, which is probably an over­
simplified case but nevertheless will give us some understanding and insight into the
phenomenon of thermal runaway.
3.1. W ave propagation in inhomogeneous media
Consider a uniform plane wave incident upon the surface o f a dielectric half-space
where z > 0 as shown in Fig. 3-1. The relative complex permittivity o f the dielectric
medium is assumed to be independent of x and y but is a function of z, i.e.
z(z) = (er( z ) - j z ' r (z))£o
(3-1)
Also, assume the only component o f the electric field is E = Ey . Since there is no field
variation in the x-y plane ( -^ - = - ~ = 0 ), the wave equation given by eq. (1-13) is
ox
dy
reduced to
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
39
E(2)
eo
E:
propagation
z=0
Fig. 3-1 Geometry for one dimensional microwave hearing
d 2£
2 - y 2(z ) £ ' = 0
y2(z) = -co2 ^ 0 e(z)
(3-2)
This equation for wave propagation in inhomogeneous media can be solved by various
approximation methods. The conventional WKB method is most convenient for finding
the solution for slowly varying media. Let
E(z) = A(z)e~i m
(3-3)
where A and 4> can be regarded as the amplitude and phase, respectively.
Substitution into the wave equation gives
A" - 2y<t>'A' - ;<(>"A + [—"y2 (z) - ((j)')2]A = 0
(3-4)
where the prime indicates a derivative with respect to z. For slowly varying media, we
make the approximation
a" = 0
Then A e ^
isan approximate solution of eq. (4) if
<j> = ± 7 t(z )
or
0 = ± yJ t z ) d z
and
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
(3-5)
40
or
A=
const
(3-6)
■ fT W f
then, the general solution is
a exp - J V
E(z) =
(z )
dz + b exp + J y ( z ) dz
'F J W )
(3-7)
o
where a and b are arbitrary constants. Let y(z) = a(z) + y‘P(z), the integrals in eq. (3-7)
representing the propagation and attenuation of two waves travelling in the positive and
negative z directions. Since the medium is infinite for z > 0, only the first wave exists.
Furthermore, additional amplitude and phase variations of this wave are given by the fac­
tor outside the brackets. The field amplitude is given by
2
\E(z)\ =
V |i(z )-y a (z )
exp - J a (z)dz
(3-8)
d$
It is interesting to note that for a low-loss dielectric, a c (J and if - j - < 0, we could have
a growing wave even in a lossy medium, at least at some short distance from the surface
where the decaying characteristic of the factor in the denominator is dominant over that
of the exponential function. On the other hand, if
dz
> 0, then we will have another
attenuation factor beyond the one caused by the loss of the dielectric. This attenuation
can be regarded as being a result of the internal reflection at every point along the z
direction caused by the growing characteristic of the permittivity. The generalized WKB
method is well described by Wait [16], who gives a more exact solution for the wave
equation. Anyway, even this modified method must be based on the assumption of slowly
varying media. For rapidly varying media, other methods are available [17], but they are
less compact and not necessarily numerically efficient A numerical approach based on
the finite difference method will be utilized in the following sections to solve practical
problems of one-dimensional microwave heating.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
41
3.2. Numerical method for coupled wave-thermal equation
Microwave heating is usually a dynamic process, in which most of the parameters
are functions of space and time. The field distribution in the dielectric medium has a
direct impact on the temperature distribution through the power dissipation given in eq.
(1-16). On the other hand, the field distribution depends on the dielectric profile, which is
a function of the temperature distribution. The time dependence of the electric field mag­
nitude and dielectric properties comes into the picture as the temperature itself is time
dependent. The transient electric field need not be considered since it practically will not
have any impact on the slow thermal process.
The one-dimensional thermal equation is given by
dT
d2T
p
3or7 = a'T9 zTz + TT"
pcp
(3-9)
where
T = T ( z , t ) : temperature [°C]
p = p (z,r) : power density [W/m3]
~
cq = ----- : thermal diffusivity [m /r]
PCp
Xc : thermal conductivity [W/(m°C)]
p : mass density [kg/m3]
cp : specific heat [J/(kg°C)]
The first term on the right hand side of eq. (9) represents the rate of change of tempera­
ture by thermalconduction and the second one represents the heatgeneration. Normally,
thermal conduction can be ignored due to the quick microwave heating process. Since eq.
(9) presents no major problem when solving numerically, thermal conduction will be
included for the sake o f generality.
Assume the temperature dependence of the relative dielectric constant and loss fac­
tor is given by
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
42
Er( z , t ) = f ( T ( z , t ) )
and
er = g ( J (z,r))
(3-10)
Combining eqs. (2),(9),(10) and (1-16), we get the coupled wave-thermal equations
dz
+ co2 Moeo[f(T) - jg (D J£ = 0
(3-11)
(3-12)
3.2.1 Boundary condition considerations
Eqs. (11) and (12) are second-order differential equations, which can be solved
with tv/o specified boundary conditions for each equation. First, consider the wave equa­
tion for the half space problem shown in Fig. 3-1. Theoretically, the wave equation can
be solved by the Runge-Kutta method [18] when the boundary conditions £ (0 ,r) and
E (0,r) are given. For a homogeneous medium (for which the exact solution is known),
E (0 ,0 and £ (0,0 can be found if £,• is given, as
£(0,/) = r£f
(3-13)
where T is the reflection coefficient at the surface.
£ (0 ,r) = - y £ ( 0 ,r)
(3-14)
The second boundary condition is valid as only a decaying wave is assumed. Of course,
there is no reason to solve the equation numerically when its analytical solution is
already known, but this discussion will outline the problems associated with the wave
equation. With the given boundary conditions above, the wave equation cannot be solved
in practice because of numerical instability. Recall that the complete solution for the
wave equation is a linear combination of decaying and growing waves. With the boun­
dary condition given in eq. (14), only a decaying wave should be selected but the inevit­
able round-off errors in the numerical solving procedure will excite the growing wave.
Eventually, this will become the dominant component causing an erroneous solution with
wild oscillations. Alternatively, one can use the finite difference method, implicit scheme
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
43
[19], to suppress the instability but the solution is of poor accuracy. At this point, it
seems that any numerical method using the above boundary conditions cannot yield a
satisfactory solution. On the other hand, as the actual medium is lossy dielectric, the field
at infinity vanishes; therefore, practically, the field is negligible at some distance z = d,
which can be estimated as follows
(3-15)
which gives
(3-16)
where 8 is some small constant and
a min= m m [a(z)]
0 < z< °°
(3-17)
Intuitively, the attenuation in any medium with a(z) > a , ^ for all z must be larger than
that in a homogeneous medium with a(z) =
for all z. Thus, eq. (16) gives the dis­
tance at which a vanishing field is ensured. With the boundary conditions E (0,r) and
E( d, t ) = 0, the wave equation (11) can be solved using the finite difference method and
the solution is guaranteed to be stable with good accuracy.
The thermal equation (12) also can be solved using the finite difference method
with the initial conditions
(3-18)
7(z, 0)=7~o(z)
and the boundary conditions for the thermal isolation case
r'(0 ,r) = 0
and
T{d, t) = 0
(3-19)
or, for the case of convection loss at the surface z = 0
T (0,t) = -^-[Tan - 7 (0 ,0 ]
and
T(d,t) = 0
(3-20)
where h is the convection heat coefficient [W/(m 2 oC)]. In both cases, the boundary
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
44
condition T (d,t) = 0 follows from the negligible field and hence power dissipation at
z = d, which means a vanishing temperature differential at that location.
3.2.2 Numerical approach
The space-time domain can be discretized into a finite difference grid as shown in
Fig. 3-2. Suppose the solution domain is restricted to (M x N ) nodes. Then
Xnun= X [ ( m - \ ) h , ( n - l ) k ]
\< m < M ,l< n< N
(3-21)
where X is any quantity varying with z and t.
Wave Equation. The finite difference formula for the second order derivative of E
at node (m ,n ), n constant, takes the form
d~E
Em+\,n
"h
Em —l.n ~
1<m <M
— —= ---------------dz2
h2
(3-22)
With substitution of eq. (22) into eq. (11), the wave equation takes the finite ditference
form at time step n
(■
•i
■>■— ■ »—
(> ■ <i ■
— <> —
— <3
f
<
(m.n) (
k
1r
►
•
•
•
.
m= 1
x
h
m==M
Fig. 3-2 Finite difference grid of space time domain
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission .
45
E m -l,n
^ m ,n + E-m+l.r. ~ 0
( . 2 +
h
\ < m < M
(3-23)
where
T ^ n = “ 2 ^ e o [ /( 7 ’m,n)-;'g(7'/n,n)]
(3-24)
Assuming £ l n and E m „ are given as boundary conditions, wehave from (23) a system
of (M - 2) equations,which
can be solved easily. Usually, only theincidentfield £,• is
given and E 1>n = E (0) must be determined from eq. (13) with
t(°) = — n(°)~ ■
Ho+ 11(0 )
(3-25)
where TJo = yPoTSo" 2nd r](0) is unknown. It should be noted that, by definition
n(z) =
E(z)
H(z)
E(z)
1 E ' {z)
J'mo
For homogeneous media, E / E = y and hence rj =
(3-26)
qj/e. For inhomogeneous media, this
is no longer valid since the solution is still unknown. Instead, we assume
E(z) - C E(z)
(3-27)
where C is an arbitrary constant. We start to solve the boundary value problem
r 32£
3z 2
te
=0
£(0) = C E (0)
and
E(d) = 0
(3-28)
Once £(z) has been found, its derivative can be determined by numerical differentiation
[20]
E(Zs) = I ^ - E ( z u) + R [ E ]
u=0
h
where
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
(3-29)
46
, if S
Lu
V_1
= U
i v
v*u
(3-30)
-1
I
I
, if 5
1 1 (Z|» —Zv)
u
v=l
V= 1
v*u
v^s.u
For five-point interpolation (/ = 5) we have
E (0) = — [-25£(0) + 48£(A) - 36E{2h) + 16£(3/r) - 3£(4/t)] + 4-A4£ (5)(3-31)
12/i
It readily follows from eqs. (26) and (27) that
tKO)=7O)Po-I^0)
(3-32)
£ (0)
E (0) is then determined from eqs. (13) and (25) and hence E (z) can be found by
E ( z) = ! ^
C
= - P ^ £ ( z)
E( 0)
(3-33)
Thermal Equation. Equation (12) involves both space and time variables for which
a widely used differencing scheme of Ciank and Nicolson can be found in [21]
T m,n + 1 ~ T 'm ,n
Ctf r
k
l
7--------- = ~ [ T
"
>•
.
1
m,n + F m,„+l] + —[Gm,n + Gr^n+i]
2
^
^
0 < m < M (3-34)
where
r m + l,n + T m - l.n “ 2 F m_„
nun
_
Gm n
coeo
pcp
2
I^-m ,n I
Eqn. (34) can be rearranged into the form
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
(3-35)
(3-36)
47
- j - T m. hn+l + (1 + r)Tmn+i - y 7 m+1.n+1 = Dn
0 <>m<M
(3‘37)
where r = -Kr and
h2
= [~2^'m+l,n ~ (1 —r )Tm,n +
+ “ [^m,n + Gm,n + 1]
(3-38)
The boundary conditions given in eqs. (19) at the time step (n+1) take the finite differ­
ence forms
T'i.n+i -T i,n + \
and
TM'H+i = r^ /_ l n+1
(3-39)
and for eq. (2 0 )
(1 + “
) r 1.„+1 = r Zn+ 1 + - ^ r am
and
r M.fl+1 = r lW_1.n+1
(3-40)
From eqs. (37) and (39) or (40), we have again a system of equations with ( M - 2)
unknowns. These are the temperature at time step n+ l expressed explicitly by constants
and quantities known from time step n as given in Dn, with one exception of GOTin+1. The
coupled wave-thermal equation can be solved iteratively as follows:
1/ At time step n with known temperature distribution Tmn, solve eq. (23) for field
distribution E ^ ni 0 <m <M.
2/ Find G ^ „ 1 <m <M from eq. (36)
3/ Assuming Gm,n+1 = Gm,„, solve eq. (37) for T„^n+X
4/ Repeat 1/ and 2/ for time step n+1; solve for Tm^n+\ with the new value of
. This procedure is repeated until Tm<n+i converges. Usually, the conver­
gence is robust and one iteration is enough.
5/ Go to 1/ for next time step.
It is noted that both system equations (23) and (37) are of tri-diagonal form. Using Gauss
elimination, the solutions can be expressed in implicit formulas and then found by back
substitution [2 2 ].
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
48
3.3. T herm al runaw ay in one-dim ensional microwave heating
Thermal runaway can be defined as a process -a which the difference between the
n.tes of temperature rise at two particular points grows with time, or mathematically
d \ T ( z a) - f ( z b)\
(3-41)
0
dt
where the dot indicates a derivative with respect to time. This definition is consistent
with the common understanding of thermal runaway as a process in which the power dis­
sipation becomes more and more localized it some particular points. The relation is obvi­
ous from eq. (9) by neglect'ng the conduction term
T(z) = E llL =:G( z )
Pcp
(3-42)
which can reasonably be done as thermal conduction is usually small in microwave heat­
ing. G (z) will be referred to as heat source density. Substitute eq. (42) into (41), we get
d \ G ( z a) - G ( z b)\
dt
(3-43)
>0
Thermal runaway can thus be illustrated graphically as in Fig. 3-3 which sketches heat
source density distributions at two different times t\ and r?, where t ] < t 2
a)
b)
a
O
z
z
Fig. 3-3 Postulated qualitative distributions of heat source density for cases of
a) thermal runaway and b) non-thermal runaway
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
49
3.3.1 Examnle 1: Cooked beef juice
The dielectric properties of cooked beef juice at / =2.45 GHz [6 ] are approximately
given as functions of temperature T
er(T) = 7 0 - j T
e” (T) = 23
in the temperature interval [5°C , 65°C]. Assuming the thermal properties of cooked beef
juice are the same as those of water, we have cp = 4180 J/(kg°C), X,
W/(m °C),
p = 1000 kg/m 3 and h = 6 W/(m2 CC) for free convection with air which causes a negli­
gible surface loss as seen bv eq. (20). Assume a sinusoidal field of magnitude
t i = 20000 Vh n incident upon the surface of the half space dielectric medium as shown
in Fig. 3-1. The distance at which the field magnitude E(d) < 10- 3£ (0 ) is determined
from eq. (16) is
^ _
InS
____________________ Int.iO-3)___________________ = 0 1 m
c min " Re(-(2rt 2.45T0 9 )24 k - 10-78.85H 0 -12(70 - j 23)}
'm
Suppose the initial temperature distribution is T (?. )) = 10 °C f o r C < z < J . With the
numerical method developed in the previous section, a simulation o f 1-D microwave
heating process is run ovei a time period of rmax = 8 .s. Tne length and time ster are
h = 0.0005 m and k = 1 s, respectively. The temperature and power density distributions
at different elapsed times are shown in Fig. 3-4a and b respectively. As the loss factor is
constant, the field distribution is directly proportional to the square root of the power
density and is not shown.
The results obtained show a weak thermal tunaway taking place at the surface of
the dielectric medium: the heat source density at the surface increases from 6.6 to 7.8°C/s
while it remains unchanged for z >0.01 m. The explanation to this observation, despite
the constant loss factor er = 23, can be sought in the decrement of er at the surface where
the heat source density is highest and hence the temperature rises fastest. First of all, the
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
so
b)
a)
65-
6-
G
4°C/s!
■r
60-
2-
0.01
0.02
0.03
7.
0.04
0.05
[m]
0
0.01
0.02
0.03
z [m]
0.05
c)
: t = 0s
: t = Is
40-Ar
T
[°C]
: t = 8s
20-
-e-
-0
0.01
0.02
0.03
0.04
0.05
z [m]
Fig. 2-* Distribution o f a) heat source
density b) dielectric constant and
c) temperature in cooked beef juice
reflection coefficient T drops for decreasing £r , i.e. the transmitted field at the surface
increases. This increasing field, however, decays quickly away from the surface because
the medium, after the initial heating, becomes inhomogeneous as shown in Fig. 3-4b.
Clearly, the increasing attenuation must be due to the growing dielectric constant profile
in the direction of propagation, which in turn causes internal reflection within the
medium as has been explained in Section 3.1.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
51
3.3.2 Example 2: Nylon 66
Nylon 66 is a typical example of a material that manifests a strong tendency for
thermal runaway [23]. Its dielectric properties are shown in Fig. 3-5. An investigation of
thermal runaway for Nylon 66 can be carried out as before assuming £,- = 120000 V/m,
/ = 2.45 GHz, cp = 4180 J/(kg°C),
= 10 W/(m °C), p = 1000 kg/m 3 and negligible con­
vection loss at the surface of the dielectric half space. Because of the low loss factor,
d = 8 m is chosen to ensure a vanishing field at that distance.The simulation time rmax is
chosen to be 8 s and the steps h =0.005 m, k — Is. The numerical results are shown in
Fig. 3-6.
Thermal runaway is obvious by the strong shift of G (z) toward the surface. This is
mainly caused by the significant increase of e r at the surface, which is brought about as a
dsr
«
result of the positive slope, ------ , of the er vs. temperature response. Since G is propordT
tional to er , tius is the positive feedback that typically causes tnermal runaway in the
3
3
3
100
UO
1G
0
160
T C O
0 5
04
0 3
0
2
01
0
T CC)
Fig. 3-5 Dielectric properties of Nylon 66
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
52
a)
b)
80 - |
156010-
20-
■0
4
3
7.
•0
1
[m]
2
3
4
3
4
z [m]
c)
d)
.15
2 .9 9 -,
t r 2.98 - i i
2 .9 7 -
0.05
2.96
0
2
4
3
-0
i
7 [m]
2
z [m]
c)
— B
: / = Os
■x- : t = 4s
100-
-A- : t = 8s
80
T
[°C]
Fig. 3-6 Distribution of a) field magnitude
b) heat source density c) dielectric
constant d) loss factor and e)
temperature in Nylon 66
60
40
-
0
1
2
3
4
z[m ]
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
53
materials that possess this property. As a secondary effect, the rising er leads to increas­
ing power absorption at the surface and the field attenuation increases causing a drop of
G for z > 0.3 m. It is noted that the transmitted field at the surface appears unchanged
despite the rapid rise of Er ; this is because er is small in comparision to e r , which almost
does not change in the actual temperature interval [40°C, 80°C] ( see Fig. 3-5 ). There is,
in fact, a factor that counteracts the runaway effect in this case; it is the decaying profile
of both £r (z) and e r (z), which causes the field to grow in the z direction. However, the
effect is too weak to balance the strong field attenuation caused by power absorption at
the surface.
3.3.3 Example 3: Thawing o f frozen water
W hen defrosted, water and water-containing substances like most food materials
undergo a phase change near 0°C, at v/hich significant changes o f dielectric as well as
thermal properties are observed. As mentioned in Section 1.1, the relatively small er and
e r of ice make the penetration depth o f electric field much larger than that in water. In a
thawing process, when the ice at the surface of the material melts, the layer of water will
absorb most of the microwave energy, preventing further penetration into the material.
Of special importance in the thawing process is the significant latent heat of water at 0°C;
that is, the he^'. required for the ice-to-water phase change. In order to incorporate this
factor into the simulation, it seems most natural to reformulate the thermal equation (9)
in terms of enthalphy and temperature
1 dH
d2T
p
— TT + _
cp dt
1 dz2
pcp
(3-44)
where H is given by
\c pT
t
H = y CpT + La
H : enthalphy [J/kg],
<0
T> 0
Lc\ latent heat [J/kg]
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
(3' 45)
54
The thermal equation given by eqs. (44) and (45) will be solved by the finite differ­
ence method as before. For water, La = 3.33T0 5 [J/kg] and unlike the cases studied
before, Xt , p and cp can no longer be considered as constants. The temperature
dependences of these quantities [5] are here approximated by piecewise linear
0
W
t
%
functions as shown in Fig. 3-7 along with those for s r and £r . In this simulation we
choose Ei = 20000 V/m, / = 2.45 GHz, r max = 8 s, h = 0.002 m, k - I s .
The
numerical results are shown in Fig. 3-8.
Looking at the field distributions, we see a smooth curve that dies out at
around z = 0.1 m for t - 0. At this starting moment, the medium is homogeneous
because the temperature everywhere is —8 °C, so the field distribution is exponential
as expected. Because of the high power density, the temperature at the surface rises
fastest and so do er and e r . For t = 1 s, the surface temperature is about —1°C and
er, e'r have risen so much that the transmitted field has been reduced significantly.
On the other hand, e'r and e r' still remain low for z > 0.02m, which gives rise to a
local maximum of the field there marking the transition from a medium with high
permittivity to another with low permittivity. For r = 2s, the temperature at the sur­
face and up to a depth of 0.01m has reached 0°C indicating a phase change is tak­
ing place. The medium in this "mushy" region is a mix of ice and liquid water hav­
ing the same temperature because all the dissipated power in the region goes to
raising the enthalphy. The widening of the mushy region continues until t = 16s,
when the surface enthalphy becomes high enough to overcome the latent heat
needed to complete the phase change from ice to liquid water. After this moment,
the temperature rises quickly in the newly formed layer of water at the surface and
the mushy region moves deeper into the medium. At t = 24s, the temperature at the
surface is over 30°C while it is still 0°C just 0.005m below the surface. It should be
noted that the rate of temperature rise dT/dt is no longer approximately propor­
tional to G in the mushy region due the latent heat, so G cannot be used as an
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
55
b)
a)
60
• 10
d)
c)
C
5,
3000
4<
at
06
2300
0.4
•10
• 10
e)
I os
a?
012
0.1
•10
40
Fig. 3-7 Variation of a) relative permittivity
u) iuss factor c) tiierniai cunducdvity
u) specific heal and c) mass density
of water on temperature around 0°C
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission .
70
56
b)
a)
10
10-
G
[°C/s]
5
5-
0
0
0.02
0.04
0.06
0.08
1
■0
0.02
0.04
0.06
z [m]
z [m]
c)
d)
0.08
1
80
60
20-
40
10-
20
■0
■0
0.02
0.04
0.06
0.08
1
-0
0.02
0.04
z[m ]
0.06
0.08
z[m ]
— B : t = Os
c)
1- : r = Is
30- \
— A- : t = 2s
— x- : t = 16s
20-
— e- : t = 24s
Fig. 3-8 Distribution of a) field magnitude
b) heat source density c) dielectric
constant d) loss factor and e)
temperature in water
-10
-0
0.02
0.04
0.06
0.08
z [m]
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
57
indicator for the thermal runaway effect as before. However, one can see that the
difference between the rates of temperature rise at the surface and at some point in
the mushy region, say at z = 0.01m, grows when t passes 15s. According to the
definition given by eq. (41), this exhibits thermal runaway. The field distribution at
t = 24s decays quickly in the water and mushy regions and slowly dies out in the
ice region to the right. From the analysis carried out above, one can easily see why
it is so difficult to heat any water-containing substance from the frozen state; the
latent heat together with the high loss factor of water above the freezing point are
in effect a barrier that allows the temperature to rise quickly on one side where
T >0°C while suppressing the temperature rise on the other side where T < 0°C. It
is noted that the problem of one-dimensional microwave heating of frozen sub­
stances has also been studied by Coleman [31] using both analytical and numerical
approaches. The results he obtained are qualitatively similar to the ones given
above. Coleman concentrated in solving only the thermal equation while ignoring
the wave equation by assuming exponential field distributions in the water and ice
regions. This assumption is acceptable but strictly incorrect because of the field
reflection at the interface of the two regions. One interesting result he obtained is
that, for some materials with relatively high thermal diffusivity, the mushy region
between the thawed and frozen regions can be totally eliminated, especially when
the microwave heating is supplemented by convective heating at the surface of the
medium.
The convergence of the numerical solutions has been tested for all three cases
investigated above. For cooked beef juice and Nylon 66 , the numerical solutions
converge well even with the length steps and time steps twice as large as the values
chosen above but the results are not shown here. Instead, we look at the case of
ice-water for which the numerical convergence could be expected to be the worst
among the three cases because o f the rapid changes o f many parameters around the
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
58
(x) : h=0.0005m k=0.5s
(o): h=0 .001 m k= 1.0s
(+) : h=0 .002 m k=2 .0 s
(*): h=0.004m k=4.0s
T[°C\
h: length step k: time step
-10
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
0.09
0.1
z [m]
Fig. 3-9 Convergence o f numerical solution for temperature
distribution in water at r = 24s
freezing point. Fig. 3-9 shows, nevertheless, a robust convergence of the solution of
T for h = 0.001 m and k = 1 s. A further reduction in length and time step yields
only a negligible difference.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission .
59
4. TRANSVERSE MICROWAVE HEATING - THERMAL RUNAWAY
The one-dimensional microwave heating considered in Chapter 3 is not just o f
illustrative value but also it applies approximately to the situations where the dielectric
bodies are very large in wavelengths. In such situations, the incident microwave has no
other way to reach the interior of the body than penetrating through the surface and dissi­
pating most of its power there. The power penetration, thermal runaway and all the
characteristics of one-dimensional microwave heating depend totally on the dielectric
and thermal properties o f the body. In practice, the physical dimensions of most dielec­
tric bodies subject to microwave heating are of several wavelengths or less. The finite
size of the bodies necessitates analysis of the field inside the Fl/ uos in three dimensions.
Here, beside the material properties, the body geometry plays an important role because
o f resonances or modes that exist in the bodies. The three-dimensional analytical treat­
ment of the electromagnetic field in a dielectric body is quite complicated even for
bodies with regular shapes like circular or rectangular cylir ders.
In this chapter, we will again use the wire grid modeling method to compute the
power dissipation for several lossy dielectric slabs. Such bodies are chosen to illustrate
some frequently observed phenomena like edge heating and thermal runaway of a cheese
slice in a microwave oven. Fortunately, _hese phenomena can be explained by surface
wave theory without involving the complication of the three-dimensional analysis o f
dielectric resonators. In contrast to the case of one-dimensional heating, the thicknesses
o f the slabs to be considered a*e just a fraction of one wavelength so that the incident
fields can easily penetrate the slabs and non-uniform heating only happens in the
transverse planes. We therefore refer to this kind o f heating as transverse microwave
heating.
4.1. Surface wave theory
The established theory of surface waves on a dielectric slab will be summarised
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
60
here because of its direct applicability to the problem at hand. Consider the infinite flat
dielectric slab in free space, shown in Fig. 4-1. Wave propagation is taken to be in the x
—v r
direction and the x - variation of the fields is assumed to take the form X ( x ) x e
,
where yx is the complex propagation constant in the x direction. As in the analysis of a
parallel-plate waveguide, the wave fields are assumed to be independent of the y coordinate, so that all y - derivatives in M axwell’s equations are zero (
dy
= 0 ). With these
assumptions, Maxwell’s equations take the following form, in which e = £i inside the
slab and £ = £o outside:
dHx
dEx
- f c - + IxH z =
-= £ -
-jJxHy = ja x E z
dHy
+
=
ju \iH y
- j y xEy =-ja>\iH2
dEy
= j coe£z
(4-1)
= j($\xHx
The permeability p is taken to be the same everywhere, with p = Po- The field inside and
outside the slab must have the same propagation constant yx in order that the field con­
tinuity conditions can be satisfied at all points on the two sides o f the slab suiface.
z
Eo
t
I/
E’
y vc.*
Propagation
----------- ►
\
\
y
-t
Fig. 4-1 Edge view o f infinite slab supporting surface waves
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
61
Maxwell’s equations (1) are readily separable in two groups, the TE group consisting of
the second and the third equation in the right column together with the first equation in
the left column, and the TM group consisting of the remaining three equations.
4.1.1. TE waves
From M axwell’s equations (1) the relevant TE field components are Hx , Ey and
Hz. As in the analysis of rectangular waveguide, we use the separation-of-variables tech­
nique, assuming that Hx(z,x) = Z(z)X(x). The solution for Hx inside the slab takes the
form:
Hx l = C {Z (z)e~v
(4-2)
where C i is an arbitrary constant and Z (z) is proportional to sinh(n \ z ) or cosh(w \ z ) or a
linear combination of the two functions. The solution is most conveniently expressed in
terms of these funtions because u t , defined by the separation relation
u \ = -(co 2p.£1+y^)
(4-3)
in general, is a complex constant.
Even Modes. Choose Z (z) = sinh(u i z). Inside the slab, the complex solution for TE-even
modes is
Hx i(z,x) = C isinh(w iz)e ^
(4-4)
£ y i ( z ,x ) = ^ ^ C 1c o s h ( u 1z )e
(4-5)
J
ti
.
Hz \(z,x) = — C x C o s h ( u i z)e YxZ
“l
Outside the slab, for
z>t
the solution is of the form
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
(4-6)
where Co is a constant and uq is defined by
uq
= -(co : 'i£o+ ^)
(4 10)
Note that Im («o) > 0 “ a Re {«o ) > 0 to assure a vanishing field for z —
At the boundary surface z = r, the condition Hx y = Hxq gives
C isin h (u !r) = Co^ u '
(4-11)
The boundary condition Ey \ = Ey o at z = r gives
Ci
ui
cosh(uiz) =
Co - Ur,i
e
uo
(4-12)
Eliminating C \ I C q gives
witanh(uir)+uo = 0
(4-13)
With u\ and u q given in eqs. (3) and (10) respectively, the transcendental equation (13)
can be solved to give the propagation constant yx for a given mode. The Muller method
[25] is most convenient for seeking a complex root of the equation.
Odd Modes. Choose Z(z) = cosh(u i). The soludon for TE-odd modes inside the slab is
Hx \ (z,x) = Ci Cosh(ui z)e' yiX
Ey i (z,x)=
] sinh(«\z )e ^ llX
u!
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
(4-14)
(4-15)
Outside the slab, the solution is given by eqs. (7)-(9) as before. Matching the tangential
fields at z = t gives the transcendental equation for TE-odd modes
u 1coth(t' it>fUo = 0
(4-17)
It turns out that, for thin slabs with thickness much smaller than the intrinsic
wavelength of the dielectric medium, only the lowest even mode TEo can propagate.
4.1.2. 7 M M/.<7v*».r
From M axwell’s equations (1), the relevant TM field components are Ex , Hy and
Ez. Following the derivation procedure as in the TE case, we get the solutions with asso­
ciated transcendental equations for TM even and odd modes.
Even M -ie s. The complex solution forTM -even modes inside the slab is
Ex \(z,x) = C \S\n\i{u\z)e v
Hy i ( z , x ) =
J (0£l
_yr
C;Cosh(u!z)e f
(4-18)
(4-19)
“1
£, i( z ,x ) = —— CiCosh(uiz)e YrJ:
Ul
(4-20)
where yx is the root or the transcendental equation
Ui
£i
and
u
\ ,
uq
tanh(u 11)H
Uq
=0
Eo
are defined as above.
Odd Modes. The complex solution for TM odd modes inside the slab is
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
(4-21)
64
£ ^ 1(z,x) = CiCOsh.''w12)e Yp:
(4-22)
(4-23)
Ez i(z,x) = — ^ -C 1sinh(u 1z )f ^
Ml
(4-24)
where yx is the root of the transcendental equation
ill
(4-25)
El
and u i , Uq are defined as above.
4.2. Edge heating in infinitely long dielectric strips
For lossless dielectric slabs, yx is imaginary and the surface wave propagates
without attenuation. For lossy slabs, yx is generally a complex constant and the surface
wave dies out with the attenuation constant a x, where cc* = Re {yx }. Consider the
infinitely long dielectric strip shown in Fig. 4-2. The strip can be thought o f as a segment
of the infinite slab shown in Fig. 4-1. Assume the field in the strip can be expressed in
terms of the surface wave modes in an infinite slab. We are interested in investigating
the field variation inside the strip when a wave is incident upon it. It is o f interest to con­
sider the two cases, TE and TM:
4.2.1. The TE case
Let the strip be illuminated by some incident wave with the electric vector field
parallel to the edges of the strip. The tuiai field inside the strip consists o f the "forced
wavv" and surface waves. By "forced wave", we mean the wave that has a variation only
in the z direction b' i no variation along the x-direction. The existence o f this wave can
be verified by modifying e -YrX with e ^ + C, where C is a constant, in eqs. (4)-(6) and
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
65
y
-a
Fig. 4-2 Geometry o f infinite dielectric strip with TE surface waves
Fig. 4-3 Normalized power density distribution in infinite cheese strip
(er = 5 4 -/2 1 ) due to TEq surface waves.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
66
still satisfying the Maxwell’s equation given in eq. (1). The other waves are TE surface
waves which are induced by diffraction of the incident wave at the edges of the strip.
(Because of the infinite length of the strip, the total E-field has no x — component and
hence no TM wave is present in the total fieldV The TE surface waves propagate either in
the positive or in the negative x direction and the propagation constants for these waves
are determined by eq. (13). Usually, only the lowest even modes TEo can propagate on a
thin dielectric strip; other non-propagating modes may contribute to the field at the edges
but they most probably die out very quickly. Richmond [26] investigated the same prob­
lem using Galerkin’s method to determine the amplitudes of different modes existing in
the strip. He found that it is accurate to express the field in thin strips with TEo waves
and the forced wave only. If the forced wave is small compared to the TEo waves, the
field variation along the x direction can be approximately expressed by
X (x) = A e ^ + B
e
(4-26)
where A and B are constants. For symmetric excitation, A = B and we have
X ( x ) = 2A c o s h ^ x )
(4-27)
The normalized power density distribution is given by
IX CO | 2
2
|cosh(yxs ) |
t = '— ; |T
|X (n ) | 2
|cosh(Yra ) |
(4-28)
Obviously, if ct* is large, that is a xa :» 1, the field attenuates quickly from the edges so
the field magnitudes at the edges will be much larger than at the center of the strip, i.e.
w*. have edge heating. Otherwise, if a xa < 1, p(x) may have a large local maximum at
the center, i.e. we have center heating.
Example. Dielectric properties of cheese at / = 2.45GHz and T —50°C are e / = 54 and
?/■’ = 21. For a cheese slab of thickness 2r = 2mm, the propagation constant for the TEo
surface wave is found from eq. (13) to be yx = 44 + j 139 m-1. A plot of p (x) for a cheese
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
67
strip of width 10cm is ^ ° w n in Fig. 4-3. It is seen that the power density dies out quickly
from the edges; at distances 1cm from the edges, the power density is reduced to less
than half of the maximum value.
It is worth noting that the nature of decaying surface waves on a thin lossy strip
agrees well with the edge behaviour of surface current on a strip conductor. It can be
shown that for a perfectly conducting strip, the attenuation constant for TE surface waves
is infinitely large, which means the surface current dies out "instantly" in the vicinity of
the edges. Since the surface current at the center of the strip is nonzero, the surface
currents at the edges are forced to infinity. The singularity of edge current is well known
and was observed by several authors [27], [28] when they used different methods to cal­
culate the surface current distribution on a perfectly conducting strip being illuminated
by plane waves with their electric field vectors parallel to the edges of the strip.
\
4.2.2. The T M case
;
Consider the same dielectric strip as before; the electric field vector is now perpen­
dicular to the edges of the strip. It is clear that only TM surface waves are induced in this
!
case. For strips with high dielectric constant or loss factor, it is found that the propaga­
tion constant fc. TM waves is very close to that o f plane wave in free space. For the
;
cheese strip in the above example, the propagation constant for TMo wave is found from
eq. (17) to be yx = 0.0038 + j 51.42 m_1. For comparision, the propagation constant of a
plane wave in free space ;s y’Po = yco^poEo = j 51.35 m_1 at / = 2.45GHz. In contrast to
the TE case, the attenuation constant for the 1 ?vl wave is very small, which does not lead
to edge heating as in TE case. For normal illumination, the two TM waves propagating in
the positive and negative x directions have the same magnitude arid the resultant field
distribution has a form similar to the standing wave pattern on an electric dipole in free
r
space.
r
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
]
i
68
4.3. Edge heating in rectangular dielectric slabs
The infinite strip considered above was useful to illustrate the edge heating
phenomenon in lossy dielectric. In practice, any dielectric material subject to microwave
heating has finite dimensions and the heating pattern for such a three-dimensional object
becomes understandably more complicated. It is most natural in the study of edge heating
to consider a rectangular dielectric slab as shown in Fig. 4-4. The slab can be thought of
as being a segment of the strip considered above. However, no straightforward derivation
can be made for the slab from the analysis of the strip because the assumption
dy
= 0 is
no longer valid. An appropriate analytical treatment of the field inside the slab would
involve the theory o f dielectric resonators which is beyond the scope of this thesis. Even
if the propagation constants for different modes in the slab cavity could be determined
analytically, the simultaneous existence of many modes with unknown magnitudes
would not help one in gaining an insight into the true field existing in the slab. In this
section, we will instead use the wire grid modeling method to study the heating pattern
for some dielectric slabs of food materials with known dielectric properties. The surface
Z
Fig. 4-4 Geometry of rectangular dielectric slab
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
1
69
wave theory can still be used approximately to choose the approptiate dimensions for the
slabs and to discuss the results obtained.
To demonstrate edge heating in a rectangular slab, we choose to model three slabs
with the same dimensions but different food materials. The dielectric properties for these
materials are estimated from the data given by Ohlsson and Bengtsson [29]. The related
propagation constants and wavelengths for these materials are listed in Table 4-1, where
_ y = ct + v'p = (-cd2 [Io£o)1/2 is the propagation constant of a plane wave in an
infinite dielectric medium
_ Is = a j +
is the propagation constant of a TEo surface wave in an infinite
dielectric slab of thickness 2 1 = 2 mm
It should be noted that the listed values of e / and er" are just estimates at / = 2.45GHz
and T = 50°C for the chosen materials because of the limited quantity of available data.
This, however, does not affect the theoretical aspect of the study. The phase-shift con­
stants P are of the same order for the three kinds of food but there are significant differ­
ences between the attenuation constants for these niaterials, with ham being the lossiest
one. Similar observations can be made for P5 and a*. The choice of the materials has in
fact been based on these criteria to achieve slabs with the same physical dimensions as
well as the same wire grid model (but different loads). The results obtained for these
slabs would then be most easily compared and interpreted.
Since the surface
wavelengths Xs of all three materials are about 0.045 m, an appropriate length for the
slab would be 2a = 0.064 m which is about 3/2 Xs. This is certainly the minimum length
that would allow us to see the edge heating phenomenon in the slab. The width 2b of the
slab can be chosen more freely but it should be much larger than the thickness o f the slab
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
70
to keep the assumption of surface waves propagating on thin slab valid. We choose
b = a/2 = 0.032 m. The thickness of the slab has been chosen to be 2 mm, which is about
1/8X or the maximum allowpble wire grid spacing. The dimensions of the full wire grid
model for the slab are then (2ax2£x2c) = (32x16x1) d, where d = 2 mm is the grid spac­
ing. Using symmetry, the size of the quarter wire grid is (16x8x1) d. The source dipole is
y-oriented and placed lm above the slab. This large distance together with the utilized
ME-symmetry guarantee a virtually uniform plane wave incident normally upon the slab.
The whole configuration is shown in Fig. 4-5. The wire grid model for each material is
loaded with corresponding lumped impedances calculated from eqs. (2-4)- (2-6). Volume
resistance Rv and capacitance Cv are listed in Table 4-1. The computed power density
distributions for the three cases are shown in Figs. 4-6. The expected edge heating is
clearly shown in Fig. 4-6a, where the power density in the cheese slab has its maxima at
the two edges (x = ±a) and decays toward a minimum at the center (x = 0). The yvariation of the power density takes the form of a standing wave with its maximum at
y = 0. It is seen here that the concept of decaying surface waves is still applicable to
ex p iJ^
; edge heating phenomenon although the attenuation and phase-shift constants
in the x direction are different from those in case of the infinite strip. The evidence for
surface waves is more obvious in the less lossy pineapple slab for which the power den­
sity is shown in Fig. 4-6b. Here, we again see two maxima at the edges but what distin­
guishes this case from the previous one is the local maximum at the center of the slab. It
should be noted from the previous section that the field distribution on an infinite strip
always has a local maximum at the center if the strip is illuminated symmetrically. How­
ever, if the strip is very lossy or wide, the surface waves will effectively die out at the
center and the local maximum becomes negligible. The same can reasonably be stated
for a rectangular slab; we see that no local maximum is present at the center of the cheese
slab because the surface wave’s attenuation constant a* for cheese is much larger than
that for pineapple (see Table 4-1). The distribudon of power density in a rectangular slab
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
71
Table 4-1 Dielectric properties for some food materials at / = 2.45 GHz and related
propagation constants, wavelengths and volume lumped loads for wire g ud model
Food
£r
Y [m-1 ]
X [m]
Yj [m l ]
Xj [m]
Rv [« ]
Cv [pF]
Cheese
54-j21
72+J384
0.0164
44+jl39
0.0450
88
1.88
Pineapple
51-j9
32+j368
0.0171
19+jl45
0.0434
204
1.80
Ham
50-j38
130+j386
0.0163
81+jl31
0.0470
48
1.76
(0 .02, 0 . 02, 1.00 )
2 cm source dipole
Fig. 4-5 Quarter grid model for (64x32x2) mm dielectric slab with source dipole
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
72
a)
b)
Fig. 4-6 Power density distribution in a) cheese slab b) pineapple slab
c) ham slab o f (64x32x2)mm
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
73
’s admittedly much more complicates than in the case of a corresponding infinite strip
because it is a combination o f many modes that exist in the slab. With ham being the lossiest dielectric considered, the power density distribution in Fig. 4-6c shows an extreme
case of edge heating where the power density highly concentrates at two 15mm wide
strips at the edges.
As a final investigation on edge heating, consider a square cheese slab of
(16x16x2) d being illuminated by a uniform plane wave. The whole configuration
shown in Fig. 4-7 is identical to the one shown in Fig. 4-5 except the length of the slab
has been reduced to 2a = 16 d = 32 mm which is only 3/4 Xs . The power density distribu­
tion is shown in Fig. 4-8. In contrast to the cases studied above, the power distribution
has its maximum at the center of the slab, which can be reasonably explained in the same
way as before. Due to the relatively small length of he slab, the surface waves travelling
from the edges have not decayed significantly before they interfere and build up the
standing wave peak at the center. One special characteristic of power density in the slab
is the peaks that are highly localized at the comers. Simulations of microwave heating in
other slabs showed that the peaks only show up for small slabs but not for those whose
dimensions are in order of wavelengths (the results are not shown). Similar observations
have been made by Ko and Mittra [28] for surface currents *n perfectly conducting
plates. There seems to be no exact explanation for this phenomenon at the present time.
4.4. Thermal runaway in transverse microwave heating
Just as in the case of one-dimensional heating, initial non-uniform heating can lead
to thermal runaway in the material being illuminated transversely. In microwave
transverse heating, the source of initial non-uniform heating can be either uneven illumi­
nation or an edge heating effect. Whatever the reason, the temperature dependence o f the
dielectric properties determines the tendency of a material to develop thermal runaway.
Thermal runaway is typically observed by heating a rectangular cheese slice in a
microwave oven. Usually, there are two parallel edges being overheated. As time passes,
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
74
(0.02 , 0 .02, 1.00)
2 cm source dipole
d = 2 mm
y
Fig. 4-7 Quarter grid modei for (32x32x2) mm cheese slab with source dipole
f
i
a
Fig. 4-8 Power density distribution in (32x32x2) m m cheese slab
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
75
the borders between these edges and the rest become more and more well defined; the
edges finally melt and degrade while the rest is bareiy warmed up. It is o f interest to con­
sider the cheese slab with the wire grid model shown in Fig. 4-5. The temperature depen­
dence of the dielectric properties of cheese can be described by the linear approximations
[291:
e / = 54 + ^ ( 7 - 5 0 )
er" = 21 + - j y ( r - 5 0 )
in the temperature range T = 50 —» 125°C.
After an initial uneven heating, the slab becomes dielectrically inhomogeneous due
to the non-uniform temperature distribution in the slab. The wire grid model for the slab
must therefore be loaded non-uniformly. The new lumped impedances for the wire grid
model can be calculated as described previously in Chapter 2 with the assumption that
the temperature and hence the dielectric properties are constant over any given elemen­
tary dielectric cube of the model. With the new loads, a new power distribution must be
computed. With the assumption that the heat exchange in the slab is negligible due to the
quick heating process, the dynamic heating pattern in the slab can be simulated as fol­
lows:
1/ Set up the wire grid model for the homogeneous slab.
2/ Compute the power density distribution
3/ Compute the temperature rise in the slab AT =
pC
Usually, pC is constant
with temperature and the time step can be scaled to give a maximum temperature
rise ( A
=
— • The temperature rise is then generally
Pc
given
A T = - 2 - ( A T ) max.
P m ax
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
by
75
the borders between these edges and the rest become more arid more well defined; the
edges finally melt and degrade while the rest is barely wanned up. It is of interest to con­
sider the cheese slab with the wire grid model shown in Fig. 4-5. The temperature depen­
dence of the dielectric propenies of cheese can be described by the linear approximations
[291:
e / = 54 4 - ^ ( 1 - 5 0 )
e / ' = 21 + - ^ - ( r - 5 0 )
in the temperature range T = 50 —» 125°C.
After an initial uneven heating, the slab becomes dielectrically inhomogeneous due
to the non-uniform temperature distribution in the slab. The wire grid model for the slab
must therefore be loaded non-uniformly. The new lumped impedances for the wire grid
model can be calculated as described previously in Chapter 2 with the assumption that
the temperature ana hence the dielectric properties are constat
over any given elemen­
tary dielectric cube of the model. With the new loads, a new power distribution must be
computed. With the assumption that the heat exchange in the slab is negligible due ro the
quick heating process, the dynamic heating pattern in the slab can be simulated as fol­
lows:
1/ Set up the wire grid model for the homogeneous slab.
2/ Compute the power density distribution
p A'
3/ Compute the temperature rise in the dab .jT --—— ' Is iallv, ->C is constant
DC
with temperature and the time step can b" 'T e d to g ve a m'ximum t mperature
CmaxAr
rise (AT)mM = ----- -— . The temper:-1 '
pC
rise is d u n
. lerall • given by
A7 = - £ - ( A r w
P nax
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission.
76
4 / Calculate the new distribution of dielectric properties and hence the new lumped
loads for the model.
5/ Go to 2/ for the next time step.
If the time step is small such that the temperature rise and hence the dielectric
change are small then the simuladon should approach the rea* evolution of power density
in the slab. For the cheese slab, a maximum temperature rise is chosen to be 25°C for
each time step. Assuming yC = 1.0-10'J/(m3oC) the results are shown in Fig. 4-9.
Thermal runaway is evident by the rising power density at the edges and the
remaining low power density it the center. Thermal runaway is expected in :his case
because of the rising loss factor er" with temperature, which means the hot regions in the
slab get lossier and become hotter and so on. Furthermore, the power carried by the <urface waves travelling from the edges to the center gets absorbed more and more n r r the
edges leaving virtually nothing for the center. However, thermal runaway is not very
strong in this case probably because the dielectric constant e / is also rising with tempera­
ture, which in effect reduces the penetrability of waves into the hottest regions. Admit­
tedly, the investigation conducted above is very' sirr~'e and is only one of many cases
that show thermal runaway; a rigorous treatment of thermal runaway would require much
more research in this virtually unexplored area.
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
77
a)
a
b)
c)
Fig. 4-9 Evolution of power density distribution in cheese slab:
a) r = Os, b) r = 76s, c) t — 146s
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
78
CONCLUSIONS
The wire grid modeling method has been shown to be successful in modeling con­
ductors and dielectric materials. Simulation of microwave heating of a dielectric in a
small, deep rectangular waveguide model shows that the numerical results obtained by
the method are comparable to know: analytical solutions. The method was then used to
investigate the microwave heating properties of various dielectric models, including edge
heating and thermal runaway phenomena in common food materials. The predicted heat­
ing patterns in these models are consistent with the practical heating patterns of food
conanonly observed in commercial microwave ovens. General aspects of one dimen­
sional microwave heating have also been studied to give some insights into the physics of
the thermal runaway effect. In particular, the finite difference method was used to solve
coupled wave-thermal equations witti some simulations of one-dimensional microwave
heating revealing the tendency of several dielectric materials to develop thermal runawav
With the ability of modeling one- two- and three-uimensionai objects for a wide
variety of materials, the wire grid modeling method is shown to be a powerful tool, not
only for the research in microwave heating but also for solving general electromagnetic
scattering problems. Despite the limitation that only small objects of several wavelengths
can be modeled with the wire grid modeling method, it is hoped that the rapid evolution
o f the modem computer technology will enaLe larger models to be run on machines with
higher speed and larger memory and this will certainly open up a new spectrum of prob­
lems that can be solved by the method.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
79
APPENDIX A
IM PLEM ENTATION O F ELEC TR IC AND M AGNETIC SYM METRIES
IN TH E TH IN -W IRE PROGRAM
In electromagnetic analysis, the theoretical perfect electric conductor is frequently
used to get the simple boundary condition nxE = 0 or E tan = 0. As a dual to the perfect
electric conductor, the hypothetical perfect magnetic conductor is the one that satisfies
the boundary condition -n x H = 0 or H tan = 0. The images of electric current elements
normal and parallel to these perfectly conducting planes are shown in Fig. A -l
The concept of electric and magnetic symmetries is extremely useful in the
perfect electric conducting plane
(E-planc)
positive image
negative image
pcncct magnetic conducting plane
(M-planc)
negative image
positive image
Fig. A -l
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited without p erm ission .
80
reduction of unknown dipole current modes when solving a problem using the thin wire
program. If the wire structure is symmetric with respect to the X=0 and Y -0 planes for
instance, electric and magnetic symmetries can be imposed on these planes to reduce the
unknown currents to the ones flowing in the first quadrant only. The currents in other
quadrants are images of the first quadrant currents. The physical interpretation o f impos­
ing a symmetry on a plane is that the wire structure is symmetrically excited or
illuminated by electromagnetic sources around that plane. Consider the simple wire
structure consisting of two dipoles shown in Fig. A-2a. System equations for the wire
dipole 1 i~
dipole 2 i—.
i i
V'i
V2
a i
i
Fig. A-2a
i
l
— >dipole
E-planc
-----------------------------------------image
Fig. A-2b
-T,
-c-------
A
dipole
V'.
M-planc — — — — — — — —
image
Fig. A-2c
c ::::::
-----I
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
81
structure are
Z i i / i + Z i 2/ 2 = Vi
(A-l)
^ 21^2 "h ^ 22^2 = ^ 2
(A-2)
We always have Z 12 = Z 2 i. And Z n = Z 22 if the dipoles are identical. V l and V 2 are
excitation voltages.
Electric symmetry consideration.
If V x = - V
2
then it readily follows 1 1 = - / 2 and the configuration corresponds to
the one shown in Fig. A-2b with
If the dipole and its image coincide, that is if the dipole is on the symmetry plane then
Z\ \ = Z 12 and I \ is undetermined. In general, therefore, when solving for a thin wire
problem using electric symmetries, the segments lying on these planes must be removed
to avoid the above problem. This is possible because the total currents on these segments
are zeros and removing them does not effect the solution for the other currents on the
wire structrure. Another advantage is that the number of unknown currents will be then
reduced.
Magnetic symmetry consideration.
U V X = V 2 then we have I x =1 2, which corresponds to the configuration with mag­
netic symmetry shown in Fig. A-2c, and
Vi
'1 = 7 — 7 Z 11 + -M2
Now, if the dipole is on the symmetry plane then Z \ \ = Z i 2 and
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
(A-4)
82
(A-5)
2Z 11
2
where I \ is the current flowing in the dipole in free space.
If the dipole is loaded with Z i ,
/
^1—
2Zn + Z l
1
(A- 6 )
and the dipole in free space is loaded with Z' \ ,
Vi
/' i =
(A-7)
Z n +Z\
Since the total current on the segment lying on the symmetry plane is 21 [ which should
be equal to l \ for all values of V i , it follows that Z t = 2Z 'i. When using symmetries,
this means the segments lying on magnetic symmetry planes should be loaded witn twice
the impedance these segments are loaded with when no symmetry is used.
For dipoles perpendicular to symmetry planes, it is worth noting that currents are
able to flow normally through an electric plane but are shorted in the case of a magnetic
plane. In practice, this means, if a dipole is placed perpendicular next to ar. electric plane,
an extra current mode must be set up at the plane. In contrast, no extra current mode is
needed if the plane is magnetic. See Fig. A-3.
M-plane
E-plane
r>
11'
II '
u
n
extra current mode1 n
II
“U r /
" ,
ly
Fig. A-3
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
83
APPENDIX B
COM BINATION OF SYM M ETRIES
The use of electric and magnetic symmetries for a symmetrical wire structure
implies that it must be symmetrically excited. However, nonsymmetrical excitation can
be achieved by combining different symmetry configurations together since electromag­
netic quantities are linear and additive. Consider the simple wire structure shown in Fig.
B-la. Dipole 1 is fed with current I and dipole 2 with a zero current. The configuration is
equivalent to the linear combination of the electric and magnetic symmetries as shown in
Fig. B -lb.
Combination of two-dimensional symmetries can also be carried out and is shown
in the following. The four basic two-dimensional symmetry configurations with x- and
y-oriented current elements are shown in Fig. B-2. Let x ,, x 2, x 2 and X 4 be the
coefficients of the basic configurations MEx, EMx, MMx and EEx respectively in the
combination to create a desired excitation configuration with x-oriented currents having
the values lclt k 2, £3 and k.4 in the 4 quadrants. Referring to Fig. B-2, we can write the
combination equations in matrix notation as follows,
/ , =1
Fig. B -la
/ 2 =0
1/2
1/2
E-plane
+
a
M-plane
Fig. B-lb
c
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
84
M
M
M Ex
MEy
EM x
EMy
M
M
M
MMx
MMy
— M
M
EEx
EEy
Fig. B-2
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
85
1
1
-1
il
1
1
1 1 -1 -11
-1
-
1
1 - 11
'
'k{
*1
ki
*2
“*
(B -l)
*3
*3
X4
k4
Solving for x x to x 4, we get
’ 1/4 -1/4 1/4 —1/4
*2
1/4 1/4 1/4 1/4
“
1/4 -1/4 -1 /4 1/4
*3
1/4 1/4 -1 /4 -1 /4
x4
'
k
f
ki
ki
(B-2)
k4
Similarly, let y i , >’2 , >3 and y 4 be the coefficients for MEy, EMy, MMy and EEy
configurations respectively. And let l \ , I 2 , Ij, and l4 be the values of y-oriented currents
in the 4 quadrants. Referring to Fig. B-2, we have the matrix equation
1
1
1 -1 -1
il
'yi
1
yi
12
1
y
h
1
ya
3
• / r
(B-3)
U
which gives
yi
yi
y3
y4
' 1/4 1/4 1/4 1/4
1/4 -1/4 1/4 -1/4
1/4 1/4 -1/4 -1/4
1/4-1/4-1/4 1/4
-li
h
h
U
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
('3-4,
86
REFERENCES
[1]
W.L. Stutzman and G.A. Thiele, "Antenna theory and design", John Wiley &
Sons, New York, 1981, pp. 356-364.
[2]
J.H. Richmond, "Radiation and scattering by thin-wire structures in the complex
frequency domain", National Technical Information Service, Springfield, VA
22151, NASA CR-2396, May 1974.
[3]
C. Aversa, "Microwave properties of wire grids", B.A.Sc., thesis, Div. of
Engineering Science, University of Toronto, Toronto, Canada, 1988.
[4]
A.C. Metaxas and R.J. Meredith, "Industrial microwave heating", Peter Peregrir.us, London, 1983, pp. 5-24.
[5]
N.E. Bengston and P.D. Risman, "Dielectric properties of foods at 2.8 GHz as
determined by a cavity pertubation technique", J. Microwave Power, vol. 6 , no. 2,
p. 107, 1971.
[6 ]
E.C. To, R.E. Mudgett, D.I.C. Wang, S.A. Goldblith and R.V. Decareau, "Dielec­
tric properties of cood materials", J. Microwave Power, vol. 9, no. 4, p. 303, 1974.
[7]
R.F. Harrington, "Field computation by moment methcds", Macmillan, NewYork,
1968.
[8]
R. Mittra, ed., "Computer techniques for electromagnetics", Pergamon Press,
Oxford, 1973.
[9]
K.S.H. Lee, L. Marin and J.P. Castillo, "Limitation of wire grid modeling o f a
closed surface", IEEE Trans. Electromag. Compat., vol. EMC-18, no. 3, pp. 123129, Aug. 1976.
[10]
J. Moore and R. Pizer, "Moment methods in electromagnetics", Research Studies
Tress, New York, 1984, secs. 1.1.2 and 6 4.
[11]
A.C. Ludwig, "Wire grid modeling of surfaces", IEEE Trans. Antennas Propagat.,
vol. AP-35, no. 9, pp. 1045-i048, Sept. 1987.
[12]
E.C. Jordan and K.G. Balmain, "Electromagnetic waves and radiating systems",
* entice-Hall, Toronto, 1968, pp. 470-473.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
87
[13]
M.A. Tilston, "Thin-wire reciprocal multiradius implementation o f the elec­
tromagnetic moment method", Ph.D, thesis, University o f Toronto, Toronto,
Canada, 1989.
[14]
P.C. Clemmow, "The plane wave spectrum representation o f the electromagnetic
field", Pergamon Press, Oxford, 1966.
[15]
D.T. Paris and F.K. Hurd, "Basic electromagnetic theory", McGraw-Hill, 1969,
New York, pp. 367-370.
[16]
J.R. Wait, "Electromagnetic waves in stratified media", Pergamon Press, New
York, 1962, pp. 88-91.
[17]
in [16], pp. 100-105.
[18]
D.M. Young and R.T. Gregory, "A survey of numerical mathematics", AddisonWesley Publishing Company, New York, 1972, pp. 473-477.
[19]
"Numerical recipes: the art of scientific computing", Cambridge University Press,
1986.
[20]
in [18], pp. 357-358.
[2 1 ] in [18], pp. 1086-1088.
[22]
in [18], pp. 586-587.
[23]
in[4], p. 55.
[24]
E.N. Ganic, J.P. Hartnett and W.M. Rohsenow, "Handbook of heat transfer funda­
mentals", McGraw-Hill, New York, 1985.
[25]
D.E. Muller, "A method for solving algebraic equations using an automatic com­
puter", Mathematical Tables and Aids to Computation, vol. 10, 1956, pp. 208-215.
[26]
J.H. Richmond, "Scattering by thin dielectric strips", IEEE Trans. Antennas Propagat., vol. AP-'J>, no. 1, op. 64-68, Jan. 1985.
[27]
S.M. Rao, D.R. V/ilton and A.W. Glisson, "Electromagnetic scattering by surfaces
o f arbitraiy shape , IEEE Trans. Antennas Propagat., vol. AP-30, no. 3, pp 416-
R ep ro d u ced with p erm ission o f th e copyright ow ner. Further reproduction prohibited w ithout perm ission.
88
417, May 1982.
[28]
W.L. Ko and R. Mittra, "A new approach based on a combination of integral
equation and asymptotic techniques for solving electromagnetic scattering prob­
lems", EEEE Trans. Antennas PropagaL, vol. AP-25, no. 2, pp. 187-197, March
1977.
[29]
Th. Ohlsson and N.E. Bengtsson, "Dielectric food data for microwave sterilization
processing", J. Microwave Power, vol. 10, no. 1, 1975.
[30]
J.T. Mayhan, "Characteristic modes and wire grid modeling", IEEE Trans. Anten­
nas PropagaL, vol. 38, no. 4, p. 459, April 199G.
[31]
C.J. Coleman, "The microwave heating of frozen substances", Appl. Math. Model­
ling, vol. 14, pp. 439-443, Aug. 1990.
R ep ro d u ced with p erm ission o f the copyright ow ner. Further reproduction prohibited w ithout p erm ission.
Документ
Категория
Без категории
Просмотров
0
Размер файла
3 314 Кб
Теги
sdewsdweddes
1/--страниц
Пожаловаться на содержимое документа