《电磁学》电子书（英文版）：Appendix-A

Appendix A Mathematical appendix A.1 The Fourier transform The Fourier transform permits us to decompose a complicated ?eld structure into elemental components. This can simplify the computation of ?elds and provide physical insight into their spatiotemporal behavior. In this section we review the properties of the transform and demonstrate its usefulness in solving ?eld equations. One-dimensional case Let f be a function of a single variable x. The Fourier transform of f(x)is the function F(k) de?ned by the integral F { f(x)}=F(k)= integraldisplay ∞ ?∞ f(x)e ?jkx dx. (A.1) Note that x and the corresponding transform variable k must have reciprocal units: if x is time in seconds, then k is a temporal frequency in radians per second; if x is a length in meters, then k is aspatialfrequency in radians per meter. We sometimes refer to F(k) as the frequencyspectrum of f(x). Not every function has a Fourier transform. The existence of (A.1) can be guaranteed by a set of su?cient conditions such as the following: 1. f is absolutely integrable: integraltext ∞ ?∞ | f(x)|dx <∞; 2. f has no in?nite discontinuities; 3. f has at most ?nitely many discontinuities and ?nitely many extrema in any ?nite interval (a,b). While such rigor is certainly of mathematical value, it may be of less ultimate use to the engineer than the following heuristic observation o?ered by Bracewell [22]: a good mathematicalmodelofaphysicalprocessshouldbeFouriertransformable. That is, if the Fourier transform of a mathematical model does not exist, the model cannot precisely describe a physical process. The usefulness of the transform hinges on our ability to recover f through the inverse transform: F ?1 {F(k)}= f(x)= 1 2π integraldisplay ∞ ?∞ F(k)e jkx dk. (A.2) When this is possible we write f(x) ? F(k) and say that f(x) and F(k) form a Fourier transform pair. TheFourierintegraltheorem states that FF ?1 { f(x)}= F ?1 F { f(x)}= f(x), except at points of discontinuity of f . At a jump discontinuity the inversion formula returns the average value of the one-sided limits f(x + ) and f(x ? ) of f(x). At points of continuity the forward and inverse transforms are unique. Transform theorems and properties. We now review some basic facts pertaining to the Fourier transform. Let f(x) ? F(k)= R(k)+ jX(k), and g(x) ? G(k). 1. Linearity. αf(x)+βg(x) ? αF(k)+βG(k) if α and β are arbitrary constants. This follows directly from the linearity of the transform integral, and makes the transform useful for solving linear di?erential equations (e.g., Maxwell’s equations). 2. Symmetry. The property F(x) ? 2πf(?k) is helpful when interpreting transform tables in which transforms are listed only in the forward direction. 3. Conjugatefunction. We have f ? (x) ? F ? (?k). 4. Realfunction. If f is real, then F(?k)= F ? (k). Also, R(k)= integraldisplay ∞ ?∞ f(x)cos kx dx, X(k)=? integraldisplay ∞ ?∞ f(x)sin kx dx, and f(x)= 1 π Re integraldisplay ∞ 0 F(k)e jkx dk. A real function is completely determined by its positive frequency spectrum. It is obviously advantageous to know this when planning to collect spectral data. 5. Realfunctionwithre?ectionsymmetry. If f is real and even, then X(k)≡ 0 and R(k)= 2 integraldisplay ∞ 0 f(x)cos kx dx, f(x)= 1 π integraldisplay ∞ 0 R(k)cos kx dk. If f is real and odd, then R(k)≡ 0 and X(k)=?2 integraldisplay ∞ 0 f(x)sin kx dx, f(x)=? 1 π integraldisplay ∞ 0 X(k)sin kx dk. (Recallthat f isevenif f(?x)= f(x)forall x. Similarly f isoddif f(?x)=?f(x) for all x.) 6. Causalfunction. Recall that f is causal if f(x)= 0 for x < 0. (a) If f is real and causal, then X(k)=? 2 π integraldisplay ∞ 0 integraldisplay ∞ 0 R(k prime )cos k prime x sin kx dk prime dx, R(k)=? 2 π integraldisplay ∞ 0 integraldisplay ∞ 0 X(k prime )sin k prime x cos kx dk prime dx. (b) If f is real and causal, and f(0) is ?nite, then R(k) and X(k) are related by the Hilberttransforms X(k)=? 1 π P.V. integraldisplay ∞ ?∞ R(k) k ? k prime dk prime , R(k)= 1 π P.V. integraldisplay ∞ ?∞ X(k) k ? k prime dk prime . (c) If f is causal and has ?nite energy, it is not possible to have F(k) = 0 for k 1 < k < k 2 . That is, the transform of a causal function cannot vanish over an interval. A causal function is completely determined by the real or imaginary part of its spectrum. As with item 4, this is helpful when performing calculations or mea- surements in the frequency domain. If the function is not band-limited however, truncation of integrals will give erroneous results. 7. Time-limitedvs.band-limitedfunctions. Assume t 2 > t 1 .Iff(t)= 0 for both t < t 1 and t > t 2 , then it is not possible to have F(k) = 0 for both k < k 1 and k > k 2 where k 2 > k 1 . That is, a time-limited signal cannot be band-limited. Similarly, a band-limited signal cannot be time-limited. 8. Null function. If the forward or inverse transform of a function is identically zero, then the function is identically zero. This important consequence of the Fourier integral theorem is useful when solving homogeneous partial di?erential equations in the frequency domain. 9. Spaceortimeshift. For any ?xed x 0 , f(x ? x 0 ) ? F(k)e ?jkx 0 . (A.3) A temporal or spatial shift a?ects only the phase of the transform, not the magni- tude. 10. Frequencyshift. For any ?xed k 0 , f(x)e jk 0 x ? F(k ? k 0 ). Note that if f ? F where f is real, then frequency-shifting F causes f to be- come complex — again, this is important if F has been obtained experimentally or through computation in the frequency domain. 11. Similarity. We have f(αx) ? 1 |α| F parenleftbigg k α parenrightbigg , where α is any real constant. “Reciprocal spreading” is exhibited by the Fourier transform pair; dilation in space or time results in compression in frequency, and vice versa. 12. Convolution. We have integraldisplay ∞ ?∞ f 1 (x prime )f 2 (x ? x prime )dx prime ? F 1 (k)F 2 (k) and f 1 (x)f 2 (x) ? 1 2π integraldisplay ∞ ?∞ F 1 (k prime )F 2 (k ? k prime )dk prime . The ?rst of these is particularly useful when a problem has been solved in the frequency domain and the solution is found to be a product of two or more functions of k. 13. Parseval’sidentity. We have integraldisplay ∞ ?∞ | f(x)| 2 dx = 1 2π integraldisplay ∞ ?∞ |F(k)| 2 dk. Computations of energy in the time and frequency domains always give the same result. 14. Di?erentiation. We have d n f(x) dx n ? (jk) n F(k) and (?jx) n f(x) ? d n F(k) dk n . The Fourier transform can convert a di?erential equation in the x domain into an algebraic equation in the k domain, and vice versa. 15. Integration. We have integraldisplay x ?∞ f(u)du ? πF(k)δ(k)+ F(k) jk where δ(k) is the Dirac delta or unit impulse. Generalized Fourier transforms and distributions. It is worth noting that many useful functions are not Fourier transformable in the sense given above. An example is the signum function sgn(x)= braceleftBigg ?1, x < 0, 1, x > 0. Although this function lacks a Fourier transform in the usual sense, for practical purposes it may still be safely associated with what is known as ageneralizedFouriertransform.A treatment of this notion would be out of place here; however, the reader should certainly be prepared to encounter an entry such as sgn(x) ? 2/jk in a standard Fourier transform table. Other functions can be regarded as possessing transforms when generalized functions are permitted into the discussion. An important example of a generalized function is the Dirac delta δ(x), which has enormous value in describing distributions that are very thin, such as the charge layers often found on conductor surfaces. We shall not delve into the intricacies of distribution theory. However, we can hardly avoid dealing with generalized functions; to see this we need look no further than the simple function cos k 0 x with its transform pair cos k 0 x ? π[δ(k + k 0 )+δ(k ? k 0 )]. The reader of this book must therefore know the standard facts about δ(x): that it acquires meaning only as part of an integrand, and that it satis?es the siftingproperty integraldisplay ∞ ?∞ δ(x ? x 0 )f(x)dx = f(x 0 ) for any continuous function f . With f(x)= 1 we obtain the familiar relation integraldisplay ∞ ?∞ δ(x)dx = 1. With f(x)= e ?jkx we obtain integraldisplay ∞ ?∞ δ(x)e ?jkx dx = 1, thus δ(x) ? 1. It follows that 1 2π integraldisplay ∞ ?∞ e jkx dk =δ(x). (A.4) Useful transform pairs. Some of the more common Fourier transforms that arise in the study of electromagnetics are given in Appendix C. These often involve the simple functions de?ned here: 1. Unitstepfunction U(x)= braceleftBigg 1, x < 0, 0, x > 0. (A.5) 2. Signumfunction sgn(x)= braceleftBigg ?1, x < 0, 1, x > 0. (A.6) 3. Rectangularpulsefunction rect(x)= braceleftBigg 1, |x|< 1, 0, |x|> 1. (A.7) 4. Triangularpulsefunction Lambda1(x)= braceleftBigg 1 ?|x|, |x|< 1, 0, |x|> 1. (A.8) 5. Sincfunction sinc(x)= sin x x . (A.9) Transforms of multi-variable functions Fourier transformations can be performed over multiple variables by successive appli- cations of (A.1). For example, the two-dimensional Fourier transform over x 1 and x 2 of the function f(x 1 ,x 2 ,x 3 ,...,x N ) is the quantity F(k x 1 ,k x 2 ,x 3 ,...,x N ) given by integraldisplay ∞ ?∞ bracketleftbiggintegraldisplay ∞ ?∞ f(x 1 ,x 2 ,x 3 ,...,x N )e ?jk x 1 x 1 dx 1 bracketrightbigg e ?jk x 2 x 2 dx 2 = integraldisplay ∞ ?∞ integraldisplay ∞ ?∞ f(x 1 ,x 2 ,x 3 ,...,x N )e ?jk x 1 x 1 e ?jk x 2 x 2 dx 1 dx 2 . The two-dimensional inverse transform is computed by multiple application of (A.2), recovering f(x 1 ,x 2 ,x 3 ,...,x N ) through the operation 1 (2π) 2 integraldisplay ∞ ?∞ integraldisplay ∞ ?∞ F(k x 1 ,k x 2 ,x 3 ,...,x N )e jk x 1 x 1 e jk x 2 x 2 dk x 1 dk x 2 . Higher-dimensional transforms and inversions are done analogously. Transforms of separable functions. If we are able to write f(x 1 ,x 2 ,x 3 ,...,x N )= f 1 (x 1 ,x 3 ,...,x N )f 2 (x 2 ,x 3 ,...,x N ), then successive transforms on the variables x 1 and x 2 result in f(x 1 ,x 2 ,x 3 ,...,x N ) ? F 1 (k x 1 ,x 3 ,...,x N )F 2 (k x 2 ,x 3 ,...,x N ). In this case a multi-variable transform can be obtained with the help of a table of one- dimensional transforms. If, for instance, f(x,y,z)=δ(x ? x prime )δ(y ? y prime )δ(z ? z prime ), then we obtain F(k x ,k y ,k z )= e ?jk x x prime e ?jk y y prime e ?jk z z prime by three applications of (A.1). A more compact notation for multi-dimensional functions and transforms makes use of the vector notation k = ?xk x + ?yk y + ?zk z and r = ?xx + ?yy + ?zz where r is the position vector. In the example above, for instance, we could have written δ(x ? x prime )δ(y ? y prime )δ(z ? z prime )=δ(r ? r prime ), and F(k)= integraldisplay ∞ ?∞ integraldisplay ∞ ?∞ integraldisplay ∞ ?∞ δ(r ? r prime )e ?jk·r dx dydz = e ?jk·r prime . Fourier–Bessel transform. If x 1 and x 2 have the same dimensions, it may be con- venient to recast the two-dimensional Fourier transform in polar coordinates. Let x 1 = ρcosφ, k x 1 = p cosθ, x 2 =ρsinφ, and k x 2 = p sinθ,where p andρ are de?ned on(0,∞) and φ and θ are de?ned on (?π,π). Then F(p,θ,x 3 ,...,x N )= integraldisplay π ?π integraldisplay ∞ 0 f(ρ,φ,x 3 ,...,x N )e ?jpρcos(φ?θ) ρdρdφ. (A.10) If f is independent of φ (due to rotational symmetry about an axis transverse to x 1 and x 2 ), then the φ integral can be computed using the identity J 0 (x)= 1 2π integraldisplay π ?π e ?jx cos(φ?θ) dφ. Thus (A.10) becomes F(p,x 3 ,...,x N )= 2π integraldisplay ∞ 0 f(ρ,x 3 ,...,x N )J 0 (ρp)ρdρ, (A.11) showing that F is independent of the angular variable θ. Expression (A.11) is termed the Fourier–Besseltransform of f . The reader can easily verify that f can be recovered from F through f(ρ,x 3 ,...,x N )= integraldisplay ∞ 0 F(p,x 3 ,...,x N )J 0 (ρp) pdp, the inverse Fourier–Bessel transform. A review of complexcontour integration Some powerful techniques for the evaluation of integrals rest on complex variable the- ory. In particular, the computation of the Fourier inversion integral is often aided by these techniques. We therefore provide a brief review of this material. For a fuller discussion the reader may refer to one of many widely available textbooks on complex analysis. We shall denote by f(z) a complex valued function of a complex variable z. That is, f(z)= u(x,y)+ jv(x,y), where the real and imaginary parts u(x,y)andv(x,y)of f are each functions of the real and imaginary parts x and y of z: z = x + jy = Re(z)+ j Im(z). Here j = √ ?1, as is mostly standard in the electrical engineering literature. Limits, di?erentiation, and analyticity. Let w = f(z), and let z 0 = x 0 + jy 0 and w 0 = u 0 + jv 0 be points in the complex z and w planes, respectively. We say that w 0 is the limit of f(z) as z approaches z 0 , and write lim z→z 0 f(z)=w 0 , if and only if both u(x,y)→ u 0 and v(x,y)→v 0 as x → x 0 and y → y 0 independently. The derivative of f(z) at a point z = z 0 is de?ned by the limit f prime (z 0 )= lim z→z 0 f(z)? f(z 0 ) z ? z 0 , ifitexists. Existencerequiresthatthederivativebeindependentofdirectionofapproach; that is, f prime (z 0 )cannot depend on the manner in which z → z 0 in the complex plane. (This turns out to be a much stronger condition than simply requiring that the functions u and v be di?erentiable with respect to the variables x and y.) We say that f(z) is analytic at z 0 if it is di?erentiable at z 0 and at all points in some neighborhood of z 0 . If f(z) is not analytic at z 0 but every neighborhood of z 0 contains a point at which f(z) is analytic, then z 0 is called a singularpointof f(z). Laurent expansions and residues. Although Taylor series can be used to expand complex functions around points of analyticity, we must often expand functions around points z 0 at or near which the functions fail to be analytic. For this we use the Laurent expansion, a generalization of the Taylor expansion involving both positive and negative powers of z ? z 0 : f(z)= ∞ summationdisplay n=?∞ a n (z ? z 0 ) n = ∞ summationdisplay n=1 a ?n (z ? z 0 ) n + ∞ summationdisplay n=0 a n (z ? z 0 ) n . The numbers a n are the coe?cients of the Laurent expansion of f(z) at point z = z 0 . The ?rst series on the right is theprincipalpart of the Laurent expansion, and the second series is theregularpart. The regular part is an ordinary power series, hence it converges in some disk |z?z 0 |< R where R ≥ 0. Puttingζ = 1/(z?z 0 ), the principal part becomes summationtext ∞ n=1 a ?n ζ n ; this power series converges for |ζ|<ρwhereρ ≥ 0, hence the principal part converges for |z ? z 0 | > 1/ρdefinesr. When r < R, the Laurent expansion converges in the annulus r <|z ? z 0 |< R; when r > R, it diverges everywhere in the complex plane. The function f(z) has an isolated singularity at point z 0 if f(z) is not analytic at z 0 but is analytic in the “punctured disk” 0 < |z ? z 0 | < R for some R > 0. Isolated singularities are classi?ed by reference to the Laurent expansion. Three types can arise: 1. Removablesingularity. Thepoint z 0 isaremovablesingularityof f(z)iftheprincipal part of the Laurent expansion of f(z) about z 0 is identically zero (i.e., if a n = 0 for n =?1,?2,?3,...). 2. Poleoforder k. The point z 0 is a pole of order k if the principal part of the Laurent expansion about z 0 contains only ?nitely many terms that form a polynomial of degree k in (z ? z 0 ) ?1 . A pole of order 1is called a simplepole. 3. Essentialsingularity. The point z 0 is an essential singularity of f(z) if the principal part of the Laurent expansion of f(z) about z 0 contains in?nitely many terms (i.e., if a ?n negationslash= 0 for in?nitely many n). The coe?cient a ?1 in the Laurent expansion of f(z) about an isolated singular point z 0 is the residueof f(z) at z 0 . It can be shown that a ?1 = 1 2πj contintegraldisplay Gamma1 f(z)dz (A.12) whereGamma1isanysimpleclosedcurveorientedcounterclockwiseandcontaininginitsinterior z 0 and no other singularity of f(z). Particularly useful to us is the formula for evaluation of residues at pole singularities. If f(z) has a pole of order k at z = z 0 , then the residue of f(z) at z 0 is given by a ?1 = 1 (k ? 1)! lim z→z 0 d k?1 dz k?1 [(z ? z 0 ) k f(z)]. (A.13) Cauchy–Goursat and residue theorems. It can be shown that if f(z) is analytic at all points on and within a simple closed contour C, then contintegraldisplay C f(z)dz = 0. This central result is known as the Cauchy–Goursattheorem. We shall not o?er a proof, but shall proceed instead to derive a useful consequence known as the residue theorem. Figure A.1: Derivation of the residue theorem. FigureA.1 depictsasimpleclosedcurve C enclosing n isolatedsingularitiesofafunction f(z). We assume that f(z) is analytic on and elsewhere within C. Around each singular point z k we have drawn a circle C k so small that it encloses no singular point other than z k ; taken together, the C k (k = 1,...,n) and C form the boundary of a region in which f(z) is everywhere analytic. By the Cauchy–Goursat theorem integraldisplay C f(z)dz+ n summationdisplay k=1 integraldisplay C k f(z)dz = 0. Hence 1 2πj integraldisplay C f(z)dz = n summationdisplay k=1 1 2πj integraldisplay C k f(z)dz, where now the integrations are all performed in a counterclockwise sense. By (A.12) integraldisplay C f(z)dz = 2πj n summationdisplay k=1 r k (A.14) where r 1 ,...,r n are the residues of f(z) at the singularities within C. Contour deformation. Suppose f is analytic in a region D and Gamma1 is a simple closed curve in D.IfGamma1 can be continuously deformed to another simple closed curveGamma1 prime without passing out of D, then integraldisplay Gamma1 prime f(z)dz = integraldisplay Gamma1 f(z)dz. (A.15) Toseethis,considerFigureA.2wherewehaveintroducedanothersetofcurves±γ; these new curves are assumed parallel and in?nitesimally close to each other. Let C be the composite curve consisting ofGamma1, +γ, ?Gamma1 prime , and ?γ, in that order. Since f is analytic on and within C, we have integraldisplay C f(z)dz = integraldisplay Gamma1 f(z)dz+ integraldisplay +γ f(z)dz+ integraldisplay ?Gamma1 prime f(z)dz+ integraldisplay ?γ f(z)dz = 0. But integraltext ?Gamma1 prime f(z)dz =? integraltext Gamma1 prime f(z)dz and integraltext ?γ f(z)dz =? integraltext +γ f(z)dz, hence (A.15) follows. The contour deformation principle often permits us to replace an integration contour by one that is more convenient. Figure A.2: Derivation of the contour deformation principle. Principal value integrals. We must occasionally carry out integrations of the form I = integraldisplay ∞ ?∞ f(x)dx where f(x)has a ?nite number of singularities x k (k = 1,...,n) along the real axis. Such singularities in the integrand force us to interpret I as an improper integral. With just one singularity present at point x 1 , for instance, we de?ne integraldisplay ∞ ?∞ f(x)dx = lim ε→0 integraldisplay x 1 ?ε ?∞ f(x)dx + lim η→0 integraldisplay ∞ x 1 +η f(x)dx provided that both limits exist. When both limits do not exist, we may still be able to obtain a well-de?ned result by computing lim ε→0 parenleftbiggintegraldisplay x 1 ?ε ?∞ f(x)dx + integraldisplay ∞ x 1 +ε f(x)dx parenrightbigg (i.e.,bytakingη=εsothatthelimitsare“symmetric”).Thisquantityiscalledthe Cauchyprincipalvalue of I and is denoted P.V. integraldisplay ∞ ?∞ f(x)dx. More generally, we have P.V. integraldisplay ∞ ?∞ f(x)dx = lim ε→0 parenleftbiggintegraldisplay x 1 ?ε ?∞ f(x)dx + integraldisplay x 2 ?ε x 1 +ε f(x)dx + +···+ integraldisplay x n ?ε x n?1 +ε f(x)dx + integraldisplay ∞ x n +ε f(x)dx parenrightbigg for n singularities x 1 <···< x n . In a large class of problems f(z) (i.e., f(x) with x replaced by the complex variable z) is analytic everywhere except for the presence of ?nitely many simple poles. Some of these may lie on the real axis (at points x 1 <···< x n ,say),andsomemaynot. Considernowtheintegrationcontour C showninFigureA.3.Wechoose R solargeand ε so small that C encloses all the poles of f that lie in the upper half of the complex Figure A.3: Complex plane technique for evaluating a principal value integral. plane. In many problems of interest the integral of f around the large semicircle tends to zero as R →∞and the integrals around the small semicircles are well-behaved as ε → 0. It may then be shown that P.V. integraldisplay ∞ ?∞ f(x)dx =πj n summationdisplay k=1 r k + 2πj summationdisplay UHP r k where r k is the residue at the kth simple pole. The ?rst sum on the right accounts for the contributions of those poles that lie on the real axis; note that it is associated with a factor πj instead of 2πj, since these terms arose from integrals over semicircles rather than over full circles. The second sum, of course, is extended only over those poles that reside in the upper half-plane. Fourier transform solution of the 1-D wave equation Successive applications of the Fourier transform can reduce a partial di?erential equa- tion to an ordinary di?erential equation, and ?nally to an algebraic equation. After the algebraic equation is solved by standard techniques, Fourier inversion can yield a solution to the original partial di?erential equation. We illustrate this by solving the one-dimensional inhomogeneous wave equation parenleftbigg ? 2 ?z 2 ? 1 c 2 ? 2 ?t 2 parenrightbigg ψ(x,y,z,t)= S(x,y,z,t), (A.16) where the ?eldψ is the desired unknown and S is the known source term. For uniqueness of solution we must specify ψ and ?ψ/?z over some z = constant plane. Assume that ψ(x,y,z,t) vextendsingle vextendsingle vextendsingle z=0 = f(x,y,t), (A.17) ? ?z ψ(x,y,z,t) vextendsingle vextendsingle vextendsingle z=0 = g(x,y,t). (A.18) We begin by positing inverse temporal Fourier transform relationships for ψ and S: ψ(x,y,z,t)= 1 2π integraldisplay ∞ ?∞ ? ψ(x,y,z,ω)e jωt dω, S(x,y,z,t)= 1 2π integraldisplay ∞ ?∞ ? S(x,y,z,ω)e jωt dω. Substituting into (A.16), passing the derivatives through the integral, calculating the derivatives, and combining the inverse transforms, we obtain 1 2π integraldisplay ∞ ?∞ bracketleftbiggparenleftbigg ? 2 ?z 2 + k 2 parenrightbigg ? ψ(x,y,z,ω)? ? S(x,y,z,ω) bracketrightbigg e jωt dω= 0 where k =ω/c. By the Fourier integral theorem parenleftbigg ? 2 ?z 2 + k 2 parenrightbigg ? ψ(x,y,z,ω)? ? S(x,y,z,ω)= 0. (A.19) We have thus converted a partial di?erential equation into an ordinary di?erential equa- tion. A spatial transform on z will now convert the ordinary di?erential equation into an algebraic equation. We write ? ψ(x,y,z,ω)= 1 2π integraldisplay ∞ ?∞ ? ψ z (x,y,k z ,ω)e jk z z dk z , ? S(x,y,z,ω)= 1 2π integraldisplay ∞ ?∞ ? S z (x,y,k z ,ω)e jk z z dk z , in (A.19), pass the derivatives through the integral sign, compute the derivatives, and set the integrand to zero to get (k 2 ? k 2 z ) ? ψ z (x,y,k z ,ω)? ? S z (x,y,k z ,ω)= 0; hence ? ψ z (x,y,k z ,ω)=? ? S z (x,y,k z ,ω) (k z ? k)(k z + k) . (A.20) The price we pay for such an easy solution is that we must now perform a two- dimensional Fourier inversion to obtain ψ(x,y,z,t) from ? ψ z (x,y,k z ,ω). It turns out to be easiest to perform the spatial inverse transform ?rst, so let us examine ? ψ(x,y,z,ω)= 1 2π integraldisplay ∞ ?∞ ? ψ z (x,y,k z ,ω)e jk z z dk z . By (A.20) we have ? ψ(x,y,z,ω)= 1 2π integraldisplay ∞ ?∞ [ ? S z (x,y,k z ,ω)] bracketleftbigg ?1 (k z ? k)(k z + k) bracketrightbigg e jk z z dk z , where the integrand involves a product of two functions. With ?g z (k z ,ω)= ?1 (k z ? k)(k z + k) , the convolution theorem gives ? ψ(x,y,z,ω)= integraldisplay ∞ ?∞ ? S(x,y,ζ,ω)?g(z ?ζ,ω)dζ (A.21) Figure A.4: Contour used to compute inverse transform in solution of the 1-D wave equation. where ?g(z,ω)= 1 2π integraldisplay ∞ ?∞ ?g z (k z ,ω)e jk z z dk z = 1 2π integraldisplay ∞ ?∞ ?1 (k z ? k)(k z + k) e jk z z dk z . To compute this integral we use complex plane techniques. The domain of integration extends along the real k z -axis in the complex k z -plane; because of the poles at k z =±k, we must treat the integral as a principal value integral. Denoting I(k z )= ?e jk z z 2π(k z ? k)(k z + k) , we have integraldisplay ∞ ?∞ I(k z )dk z = lim integraldisplay Gamma1 r I(k z )dk z = lim integraldisplay ?k?δ ?Delta1 I(k z )dk z + lim integraldisplay k?δ ?k+δ I(k z )dk z + lim integraldisplay Delta1 k+δ I(k z )dk z where the limits takeδ → 0 andDelta1→∞. Our k z -plane contour takes detours around the poles using semicircles of radius δ, and is closed using a semicircle of radius Delta1 (Figure A.4).Notethatif z> 0,wemustclosethecontourintheupperhalf-plane. By Cauchy’s integral theorem integraldisplay Gamma1 r I(k z )dk z + integraldisplay Gamma1 1 I(k z )dk z + integraldisplay Gamma1 2 I(k z )dk z + integraldisplay Gamma1 Delta1 I(k z )dk z = 0. Thus integraldisplay ∞ ?∞ I(k z )dk z =?lim δ→0 integraldisplay Gamma1 1 I(k z )dk z ? lim δ→0 integraldisplay Gamma1 2 I(k z )dk z ? lim Delta1→∞ integraldisplay Gamma1 Delta1 I(k z )dk z . The contribution from the semicircle of radius Delta1 can be computed by writing k z in polar coordinates as k z =Delta1e jθ : lim Delta1→∞ integraldisplay Gamma1 Delta1 I(k z )dk z = 1 2π lim Delta1→∞ integraldisplay π 0 ?e jzDelta1e jθ (Delta1e jθ ? k)(Delta1e jθ + k) jDelta1e jθ dθ. Using Euler’s identity we can write lim Delta1→∞ integraldisplay Gamma1 Delta1 I(k z )dk z = 1 2π lim Delta1→∞ integraldisplay π 0 ?e ?Delta1z sinθ e jDelta1z cosθ Delta1 2 e 2 jθ jDelta1e jθ dθ. Thus, as long as z > 0 the integrand will decay exponentially as Delta1→∞, and lim Delta1→∞ integraldisplay Gamma1 Delta1 I(k z )dk z → 0. Similarly, integraltext Gamma1 Delta1 I(k z )dk z → 0 when z < 0 if we close the semicircle in the lower half-plane. Thus, integraldisplay ∞ ?∞ I(k z )dk z =?lim δ→0 integraldisplay Gamma1 1 I(k z )dk z ? lim δ→0 integraldisplay Gamma1 2 I(k z )dk z . (A.22) The integrals around the poles can also be computed by writing k z in polar coordinates. Writing k z =?k +δe jθ we ?nd lim δ→0 integraldisplay Gamma1 1 I(k z )dk z = 1 2π lim δ→0 integraldisplay 0 π ?e jz(?k+δe jθ )jδe jθ (?k +δe jθ ? k)(?k +δe jθ + k) dθ = 1 2π integraldisplay π 0 e ?jkz ?2k jdθ =? j 4k e ?jkz . Similarly, using k z = k +δe jθ , we obtain lim δ→0 integraldisplay Gamma1 2 I(k z )dk z = j 4k e jkz . Substituting these into (A.22) we have ?g(z,ω)= j 4k e ?jkz ? j 4k e jkz = 1 2k sin kz, (A.23) valid for z > 0.Forz < 0, we close in the lower half-plane instead and get ?g(z,ω)=? 1 2k sin kz. (A.24) Substituting (A.23) and (A.24) into (A.21) we obtain ? ψ(x,y,z,ω)= integraldisplay z ?∞ ? S(x,y,ζ,ω) sin k(z ?ζ) 2k dζ ? 1 2k integraldisplay ∞ z ? S(x,y,ζ,ω) sin k(z ?ζ) 2k dζ where we have been careful to separate the two cases considered above. To make things a bit easier when we apply the boundary conditions, let us rewrite the above expression. Splitting the domain of integration we write ? ψ(x,y,z,ω)= integraldisplay 0 ?∞ ? S(x,y,ζ,ω) sin k(z ?ζ) 2k dζ + integraldisplay z 0 ? S(x,y,ζ,ω) sin k(z ?ζ) k dζ ? ? integraldisplay ∞ 0 ? S(x,y,ζ,ω) sin k(z ?ζ) 2k dζ. Expansion of the trigonometric functions then gives ? ψ(x,y,z,ω)= integraldisplay z 0 ? S(x,y,ζ,ω) sin k(z ?ζ) k dζ + + sin kz 2k integraldisplay 0 ?∞ ? S(x,y,ζ,ω)cos kζ dζ ? cos kz 2k integraldisplay 0 ?∞ ? S(x,y,ζ,ω)sin kζ dζ ? ? sin kz 2k integraldisplay ∞ 0 ? S(x,y,ζ,ω)cos kζ dζ + cos kz 2k integraldisplay ∞ 0 ? S(x,y,ζ,ω)sin kζ dζ. The last four integrals are independent of z, so we can represent them with functions constant in z. Finally, rewriting the trigonometric functions as exponentials we have ? ψ(x,y,z,ω)= integraldisplay z 0 ? S(x,y,ζ,ω) sin k(z ?ζ) k dζ + ? A(x,y,ω)e ?jkz + ? B(x,y,ω)e jkz . (A.25) This formula for ? ψ was found as a solution to the inhomogeneous ordinary di?erential equation (A.19). Hence, to obtain the complete solution we should add any possible solutions of the homogeneous di?erential equation. Since these are exponentials, (A.25) in fact represents the complete solution, where ? A and ? B are considered unknown and can be found using the boundary conditions. If we are interested in the frequency-domain solution to the wave equation, then we are done. However, since our boundary conditions (A.17) and (A.18) pertain to the time domain, we must temporally inverse transform before we can apply them. Writing the sine function in (A.25) in terms of exponentials, we can express the time-domain solution as ? ψ(x,y,z,t)= integraldisplay z 0 F ?1 braceleftbigg c 2 ? S(x,y,ζ,ω) jω e j ω c (z?ζ) ? c 2 ? S(x,y,ζ,ω) jω e ?j ω c (z?ζ) bracerightbigg dζ + + F ?1 braceleftbig ? A(x,y,ω)e ?j ω c z bracerightbig + F ?1 braceleftbig ? B(x,y,ω)e j ω c z bracerightbig . (A.26) A combination of the Fourier integration and time-shifting theorems gives the general identity F ?1 braceleftbigg ? S(x,y,ζ,ω) jω e ?jωt 0 bracerightbigg = integraldisplay t?t 0 ?∞ S(x,y,ζ,τ)dτ, (A.27) where we have assumed that ? S(x,y,ζ,0)= 0. Using this in (A.26) along with the time- shifting theorem we obtain ψ(x,y,z,t)= c 2 integraldisplay z 0 braceleftBigg integraldisplay t? ζ?z c ?∞ S(x,y,ζ,τ)dτ ? integraldisplay t? z?ζ c ?∞ S(x,y,ζ,τ)dτ bracerightBigg dζ + + a parenleftBig x,y,t ? z c parenrightBig + b parenleftBig x,y,t + z c parenrightBig , or ψ(x,y,z,t)= c 2 integraldisplay z 0 integraldisplay t+ z?ζ c t? z?ζ c S(x,y,ζ,τ)dτ dζ + a parenleftBig x,y,t ? z c parenrightBig + b parenleftBig x,y,t + z c parenrightBig (A.28) where a(x,y,t)= F ?1 [ ? A(x,y,ω)], b(x,y,t)= F ?1 [ ? B(x,y,ω)]. To calculate a(x,y,t) and b(x,y,t), we must use the boundary conditions (A.17) and (A.18). To apply (A.17), we put z = 0 into (A.28) to give a(x,y,t)+ b(x,y,t)= f(x,y,t). (A.29) Using (A.18) is a bit more complicated since we must compute ?ψ/?z, and z is a pa- rameter in the limits of the integral describing ψ. To compute the derivative we apply Leibnitz’ rule for di?erentiation: d dα integraldisplay θ(α) φ(α) f(x,α)dx = parenleftbigg dθ dα parenrightbigg f (θ(α),α)? parenleftbigg dφ dα parenrightbigg f (φ(α),α)+ integraldisplay θ(α) φ(α) ?f ?α dx. (A.30) Using this on the integral term in (A.28) we have ? ?z bracketleftBigg c 2 integraldisplay z 0 parenleftBigg integraldisplay t+ z?ζ c t? z?ζ c S(x,y,ζ,τ)dτ parenrightBigg dζ bracketrightBigg = c 2 integraldisplay z 0 ? ?z parenleftBigg integraldisplay t+ z?ζ c t? z?ζ c S(x,y,ζ,τ)dτ parenrightBigg dζ, which is zero at z = 0.Thus ?ψ ?z vextendsingle vextendsingle vextendsingle z=0 = g(x,y,t)=? 1 c a prime (x,y,t)+ 1 c b prime (x,y,t) where a prime =?a/?t and b prime =?b/?t. Integration gives ? a(x,y,t)+ b(x,y,t)= c integraldisplay t ?∞ g(x,y,τ)dτ. (A.31) Equations (A.29) and (A.31) represent two algebraic equations in the two unknown functions a and b. The solutions are 2a(x,y,t)= f(x,y,t)? c integraldisplay t ?∞ g(x,y,τ)dτ, 2b(x,y,t)= f(x,y,t)+ c integraldisplay t ?∞ g(x,y,τ)dτ. Finally, substitution of these into (A.28) gives us the solution to the inhomogeneous wave equation ψ(x,y,z,t)= c 2 integraldisplay z 0 integraldisplay t+ z?ζ c t? z?ζ c S(x,y,ζ,τ)dτ dζ + 1 2 bracketleftBig f parenleftBig x,y,t ? z c parenrightBig + f parenleftBig x,y,t + z c parenrightBigbracketrightBig + + c 2 integraldisplay t+ z c t? z c g(x,y,τ)dτ. (A.32) This is known as the D’Alembert solution. The terms f(x,y,t ? z/c) contribute to ψ as waves propagating away from the plane z = 0 in the ±z-directions, respectively. The integral over the forcing term S is seen to accumulate values of S over a time interval determined by z ?ζ. The boundary conditions could have been applied while still in the temporal frequency domain (but not the spatial frequency domain, since the spatial position z is lost). But to do this, we would need the boundary conditions to be in the temporal frequency domain. This is easily accomplished by transforming them to give ? ψ(x,y,z,ω) vextendsingle vextendsingle vextendsingle z=0 = ? f(x,y,ω), ? ?z ? ψ(x,y,z,ω) vextendsingle vextendsingle vextendsingle z=0 = ?g(x,y,ω). Applying these to (A.25) (and again using Leibnitz’ rule) we have ? A(x,y,ω)+ ? B(x,y,ω)= ? f(x,y,ω), ?jk ? A(x,y,ω)+ jk ? B(x,y,ω)= ?g(x,y,ω), hence 2 ? A(x,y,ω)= ? f(x,y,ω)? c ?g(x,y,ω) jω , 2 ? B(x,y,ω)= ? f(x,y,ω)+ c ?g(x,y,ω) jω . Finally, substituting these back into (A.25) and expanding the sine function we obtain the frequency-domain solution that obeys the given boundary conditions: ? ψ(x,y,z,ω)= c 2 integraldisplay z 0 bracketleftbigg ? S(x,y,ζ,ω)e j ω c (z?ζ) jω ? ? S(x,y,ζ,ω)e ?j ω c (z?ζ) jω bracketrightbigg dζ + + 1 2 bracketleftbig ? f(x,y,ω)e j ω c z + ? f(x,y,ω)e ?j ω c z bracketrightbig + + c 2 bracketleftbigg ?g(x,y,ω)e j ω c z jω ? ?g(x,y,ω)e ?j ω c z jω bracketrightbigg . This is easily inverted using (A.27) to give (A.32). Fourier transform solution of the 1-D homogeneous wave equation for dissipative media Wave propagation in dissipative media can be studied using the one-dimensional wave equation parenleftbigg ? 2 ?z 2 ? 2Omega1 v 2 ? ?t ? 1 v 2 ? 2 ?t 2 parenrightbigg ψ(x,y,z,t)= S(x,y,z,t). (A.33) This equation is nearly identical to the wave equation for lossless media studied in the previous section, except for the addition of the ?ψ/?t term. This extra term will lead to important physical consequences regarding the behavior of the wave solutions. We shall solve (A.33) using the Fourier transform approach of the previous section, but to keep the solution simple we shall only consider the homogeneous problem. We begin by writing ψ in terms of its inverse temporal Fourier transform: ψ(x,y,z,t)= 1 2π integraldisplay ∞ ?∞ ? ψ(x,y,z,ω)e jωt dω. Substituting this into the homogeneous version of (A.33) and taking the time derivatives, we obtain 1 2π integraldisplay ∞ ?∞ bracketleftbigg (jω) 2 + 2Omega1(jω)?v 2 ? 2 ?z 2 bracketrightbigg ? ψ(x,y,z,ω)e jωt dω= 0. The Fourier integral theorem leads to ? 2 ? ψ(x,y,z,ω) ?z 2 ?κ 2 ? ψ(x,y,z,ω)= 0 (A.34) where κ = 1 v radicalbig p 2 + 2Omega1p with p = jω. We can solve the homogeneous ordinary di?erential equation (A.34) by inspection: ? ψ(x,y,z,ω)= ? A(x,y,ω)e ?κz + ? B(x,y,ω)e κz . (A.35) Here ? A and ? B are frequency-domain coe?cients to be determined. We can either specify these coe?cients directly, or solve for them by applying speci?c boundary conditions. We examine each possibility below. Solution to wave equation by direct application of boundary conditions. The solution to the wave equation (A.33) will be unique if we specify functions f(x,y,t) and g(x,y,t) such that ψ(x,y,z,t) vextendsingle vextendsingle vextendsingle z=0 = f(x,y,t), ? ?z ψ(x,y,z,t) vextendsingle vextendsingle vextendsingle z=0 = g(x,y,t). (A.36) AssumingtheFouriertransformpairs f(x,y,t) ? ? f(x,y,ω)and g(x,y,t) ? ?g(x,y,ω), we can apply the boundary conditions (A.36) in the frequency domain: ? ψ(x,y,z,ω) vextendsingle vextendsingle vextendsingle z=0 = ? f(x,y,ω), ? ?z ? ψ(x,y,z,ω) vextendsingle vextendsingle vextendsingle z=0 = ?g(x,y,ω). From these we ?nd ? A + ? B = ? f, ?κ ? A +κ ? B = ?g,v or ? A = 1 2 bracketleftbigg ? f ? ?g κ bracketrightbigg , ? B = 1 2 bracketleftbigg ? f + ?g κ bracketrightbigg . Substitution into (A.35) gives ? ψ(x,y,z,ω)= ? f(x,y,ω)coshκz + ?g(x,y,ω) sinhκz κ = ? f(x,y,ω) ? ?z ? Q(x,y,z,ω)+ ?g(x,y,ω) ? Q(x,y,z,ω) = ? ψ 1 (x,y,z,ω)+ ? ψ 2 (x,y,z,ω) where ? Q = sinhκz/κ. Assuming that Q(x,y,z,t) ? ? Q(x,y,z,ω), we can employ the convolution theorem to immediately write down ψ(x,y,z,t): ψ(x,y,z,t)= f(x,y,t)? ? ?z Q(x,y,z,t)+ g(x,y,z,t)? Q(x,y,z,t) =ψ 1 (x,y,z,t)+ψ 2 (x,y,z,t). (A.37) To ?nd ψ we must ?rst compute the inverse transform of ? Q. Here we resort to a tabulated result [26]: sinh bracketleftbig a √ p +λ √ p +μ bracketrightbig √ p +λ √ p +μ ? 1 2 e ? 1 2 (μ+λ)t J 0 parenleftbigg 1 2 (λ?μ) radicalbig a 2 ? t 2 parenrightbigg , ?a < t < a. Here a is a positive, ?nite real quantity, and λ and μ are ?nite complex quantities. Outside the range |t|< a the time-domain function is zero. Letting a = z/v, μ= 0, and λ= 2Omega1 in the above expression, we ?nd Q(x,y,z,t)= v 2 e ?Omega1t J 0 parenleftbigg Omega1 v radicalbig z 2 ?v 2 t 2 parenrightbigg [U(t + z/v)? U(t ? z/v)] (A.38) where U(x) is the unit step function (A.5). From (A.37) we see that ψ 2 (x,y,z,t)= integraldisplay ∞ ?∞ g(x,y,t ?τ)Q(x,y,z,τ)dτ = integraldisplay z/v ?z/v g(x,y,t ?τ)Q(x,y,z,τ)dτ. Using the change of variables u = t ?τ and substituting (A.38), we then have ψ 2 (x,y,z,t)= v 2 e ?Omega1t integraldisplay t+ z v t? z v g(x,y,u)e Omega1u J 0 parenleftbigg Omega1 v radicalbig z 2 ?(t ? u) 2 v 2 parenrightbigg du. (A.39) To ?nd ψ 1 we must compute ?Q/?z. Using the product rule we have ?Q(x,y,z,t) ?z = v 2 e ?Omega1t J 0 parenleftbigg Omega1 v radicalbig z 2 ?v 2 t 2 parenrightbigg ? ?z [U(t + z/v)? U(t ? z/v)] + + v 2 e ?Omega1t [U(t + z/v)? U(t ? z/v)] ? ?z J 0 parenleftbigg Omega1 v radicalbig z 2 ?v 2 t 2 parenrightbigg . Next, using dU(x)/dx = δ(x) and remembering that J prime 0 (x) =?J 1 (x) and J 0 (0) = 1,we can write ?Q(x,y,z,t) ?z = 1 2 e ?Omega1t [δ(t + z/v)+δ(t ? z/v)] ? ? zOmega1 2 2v e ?Omega1t J 1 parenleftBig Omega1 v √ z 2 ?v 2 t 2 parenrightBig Omega1 v √ z 2 ?v 2 t 2 [U(t + z/v)? U(t ? z/v)]. Convolving this expression with f(x,y,t) we obtain ψ 1 (x,y,z,t)= 1 2 e ? Omega1 v z f parenleftBig x,y,t ? z v parenrightBig + 1 2 e Omega1 v z f parenleftBig x,y,t + z v parenrightBig ? ? zOmega1 2 2v e ?Omega1t integraldisplay t+ z v t? z v f(x,y,u)e Omega1u J 1 parenleftBig Omega1 v radicalbig z 2 ?(t ? u) 2 v 2 parenrightBig Omega1 v radicalbig z 2 ?(t ? u) 2 v 2 du. (A.40) Finally, adding (A.40) and (A.39), we obtain ψ(x,y,z,t)= 1 2 e ? Omega1 v z f parenleftBig x,y,t ? z v parenrightBig + 1 2 e Omega1 v z f parenleftBig x,y,t + z v parenrightBig ? ? zOmega1 2 2v e ?Omega1t integraldisplay t+ z v t? z v f(x,y,u)e Omega1u J 1 parenleftBig Omega1 v radicalbig z 2 ?(t ? u) 2 v 2 parenrightBig Omega1 v radicalbig z 2 ?(t ? u) 2 v 2 du + + v 2 e ?Omega1t integraldisplay t+ z v t? z v g(x,y,u)e Omega1u J 0 parenleftbigg Omega1 v radicalbig z 2 ?(t ? u) 2 v 2 parenrightbigg du. (A.41) Note that when Omega1= 0 this reduces to ψ(x,y,z,t)= 1 2 f parenleftBig x,y,t ? z v parenrightBig + 1 2 f parenleftBig x,y,t + z v parenrightBig + v 2 integraldisplay t+ z v t? z v g(x,y,u)du, which matches (A.32) for the homogeneous case where S = 0. Solution to wave equation by speci?cation of wave amplitudes. An alternative to direct speci?cation of boundary conditions is speci?cation of the amplitude functions ? A(x,y,ω)and ? B(x,y,ω)ortheirinversetransforms A(x,y,t)and B(x,y,t). Ifwespecify the time-domain functions we can write ψ(x,y,z,t) as the inverse transform of (A.35). For example, a wave traveling in the +z-direction behaves as ψ(x,y,z,t)= A(x,y,t)? F + (x,y,z,t) (A.42) where F + (x,y,z,t) ? e ?κz = e ? z v √ p 2 +2Omega1p . We can ?nd F + using the following Fourier transform pair [26]): e ? x v √ (p+ρ) 2 ?σ 2 ? e ? ρ v x δ(t ? x/v)+ σx v e ?ρt I 1 parenleftBig σ radicalbig t 2 ?(x/v) 2 parenrightBig radicalbig t 2 ?(x/v) 2 , x v < t. (A.43) Here x is real and positive and I 1 (x) is the modi?ed Bessel function of the ?rst kind and order 1. Outside the range x/v< t the time-domain function is zero. Letting ρ =Omega1 and σ =Omega1 we ?nd F + (x,y,z,t)= Omega1 2 z v e ?Omega1t I 1 (Omega1 radicalbig t 2 ?(z/v) 2 ) Omega1 radicalbig t 2 ?(z/v) 2 U(t ? z/v)+ e ? Omega1 v z δ(t ? z/v). (A.44) Note that F + is a real functions of time, as expected. Substituting (A.44) into (A.42) and writing the convolution in integral form we have ψ(x,y,z,t)= integraldisplay ∞ z/v A(x,y,t ?τ) bracketleftBigg Omega1 2 z v e ?Omega1τ I 1 (Omega1 radicalbig τ 2 ?(z/v) 2 ) Omega1 radicalbig τ 2 ?(z/v) 2 bracketrightBigg dτ + + e ? Omega1 v z A parenleftBig x,y,t ? z v parenrightBig , z > 0. (A.45) The 3-D Green’s function for waves in dissipative media To understand the ?elds produced by bounded sources within a dissipative medium we may wish to investigate solutions to the wave equation in three dimensions. The Green’s function approach requires the solution to parenleftbigg ? 2 ? 2Omega1 v 2 ? ?t ? 1 v 2 ? 2 ?t 2 parenrightbigg G(r|r prime ;t)=?δ(t)δ(r ? r prime ) =?δ(t)δ(x ? x prime )δ(y ? y prime )δ(z ? z prime ). That is, we are interested in the impulse response of a point source located at r = r prime . We begin by substituting the inverse temporal Fourier transform relations G(r|r prime ;t)= 1 2π integraldisplay ∞ ?∞ ? G(r|r prime ;ω)e jωt dω, δ(t)= 1 2π integraldisplay ∞ ?∞ e jωt dω, obtaining 1 2π integraldisplay ∞ ?∞ bracketleftbiggparenleftbigg ? 2 ? jω 2Omega1 v 2 ? 1 v 2 (jω) 2 parenrightbigg ? G(r|r prime ;ω)+δ(r ? r prime ) bracketrightbigg e jωt dω= 0. By the Fourier integral theorem we have (? 2 + k 2 ) ? G(r|r prime ;ω)=?δ(r ? r prime ). (A.46) This is known as the Helmholtzequation. Here k = 1 v radicalbig ω 2 ? j2ωOmega1 (A.47) is called the wavenumber. To solve the Helmholtz equation we write ? G in terms of a 3-dimensional inverse Fourier transform. Substitution of ? G(r|r prime ;ω)= 1 (2π) 3 integraldisplay ∞ ?∞ ? G r (k|r prime ;ω)e jk·r d 3 k, δ(r ? r prime )= 1 (2π) 3 integraldisplay ∞ ?∞ e jk·(r?r prime ) d 3 k, into (A.46) gives 1 (2π) 3 integraldisplay ∞ ?∞ bracketleftBig ? 2 parenleftbig ? G r (k|r prime ;ω)e jk·r parenrightbig + k 2 ? G r (k|r prime ;ω)e jk·r + e jk·(r?r prime ) bracketrightBig d 3 k = 0. Here k = ?xk x + ?yk y + ?zk z with |k| 2 = k 2 x + k 2 y + k 2 z = K 2 . Carrying out the derivatives and invoking the Fourier integral theorem we have (K 2 ? k 2 ) ? G r (k|r prime ;ω)= e ?jk·r prime . Solving for ? G and substituting it into the inverse transform relation we have ? G(r|r prime ;ω)= 1 (2π) 3 integraldisplay ∞ ?∞ e jk·(r?r prime ) (K ? k)(K + k) d 3 k. (A.48) Tocomputetheinversetransformintegralin(A.48)wewritethe3-Dtransformvariable in spherical coordinates: k ·(r ? r prime )= KRcosθ, d 3 k = K 2 sinθ dK dθ dφ, where R =|r ? r prime | and θ is the angle between k and r ? r prime . Hence (A.48) becomes ? G(r|r prime ;ω)= 1 (2π) 3 integraldisplay ∞ 0 K 2 dK (K ? k)(K + k) integraldisplay 2π 0 dφ integraldisplay π 0 e jKRcosθ sinθ dθ = 2 (2π) 2 R integraldisplay ∞ 0 K sin(KR) (K ? k)(K + k) dK, or, equivalently, ? G(r|r prime ;ω)= 1 2 jR(2π) 2 integraldisplay ∞ ?∞ e jKR (K ? k)(K + k) KdK? ? 1 2 jR(2π) 2 integraldisplay ∞ ?∞ e ?jkR (K ? k)(K + k) KdK. We can compute the integrals over K using the complex plane technique. We consider K to be a complex variable, and note that for dissipative media we have k = k r + jk i , where k r > 0 and k i < 0. Thus the integrand has poles at K =±k. For the integral involving e +jKR we close the contour in the upper half-plane using a semicircle of radius Delta1 and use Cauchy’s residue theorem. Then at all points on the semicircle the integrand decays exponentially as Delta1 →∞, and there is no contribution to the integral from this part of the contour. The real-line integral is thus equal to 2πj times the residue at K =?k: integraldisplay ∞ ?∞ e jKR (K ? k)(K + k) KdK= 2πj e ?jkR ?2k (?k). For the term involving e ?jKR we close in the lower half-plane and again the contribution from the in?nite semicircle vanishes. In this case our contour is clockwise and so the real line integral is ?2πj times the residue at K = k: integraldisplay ∞ ?∞ e ?jKR (K ? k)(K + k) KdK=?2πj e ?jkR 2k k. Thus ? G(r|r prime ;ω)= e ?jkR 4πR . (A.49) Note that if Omega1= 0 then this reduces to ? G(r|r prime ;ω)= e ?jωR/v 4πR . (A.50) Our last step is to ?nd the temporal Green’s function. Let p = jω. Then we can write ? G(r|r prime ;ω)= e κR 4πR where κ =?jk = 1 v radicalbig p 2 + 2Omega1p. We may ?nd the inverse transform using (A.43). Letting x = R, ρ = Omega1, and σ = Omega1 we ?nd G(r|r prime ;t)= e ? Omega1 v R δ(t ? R/v) 4πR + Omega1 2 4πv e ?Omega1t I 1 parenleftBig Omega1 radicalbig t 2 ?(R/v) 2 parenrightBig Omega1 radicalbig t 2 ?(R/v) 2 U parenleftbigg t ? R v parenrightbigg . We note that in the case of no dissipation where Omega1= 0 this reduces to G(r|r prime ;t)= δ(t ? R/v) 4πR which is the inverse transform of (A.50). Fourier transform representation of the static Green’s function In the study of static ?elds, we shall be interested in the solution to the partial di?er- ential equation ? 2 G(r|r prime )=?δ(r ? r prime )=?δ(x ? x prime )δ(y ? y prime )δ(z ? z prime ). (A.51) Here G(r|r prime ), called the “static Green’s function,” represents the potential at point r produced by a unit point source at point r prime . In Chapter 3 we ?nd that G(r|r prime ) = 1/4π|r ? r prime |. In a variety of problems it is also useful to have G written in terms of an inverse Fourier transform over the variables x and y. Letting G r form a three-dimensional Fourier transform pair with G, we can write G(r|r prime )= 1 (2π) 3 integraldisplay ∞ ?∞ G r (k x ,k y ,k z |r prime )e jk x x e jk y y e jk z z dk x dk y dk z . Substitution into (A.51) along with the inverse transformation representation for the delta function (A.4) gives 1 (2π) 3 ? 2 integraldisplay ∞ ?∞ G r (k x ,k y ,k z |r prime )e jk x x e jk y y e jk z z dk x dk y dk z =? 1 (2π) 3 integraldisplay ∞ ?∞ e jk x (x?x prime ) e jk y (y?y prime ) e jk z (z?z prime ) dk x dk y dk z . We then combine the integrands and move the Laplacian operator through the integral to obtain 1 (2π) 3 integraldisplay ∞ ?∞ bracketleftBig ? 2 parenleftbig G r (k|r prime )e jk·r parenrightbig + e jk·(r?r prime ) bracketrightBig d 3 k = 0, where k = ?xk x + ?yk y + ?zk z . Carrying out the derivatives, 1 (2π) 3 integraldisplay ∞ ?∞ bracketleftBig parenleftbig ?k 2 x ? k 2 y ? k 2 z parenrightbig G r (k|r prime )+ e ?jk·r prime bracketrightBig e jk·r d 3 k = 0. Letting k 2 x + k 2 y = k 2 ρ and invoking the Fourier integral theorem we get the algebraic equation parenleftbig ?k 2 ρ ? k 2 z parenrightbig G r (k|r prime )+ e ?jk·r prime = 0, which we can easily solve for G r : G r (k|r prime )= e ?jk·r prime k 2 ρ + k 2 z . (A.52) Equation (A.52) gives us a 3-D transform representation for the Green’s function. Since we desire the 2-D representation, we shall have to perform the inverse transform over k z . Writing G xy (k x ,k y ,z|r prime )= 1 2π integraldisplay ∞ ?∞ G r (k x ,k y ,k z |r prime )e jk z z dk z we have G xy (k x ,k y ,z|r prime )= 1 2π integraldisplay ∞ ?∞ e ?jk x x prime e ?jk y y prime e jk z (z?z prime ) k 2 ρ + k 2 z dk z . (A.53) To compute this integral, we let k z be a complex variable and consider a closed contour in the complex plane, consisting of a semicircle and the real axis. As previously discussed, we compute the principal value integral as the semicircle radius Delta1 →∞, and ?nd that the contribution along the semicircle reduces to zero. Hence we can use Cauchy’s residue theorem (A.14) to obtain the real-line integral: G xy (k x ,k y ,z|r prime )= 2πj res braceleftBigg 1 2π e ?jk x x prime e ?jk y y prime e jk z (z?z prime ) k 2 ρ + k 2 z bracerightBigg . Here res{ f(k z )} denotes the residues of the function f(k z ). The integrand in (A.53) has poles of order 1at k z =±jk ρ , k ρ ≥ 0.Ifz ? z prime > 0 we close in the upper half-plane and enclose only the pole at k z = jk ρ . Computing the residue using (A.13), we obtain G xy (k x ,k y ,z|r prime )= j e ?jk x x prime e ?jk y y prime e ?k ρ (z?z prime ) 2 jk ρ , z > z prime . Since z > z prime this function decays for increasing z, as expected physically. For z?z prime < 0 we close in the lower half-plane, enclosing the pole at k z =?jk ρ and incurring an additional negative sign since our contour is now clockwise. Evaluating the residue we have G xy (k x ,k y ,z|r prime )=?j e ?jk x x prime e ?jk y y prime e k ρ (z?z prime ) ?2 jk ρ , z < z prime . We can combine both cases z > z prime and z < z prime by using the absolute value function: G xy (k x ,k y ,z|r prime )= e ?jk x x prime e ?jk y y prime e ?k ρ |z?z prime | 2k ρ . (A.54) Finally, wesubstitute(A.54)intotheinversetransformformula. ThisgivestheGreen’s function representation G(r|r prime )= 1 4π|r ? r prime | = 1 (2π) 2 integraldisplay ∞ ?∞ e ?k ρ |z?z prime | 2k ρ e jk ρ ·(r?r prime ) d 2 k ρ , (A.55) where k ρ = ?xk x + ?yk y , k ρ =|k ρ |, and d 2 k ρ = dk x dk y . On occasion we may wish to represent the solution of the homogeneous (Laplace) equation ? 2 ψ(r)= 0 intermsofa2-DFouriertransform. Inthiscasewerepresentψ asa2-Dinversetransform and substitute to obtain 1 (2π) 2 integraldisplay ∞ ?∞ ? 2 parenleftbig ψ xy (k x ,k y ,z)e jk x x e jk y y parenrightbig dk x dk y = 0. Carrying out the derivatives and invoking the Fourier integral theorem we ?nd that parenleftbigg ? 2 ?z 2 ? k 2 ρ parenrightbigg ψ xy (k x ,k y ,z)= 0. Hence ψ xy (k x ,k y ,z)= Ae k ρ z + Be ?k ρ z where A and B are constants with respect to z. Inverse transformation gives ψ(r)= 1 (2π) 2 integraldisplay ∞ ?∞ bracketleftbig A(k ρ )e k ρ z + B(k ρ )e ?k ρ z bracketrightbig e jk ρ ·r d 2 k ρ . (A.56) A.2 Vector transport theorems We are often interested in the time rate of change of some ?eld integrated over a moving volume or surface. Such a derivative may be used to describe the transport of a physical quantity (e.g., charge, momentum, energy) through space. Many of the relevant theorems are derived in this section. The results ?nd application in the development of the large-scale forms of Maxwell equations, the continuity equation, and the Poynting theorem. Partial, total, and material derivatives The key to understanding transport theorems lies in the di?erence between the various meansoftime-di?erentiatinga?eld. Considerascalar?eld T(r,t)(whichcouldrepresent one component of a vector or dyadic ?eld). If we ?x our position within the ?eld and examine how the ?eld varies with time, we describe thepartialderivative of T. However, this may not be the most useful means of measuring the time rate of change of a ?eld. For instance, in mechanics we might be interested in the rate at which water cools as it sinks to the bottom of a container. In this case, T could represent temperature. We could create a “depth pro?le” at any given time (i.e., measure T(r,t 0 ) for some ?xed t 0 ) by taking simultaneous data from a series of temperature probes at varying depths. We could also create a temporal pro?le at any given depth (i.e., measure T(r 0 ,t) for some ?xed r 0 ) by taking continuous data from a probe ?xed at that depth. But neither of these would describe how an individual sinking water particle “experiences” a change in temperature over time. Instead, we could use a probe that descends along with a particular water packet (i.e., volume element), measuring the time rate of temperature change of that element. This rate of change is called the convective or material derivative, since it corresponds to a situation in which a physical material quantity is followed as the derivative is calculated. We anticipate that this quantity will depend on (1) the time rate of change of T at each ?xed point that the particle passes, and (2) the spatial rate of change of T as well as the rapidity with which the packet of interest is swept through that space gradient. The faster the packet descends, or the faster the temperature cools with depth, the larger the material derivative should be. To compute the material derivative we describe the position of a water packet by the vector r(t)= ?xx(t)+ ?yy(t)+ ?zz(t). Because no two packets can occupy the same place at the same time, the speci?cation of r(0)= r 0 uniquely describes (or “tags”) a particular packet. The time rate of change of r with r 0 held constant (the material derivative of the position vector) is thus the velocity ?eld u(r,t) of the ?uid: parenleftbigg dr dt parenrightbigg r 0 = Dr Dt = u. (A.57) Here we use the “big D” notation to denote the material derivative, thereby avoiding confusion with the partial and total derivatives described below. To describe the time rate of change of the temperature of a particular water packet, we only need to hold r 0 constant while we examine the change. If we write the temperature as T(r,t)= T(r(r 0 ,t),t)= T [x(r 0 ,t),y(r 0 ,t),z(r 0 ,t),t], then we can use the chain rule to ?nd the time rate of change of T with r 0 held constant: DT Dt = parenleftbigg dT dt parenrightbigg r 0 = parenleftbigg ?T ?x parenrightbiggparenleftbigg dx dt parenrightbigg r 0 + parenleftbigg ?T ?y parenrightbiggparenleftbigg dy dt parenrightbigg r 0 + parenleftbigg ?T ?z parenrightbiggparenleftbigg dz dt parenrightbigg r 0 + ?T ?t . We recognize the partial derivatives of the coordinates as the components of the material velocity (A.57), and thus can write DT Dt = ?T ?t + u x ?T ?x + u y ?T ?y + u z ?T ?z = ?T ?t + u ·?T. As expected, the material derivative depends on both the local time rate of change and the spatial rate of change of temperature. Suppose next that our probe is motorized and can travel about in the sinking water. If the probe sinks faster than the surrounding water, the time rate of change (measured by the probe) should exceed the material derivative. Let the probe position and velocity be r(t)= ?xx(t)+ ?yy(t)+ ?zz(t), v(r,t)= ?x dx(t) dt + ?y dy(t) dt + ?z dz(t) dt . We can use the chain rule to determine the time rate of change of the temperature observed by the probe, but in this case we do not constrain the velocity components to represent the moving ?uid. Thus, we merely obtain dT dt = ?T ?x dx dt + ?T ?y dy dt + ?T ?z dz dt + ?T ?t = ?T ?t + v ·?T. This is called the totalderivative of the temperature ?eld. In summary, the time rate of change of a scalar ?eld T seen by an observer moving with arbitrary velocity v is given by the total derivative dT dt = ?T ?t + v ·?T. (A.58) If the velocity of the observer happens to match the velocity u of a moving substance, the time rate of change is the material derivative DT Dt = ?T ?t + u ·?T. (A.59) Wecanobtainthematerialderivativeofavector?eld F bycomponent-wiseapplication of (A.59): DF Dt = D Dt bracketleftbig ?xF x + ?yF y + ?zF z bracketrightbig = ?x ?F x ?t + ?y ?F y ?t + ?z ?F z ?t + ?x [u ·(?F x )] + ?y bracketleftbig u ·(?F y ) bracketrightbig + ?z [u ·(?F z )]. Figure A.5: Derivation of the Helmholtz transport theorem. Using the notation u ·?=u x ? ?x + u y ? ?y + u z ? ?z we can write DF Dt = ?F ?t +(u ·?)F. (A.60) This is the material derivative of a vector ?eld F when u describes the motion of a physical material. Similarly, the total derivative of a vector ?eld is dF dt = ?F ?t +(v ·?)F where v is arbitrary. The Helmholtz and Reynolds transport theorems We choose the intuitive approach taken by Tai [190] and Whitaker [214]. Consider an open surface S(t) moving through space and possibly deforming as it moves. The velocity of the points comprising the surface is given by the vector ?eld v(r,t). We are interested in computing the time derivative of the ?ux of a vector ?eld F(r,t) through S(t): ψ(t)= d dt integraldisplay S(t) F(r,t)· dS = lim Delta1t→0 integraltext S(t+Delta1t) F(r,t +Delta1t)· dS ? integraltext S(t) F(r,t)· dS Delta1t . (A.61) Here S(t+Delta1t)= S 2 isfoundbyextendingeachpointon S(t)= S 1 throughadisplacement vDelta1t,asshowninFigureA.5.SubstitutingtheTaylorexpansion F(r,t +Delta1t)= F(r,t)+ ?F(r,t) ?t Delta1t +··· into (A.61), we ?nd that only the ?rst two terms give non-zero contributions to the integral and ψ(t)= integraldisplay S(t) ?F(r,t) ?t · dS + lim Delta1t→0 integraltext S 2 F(r,t)· dS ? integraltext S 1 F(r,t)· dS Delta1t . (A.62) ThesecondtermontherightcanbeevaluatedwiththehelpofFigureA.5.Asthesurface moves through a displacement vDelta1t it sweeps out a volume region Delta1V that is bounded on the back by S 1 , on the front by S 2 , and on the side by a surface S 3 = Delta1S. We can thus compute the two surface integrals in (A.62) as the di?erence between contributions from the surface enclosing Delta1V and the side surface Delta1S (remembering that the normal to S 1 in (A.62) points into Delta1V). Thus ψ(t)= integraldisplay S(t) ?F(r,t) ?t · dS + lim Delta1t→0 contintegraltext S 1 +S 2 +Delta1S F(r,t)· dS ? integraltext Delta1S F(r,t)· dS 3 Delta1t = integraldisplay S(t) ?F(r,t) ?t · dS + lim Delta1t→0 integraltext Delta1V ?·F(r,t)dV 3 ? integraltext Delta1S F(r,t)· dS 3 Delta1t by the divergence theorem. To compute the integrals over Delta1S and Delta1V we note from FigureA.5thattheincrementalsurfaceandvolumeelementsarejust dS 3 = dl ×(vDelta1t), dV 3 =(vDelta1t)· dS. Then, since F · [dl ×(vDelta1t)] =Delta1t(v × F)· dl, we have ψ(t)= integraldisplay S(t) ?F(r,t) ?t · dS + lim Delta1t→0 Delta1t integraltext S(t) [v?·F(r,t)] · dS Delta1t ? lim Delta1t→0 Delta1t contintegraltext Gamma1 [v × F(r,t)] · dl Delta1t . Taking the limit and using Stokes’s theorem on the last integral we have ?nally d dt integraldisplay S(t) F · dS = integraldisplay S(t) bracketleftbigg ?F ?t + v?·F ??×(v × F) bracketrightbigg · dS, (A.63) which is the Helmholtztransporttheorem [190, 43]. In case the surface corresponds to a moving physical material, we may wish to write the Helmholtz transport theorem in terms of the material derivative. We can set v = u and use ?×(u × F)= u(?·F)? F(?·u)+(F ·?)u ?(u ·?)F and (A.60) to obtain d dt integraldisplay S(t) F · dS = integraldisplay S(t) bracketleftbigg DF Dt + F(?·u)?(F ·?)u bracketrightbigg · dS. If S(t) in (A.63) is closed, enclosing a volume region V(t), then contintegraldisplay S(t) [?×(v × F)] · dS = integraldisplay V(t) ?·[?×(v × F)] dV = 0 by the divergence theorem and (B.49). In this case the Helmholtz transport theorem becomes d dt contintegraldisplay S(t) F · dS = contintegraldisplay S(t) bracketleftbigg ?F ?t + v?·F bracketrightbigg · dS. (A.64) We now come to an essential tool that we employ throughout the book. Using the divergence theorem we can rewrite (A.64) as d dt integraldisplay V(t) ?·F dV = integraldisplay V(t) ?· ?F ?t dV + contintegraldisplay S(t) (?·F)v · dS. Replacing ?·F by the scalar ?eld ρ we have d dt integraldisplay V(t) ρdV = integraldisplay V(t) ?ρ ?t dV + contintegraldisplay S(t) ρv · dS. (A.65) In this general form of the transport theorem v is an arbitrary velocity. In most appli- cations v = u describes the motion of a material substance; then D Dt integraldisplay V(t) ρdV = integraldisplay V(t) ?ρ ?t dV + contintegraldisplay S(t) ρu · dS, (A.66) which is the Reynolds transport theorem [214]. The D/Dt notation implies that V(t) retains exactly the same material elements as it moves and deforms to follow the material substance. We may rewrite the Reynolds transport theorem in various forms. By the divergence theorem we have d dt integraldisplay V(t) ρdV = integraldisplay V(t) bracketleftbigg ?ρ ?t +?·(ρv) bracketrightbigg dV. Setting v = u, using (B.42), and using (A.59) for the material derivative of ρ, we obtain D Dt integraldisplay V(t) ρdV = integraldisplay V(t) bracketleftbigg Dρ Dt +ρ?·u bracketrightbigg dV. (A.67) We may also generate a vector form of the general transport theorem by taking ρ in (A.65) to be a component of a vector. Assembling all of the components we have d dt integraldisplay V(t) A dV = integraldisplay V(t) ?A ?t dV + contintegraldisplay S(t) A(v · ?n)dS. (A.68) A.3 Dyadic analysis Dyadic analysis was introduced in the late nineteenth century by Gibbs to generalize vector analysis to problems in which the components of vectors are related in a linear manner. It has now been widely supplanted by tensor theory, but maintains a foothold in engineering where the transformation properties of tensors are not paramount (except, of course, in considerations such as those involving special relativity). Terms such as “tensor permittivity” and “dyadic permittivity” are often used interchangeably. Component form representation. We wish to write one vector ?eld A(r,t) as a linear function of another vector ?eld B(r,t): A = f(B). By this we mean that each component of A is a linear combination of the components of B: A 1 (r,t)= a 11 prime B 1 prime(r,t)+ a 12 prime B 2 prime(r,t)+ a 13 prime B 3 prime(r,t), A 2 (r,t)= a 21 prime B 1 prime(r,t)+ a 22 prime B 2 prime(r,t)+ a 23 prime B 3 prime(r,t), A 3 (r,t)= a 31 prime B 1 prime(r,t)+ a 32 prime B 2 prime(r,t)+ a 33 prime B 3 prime(r,t). Here the a ij prime may depend on space and time (or frequency). The prime on the second index indicates that A and B may be expressed in distinct coordinate frames ( ? i 1 , ? i 2 , ? i 3 ) and ( ? i 1 prime, ? i 2 prime, ? i 3 prime), respectively. We have A 1 = parenleftbig a 11 prime ? i 1 prime + a 12 prime ? i 2 prime + a 13 prime ? i 3 prime parenrightbig · parenleftbig ? i 1 prime B 1 prime + ? i 2 prime B 2 prime + ? i 3 prime B 3 prime parenrightbig , A 2 = parenleftbig a 21 prime ? i 1 prime + a 22 prime ? i 2 prime + a 23 prime ? i 3 prime parenrightbig · parenleftbig ? i 1 prime B 1 prime + ? i 2 prime B 2 prime + ? i 3 prime B 3 prime parenrightbig , A 3 = parenleftbig a 31 prime ? i 1 prime + a 32 prime ? i 2 prime + a 33 prime ? i 3 prime parenrightbig · parenleftbig ? i 1 prime B 1 prime + ? i 2 prime B 2 prime + ? i 3 prime B 3 prime parenrightbig , and since B = ? i 1 prime B 1 prime + ? i 2 prime B 2 prime + ? i 3 prime B 3 prime we can write A = ? i 1 (a prime 1 · B)+ ? i 2 (a prime 2 · B)+ ? i 3 (a prime 3 · B) where a prime 1 = a 11 prime ? i 1 prime + a 12 prime ? i 2 prime + a 13 prime ? i 3 prime, a prime 2 = a 21 prime ? i 1 prime + a 22 prime ? i 2 prime + a 23 prime ? i 3 prime, a prime 3 = a 31 prime ? i 1 prime + a 32 prime ? i 2 prime + a 33 prime ? i 3 prime. In shorthand notation A = ˉa · B (A.69) where ˉa = ? i 1 a prime 1 + ? i 2 a prime 2 + ? i 3 a prime 3 . (A.70) Written out, the quantity ˉa looks like ˉa = a 11 prime( ? i 1 ? i 1 prime)+ a 12 prime( ? i 1 ? i 2 prime)+ a 13 prime( ? i 1 ? i 3 prime)+ + a 21 prime( ? i 2 ? i 1 prime)+ a 22 prime( ? i 2 ? i 2 prime)+ a 23 prime( ? i 2 ? i 3 prime)+ + a 31 prime( ? i 3 ? i 1 prime)+ a 32 prime( ? i 3 ? i 2 prime)+ a 33 prime( ? i 3 ? i 3 prime). Terms such as ? i 1 ? i 1 prime are called dyads, while sums of dyads such as ˉa are called dyadics. The components a ij prime of ˉa may be conveniently placed into an array: [ˉa] = ? ? a 11 prime a 12 prime a 13 prime a 21 prime a 22 prime a 23 prime a 31 prime a 32 prime a 33 prime ? ? . Writing [A] = ? ? A 1 A 2 A 3 ? ? , [B] = ? ? B 1 prime B 2 prime B 3 prime ? ? , we see that A = ˉa · B can be written as [A] = [ˉa][B] = ? ? a 11 prime a 12 prime a 13 prime a 21 prime a 22 prime a 23 prime a 31 prime a 32 prime a 33 prime ? ? ? ? B 1 prime B 2 prime B 3 prime ? ? . Note carefully that in (A.69) ˉa operates on B from the left. A reorganization of the components of ˉa allows us to write ˉa = a 1 ? i 1 prime + a 2 ? i 2 prime + a 3 ? i 3 prime (A.71) where a 1 = a 11 prime ? i 1 + a 21 prime ? i 2 + a 31 prime ? i 3 , a 2 = a 12 prime ? i 1 + a 22 prime ? i 2 + a 32 prime ? i 3 , a 3 = a 13 prime ? i 1 + a 23 prime ? i 2 + a 33 prime ? i 3 . We may now consider using ˉa to operate on a vector C = ? i 1 C 1 + ? i 2 C 2 + ? i 3 C 3 from the right: C · ˉa =(C · a 1 ) ? i 1 prime +(C · a 2 ) ? i 2 prime +(C · a 3 ) ? i 3 prime. In matrix form C · ˉa is [ˉa] T [C] = ? ? a 11 prime a 21 prime a 31 prime a 12 prime a 22 prime a 32 prime a 13 prime a 23 prime a 33 prime ? ? ? ? C 1 C 2 C 3 ? ? where the superscript “T” denotes the matrix transpose operation. That is, C · ˉa = ˉa T · C where ˉa T is the transpose of ˉa. If the primed and unprimed frames coincide, then ˉa = a 11 ( ? i 1 ? i 1 )+ a 12 ( ? i 1 ? i 2 )+ a 13 ( ? i 1 ? i 3 )+ + a 21 ( ? i 2 ? i 1 )+ a 22 ( ? i 2 ? i 2 )+ a 23 ( ? i 2 ? i 3 )+ + a 31 ( ? i 3 ? i 1 )+ a 32 ( ? i 3 ? i 2 )+ a 33 ( ? i 3 ? i 3 ). In this case we may compare the results of ˉa · B and B · ˉa for a given vector B = ? i 1 B 1 + ? i 2 B 2 + ? i 3 B 3 . We leave it to the reader to verify that in general B · ˉa negationslash= ˉa · B. Vector form representation. We can express dyadics in coordinate-free fashion if we expand the concept of a dyad to permit entities such as AB. Here A and B are called the antecedent and consequent, respectively. The operation rules (AB)· C = A(B · C), C ·(AB)=(C · A)B, de?ne the anterior and posterior products of AB with a vector C, and give results consistent with our prior component notation. Sums of dyads such as AB + CD are called dyadicpolynomials, or dyadics. The simple dyadic AB =(A 1 ? i 1 + A 2 ? i 2 + A 3 ? i 3 )(B 1 prime ? i 1 prime + B 2 prime ? i 2 prime + B 3 prime ? i 3 prime) can be represented in component form using AB = ? i 1 a prime 1 + ? i 2 a prime 2 + ? i 3 a prime 3 where a prime 1 = A 1 B 1 prime ? i 1 prime + A 1 B 2 prime ? i 2 prime + A 1 B 3 prime ? i 3 prime, a prime 2 = A 2 B 1 prime ? i 1 prime + A 2 B 2 prime ? i 2 prime + A 2 B 3 prime ? i 3 prime, a prime 3 = A 3 B 1 prime ? i 1 prime + A 3 B 2 prime ? i 2 prime + A 3 B 3 prime ? i 3 prime, or using AB = a 1 ? i 1 prime + a 2 ? i 2 prime + a 3 ? i 3 prime where a 1 = ? i 1 A 1 B 1 prime + ? i 2 A 2 B 1 prime + ? i 3 A 3 B 1 prime, a 2 = ? i 1 A 1 B 2 prime + ? i 2 A 2 B 2 prime + ? i 3 A 3 B 2 prime, a 3 = ? i 1 A 1 B 3 prime + ? i 2 A 2 B 3 prime + ? i 3 A 3 B 3 prime. Note that if we write ˉa = AB then a ij = A i B j prime. A simple dyad AB by itself cannot represent a general dyadic ˉa; only six independent quantities are available in AB (the three components of A and the three components of B), while an arbitrary dyadic has nine independent components. However, it can be shown that any dyadic can be written as a sum of three dyads: ˉa = AB + CD + EF. This is called a vectorrepresentation of ˉa.IfV is a vector, the distributive laws ˉa · V =(AB + CD + EF)· V = A(B · V)+ C(D · V)+ E(F · V), V · ˉa = V ·(AB + CD + EF)=(V · A)B +(V · C)D +(V · E)F, apply. Dyadic algebra and calculus. The cross product of a vector with a dyadic produces another dyadic. If ˉa = AB + CD + EF then by de?nition ˉa × V = A(B × V)+ C(D × V)+ E(F × V), V × ˉa =(V × A)B +(V × C)D +(V × E)F. The corresponding component forms are ˉa × V = ? i 1 (a prime 1 × V)+ ? i 2 (a prime 2 × V)+ ? i 3 (a prime 3 × V), V × ˉa =(V × a 1 ) ? i 1 prime +(V × a 2 ) ? i 2 prime +(V × a 3 ) ? i 3 prime, where we have used (A.70) and (A.71), respectively. Interactions between dyads or dyadics may also be de?ned. The dot product of two dyads AB and CD is a dyad given by (AB)·(CD)= A(B · C)D =(B · C)(AD). The dot product of two dyadics can be found by applying the distributive property. If α is a scalar, then the product αˉa is a dyadic with components equal to α times the components of ˉa. Dyadic addition may be accomplished by adding individual dyadic components as long as the dyadics are expressed in the same coordinate system. Sub- traction is accomplished by adding the negative of a dyadic, which is de?ned through scalar multiplication by ?1. Some useful dyadic identities appear in Appendix B. Many more can be found in Van Bladel [202]. The various vector derivatives may also be extended to dyadics. Computations are easiest in rectangular coordinates, since ? i 1 = ?x, ? i 2 = ?y, and ? i 3 = ?z are constant with position. The dyadic ˉa = a x ?x + a y ?y + a z ?z has divergence ?·ˉa =(?·a x )?x +(?·a y )?y +(?·a z )?z, and curl ?×ˉa =(?×a x )?x +(?×a y )?y +(?×a z )?z. Note that the divergence of a dyadic is a vector while the curl of a dyadic is a dyadic. The gradient of a vector a = a x ?x + a y ?y + a z ?z is ?a =(?a x )?x +(?a y )?y +(?a z )?z, a dyadic quantity. Thedyadicderivativesmaybeexpressedincoordinate-freenotationbyusingthevector representation. The dyadic AB has divergence ?·(AB)=(?·A)B + A ·(?B) and curl ?×(AB)=(?×A)B ? A ×(?B). The Laplacian of a dyadic is a dyadic given by ? 2 ˉa =?(?·ˉa)??×(?×ˉa). The divergence theorem for dyadics is integraldisplay V ?·ˉa dV = contintegraldisplay S ?n · ˉa dS. Some of the other common di?erential and integral identities for dyadics can be found in Van Bladel [202] and Tai [192]. Special dyadics. We say that ˉa is symmetric if B · ˉa = ˉa · B for any vector B. This requires ˉa T = ˉa, i.e., a ij prime = a ji prime. We say that ˉa is antisymmetric if B · ˉa =?ˉa · B for any B. In this case ˉa T =?ˉa. That is, a ij prime =?a ji prime and a ii prime = 0. A symmetric dyadic has only six independent components while an antisymmetric dyadic has only three. The reader can verify that any dyadic can be decomposed into symmetric and antisymmetric parts as ˉa = 1 2 parenleftbig ˉa + ˉa T parenrightbig + 1 2 parenleftbig ˉa ? ˉa T parenrightbig . A simple example of a symmetric dyadic is the unitdyadic ˉ I de?ned by ˉ I = ? i 1 ? i 1 + ? i 2 ? i 2 + ? i 3 ? i 3 . This quantity often arises in the manipulation of dyadic equations, and satis?es A · ˉ I = ˉ I · A = A for any vector A. In matrix form ˉ I is the identity matrix: [ ˉ I] = ? ? 100 010 001 ? ? . The components of a dyadic may be complex. We say that ˉa is hermitian if B · ˉa = ˉa ? · B (A.72) holds for any B. This requires that ˉa ? = ˉa T . Taking the transpose we can write ˉa =(ˉa ? ) T = ˉa ? where “?” stands for the conjugate-transpose operation. We say that ˉa isanti-hermitian if B · ˉa =?ˉa ? · B (A.73) for arbitrary B. In this case ˉa ? =?ˉa T . Any complex dyadic can be decomposed into hermitian and anti-hermitian parts: ˉa = 1 2 parenleftbig ˉa H + ˉa A parenrightbig (A.74) where ˉa H = ˉa + ˉa ? , ˉa A = ˉa ? ˉa ? . (A.75) A dyadic identity important in the study of material parameters is B · ˉa ? · B ? = B ? · ˉa ? · B. (A.76) We show this by decomposing ˉa according to (A.74), giving B · ˉa ? · B ? = 1 2 parenleftBig bracketleftbig B ? · ˉa H bracketrightbig ? + bracketleftbig B ? · ˉa A bracketrightbig ? parenrightBig · B ? where we have used (B · ˉa) ? =(B ? · ˉa ? ). Applying (A.72) and (A.73) we obtain B · ˉa ? · B ? = 1 2 parenleftBig bracketleftbig ˉa H? · B ? bracketrightbig ? ? bracketleftbig ˉa A? · B ? bracketrightbig ? parenrightBig · B ? = B ? · 1 2 parenleftbigbracketleftbig ˉa H · B bracketrightbig ? bracketleftbig ˉa A · B bracketrightbigparenrightbig = B ? · parenleftbigg 1 2 bracketleftbig ˉa H ? ˉa A bracketrightbig · B parenrightbigg . Since the term in brackets is ˉa H ? ˉa A = 2ˉa ? by (A.75), the identity is proved. A.4 Boundary value problems Many physical phenomena may be described mathematically as the solutions tobound- ary value problems. The desired physical quantity (usually called a “?eld”) in a certain region of space is found by solving one or more partial di?erential equations subject to certain conditions over the boundary surface. The boundary conditions may specify the values of the ?eld, some manipulated version of the ?eld (such as the normal derivative), or a relationship between ?elds in adjoining regions. If the ?eld varies with time as well as space, initial or ?nal values of the ?eld must also be speci?ed. Particularly important is whether a boundary value problem is well-posed and therefore has a unique solution which depends continuously on the data supplied. This depends on the forms of the dif- ferential equation and boundary conditions. The well-posedness of Maxwell’s equations is discussed in § 2.2. The importance of boundary value problems has led to an array of techniques, both analytical and numerical, for solving them. Many problems (such as boundary value problems involving Laplace’s equation) may be solved in several di?erent ways. Unique- nesspermitsanengineertofocusattentiononwhichtechniquewillyieldthemoste?cient solution. In this section we concentrate on the separation of variables technique, which is widely applied in the solution of Maxwell’s equations. We ?rst discuss eigenvalue prob- lems and then give an overview of separation of variables. Finally we consider a number of example problems in each of the three common coordinate systems. Sturm–Liouville problems and eigenvalues The partial di?erential equations of electromagnetics can often be reduced to ordinary di?erential equations. In some cases symmetry permits us to reduce the number of dimensions by inspection; in other cases, we may employ an integral transform (e.g., the Fourier transform) or separation of variables. The resulting ordinary di?erential equations may be viewed as particular cases of the Sturm–Liouvilledi?erentialequation d dx bracketleftbigg p(x) dψ(x) dx bracketrightbigg + q(x)ψ(x)+λσ(x)ψ(x)= 0, x ∈ [a,b]. (A.77) In linear operator notation L[ ψ(x)] =?λσ(x)ψ(x), (A.78) where L is the linear Sturm–Liouville operator L = parenleftbigg d dx bracketleftbigg p(x) d dx bracketrightbigg + q(x) parenrightbigg . Obviously ψ(x) = 0 satis?es (A.78). However, for certain values of λ dependent on p, q, σ, and the boundary conditions we impose, (A.78) has non-trivial solutions. Each λ that satis?es (A.78) is an eigenvalue of L , and any non-trivial solution associated with that eigenvalue is an eigenfunction. Taken together, the eigenvalues of an operator form its eigenvaluespectrum. We shall restrict ourselves to the case in which L is self-adjoint. Assume p, q, and σ are real and continuous on [a,b]. It is straightforward to show that for any two functions u(x) and v(x) Lagrange’s identity u L [v] ?v L [u] = d dx bracketleftbigg p parenleftbigg u dv dx ?v du dx parenrightbiggbracketrightbigg (A.79) holds. Integration gives Green’s formula integraldisplay b a (u L[v] ?vL[u]) dx = p parenleftbigg u dv dx ?v du dx parenrightbigg vextendsingle vextendsingle vextendsingle b a . The operator L is self-adjoint if its associated boundary conditions are such that p parenleftbigg u dv dx ?v du dx parenrightbigg vextendsingle vextendsingle vextendsingle b a = 0. (A.80) Possible sets of conditions include the homogeneous boundary conditions α 1 ψ(a)+β 1 ψ prime (a)= 0,α 2 ψ(b)+β 2 ψ prime (b)= 0, (A.81) and the periodic boundary conditions ψ(a)=ψ(b), p(a)ψ prime (a)= p(b)ψ prime (b). (A.82) By imposing one of these sets on (A.78) we obtain a Sturm–Liouvilleproblem. The self-adjoint Sturm–Liouville operator has some nice properties. Each eigenvalue is real, and the eigenvalues form a denumerable set with no cluster point. Moreover, eigen- functions corresponding to distinct eigenvalues are orthogonal, and the eigenfunctions form a complete set. Hence we can expand any su?ciently smooth function in terms of the eigenfunctions of a problem. We discuss this further below. A regular Sturm–Liouville problem involves a self-adjoint operator L with p(x)>0 and σ(x)>0 everywhere, and the homogeneous boundary conditions (A.81). If p or σ vanishes at an endpoint of [a,b], or an endpoint is at in?nity, the problem is singular. The harmonic di?erential equation can form the basis of regular problems, while prob- lems involving Bessel’s and Legendre’s equations are singular. Regular Sturm–Liouville problems have additional properties. There are in?nitely many eigenvalues. There is a smallest eigenvalue but no largest eigenvalue, and the eigenvalues can be ordered as λ 0 <λ 1 <···<λ n ···. Associatedwitheachλ n isaunique(toanarbitrarymultiplicative constant) eigenfunction ψ n that has exactly n zeros in (a,b). If a problem is singular because p = 0 at an endpoint, we can also satisfy (A.80) by demanding that ψ be bounded at that endpoint (a singularity condition) and that any regular Sturm–Liouville boundary condition hold at the other endpoint. This is the case for Bessel’s and Legendre’s equations discussed below. Orthogonality of the eigenfunctions. Let L be self-adjoint, and let ψ m and ψ n be eigenfunctions associated with λ m and λ n , respectively. Then by (A.80) we have integraldisplay b a (ψ m (x) L [ψ n (x)] ?ψ n (x) L [ψ m (x)]) dx = 0. But L [ψ n (x)] =?λ n σ(x)ψ n (x) and L [ψ m (x)] =?λ m σ(x)ψ m (x). Hence (λ m ?λ n ) integraldisplay b a ψ m (x)ψ n (x)σ(x)dx = 0, and λ m negationslash=λ n implies that integraldisplay b a ψ m (x)ψ n (x)σ(x)dx = 0. (A.83) We say that ψ m and ψ n are orthogonal with respect to the weight function σ(x). Eigenfunction expansion of an arbitrary function. If L is self-adjoint, then its eigenfunctions form acompleteset. This means that any piecewise smooth function may be represented as a weighted series of eigenfunctions. Speci?cally, if f and f prime are piece- wise continuous on [a,b], then f may be represented as the generalizedFourierseries f(x)= ∞ summationdisplay n=0 c n ψ n (x). (A.84) Convergence of the series is uniform and gives, at any point of (a,b), the average value [ f(x + )+ f(x ? )]/2 of the one-sided limits f(x + )and f(x ? )of f(x). The c n can be found using orthogonality condition (A.83): multiply (A.84) by ψ m σ and integrate to obtain integraldisplay b a f(x)ψ m (x)σ(x)dx = ∞ summationdisplay n=0 c n integraldisplay b a ψ n (x)ψ m (x)σ(x)dx, hence c n = integraltext b a f(x)ψ n (x)σ(x)dx integraltext b a ψ 2 n (x)σ(x)dx . (A.85) These coe?cients ensure that the series converges in mean to f ; i.e., the mean-square error integraldisplay b a vextendsingle vextendsingle vextendsingle vextendsingle vextendsingle f(x)? ∞ summationdisplay n=0 c n ψ n (x) vextendsingle vextendsingle vextendsingle vextendsingle vextendsingle 2 σ(x)dx is minimized. Truncation to ?nitely-many terms generally results in oscillations (Gibb’s phenomena) near points of discontinuity of f . The c n are easier to compute if the ψ n are orthonormal with integraldisplay b a ψ 2 n (x)σ(x)dx = 1 for each n. Uniqueness of the eigenfunctions. If both ψ 1 and ψ 2 are associated with the same eigenvalue λ, then L [ψ 1 (x)] +λσ(x)ψ 1 (x)= 0, L [ψ 2 (x)] +λσ(x)ψ 2 (x)= 0, hence ψ 1 (x) L [ψ 2 (x)] ?ψ 2 (x) L [ψ 1 (x)] = 0. By (A.79) we have d dx bracketleftbigg p(x) parenleftbigg ψ 1 (x) dψ 2 (x) dx ?ψ 2 (x) dψ 1 (x) dx parenrightbiggbracketrightbigg = 0 or p(x) parenleftbigg ψ 1 (x) dψ 2 (x) dx ?ψ 2 (x) dψ 1 (x) dx parenrightbigg = C where C is constant. Either of (A.81) implies C = 0, hence d dx parenleftbigg ψ 2 (x) ψ 1 (x) parenrightbigg = 0 so that ψ 1 (x) = Kψ 2 (x) for some constant K. So under homogeneous boundary condi- tions, every eigenvalue is associated with a unique eigenfunction. Thisisfalsefortheperiodicboundaryconditions(A.82). Eigenfunctionexpansionthen becomes di?cult, as we can no longer assume eigenfunction orthogonality. However, the Gram–Schmidt algorithm may be used to construct orthogonal eigenfunctions. We refer the interested reader to Haberman [79]. The harmonic di?erential equation. The ordinary di?erential equation d 2 ψ(x) dx 2 =?k 2 ψ(x) (A.86) is Sturm–Liouville with p ≡ 1, q ≡ 0, σ ≡ 1, andλ= k 2 . Suppose we take [a,b] = [0,L] and adopt the homogeneous boundary conditions ψ(0)= 0 and ψ(L)= 0. (A.87) Since p(x)>0 andσ(x)>0 on [0,L], equations(A.86)and(A.87)formaregularSturm– Liouville problem. Thus we should have an in?nite number of discrete eigenvalues. A power series technique yields the two independent solutions ψ a (x)= A a sin kx,ψ b (x)= A b cos kx, to (A.86); hence by linearity the most general solution is ψ(x)= A a sin kx + A b cos kx. (A.88) The condition at x = 0 gives A a sin 0 + A b cos 0 = 0, hence A b = 0. The other condition then requires A a sin kL = 0. (A.89) Since A a = 0 would give ψ ≡ 0, we satisfy (A.89) by choosing k = k n = nπ/L for n = 1,2,....Because λ= k 2 , the eigenvalues are λ n =(nπ/L) 2 with corresponding eigenfunctions ψ n (x)= sin k n x. Note that λ = 0 is not an eigenvalue; eigenfunctions are nontrivial by de?nition, and sin(0πx/L)≡ 0. Likewise, the di?erential equation associated with λ= 0 can be solved easily, but only its trivial solution can ?t homogeneous boundary conditions: with k = 0, (A.86) becomes d 2 ψ(x)/dx 2 = 0, giving ψ(x)= ax +b; this can satisfy (A.87) only with a = b = 0. These “eigensolutions” obey the properties outlined earlier. In particular the ψ n are orthogonal, integraldisplay L 0 sin parenleftBig nπx L parenrightBig sin parenleftBig mπx L parenrightBig dx = L 2 δ mn , and the eigenfunction expansion of a piecewise continuous function f is given by f(x)= ∞ summationdisplay n=1 c n sin parenleftBig nπx L parenrightBig where, with σ(x)= 1 in (A.85), we have c n = integraltext L 0 f(x)sin parenleftbig nπx L parenrightbig dx integraltext L 0 sin 2 parenleftbig nπx L parenrightbig dx = 2 L integraldisplay L 0 f(x)sin parenleftBig nπx L parenrightBig dx. Hence we recover the standard Fourier sine series for f(x). With little extra e?ort we can examine the eigenfunctions resulting from enforcement of the periodic boundary conditions ψ(0)=ψ(L) and ψ prime (0)=ψ prime (L). The general solution (A.88) still holds, so we have the choices ψ(x)= sin kx and ψ(x)= cos kx. Evidently both ψ(x)= sin parenleftbigg 2nπx L parenrightbigg and ψ(x)= cos parenleftbigg 2nπx L parenrightbigg satisfy the boundary conditions for n = 1,2,.... Thus each eigenvalue (2nπ/L) 2 is associated with two eigenfunctions. Bessel’s di?erential equation. Bessel’s equation d dx parenleftbigg x dψ(x) dx parenrightbigg + parenleftbigg k 2 x ? ν 2 x parenrightbigg ψ(x)= 0 (A.90) occurs when problems are solved in circular-cylindrical coordinates. Comparison with (A.77) shows thatλ= k 2 , p(x)= x, q(x)=?ν 2 /x, andσ(x)= x. We take [a,b] = [0,L] along with the boundary conditions ψ(L)= 0 and |ψ(0)|<∞. (A.91) Although the resulting Sturm–Liouville problem is singular, the speci?ed conditions (A.91) maintain satisfaction of (A.80). The eigenfunctions are orthogonal because (A.80) is satis?ed by having ψ(L)= 0 and p(x)dψ(x)/dx → 0 as x → 0. As a second-order ordinary di?erential equation, (A.90) has two solutions denoted by J ν (kx) and N ν (kx), and termed Bessel functions. Their properties are summarized in Appendix E.1. The function J ν (x), the Bessel function of the ?rst kind and orderν, is well-behaved in [0,L]. The function N ν (x), the Bessel function of the second kind and order ν, is unbounded at x = 0; hence it is excluded as an eigenfunction of the Sturm–Liouville problem. The condition at x = L shows that the eigenvalues are de?ned by J ν (kL)= 0. We denote the mth root of J ν (x)= 0 by p νm . Then k νm = radicalbig λ νm = p νm /L. The in?nitely many eigenvalues are ordered as λ ν1 <λ ν2 <.... Associated with eigen- value λ νm is a single eigenfunction J ν ( √ λ νm x). The orthogonality relation is integraldisplay L 0 J ν parenleftBig p νm L x parenrightBig J ν parenleftBig p νn L x parenrightBig xdx= 0, m negationslash= n. Since the eigenfunctions are also complete, we can expand any piecewise continuous function f in a Fourier–Besselseries f(x)= ∞ summationdisplay m=1 c m J ν parenleftBig p νm x L parenrightBig , 0 ≤ x ≤ L,ν>?1. By (A.85) and (E.22) we have c m = 2 L 2 J 2 ν+1 (p νm ) integraldisplay L 0 f(x)J ν parenleftBig p νm x L parenrightBig xdx. The associated Legendre equation. Legendre’s equation occurs when problems are solved in spherical coordinates. It is often written in one of two forms. Letting θ be the polar angle of spherical coordinates (0 ≤θ ≤π), the equation is d dθ parenleftbigg sinθ dψ(θ) dθ parenrightbigg + parenleftbigg λsinθ ? m 2 sinθ parenrightbigg ψ(θ)= 0. This is Sturm–Liouville with p(θ) = sinθ, σ(θ) = sinθ, and q(θ) =?m 2 /sinθ. The boundary conditions |ψ(0)|<∞ and |ψ(π)|<∞ de?ne a singular problem: the conditions are not homogeneous, p(θ) = 0 at both end- points, and q(θ) < 0. Despite this, the Legendre problem does share properties of a regular Sturm–Liouville problem — including eigenfunction orthogonality and complete- ness. Using x = cosθ, we can put Legendre’s equation into its other common form d dx parenleftbigg [1 ? x 2 ] dψ(x) dx parenrightbigg + parenleftbigg λ? m 2 1 ? x 2 parenrightbigg ψ(x)= 0, (A.92) where ?1 ≤ x ≤ 1. It is found that ψ is bounded at x =±1 only if λ= n(n + 1) where n ≥ m is an integer. These λ are the eigenvalues of the Sturm–Liouville problem, and the corresponding ψ n (x) are the eigenfunctions. As a second-order partial di?erential equation, (A.92) has two solutions known as associated Legendre functions. The solution bounded at both x =±1 is the associated Legendre function of the ?rst kind, denoted P m n (x). The second solution, unbounded at x =±1, is the associated Legendre function of the second kind Q m n (x). Appendix E.2 tabulates some properties of these functions. For ?xed m, each λ mn is associated with a single eigenfunction P m n (x). Since P m n (x) is bounded at x =±1, and since p(±1) = 0, the eigenfunctions obey Lagrange’s identity (A.79), hence are orthogonal on [?1,1] with respect to the weight function σ(x) = 1. Evaluation of the orthogonality integral leads to integraldisplay 1 ?1 P m l (x)P m n (x)dx =δ ln 2 2n + 1 (n + m)! (n ? m)! (A.93) or equivalently integraldisplay π 0 P m l (cosθ)P m n (cosθ)sinθ dθ =δ ln 2 2n + 1 (n + m)! (n ? m)! . For m = 0, P m n (x)is a polynomial of degree n. Each suchLegendrepolynomial, denoted P n (x), is given by P n (x)= 1 2 n n! d n (x 2 ? 1) n dx n . It turns out that P m n (x)=(?1) m (1 ? x 2 ) m/2 d m P n (x) dx m , giving P m n (x)= 0 for m > n. Because the Legendre polynomials form a complete set in the interval [?1,1],wemay expand any su?ciently smooth function in a Fourier–Legendreseries f(x)= ∞ summationdisplay n=0 c n P n (x). Convergence in mean is guaranteed if c n = 2n + 1 2 integraldisplay 1 ?1 f(x)P n (x)dx, found using (A.85) along with (A.93). In practice, the associated Legendre functions appear along with exponential functions inthesolutionstosphericalboundaryvalueproblems. Thecombinedfunctionsareknown assphericalharmonics, andformsolutionstotwo-dimensionalSturm–Liouvilleproblems. We consider these next. Higher-dimensional SL problems: Helmholtz’s equation. Replacing d/dx by ?, we generalize the Sturm–Liouville equation to higher dimensions: ?·[p(r)?ψ(r)] + q(r)ψ(r)+λσ(r)ψ(r)= 0, where q, p, σ, ψ are real functions. Of particular interest is the case q(r) = 0, p(r) = σ(r)= 1, giving the Helmholtz equation ? 2 ψ(r)+λψ(r)= 0. (A.94) In most boundary value problems, ψ or its normal derivative is speci?ed on the surface of a bounded region. We obtain a three-dimensional analogue to the regular Sturm– Liouville problem by assuming the homogeneous boundary conditions αψ(r)+β?n ·?ψ(r)= 0 (A.95) on the closed surface, where ?n is the outward unit normal. The problem consisting of (A.94) and (A.95) has properties analogous to those of the regular one-dimensional Sturm–Liouville problem. All eigenvalues are real. There are in?nitely many eigenvalues. There is a smallest eigenvalue but no largest eigenvalue. However, associated with an eigenvalue there may be many eigenfunctions ψ λ (r). The eigenfunctions are orthogonal with integraldisplay V ψ λ 1 (r)ψ λ 2 (r)dV = 0,λ 1 negationslash=λ 2 . They are also complete and can be used to represent any piecewise smooth function f(r) according to f(r)= summationdisplay λ a λ ψ λ (r), which converges in mean when a λ m = integraltext V f(r)ψ λ m (r)dV integraltext V ψ 2 λ m (r)dV . Thesepropertiesaresharedbythetwo-dimensionaleigenvalueprobleminvolvinganopen surface S with boundary contour Gamma1. Spherical harmonics. We now inspect solutions to the two-dimensional eigenvalue problem ? 2 Y(θ,φ)+ λ a 2 Y(θ,φ)= 0 over the surface of a sphere of radius a. Since the sphere has no boundary contour, we demand that Y(θ,φ) be bounded in θ and periodic in φ. In the next section we shall apply separation of variables and show that Y nm (θ,φ)= radicalBigg 2n + 1 4π (n ? m)! (n + m)! P m n (cosθ)e jmφ whereλ= n(n+1). Note that Q m n does not appear as it is not bounded atθ = 0,π. The functions Y nm are called spherical harmonics (sometimes zonal or tesseral harmonics, depending on the values of n and m). As expressed above they are in orthonormal form, becausetheorthogonalityrelationshipsfortheexponentialandassociatedLegendre functions yield integraldisplay π ?π integraldisplay π 0 Y ? n prime m prime(θ,φ)Ynm(θ,φ)sinθ dθ dφ =δn prime n δ m prime m . (A.96) As solutions to the Sturm–Liouville problem, these functions form a complete set on the surface of a sphere. Hence they can be used to represent any piecewise smooth function f(θ,φ) as f(θ,φ)= ∞ summationdisplay n=0 n summationdisplay m=?n a nm Y nm (θ,φ), where a nm = integraldisplay π ?π integraldisplay π 0 f(θ,φ)Y ? nm (θ,φ)sinθ dθ dφ by (A.96). The summation index m ranges from ?n to n because P m n = 0 for m > n.For negative index we can use Y n,?m (θ,φ)=(?1) m Y ? nm (θ,φ). Some properties of the spherical harmonics are tabulated in Appendix E.3. Separation of variables We now consider a technique that ?nds widespread application in solving boundary value problems, applying as it does to many important partial di?erential equations such as Laplace’s equation, the di?usion equation, and the scalar and vector wave equations. These equations are related to the scalar Helmholtz equation ? 2 ψ(r)+ k 2 ψ(r)= 0 (A.97) where k is a complex constant. If k is real and we supply the appropriate boundary conditions, we have the higher-dimensional Sturm–Liouville problem with λ = k 2 .We shall not pursue the extension of Sturm–Liouville theory to complex values of k. Laplace’s equation is Helmholtz’s equation with k = 0. With λ = k 2 = 0 it might appear that Laplace’s equation does not involve eigenvalues; however, separation of vari- ables does lead us to lower-dimensional eigenvalue problems to which our previous meth- ods apply. Solutions to the scalar or vector wave equations usually begin with Fourier transformation on the time variable, or with an initial separation of the time variable to reach a Helmholtz form. The separation of variables idea is simple. We seek a solution to (A.97) in the form of a product of functions each of a single variable. If ψ depends on all three spatial dimensions, then we seek a solution of the type ψ(u,v,w)= U(u)V(v)W(w), where u, v, and w are the coordinate variables used to describe the problem. If ψ depends on only two coordinates, we may seek a product solution involving two functions each dependent on a single coordinate; alternatively, may use the three-variable solution and choose constants so that the result shows no variation with one coordinate. The Helmholtz equation is considered separable if it can be reduced to a set of independent ordinary di?erential equations, each involving a single coordinate variable. The ordinary di?erential equations, generally of second order, can be solved by conventional techniques resulting in solutions of the form U(u)= A u U A (u,k u ,k v ,k w )+ B u U B (u,k u ,k v ,k w ), V(v)= A v V A (v,k u ,k v ,k w )+ B v V B (v,k u ,k v ,k w ), W(w)= A w W A (w,k u ,k v ,k w )+ B w W B (w,k u ,k v ,k w ). The constants k u ,k v ,k w are calledseparationconstants and are found, along with the am- plitudeconstants A, B,byapplyingboundaryconditionsappropriateforagivenproblem. At least one separation constant depends on (or equals) k, so only two are independent. In many cases k u , k v , and k w become the discrete eigenvalues of the respective di?er- ential equations, and correspond to eigenfunctions U(u,k u ,k v ,k w ), V(v,k u ,k v ,k w ), and W(w,k u ,k v ,k w ). In other cases the separation constants form a continuous spectrum of values, often when a Fourier transform solution is employed. The Helmholtz equation can be separated in eleven di?erent orthogonal coordinate systems [134]. Undoubtedly the most important of these are the rectangular, circular- cylindrical, and spherical systems, and we shall consider each in detail. We do note, however, that separability in a certain coordinate system does not imply that all prob- lems expressed in that coordinate system can be easily handled using the resulting solu- tions. Only when the geometry and boundary conditions are simple do the solutions lend themselves to easy application; often other solution techniques are more appropriate. Although rigorous conditions can be set forth to guarantee solvability by separation of variables [119], we prefer the following, more heuristic list: 1. Use a coordinate system that allows the given partial di?erential equation to sep- arate into ordinary di?erential equations. 2. The problem’s boundaries must be such that those boundaries not at in?nity co- incide with a single level surface of the coordinate system. 3. Use superposition to reduce the problem to one involving a single nonhomogeneous boundary condition. Then: (a) Solve the resulting Sturm–Liouville problem in one or two dimensions, with homogeneous boundary conditions on all boundaries. Then use a discrete eigenvalue expansion (Fourier series) and eigenfunction orthogonality to sat- isfy the remaining nonhomogeneous condition. (b) If a Sturm–Liouville problem cannot be formulated with the homogeneous boundary conditions (because, for instance, one boundary is at in?nity), use a Fourier integral (continuous expansion) to satisfy the remaining nonhomo- geneous condition. If a Sturm–Liouville problem cannot be formulated, discovering the form of the integral transform to use can be di?cult. In these cases other approaches, such as conformal mapping, may prove easier. Solutions in rectangular coordinates. In rectangular coordinates the Helmholtz equation is ? 2 ψ(x,y,z) ?x 2 + ? 2 ψ(x,y,z) ?y 2 + ? 2 ψ(x,y,z) ?z 2 + k 2 ψ(x,y,z)= 0. (A.98) We seek a solution of the form ψ(x,y,z) = X(x)Y(y)Z(z); substitution into (A.98) followed by division through by X(x)Y(y)Z(z) gives 1 X(x) d 2 X(x) dx 2 + 1 Y(y) d 2 Y(y) dy 2 + 1 Z(z) d 2 Z(z) dz 2 =?k 2 . (A.99) At this point we require the separation argument. The left-hand side of (A.99) is a sum of three functions, each involving a single independent variable, whereas the right-hand side is constant. But the only functions of independent variables that always sum to a constant are themselves constants. Thus we may equate each term on the left to a di?erent constant: 1 X(x) d 2 X(x) dx 2 =?k 2 x , 1 Y(y) d 2 Y(y) dy 2 =?k 2 y , (A.100) 1 Z(z) d 2 Z(z) dz 2 =?k 2 z , provided that k 2 x + k 2 y + k 2 z = k 2 . The negative signs in (A.100) have been introduced for convenience. Let us discuss the general solutions of equations (A.100). If k x = 0, the two indepen- dent solutions for X(x) are X(x)= a x x and X(x)= b x where a x and b x are constants. If k x negationslash= 0, solutions may be chosen from the list of functions e ?jk x x , e jk x x , sin k x x, cos k x x, any two of which are independent. Because sin x =(e jx ? e ?jx )/2 j and cos x =(e jx + e ?jx )/2, (A.101) the six possible solutions for k x negationslash= 0 are X(x)= ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? A x e jk x x + B x e ?jk x x , A x sin k x x + B x cos k x x, A x sin k x x + B x e ?jk x x , A x e jk x x + B x sin k x x, A x e jk x x + B x cos k x x, A x e ?jk x x + B x cos k x x. (A.102) We may base our choice on convenience (e.g., the boundary conditions may be amenable to one particular form) or on the desired behavior of the solution (e.g., standing waves vs. traveling waves). If k is complex, then so may be k x , k y ,ork z ; observe that with imaginary arguments the complex exponentials are actually real exponentials, and the trigonometric functions are actually hyperbolic functions. The solutions for Y(y) and Z(z) are identical to those for X(x). We can write, for instance, X(x)= braceleftBigg A x e jk x x + B x e ?jk x x , k x negationslash= 0, a x x + b x , k x = 0, (A.103) Y(y)= braceleftBigg A y e jk y y + B y e ?jk y y , k y negationslash= 0, a y y + b y , k y = 0, (A.104) Z(z)= braceleftBigg A z e jk z z + B z e ?jk z z , k z negationslash= 0, a z z + b z , k z = 0. (A.105) Examples. Let us begin by solving the simple equation ? 2 V(x)= 0. Since V depends only on x we can use (A.103)–(A.105) with k y = k z = 0 and a y = a z = 0. Moreover k x = 0 because k 2 x +k 2 y +k 2 z = k 2 = 0 forLaplace’sequation. Thegeneralsolution is therefore V(x)= a x x + b x . Boundary conditions must be speci?ed to determine a x and b x ; for instance, the condi- tions V(0)= 0 and V(L)= V 0 yield V(x)= V 0 x/L. Next let us solve ? 2 ψ(x,y)= 0. We produce a lack of z-dependence inψ by letting k z = 0 and choosing a z = 0. Moreover, k 2 x =?k 2 y since Laplace’s equation requires k = 0. This leads to three possibilities. If k x = k y = 0, we have the product solution ψ(x,y)=(a x x + b x )(a y y + b y ). (A.106) If k y is real and nonzero, then ψ(x,y)=(A x e ?k y x + B x e k y x )(A y e jk y y + B y e ?jk y y ). (A.107) Using the relations sinh u =(e u ? e ?u )/2 and cosh u =(e u + e ?u )/2 (A.108) along with (A.101), we can rewrite (A.107) as ψ(x,y)=(A x sinh k y x + B x cosh k y x)(A y sin k y y + B y cos k y y). (A.109) (We can reuse the constant names A x , B x , A y , B y , since the constants are unknown at this point.) If k x is real and nonzero we have ψ(x,y)=(A x sin k x x + B x cos k x x)(A y sinh k x y + B y cosh k x y). (A.110) Consider the problem consisting of Laplace’s equation ? 2 V(x,y)= 0 (A.111) holding in the region 0 < x < L 1 , 0 < y < L 2 , ?∞< z <∞, together with the boundary conditions V(0,y)= V 1 , V(L 1 ,y)= V 2 , V(x,0)= V 3 , V(x,L 2 )= V 4 . The solution V(x,y) represents the potential within a conducting tube with each wall held at a di?erent potential. Superposition applies: since Laplace’s equation is linear we can write the solution as the sum of solutions to four di?erent sub-problems. Each sub-problem has homogeneous boundary conditions on one independent variable and inhomogeneous conditions on the other, giving a Sturm–Liouville problem in one of the variables. For instance, let us examine the solutions found above in relation to the sub- problem consisting of Laplace’s equation (A.111) in the region 0 < x < L 1 , 0 < y < L 2 , ?∞< z <∞, subject to the conditions V(0,y)= V(L 1 ,y)= V(x,0)= 0, V(x,L 2 )= V 4 negationslash= 0. First we try (A.106). The boundary condition at x = 0 gives V(0,y)=(a x (0)+ b x )(a y y + b y )= 0, which holds for all y ∈(0,L 2 ) only if b x = 0. The condition at x = L 1 , V(L 1 ,y)= a x L 1 (a y y + b y )= 0, then requires a x = 0. But a x = b x = 0 gives V(x,y) = 0, and the condition at y = L 2 cannot be satis?ed; clearly (A.106) was inappropriate. Next we examine (A.109). The condition at x = 0 gives V(0,y)=(A x sinh 0 + B x cosh 0)(A y sin k y y + B y cos k y y)= 0, hence B x = 0. The condition at x = L 1 implies V(L 1 ,y)= [A x sinh(k y L 1 )](A y sin k y y + B y cos k y y)= 0. This can hold if either A x = 0 or k y = 0, but the case k y = 0 (= k x ) was already considered. Thus A x = 0 and the trivial solution reappears. Our last candidate is (A.110). The condition at x = 0 requires V(0,y)=(A x sin 0 + B x cos 0)(A y sinh k x y + B y cosh k x y)= 0, which implies B x = 0. Next we have V(L 1 ,y)= [A x sin(k x L 1 )](A y sinh k y y + B y cosh k y y)= 0. We avoid A x = 0 by setting sin(k x L 1 ) = 0 so that k x n = nπ/L 1 for n = 1,2,....(Here n = 0 is omitted because it would produce a trivial solution.) These are eigenvalues corresponding to the eigenfunctions X n (x) = sin(k x n x), and were found in § A.4 for the harmonic equation. At this point we have a family of solutions V n (x,y)= sin(k x n x)[A y n sinh(k x n y)+ B y n cosh(k x n y)], n = 1,2,.... The subscript n on the left identi?es V n as the eigensolution associated with eigenvalue k x n . It remains to satisfy boundary conditions at y = 0,L 2 .Aty = 0 we have V n (x,0)= sin(k x n x)[A y n sinh 0 + B y n cosh 0] = 0, hence B y n = 0 and V n (x,y)= A y n sin(k x n x)sinh(k x n y), n = 1,2,.... (A.112) It is clear that no single eigensolution (A.112) can satisfy the one remaining boundary condition. However, we are guaranteed that a series of solutions can represent the con- stant potential on y = L 2 ; recall that as a solution to a regular Sturm–Liouville problem, the trigonometric functions are complete (hence they could represent any well-behaved function on the interval 0 ≤ x ≤ L 1 ). In fact, the resulting series is a Fourier sine series for the constant potential at y = L 2 . So let V(x,y)= ∞ summationdisplay n=1 V n (x,y)= ∞ summationdisplay n=1 A y n sin(k x n x)sinh(k x n y). The remaining boundary condition requires V(x,L 2 )= ∞ summationdisplay n=1 A y n sin(k x n x)sinh(k x n L 2 )= V 4 . The constants A y n can be found using orthogonality; multiplying through by sin(k x m x) and integrating, we have ∞ summationdisplay n=1 A y n sinh(k x n L 2 ) integraldisplay L 1 0 sin parenleftbigg mπx L 1 parenrightbigg sin parenleftbigg nπx L 1 parenrightbigg dx = V 4 integraldisplay L 1 0 sin parenleftbigg mπx L 1 parenrightbigg dx. The integral on the left equals δ mn L 1 /2 where δ mn is the Kronecker delta given by δ mn = braceleftBigg 1, m = n, 0, n negationslash= m. After evaluating the integral on the right we obtain ∞ summationdisplay n=1 A y n δ mn sinh(k x n L 2 )= 2V 4 (1 ? cos mπ) mπ , hence A y m = 2V 4 (1 ? cos mπ) mπ sinh(k x m L 2 ) . The ?nal solution for this sub-problem is therefore V(x,y)= ∞ summationdisplay n=1 2V 4 (1 ? cos nπ) nπ sinh parenleftBig nπL 2 L 1 parenrightBig sin parenleftbigg nπx L 1 parenrightbigg sinh parenleftbigg nπy L 1 parenrightbigg . The remaining three sub-problems are left for the reader. Let us again consider (A.111), this time for 0 ≤ x ≤ L 1 , 0 ≤ y <∞, ?∞< z <∞, and subject to V(0,y)= V(L 1 ,y)= 0, V(x,0)= V 0 . Let us try the solution form that worked in the previous example: V(x,y)= [A x sin(k x x)+ B x cos(k x x)][A y sinh(k x y)+ B y cosh(k x y)]. The boundary conditions at x = 0,L 1 are the same as before so we have V n (x,y)= sin(k x n x)[A y n sinh(k x n y)+ B y n cosh(k x n y)], n = 1,2,.... To ?nd A y n and B y n we note that V cannot grow without bound as y →∞. Individually the hyperbolic functions grow exponentially. However, using (A.108) we see that B y n = ?A y n gives V n (x,y)= A y n sin(k x n x)e ?k xn y where A y n is a new unknown constant. (Of course, we could have chosen this exponential dependence at the beginning.) Lastly, we can impose the boundary condition at y = 0 on the in?nite series of eigenfunctions V(x,y)= ∞ summationdisplay n=1 A y n sin(k x n x)e ?k xn y to ?nd A y n . The result is V(x,y)= ∞ summationdisplay n=1 2V 0 πn (1 ? cos nπ)sin(k x n x)e ?k xn y . As in the previous example, the solution is a discrete superposition of eigenfunctions. The problem consisting of (A.111) holding for 0 ≤ x ≤ L 1 , 0 ≤ y <∞, ?∞< z <∞, along with V(0,y)= 0, V(L 1 ,y)= V 0 e ?ay , V(x,0)= 0, requires a continuous superposition of eigenfunctions to satisfy the boundary conditions. Let us try V(x,y)= [A x sinh k y x + B x cosh k y x][A y sin k y y + B y cos k y y]. The conditions at x = 0 and y = 0 require that B x = B y = 0.Thus V k y (x,y)= A sinh k y x sin k y y. A single function of this form cannot satisfy the remaining condition at x = L 1 .Sowe form a continuous superposition V(x,y)= integraldisplay ∞ 0 A(k y )sinh k y x sin k y ydk y . (A.113) By the condition at x = L 1 integraldisplay ∞ 0 A(k y )sinh(k y L 1 )sin k y ydk y = V 0 e ?ay . (A.114) We can ?nd the amplitude function A(k y ) by using the orthogonality property δ(y ? y prime )= 2 π integraldisplay ∞ 0 sin xysin xy prime dx. (A.115) Multiplying both sides of (A.114) by sin k prime y y and integrating, we have integraldisplay ∞ 0 A(k y )sinh(k y L 1 ) bracketleftbiggintegraldisplay ∞ 0 sin k y y sin k prime y ydy bracketrightbigg dk y = integraldisplay ∞ 0 V 0 e ?ay sin k prime y ydy. We can evaluate the term in brackets using (A.115) to obtain integraldisplay ∞ 0 A(k y )sinh(k y L 1 ) π 2 δ(k y ? k prime y )dk y = integraldisplay ∞ 0 V 0 e ?ay sin k prime y ydy, hence π 2 A(k prime y )sinh(k prime y L 1 )= V 0 integraldisplay ∞ 0 e ?ay sin k prime y ydy. We then evaluate the integral on the right, solve for A(k y ), and substitute into (A.113) to obtain V(x,y)= 2V 0 π integraldisplay ∞ 0 k y a 2 + k 2 y sinh(k y x) sinh(k y L 1 ) sin k y ydk y . Note that our application of the orthogonality property is merely a calculation of the inverse Fourier sine transform. Thus we could have found the amplitude coe?cient by reference to a table of transforms. We can use the Fourier transform solution even when the domain is in?nite in more than one dimension. Suppose we solve (A.111) in the region 0 ≤ x <∞, 0 ≤ y <∞, ?∞< z <∞, subject to V(0,y)= V 0 e ?ay , V(x,0)= 0. Because of the condition at y = 0 let us use V(x,y)=(A x e ?k y x + B x e k y x )(A y sin k y y + B y cos k y y). The solution form V k y (x,y)= B(k y )e ?k y x sin k y y satis?es the ?niteness condition and the homogeneous condition at y = 0. The remaining condition can be satis?ed by a continuous superposition of solutions: V(x,y)= integraldisplay ∞ 0 B(k y )e ?k y x sin k y ydk y . We must have V 0 e ?ay = integraldisplay ∞ 0 B(k y )sin k y ydk y . Use of the orthogonality relationship (A.115) yields the amplitude spectrum B(k y ), and we ?nd that V(x,y)= 2 π integraldisplay ∞ 0 e ?k y x k y a 2 + k 2 y sin k y ydk y . (A.116) As a ?nal example in rectangular coordinates let us consider a problem in which ψ depends on all three variables: ? 2 ψ(x,y,z)+ k 2 ψ(x,y,z)= 0 for 0 ≤ x ≤ L 1 , 0 ≤ y ≤ L 2 , 0 ≤ z ≤ L 3 , subject to ψ(0,y,z)=ψ(L 1 ,y,z)= 0, ψ(x,0,z)=ψ(x,L 2 ,z)= 0, ψ(x,y,0)=ψ(x,y,L 3 )= 0. Here k negationslash= 0 is a constant. This is a three-dimensional eigenvalue problem as described in § A.4, where λ = k 2 are the eigenvalues and the closed surface is a rectangular box. Physically, the wave function ψ represents the so-calledeigenvalue ornormalmode solutions for the “TM modes” of a rectangular cavity. Since k 2 x + k 2 y + k 2 z = k 2 ,we might have one or two separation constants equal to zero, but not all three. We ?nd, however, that the only solution with a zero separation constant that can ?t the boundary conditions is the trivial solution. In light of the boundary conditions and because we expect standing waves in the box, we take ψ(x,y,z)= [A x sin(k x x)+ B x cos(k x x)] · · [A y sin(k y y)+ B y cos(k y y)] · · [A z sin(k z z)+ B z cos(k z z)]. The conditions ψ(0,y,z) = ψ(x,0,z) = ψ(x,y,0) = 0 give B x = B y = B z = 0. The conditions at x = L 1 , y = L 2 , and z = L 3 require the separation constants to assume the discrete values k x = k x m = mπ/L 1 , k y = k y n = nπ/L 2 , and k z = k z p = pπ/L 3 , where k 2 x m + k 2 y n + k 2 z p = k 2 mnp and m,n, p = 1,2,.... Associated with each of these eigenvalues is an eigenfunction of a one-dimensional Sturm–Liouville problem. For the three-dimensional problem, an eigenfunction ψ mnp (x,y,z)= A mnp sin(k x m x)sin(k y n y)sin(k z p z) is associated with each three-dimensional eigenvalue k mnp . Each choice of m,n, p pro- duces a discrete cavity resonance frequency at which the boundary conditions can be satis?ed. Depending on the values of L 1,2,3 , we may have more than one eigenfunction associated with an eigenvalue. For example, if L 1 = L 2 = L 3 = L then k 121 = k 211 = k 112 = √ 6π/L. However, the eigenfunctions associated with this single eigenvalue are all di?erent: ψ 121 = sin(k x 1 x)sin(k y 2 y)sin(k z 1 z), ψ 211 = sin(k x 2 x)sin(k y 1 y)sin(k z 1 z), ψ 112 = sin(k x 1 x)sin(k y 1 y)sin(k z 2 z). When more than one cavity mode corresponds to a given resonant frequency, we call the modes degenerate. By completeness, we can represent any well-behaved function as f(x,y,z)= summationdisplay m,n,p A mnp sin(k x m x)sin(k y n y)sin(k z p z). The A mnp are found using orthogonality. When such expansions are used to solve prob- lems involving objects (such as excitation probes) inside the cavity, they are termed normalmodeexpansions of the cavity ?eld. Solutionsincylindricalcoordinates. IncylindricalcoordinatestheHelmholtzequa- tion is 1 ρ ? ?ρ parenleftbigg ρ ?ψ(ρ,φ,z) ?ρ parenrightbigg + 1 ρ 2 ? 2 ψ(ρ,φ,z) ?φ 2 + ? 2 ψ(ρ,φ,z) ?z 2 + k 2 ψ(ρ,φ,z)= 0. (A.117) With ψ(ρ,φ,z)= P(ρ)Phi1(φ)Z(z) we obtain 1 ρ ? ?ρ parenleftbigg ρ ?(PPhi1Z) ?ρ parenrightbigg + 1 ρ 2 ? 2 (PPhi1Z) ?φ 2 + ? 2 (PPhi1Z) ?z 2 + k 2 (PPhi1Z)= 0; carrying out the ρ derivatives and dividing through by PPhi1Z we have ? 1 Z d 2 Z dz 2 = k 2 + 1 ρ 2 Phi1 d 2 Phi1 dφ 2 + 1 ρP dP dρ + 1 P d 2 P dρ 2 . The left side depends on z while the right side depends on ρ and φ, hence both must equal the same constant k 2 z : ? 1 Z d 2 Z dz 2 = k 2 z , (A.118) k 2 + 1 ρ 2 Phi1 d 2 Phi1 dφ 2 + 1 ρP dP dρ + 1 P d 2 P dρ 2 = k 2 z . (A.119) We have separated the z-dependence from the dependence on the other variables. For the harmonic equation (A.118), Z(z)= braceleftBigg A z sin k z z + B z cos k z z, k z negationslash= 0, a z z + b z , k z = 0. (A.120) Of course we could use exponentials or a combination of exponentials and trigonometric functions instead. Rearranging (A.119) and multiplying through by ρ 2 , we obtain ? 1 Phi1 d 2 Phi1 dφ 2 = parenleftbig k 2 ? k 2 z parenrightbig ρ 2 + ρ P dP dρ + ρ 2 P d 2 P dρ 2 . The left and right sides depend only on φ and ρ, respectively; both must equal some constant k 2 φ : ? 1 Phi1 d 2 Phi1 dφ 2 = k 2 φ , (A.121) parenleftbig k 2 ? k 2 z parenrightbig ρ 2 + ρ P dP dρ + ρ 2 P d 2 P dρ 2 = k 2 φ . (A.122) The variables ρ and φ are thus separated, and harmonic equation (A.121) has solutions Phi1(φ)= braceleftBigg A φ sin k φ φ+ B φ cos k φ φ, k φ negationslash= 0, a φ φ+ b φ , k φ = 0. (A.123) Equation (A.122) is a bit more involved. In rearranged form it is d 2 P dρ 2 + 1 ρ dP dρ + parenleftBigg k 2 c ? k 2 φ ρ 2 parenrightBigg P = 0 (A.124) where k 2 c = k 2 ? k 2 z . The solution depends on whether any of k z , k φ ,ork c are zero. If k c = k φ = 0, then d 2 P dρ 2 + 1 ρ dP dρ = 0 so that P(ρ)= a ρ lnρ+ b ρ . If k c = 0 but k φ negationslash= 0, we have d 2 P dρ 2 + 1 ρ dP dρ ? k 2 φ ρ 2 P = 0 so that P(ρ)= a ρ ρ ?k φ + b ρ ρ k φ . (A.125) This includes the case k = k z = 0 (Laplace’s equation). If k c negationslash= 0 then (A.124) is Bessel’s di?erential equation. For noninteger k φ the two independent solutions are denoted J k φ (z) and J ?k φ (z), where J ν (z) is the ordinary Bessel function of the ?rst kind of order ν.For k φ an integer n, J n (z) and J ?n (z) are not independent and a second independent solution denoted N n (z) must be introduced. This is the ordinary Bessel function of the second kind, order n. As it is also independent when the order is noninteger, J ν (z) and N ν (z) are often chosen as solutions whether ν is integer or not. Linear combinations of these independent solutions may be used to produce new independent solutions. The functions H (1) ν (z) and H (2) ν (z) are the Hankelfunctions of the ?rst and second kind of order ν, and are related to the Bessel functions by H (1) ν (z)= J ν (z)+ jN ν (z), H (2) ν (z)= J ν (z)? jN ν (z). The argument z can be complex (as can ν, but this shall not concern us). When z is imaginary we introduce two new functions I ν (z) and K ν (z), de?ned for integer order by I n (z)= j ?n J n (jz), K n (z)= π 2 j n+1 H (1) n (jz). Expressions for noninteger order are given in Appendix E.1. Bessel functions cannot be expressed in terms of simple, standard functions. However, a series solution to (A.124) produces many useful relationships between Bessel functions of di?ering order and argument. The recursion relations for Bessel functions serve to connect functions of various orders and their derivatives. See Appendix E.1. Of the six possible solutions to (A.124), R(ρ)= ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? A ρ J ν (k c ρ)+ B ρ N ν (k c ρ), A ρ J ν (k c ρ)+ B ρ H (1) ν (k c ρ), A ρ J ν (k c ρ)+ B ρ H (2) ν (k c ρ), A ρ N ν (k c ρ)+ B ρ H (1) ν (k c ρ), A ρ N ν (k c ρ)+ B ρ H (2) ν (k c ρ), A ρ H (1) ν (k c ρ)+ B ρ H (2) ν (k c ρ), which do we choose? Again, we are motivated by convenience and the physical nature of the problem. If the argument is real or imaginary, we often consider large or small argument behavior. For x real and large, J ν (x)→ radicalbigg 2 πx cos parenleftBig x ? π 4 ?ν π 2 parenrightBig , N ν (x)→ radicalbigg 2 πx sin parenleftBig x ? π 4 ?ν π 2 parenrightBig , H (1) ν (x)→ radicalbigg 2 πx e j(x? π 4 ?ν π 2 ) , H (2) ν (x)→ radicalbigg 2 πx e ?j(x? π 4 ?ν π 2 ) , I ν (x)→ radicalbigg 1 2πx e x , K ν (x)→ radicalbigg π 2x e ?x , while for x real and small, J 0 (x)→ 1, N 0 (x)→ 2 π (ln x + 0.5772157 ? ln 2), J ν (x)→ 1 ν! parenleftBig x 2 parenrightBig ν , N ν (x)→? (ν? 1)! π parenleftbigg 2 x parenrightbigg ν . Because J ν (x) and N ν (x) oscillate for large argument, they can represent standing waves along the radial direction. However, N ν (x)is unbounded for small x and is inappropriate for regions containing the origin. The Hankel functions become complex exponentials for large argument, hence represent traveling waves. Finally, K ν (x)is unbounded for small x and cannot be used for regions containing the origin, while I ν (x) increases exponentially for large x and cannot be used for unbounded regions. Examples. Consider the boundary value problem for Laplace’s equation ? 2 V(ρ,φ)= 0 (A.126) in the region 0 ≤ρ ≤∞, 0 ≤φ ≤φ 0 , ?∞< z <∞, where the boundary conditions are V(ρ,0)= 0, V(ρ,φ 0 )= V 0 . Since there is no z-dependence we let k z = 0 in (A.120) and choose a z = 0. Then k 2 c = k 2 ? k 2 z = 0 since k = 0. There are two possible solutions, depending on whether k φ is zero. First let us try k φ negationslash= 0. Using (A.123) and (A.125) we have V(ρ,φ)= [A φ sin(k φ φ)+ B φ cos(k φ φ)][a ρ ρ ?k φ + b ρ ρ k φ ]. (A.127) Assuming k φ > 0 we must have b ρ = 0 to keep the solution ?nite. The condition V(ρ,0)= 0 requires B φ = 0.Thus V(ρ,φ)= A φ sin(k φ φ)ρ ?k φ . Our ?nal boundary condition requires V(ρ,φ 0 )= V 0 = A φ sin(k φ φ 0 )ρ ?k φ . Because this cannot hold for all ρ, we must resort to k φ = 0 and V(ρ,φ)=(a φ φ+ b φ )(a ρ lnρ+ b ρ ). (A.128) Proper behavior as ρ →∞dictates that a ρ = 0. V(ρ,0) = 0 requires b φ = 0.Thus V(ρ,φ)= V(φ)= b φ φ.The constant b φ is found from the remaining boundary condition: V(φ 0 )= V 0 = b φ φ 0 so that b φ = V 0 /φ 0 . The ?nal solution is V(φ)= V 0 φ/φ 0 . It is worthwhile to specialize this toφ 0 =π/2 and compare with the solution to the same problem found earlier using rectangular coordinates. With a = 0 in (A.116) we have V(x,y)= 2 π integraldisplay ∞ 0 e ?k y x sin k y y k y dk y . Despite its much more complicated form, this must be the same solution by uniqueness. Next let us solve (A.126) subject to the “split cylinder” conditions V(a,φ)= braceleftBigg V 0 , 0 <φ<π, 0, ?π<φ<0. Because there is no z-dependence we choose k z = a z = 0 and have k 2 c = k 2 ?k 2 z = 0. Since k φ = 0 would violate the boundary conditions at ρ = a, we use V(ρ,φ)=(a ρ ρ ?k φ + b ρ ρ k φ )(A φ sin k φ φ+ B φ cos k φ φ). The potential must be single-valued inφ: V(ρ,φ+2nπ)= V(ρ,φ). This is only possible if k φ is an integer, say k φ = m. Then V m (ρ,φ)= braceleftBigg (A m sin mφ+ B m cos mφ)ρ m ,ρ<a, (C m sin mφ+ D m cos mφ)ρ ?m ,ρ>a. On physical grounds we have discarded ρ ?m for ρ<a and ρ m for ρ>a. To satisfy the boundary conditions at ρ = a we must use an in?nite series of the complete set of eigensolutions. For ρ<a the boundary condition requires B 0 + ∞ summationdisplay m=1 (A m sin mφ+ B m cos mφ)a m = braceleftBigg V 0 , 0 <φ<π, 0, ?π<φ<0. Application of the orthogonality relations integraldisplay π ?π cos mφcos nφdφ = 2π epsilon1 n δ mn , m,n = 0,1,2,..., (A.129) integraldisplay π ?π sin mφsin nφdφ =πδ mn , m,n = 1,2,..., (A.130) integraldisplay π ?π cos mφsin nφdφ = 0, m,n = 0,1,2,..., (A.131) where epsilon1 n = braceleftBigg 1, n = 0, 2, n > 0, (A.132) is Neumann’s number, produces appropriate values for the constants A m and B m . The full solution is V(ρ,φ)= ? ? ? ? ? ? ? ? ? V 0 2 + V 0 π ∞ summationdisplay n=1 [1 ?(?1) n ] n parenleftBig ρ a parenrightBig n sin nφ, ρ<a, V 0 2 + V 0 π ∞ summationdisplay n=1 [1 ?(?1) n ] n parenleftbigg a ρ parenrightbigg n sin nφ, ρ>a. The boundary value problem ? 2 V(ρ,φ,z)= 0, 0 ≤ρ ≤ a, ?π ≤φ ≤π, 0 ≤ z ≤ h, V(ρ,φ,0)= 0, 0 ≤ρ ≤ a, ?π ≤φ ≤π, V(a,φ,z)= 0, ?π ≤φ ≤π, 0 ≤ z ≤ h, V(ρ,φ,h)= V 0 , 0 ≤ρ ≤ a, ?π ≤φ ≤π, describes the potential within a grounded “canister” with top at potential V 0 . Symmetry precludes φ-dependence, hence k φ = a φ = 0. Since k = 0 (Laplace’s equation) we also have k 2 c = k 2 ? k 2 z =?k 2 z . Thus we have either k z real and k c = jk z ,ork c real and k z = jk c . With k z real we have V(ρ,z)= [A z sin k z z + B z cos k z z][A ρ K 0 (k z ρ)+ B ρ I 0 (k z ρ)]; (A.133) with k c real we have V(ρ,z)= [A z sinh k c z + B z cosh k c z][A ρ J 0 (k c ρ)+ B ρ N 0 (k c ρ)]. (A.134) The functions K 0 and I 0 are inappropriate for use in this problem, and we proceed to (A.134). Since N 0 is unbounded for small argument, we need B ρ = 0. The condition V(ρ,0)= 0 gives B z = 0,thus V(ρ,z)= A z sinh(k c z)J 0 (k c ρ). The oscillatory nature of J 0 means that we can satisfy the condition at ρ = a: V(a,z)= A z sinh(k c z)J 0 (k c a)= 0 for 0 ≤ z < h if J 0 (k c a) = 0. Letting p 0m denote the mth root of J 0 (x) = 0 for m = 1,2,..., we have k c m = p 0m /a. Because we cannot satisfy the boundary condition at z = h with a single eigensolution, we use the superposition V(ρ,z)= ∞ summationdisplay m=1 A m sinh parenleftBig p 0m z a parenrightBig J 0 parenleftBig p 0m ρ a parenrightBig . We require V(ρ,h)= ∞ summationdisplay m=1 A m sinh parenleftbigg p 0m h a parenrightbigg J 0 parenleftBig p 0m ρ a parenrightBig = V 0 , (A.135) where the A m can be evaluated by orthogonality of the functions J 0 (p 0m ρ/a).Ifp νm is the mth root of J ν (x)= 0, then integraldisplay a 0 J ν parenleftBig p νm ρ a parenrightBig J ν parenleftBig p νn ρ a parenrightBig ρdρ =δ mn a 2 2 J prime2 ν (p νn )=δ mn a 2 2 J 2 ν+1 (p νn ) (A.136) where J prime ν (x)= dJ ν (x)/dx. Multiplying (A.135) by ρJ 0 (p 0n ρ/a) and integrating, we have A n sinh parenleftbigg p 0n h a parenrightbigg a 2 2 J prime2 0 (p 0n a)= integraldisplay a 0 V 0 J 0 parenleftBig p 0n ρ a parenrightBig ρdρ. Use of (E.105), integraldisplay x n+1 J n (x)dx = x n+1 J n+1 (x)+ C, allows us to evaluate integraldisplay a 0 J 0 parenleftBig p 0n ρ a parenrightBig ρdρ = a 2 p 0n J 1 (p 0n ). With this we ?nish calculating A m and have V(ρ,z)= 2V 0 ∞ summationdisplay m=1 sinh( p 0m a z)J 0 ( p 0m a ρ) p 0m sinh( p 0m a h)J 1 (p 0m ) as the desired solution. Finally, let us assume that k > 0 and solve ? 2 ψ(ρ,φ,z)+ k 2 ψ(ρ,φ,z)= 0 where 0 ≤ρ ≤ a, ?π ≤φ ≤π, and ?∞< z <∞, subject to the condition ?n ·?ψ(ρ,φ,z) vextendsingle vextendsingle vextendsingle ρ=a = ?ψ(ρ,φ,z) ?ρ vextendsingle vextendsingle vextendsingle ρ=a = 0 for ?π ≤φ ≤π and ?∞< z <∞.The solution to this problem leads to the transverse- electric (TE z ) ?elds in a lossless circular waveguide, whereψ represents the z-component ofthemagnetic?eld. Althoughthereissymmetrywithrespecttoφ, weseekφ-dependent solutions; the resulting complete eigenmode solution will permit us to expand any well- behaved function within the waveguide in terms of a normal mode (eigenfunction) series. In this problem none of the constants k, k z ,ork φ equal zero, except as a special case. However, the ?eld must be single-valued in φ and thus k φ must be an integer m.We consider our possible choices for P(ρ), Z(z), and Phi1(φ). Since k 2 c = k 2 ? k 2 z and k 2 > 0 is arbitrary, we must consider various possibilities for the signs of k 2 c and k 2 z . We can rule out k 2 c < 0 based on consideration of the behavior of the functions I m and K m . We also need not consider k c < 0, since this gives the same solution as k c > 0. We are then left with two possible cases. Writing k 2 z = k 2 ? k 2 c , we see that either k > k c and k 2 z > 0,or k < k c and k 2 z < 0.Fork 2 z > 0 we write ψ(ρ,φ,z)= [A z e ?jk z z + B z e jk z z ][A φ sin mφ+ B φ cos mφ]J m (k c ρ). Here the terms involving e ?jk z z represent waves propagating in the ±z directions. The boundary condition at ρ = a requires J prime m (k c a)= 0 where J prime m (x) = dJ m (x)/dx. Denoting the nth zero of J prime m (x) by p prime mn we have k c = k cm = p prime mn /a. This gives the eigensolutions ψ m = [A zm e ?jk z z + B zm e jk z z ][A φm sin mφ+ B φm cos mφ]k c J m parenleftbigg p prime mn ρ a parenrightbigg . The undetermined constants A zm , B zm , A ρm , B ρm could be evaluated when the individual eigensolutions are used to represent a function in terms of a modal expansion. For the case k 2 z < 0 we again choose complex exponentials in z; however, k z =?jα gives e ?jk z z = e ?αz and attenuation along z. The reader can verify that the eigensolutions are ψ m = [A zm e ?αz + B zm e αz ][A φm sin mφ+ B φm cos mφ]k c J m parenleftbigg p prime mn ρ a parenrightbigg where now k 2 c = k 2 +α 2 . We have used Bessel function completeness in the examples above. This property is a consequence of the Sturm–Liouville problem ?rst studied in § A.4. We often use Fourier– Bessel series to express functions over ?nite intervals. Over in?nite intervals we use the Fourier–Bessel transform. The Fourier–Bessel series can be generalized to Bessel functions of noninteger order, and to the derivatives of Bessel functions. Let f(ρ) be well-behaved over the interval [0,a]. Then the series f(ρ)= ∞ summationdisplay m=1 C m J ν parenleftBig p νm ρ a parenrightBig , 0 ≤ρ ≤ a,ν>?1 converges, and the constants are C m = 2 a 2 J 2 ν+1 (p νm ) integraldisplay a 0 f(ρ)J ν parenleftBig p νm ρ a parenrightBig ρdρ by (A.136). Here p νm is the mth root of J ν (x). An alternative form of the series uses p prime νm , the roots of J prime ν (x), and is given by f(ρ)= ∞ summationdisplay m=1 D m J ν parenleftBig p prime νm ρ a parenrightBig , 0 ≤ρ ≤ a,ν>?1. In this case the expansion coe?cients are found using the orthogonality relationship integraldisplay a 0 J ν parenleftbigg p prime νm a ρ parenrightbigg J ν parenleftbigg p prime νn a ρ parenrightbigg ρdρ =δ mn a 2 2 parenleftbigg 1 ? ν 2 p prime2 νm parenrightbigg J 2 ν (p prime νm ), and are D m = 2 a 2 parenleftBig 1 ? ν 2 p prime2 νm J 2 ν (p prime νm ) parenrightBig integraldisplay a 0 f(ρ)J ν parenleftbigg p prime νm a ρ parenrightbigg ρdρ. Solutions in spherical coordinates. If into Helmholtz’s equation 1 r 2 ? ?r parenleftbigg r 2 ?ψ(r,θ,φ) ?r parenrightbigg + 1 r 2 sinθ ? ?θ parenleftbigg sinθ ?ψ(r,θ,φ) ?θ parenrightbigg + + 1 r 2 sin 2 θ ? 2 ψ(r,θ,φ) ?φ 2 + k 2 ψ(r,θ,φ)= 0 we putψ(r,θ,φ)= R(r)Theta1(θ)Phi1(φ)and multiply through by r 2 sin 2 θ/ψ(r,θ,φ), we obtain sin 2 θ R(r) d dr parenleftbigg r 2 dR(r) dr parenrightbigg + sinθ Theta1(θ) d dθ parenleftbigg sinθ dTheta1(θ) dθ parenrightbigg + k 2 r 2 sin 2 θ =? 1 Phi1(φ) d 2 Phi1(φ) dφ 2 . Since the right side depends only on φ while the left side depends only on r and θ, both sides must equal some constant μ 2 : sin 2 θ R(r) d dr parenleftbigg r 2 dR(r) dr parenrightbigg + sinθ Theta1(θ) d dθ parenleftbigg sinθ dTheta1(θ) dθ parenrightbigg + k 2 r 2 sin 2 θ =μ 2 , (A.137) d 2 Phi1(φ) dφ 2 +μ 2 Phi1(φ)= 0. (A.138) We have thus separated o? the φ-dependence. Harmonic ordinary di?erential equation (A.138) has solutions Phi1(φ)= braceleftBigg A φ sinμφ+ B φ cosμφ, μnegationslash= 0, a φ φ+ b φ ,μ= 0. (We could have used complex exponentials to describe Phi1(φ), or some combination of exponentials and trigonometric functions, but it is conventional to use only trigonometric functions.) Rearranging (A.137) and dividing through by sin 2 θ we have 1 R(r) d dr parenleftbigg r 2 dR(r) dr parenrightbigg + k 2 r 2 =? 1 sinθTheta1(θ) d dθ parenleftbigg sinθ dTheta1(θ) dθ parenrightbigg + μ 2 sin 2 θ . We introduce a new constant k 2 θ to separate r from θ: 1 R(r) d dr parenleftbigg r 2 dR(r) dr parenrightbigg + k 2 r 2 = k 2 θ , (A.139) ? 1 sinθTheta1(θ) d dθ parenleftbigg sinθ dTheta1(θ) dθ parenrightbigg + μ 2 sin 2 θ = k 2 θ . (A.140) Equation (A.140), 1 sinθ d dθ parenleftbigg sinθ dTheta1(θ) dθ parenrightbigg + parenleftbigg k 2 θ ? μ 2 sin 2 θ parenrightbigg Theta1(θ)= 0, can be put into a standard form by letting η= cosθ (A.141) and k 2 θ =ν(ν+ 1) where ν is a parameter: (1 ?η 2 ) d 2 Theta1(η) dη 2 ? 2η dTheta1(η) dη + bracketleftbigg ν(ν+ 1)? μ 2 1 ?η 2 bracketrightbigg Theta1(η)= 0, ?1 ≤η≤ 1. This is the associated Legendre equation. It has two independent solutions called as- sociated Legendre functions of the ?rst and second kinds, denoted P μ ν (η) and Q μ ν (η), respectively. In these functions, all three quantities μ,ν,ηmay be arbitrary complex constants as long as ν+μnegationslash= ?1,?2,.... But (A.141) shows that η is real in our discus- sion; μ will generally be real also, and will be an integer whenever Phi1(φ)is single-valued. The choice of ν is somewhat more complicated. The function P μ ν (η) diverges at η =±1 unless ν is an integer, while Q μ ν (η) diverges at η=±1 regardless of whether ν is an inte- ger. In § A.4 we required that P μ ν (η) be bounded on [?1,1] to have a Sturm–Liouville problem with suitable orthogonality properties. By (A.141) we must exclude Q μ ν (η) for problems containing the z-axis, and restrict ν to be an integer n in P μ ν (η) for such prob- lems. In case the z-axis is excluded, we chooseν = n whenever possible, because the ?nite sums P m n (η) and Q m n (η) are much easier to manipulate than P μ ν (η) and Q μ ν (η). In many problems we must count on completeness of the Legendre polynomials P n (η)= P 0 n (η) or spherical harmonics Y mn (θ,φ) in order to satisfy the boundary conditions. In this book we shall consider only those boundary value problems that can be solved using integer values of ν and μ, hence choose Theta1(θ)= A θ P m n (cosθ)+ B θ Q m n (cosθ). (A.142) Single-valuedness inPhi1(φ)is a consequence of havingμ= m, andφ = constant boundary surfaces are thereby disallowed. The associated Legendre functions have many important properties. For instance, P m n (η)= ? ? ? 0, m > n, (?1) m (1 ?η 2 ) m/2 2 n n! d n+m (η 2 ? 1) n dη n+m , m ≤ n. (A.143) The case m = 0 receives particular attention because it corresponds to azimuthal invari- ance (φ-independence). We de?ne P 0 n (η)= P n (η)where P n (η)is the Legendre polynomial of order n. From (A.143), we see that 4 P n (η)= 1 2 n n! d n (η 2 ? 1) n dη n is a polynomial of degree n, and that P m n (η)=(?1) m (1 ?η 2 ) m/2 d m dη m P n (η). BoththeassociatedLegendrefunctionsandtheLegendrepolynomialsobeyorthogonality relations and many recursion formulas. In problems where the z-axis is included, the product Theta1(θ)Phi1(φ) is sometimes de?ned as the spherical harmonic Y nm (θ,φ)= radicalBigg 2n + 1 4π (n ? m)! (n + m)! P m n (cosθ)e jmθ . These functions, which are complete over the surface of a sphere, were treated earlier in this section. Remembering that k 2 r =ν(ν+ 1), the r-dependent equation (A.139) becomes 1 r 2 d dr parenleftbigg r 2 dR(r) dr parenrightbigg + parenleftbigg k 2 + n(n + 1) r 2 parenrightbigg R(r)= 0. (A.144) When k = 0 we have d 2 R(r) dr 2 + 2 r dR(r) dr ? n(n + 1) r 2 R(r)= 0 so that R(r)= A r r n + B r r ?(n+1) . When k negationslash= 0, the substitution ˉ R(r)= √ kr R(r) puts (A.144) into the form r 2 d 2 ˉ R(r) dr 2 + r d ˉ R(r) dr + bracketleftBigg k 2 r 2 ? parenleftbigg n + 1 2 parenrightbigg 2 bracketrightBigg ˉ R(r)= 0, which we recognize as Bessel’s equation of half-integer order. Thus R(r)= ˉ R(r) √ kr = Z n+ 1 2 (kr) √ kr . For convenience we de?ne the sphericalBesselfunctions j n (z)= radicalbigg π 2z J n+ 1 2 (z), n n (z)= radicalbigg π 2z N n+ 1 2 (z)=(?1) n+1 radicalbigg π 2z J ?(n+ 1 2 ) (z), h (1) n (z)= radicalbigg π 2z H (1) n+ 1 2 (z)= j n (z)+ jn n (z), h (2) n (z)= radicalbigg π 2z H (2) n+ 1 2 (z)= j n (z)? jn n (z). 4 Care must be taken when consulting tables of Legendre functions and their properties. In particular, one must be on the lookout for possible disparities regarding the factor (?1) m (cf., [76, 1, 109, 8] vs. [5, 187]). Similar care is needed with Q m n (x). These can be written as ?nite sums involving trigonometric functions and inverse powers of z. We have, for instance, j 0 (z)= sin z z , n 0 (z)=? cos z z , j 1 (z)= sin z z 2 ? cos z z , n 1 (z)=? cos z z 2 ? sin z z . We can now write R(r) as a linear combination of any two of the spherical Bessel functions j n ,n n ,h (1) n ,h (2) n : R(r)= ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? A r j n (kr)+ B r n n (kr), A r j n (kr)+ B r h (1) n (kr), A r j n (kr)+ B r h (2) n (kr), A r n n (kr)+ B r h (1) n (kr), A r n n (kr)+ B r h (2) n (kr), A r h (1) n (kr)+ B r h (2) n (kr). (A.145) Imaginary arguments produce modi?ed spherical Bessel functions; the interested reader is referred to Gradsteyn [76] or Abramowitz [1]. Examples. The problem ? 2 V(r,θ,φ)= 0,θ 0 ≤θ ≤π/2, 0 ≤ r <∞, ?π ≤φ ≤π, V(r,θ 0 ,φ)= V 0 , ?π ≤φ ≤π, 0 ≤ r <∞, V(r,π/2,φ)= 0, ?π ≤φ ≤π, 0 ≤ r <∞, gives the potential ?eld between a cone and the z = 0 plane. Azimuthal symmetry prompts us to choose μ= a φ = 0. Since k = 0 we have R(r)= A r r n + B r r ?(n+1) . (A.146) Noting that positive and negative powers of r are unbounded for large and small r, respectively, we take n = B r = 0. Hence the solution depends only on θ: V(r,θ,φ)= V(θ)= A θ P 0 0 (cosθ)+ B θ Q 0 0 (cosθ). We must retain Q 0 0 since the solution region does not contain the z-axis. Using P 0 0 (cosθ)= 1 and Q 0 0 (cosθ)= ln cot(θ/2) (cf., Appendix E.2), we have V(θ)= A θ + B θ ln cot(θ/2). A straightforward application of the boundary conditions gives A θ = 0 and B θ = V 0 /ln cot(θ 0 /2), hence V(θ)= V 0 ln cot(θ/2) ln cot(θ 0 /2) . Next we solve the boundary value problem ? 2 V(r,θ,φ)= 0, V(a,θ,φ)=?V 0 ,π/2 ≤θ<π,?π ≤φ ≤π, V(a,θ,φ)=+V 0 , 0 <θ≤π/2, ?π ≤φ ≤π, for both r > a and r < a. This yields the potential ?eld of a conducting sphere split into top and bottom hemispheres and held at a potential di?erence of 2V 0 . Azimuthal symmetry gives μ= 0. The two possible solutions for Theta1(θ) are Theta1(θ)= braceleftBigg A θ + B θ ln cot(θ/2), n = 0, A θ P n (cosθ), n negationslash= 0, where we have discarded Q 0 0 (cosθ)because the region of interest contains the z-axis. The n = 0 solution cannot match the boundary conditions; neither can a single term of the type A θ P n (cosθ), but a series of these latter terms can. We use V(r,θ)= ∞ summationdisplay n=0 V n (r,θ)= ∞ summationdisplay n=0 [A r r n + B r r ?(n+1) ]P n (cosθ). (A.147) The terms r ?(n+1) and r n are not allowed, respectively, for r < a and r > a.Forr < a then, V(r,θ)= ∞ summationdisplay n=0 A n r n P n (cosθ). Letting V 0 (θ)be the potential on the surface of the split sphere, we impose the boundary condition: V(a,θ)= V 0 (θ)= ∞ summationdisplay n=0 A n a n P n (cosθ), 0 ≤θ ≤π. This is a Fourier–Legendre expansion of V 0 (θ). The A n are evaluated by orthogonality. Multiplying by P m (cosθ)sinθ and integrating from θ = 0 to π, we obtain ∞ summationdisplay n=0 A n a n integraldisplay π 0 P n (cosθ)P m (cosθ)sinθ dθ = integraldisplay π 0 V 0 (θ)P m (cosθ)sinθ dθ. Using orthogonality relationship (A.93) and the given V 0 (θ) we have A m a m 2 2m + 1 = V 0 integraldisplay π/2 0 P m (cosθ)sinθ dθ ? V 0 integraldisplay π π/2 P m (cosθ)sinθ dθ. The substitution η= cosθ gives A m a m 2 2m + 1 = V 0 integraldisplay 1 0 P m (η)dη? V 0 integraldisplay 0 ?1 P m (η)dη = V 0 integraldisplay 1 0 P m (η)dη? V 0 integraldisplay 1 0 P m (?η)dη; then P m (?η)=(?1) m P m (η) gives A m = a ?m 2m + 1 2 V 0 [1 ?(?1) m ] integraldisplay 1 0 P m (η)dη. Because A m = 0 for m even, we can put m = 2n + 1 (n = 0,1,2,...) and have A 2n+1 = (4n + 3)V 0 a 2n+1 integraldisplay 1 0 P 2n+1 (η)dη= V 0 (?1) n a 2n+1 4n + 3 2n + 2 (2n!) (2 n n!) 2 by (E.176). Hence V(r,θ)= ∞ summationdisplay n=0 V 0 (?1) n 4n + 3 2n + 2 (2n!) (2 n n!) 2 parenleftBig r a parenrightBig 2n+1 P 2n+1 (cosθ) for r < a. The case r > a is left to the reader. Finally, consider ? 2 ψ(x,y,z)+ k 2 ψ(x,y,z)= 0, 0 ≤ r ≤ a, 0 ≤θ ≤π,?π ≤φ ≤π, ψ(a,θ,φ)= 0, 0 ≤θ ≤π, ?π ≤φ ≤π, where k negationslash= 0 is constant. This is a three-dimensional eigenvalue problem. Wave function ψ representsthesolutionsfortheelectromagnetic?eldwithinasphericalcavityformodes TE to r. Despite the prevailing symmetry, we choose solutions that vary with both θ and φ. We are motivated by a desire to solve problems involving cavity excitation, and eigenmode completeness will enable us to represent any piecewise continuous function within the cavity. We employ spherical harmonics because the boundary surface is a sphere. These exclude Q n m (cosθ), which is appropriate since our problem contains the z-axis. Since k negationslash= 0 we must choose a radial dependence from (A.145). Small-argument behavior rules out n n , h (1) n , and h (2) n , leaving us with ψ(r,θ,φ)= A mn j n (kr)Y nm (θ,φ) or, equivalently, ψ(r,θ,φ)= A mn j n (kr)P m n (cosθ)e jmφ . The eigenvalues λ= k 2 are found by applying the condition at r = a: ψ(a,θ,φ)= A mn j n (ka)Y nm (θ,φ)= 0, requiring j n (ka) = 0. Denoting the qth root of j n (x) = 0 by α nq , we have k nq = α nq /a and corresponding eigenfunctions ψ mnq (r,θ,φ)= A mnq j n (k nq r)Y nm (θ,φ). The eigenvalues are proportional to the resonant frequencies of the cavity and the eigen- functions can be used to ?nd the modal ?eld distributions. Since the eigenvalues are independent of m, we may have several eigenfunctions ψ mnq associated with each k mnq . The only limitation is that we must keep m ≤ n to have P n m (cosθ) nonzero. This is another instance of mode degeneracy. There are 2n degenerate modes associated with each resonant frequency (one for each of e ±jnφ ). By completeness we can expand any piecewise continuous function within or on the sphere as a series f(r,θ,φ)= summationdisplay m,n,q A mnq j n (k nq r)Y nm (θ,φ).