{'fraction_non_alphanumeric': 0.04924627762878017, 'fraction_numerical': 0.01780264496439471, 'mean_word_length': 4.009728978457262, 'pattern_counts': {'":': 0, '<': 14, '<?xml version=': 0, '>': 12, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Radiation from a semi-infinite unflanged planar dielectric waveguide

Radiative emission from an open-ended unflanged planar dielectric waveguide is studied for the case of TM polarization on the basis of an iterative scheme. The first step of the scheme leads to approximate values for the reflection coefficients and electromagnetic fields inside and outside the waveguide. It is shown that for the related problems of reflection from a step potential in one-dimensional quantum mechanics and of Fresnel reflection of an electromagnetic plane wave from a half-space the iterative scheme is in accordance with the exact solution.


INTRODUCTION
In a classic paper Levine and Schwinger [1] studied the radiation of sound from an unflanged circular pipe. Later they extended their theory to electromagnetic radiation [2]. Their work constituted the first major advance in the theory of diffraction after Sommerfeld's exact solution of the problem of plane wave diffraction by an ideally conducting half-plane [3]. In the theory of diffraction of sound, microwaves, or light, the radiation is assumed to propagate in uniform space with reflection by rigid objects or idealized boundaries. For a lucid introduction to the theory of diffraction we refer to Sommerfeld's lecture notes [4]. The early theory of diffraction was reviewed by Bouwkamp [5]. Later developments are discussed by Born and Wolf [6]. A brief review of the principles and applications of open-ended waveguides with idealized walls was presented by Gardiol [7].
The invention of the dielectric waveguide by Hondros and Debye [8] has led to the development of optical fibers and the subsequent advances in telecommunication. In the theoretical determination of running wave solutions to Maxwell's equations in a spatially inhomogeneous medium the radiation is assumed to propagate in a guiding structure of infinite length. The problem of emergence of radiation from a semi-infinite waveguide into a half-space is of obvious technical interest. In the case of sound the analysis is based on the exact solution of Levine and Schwinger [1], [9][10][11][12]. For a dielectric waveguide the observation of the emerging radiation can be used as a tool to study the nature of the driving incident wave [13].
It is advantageous to simplify the theoretical analysis by the use of planar symmetry. The theory of Levine and Schwinger was extended to planar geometry by Heins [14,15]. In the following we study radiation emerging from a semi-infinite planar dielectric waveguide. As an intermediate step the wavefunction in the exit plane must be calculated. In this case the integral equation technique of Schwinger [16] cannot be implemented, because of the complicated nature of the integral kernel. We evaluate the emitted radiation and the coefficients of reflection back into the waveguide approximately in a first step of an iterative scheme.
We show in two appendices that the iterative scheme converges to the exact solution in the related problems of reflection by a step potential in one-dimensional quantum mechanics and of Fresnel reflection of an electromagnetic plane wave by a half-space. For the planar dielectric waveguide it does not seem practically possible to go beyond the first step of the iterative scheme.
In a numerical example we study a planar waveguide consisting of a slab of uniform dielectric constant, bounded on both sides by a medium with a smaller dielectric constant [17][18][19]. In the case studied the first step of the iterative scheme leads to a modification of the wavefunction at the exit plane which is relatively small in comparison with the incident wave. Hence we may expect that the calculation provides a reasonable approximation to the exact solution.

II. PLANAR OPEN END GEOMETRY
We employ Cartesian coordinates (x, y, z) and consider a planar waveguide in the halfspace z < 0 with stratified dielectric constant ε(x) and uniform magnetic permeability µ 1 . In the half-space z > 0 the dielectric constant is uniform with value ε ′ and the magnetic permeability is µ 1 . We consider solutions of Maxwell's equations which do not depend on the coordinate y and depend on time t through a factor exp(−iωt). Waves traveling to the right in the left half-space will be partly reflected at the plane z = 0, and partly transmitted into the right half-space. The solutions of Maxwell's equations may be decomposed according to two polarizations . For TE-polarization the components E x , E z , and H y vanish, and the equations may be combined into the single equation We have used gaussian units, and k = ω/c is the vacuum wavenumber. For TM-polarization the components E y , H x , and H z vanish, and the equations may be combined into the single equation We assume that the profile ε(x) is symmetric, ε(−x) = ε(x), and has a simple form with ε(x) increasing monotonically for x < 0 from value ε 1 to a maximum value ε 2 at x = 0. An example of the geometry under consideration is shown in Fig. 1. In the example the dielectric constant in the left half-space equals a constant ε 2 for −d < x < d and a constant ε 1 < ε 2 for x < −d and x > d.
For definiteness we consider only TM-polarization. It is convenient to denote the magnetic field component H y (x, z) for z < 0 as u(x, z) and for z > 0 as v(x, z). The continuity conditions at z = 0 are We consider a solution u 0n (x, z) of Eq. (2.2) given by a guided mode solution u 0n (x, z) = ψ n (x) exp(ip n z), (2.4) where ψ n (x) is the guided mode wavefunction, and p n the guided mode wavenumber. We assume p n > 0, so that the wave u 0n (x, z) exp(−iωt) is traveling to the right. The complete solution takes the form where u 1n (x, z) and v n (x, z) must be determined such that the continuity conditions Eq. (2.3) are satisfied. The function u 1n (x, z) describes the reflected wave, and v n (x, z) describes the wave radiated into the right-hand half-space. Since the right-hand half-space is uniform the solution v n (x, z) takes a simple form, and can be expressed as (2.6) The contribution from the interval − √ ε ′ µ 1 |k| < q < √ ε ′ µ 1 |k| corresponds to waves traveling to the right, the contribution from |q| > √ ε ′ µ 1 |k| corresponds to evanescent waves.
Similarly the solution u 1n (x, z) in the left half-space can be expressed as where the sum corresponds to guided waves traveling to the left, with n m the number of such guided modes possible at the given frequency ω, and the integral corresponds to waves radiating towards the left. We require that the mode solutions are normalized such that [20] The guided mode solutions {ψ m (x)} can be taken to be real. Orthogonality follows from Eq. (2.2). We show in the next section how the functions u 1n (x, z) and v n (x, z) may in principle be evaluated from an iterative scheme. The coefficients {R mn } and the amplitude function R n (q) also follow from the scheme.

III. ITERATIVE SCHEME
The iterative scheme is based on successive approximations to the scattering solution. Thus we write the exact solution as infinite sums n (x, z). In zeroth approximation we identify u (0) n (x, z) with the incident wave, The corresponding v (0) n (x, z) will be determined by continuity at the exit plane z = 0. From Eq. (3.2) we have u (0) n (x, 0−) = ψ n (x). This has the Fourier transform Using continuity of the wavefunction at z = 0 and the expression (2.6) we find correspondingly so that in zeroth approximation F (0) n (q) = φ n (q). Clearly the zeroth approximation does not satisfy the second continuity equation in Eq. (2.3), and we must take care of this in the next approximation.
For the difference of terms in Eq. (2.3) we find By symmetry ρ n (x) is symmetric in x for n even, antisymmetric in x for n odd.
The next approximation u (1) n (x, z) can be found by comparison with the solution of the problem where the profile ε(x) extends over all space and radiation is generated by a source ε(x)ρ(x)δ(z) with a Sommerfeld radiation condition, so that radiation is emitted to the right for z > 0 and to the left for z < 0. This antenna solution can be expressed as with kernel K(x, x ′ , z). The latter can be calculated from the Fourier decomposition in terms of the integral with the prescription that the path of integration in the complex p plane runs just above the negative real axis and just below the positive real axis. The Green function G(x, x ′ , p) can be found from the solution of the one-dimensional wave equation, The solution takes the form [20] G is the smaller (larger) of x and x ′ , and the remaining quantities will be specified in the next section. The Green function satisfies the symmetry properties 11) and the reciprocity relation Consequently the kernel K(x, x ′ , z) has the properties 14) The function u (1) n (x, z) is now identified as The minus sign is needed to provide near cancellation of the source density between the zeroth and first order solutions, ρ n (x) ≈ 0. We find the first order function F (1) n (q) by Fourier transform from the value at z = 0 in the form The corresponding function v (1) n (x, z) is found from Eq (2.6). The first order source density ρ (1) n (x) is found to be In principle the first order function u (1) n (x, 0) in the exit plane z = 0 may be regarded as the result of a linear operator R (1) acting on the state ψ n (x) given by the incident wave. The iterated solution then corresponds to the action with the operator R = R (1) where I is the identity operator. In order j of the geometric series corresponding to the operator R the wavefunctions u n (x, z) in the left and right half-space can be found by completing the function u 1n (x, 0) in the exit plane by left and right running waves respectively.
Assuming that the scheme has been extended to all orders we obtain the solutions n (x, 0). In the limit we must have ∞ j=0 ρ (j) n (x) = 0, (3.18) so that the continuity conditions Eq. (2.3) are exactly satisfied. In Appendix A we show how the iterative scheme reproduces the exact solution for reflection from a step potential in one-dimensional quantum mechanics. In Appendix B we show the same for Fresnel reflection.
In the integral in Eq. (3.15) it is convenient to perform the integral over p first, since ρ The second term in Eq. (2.7) corresponds to the remainder of the integral, after subtraction of the simple pole contributions. The function R n (q) will be discussed in the next section. In the calculation of F where δF (1) n (q) is the contribution from the remainder of the integral over p, after subtraction of the simple pole contributions.
Formally, in the complete solution Eq. (2.7) the reflection coefficients R mn and the amplitude function R n (q) are found as R mn = (ψ m , u 1n (0)), R n (q) = (ψ(q), u 1n (0)). (3.21) with the scalar product as given by Eq. (2.8). The first step of the iterative scheme yields The continuum states ψ(q) can be discretized in the usual way, so that the expressions in Eq. (3.22) can be regarded as elements of a matrix R (1) . As indicated above, the iterative scheme corresponds to a geometric series, so that the reflection coefficients in Eq. (3.21) can be found as elements of the matrix where I is the identity matrix. Finally the function F n (q) can be found from the corresponding state u n (0) as in Eq. (3.16). The function F n (q) can be related to the radiation scattered into the right halfspace. The scattering angle θ is related to the component q by (3.24) Defining the scattering cross section σ n (θ) by σ n (θ) sin θdθ = |F n (q)| 2 qdq, (3.25) we find the relation σ n (θ) = ε ′ µ 1 k 2 ε ′ µ 1 k 2 − q 2 |F n (q)| 2 . (3.26) In lowest approximation the cross section is proportional to the absolute square of the Fourier transform of the guided mode ψ n (x). To higher order the cross section is affected by the reflection into other modes.

IV. CONTINUOUS SPECTRUM
The calculation of the function R n (q) corresponding to the contribution from the continuous spectrum requires a separate discussion. The wave equation (3.9) is related to a quantummechanical Schrödinger equation for a particle in a potential. The bound states of the Schrödinger problem correspond to the guided modes, and the scattering states correspond to a continuous spectrum of radiation modes. The eigenstates of the Hamilton operator of the Schrödinger problem satisfy a completeness relation which can be usefully employed in the waveguide problem.
Explicitly the homogeneous one-dimensional Schrödinger equation corresponding to Eq. (3.9) via the relation ψ(x) = ε(x)f (x) reads [20] where the function V (x) is given by By comparison with the quantummechanical Schrödinger equation we see that where k 1 = √ ε 1 µ 1 k, may be identified as the potential. The bound state energies correspond to {k 2 1 − p 2 n }. It is convenient to assume that the dielectric profile ε(x) equals ε 1 for x < −x 1 and x > x 1 , so that the potential U(x) vanishes for |x| > x 1 . We define three independent solutions of the Schrödinger equation (4.1) with specified behavior for |x| > x 1 . The behavior of the function f 1 (x, p) is specified as with wavenumber The behavior of the function f 2 (x, p) is specified as Similarly, the behavior of the function f 3 (x, p) is specified as The coefficients W 12 and W 22 are elements of the transfer matrix of the planar structure, Because of the assumed symmetry of the dielectric profile we have in the present case R ′ = R and T ′ = T , so that W 21 = −W 12 . Moreover The functions f 2 (x, p) and f 3 (x, p) were used in the calculation of the Green function in Eq. (3.10). The denominator in that expression is given by From the solution of the inhomogeneous Schrödinger equation it follows that the completeness relation of the normal mode solutions may be expressed as [21] nm−1 Correspondingly the Green function may be decomposed as Hence we find for the antenna kernel K(x, x ′ , z) from Eq. (3.8) by use of the integration prescription (4.13) The first order left-hand wavefunction is therefore found from Eq. (3.15) as (4.14) with scalar product as given by Eq. (2.8). The wavefunction is the sum of guided modes running to the left and of radiation into the left-hand half-space. The first term agrees with the reflection coefficient given by Eq. (3.19). By symmetry the matrix element (ψ m , ερ (0) n ) vanishes unless m and n are both even or both odd. The integral provides an alternative expression for the remainder δu (1) n (x, z). From Eq. (4.11) we may identify With this definition the function R (1) n (q) is given by The second line has the same structure as Eq. (3.19). Although the decomposition in Eq. (4.14) is of theoretical interest, the calculation of the function u (1) n (x, z) is performed more conveniently as indicated in Eq. (3.15), with the integral over p performed first. V.

NUMERICAL EXAMPLE
We demonstrate the effectiveness of the scheme on a numerical example. We consider a flat dielectric profile defined by ε(x) = ε 2 for −d < x < d and ε(x) = ε 1 for |x| > d with values ε 2 = 2.25 and ε 1 = 2.13. In the right half-space we put ε ′ = 1, and we put µ 1 = 1 everywhere. The geometry is shown in Fig. 1.
By symmetry the guided modes in infinite space are either symmetric or antisymmetric in x. The explicit expressions for the mode wavefunctions can be found by the transfer matrix method [20]. At each of the two discontinuities the coefficients of the plane waves exp(iq i x) and exp(−iq i x) are transformed into coefficients of the plane waves exp(iq j x) and exp(−iq j x) by a matrix involving Fresnel coefficients given by The wavenumbers p n (k) of the guided modes are found as zeros of the transfer matrix element W 22 (p, k), which takes the explicit form The guided mode wavefunctions {ψ n (x)}, their Fourier transforms {φ n (q)}, and the Green function G(x, x ′ , p) can be found in explicit form. In Fig. 2 we show the ratio of wavenumbers p n (k)/k as a function of kd for the first few guided modes. We choose the frequency corresponding to kd = 12. In that case there are two symmetric modes, denoted as TM0 and TM2, and one antisymmetric mode, denoted as TM1. We assume the incident wave to be symmetric in x. Then it is not necessary to consider the antisymmetric mode. In Fig. 3 we show the corresponding normalized wavefunctions ψ 0 (x) and ψ 2 (x). In Fig. 4 we show their Fourier transforms φ 0 (q) and φ 2 (q). The wavenumbers at kd = 12 are p 0 = 17.955/d and p 2 = 17.624/d. The edge of the continuum is given by k 1 d = 17.513, and the corresponding value for ε 2 is k 2 d = 18.
In Fig. 5 we show the source density −iρ n (x) of the zeroth approximation for n = 0, 2 as a function of x, as given by Eq. (3.5). The coefficients of the simple pole contributions can be calculated from Eq. (3.19). We find numerically for the discrete part R For the first correction to the emitted radiation we need to calculate the function u (1) n (x, 0). The kernel K(x, x ′ , 0) in Eq. (3.15) can be evaluated numerically. On account of the symmetry in ±p it is sufficient to calculate twice the integral along the positive real p axis, with path of integration just below the axis. In the numerical integration over p in Eq. (3.8) the simple poles at {p m } cause problems. In order to avoid the simple poles we therefore integrate instead along a contour consisting of the line from 0 to k 1 just below the axis, a semi-circle in the lower half of the complex p plane centered at (k 1 + k 2 )/2 and of radius (k 2 − k 1 )/2, and the line just below the real axis from k 2 to +∞. In Fig. 6 we plot as an example the real part of K(x, 0, 0) as a function of x. The plot of the imaginary part is similar.
In Fig. 7 we show the imaginary part of the function u (1) 0 (x, 0), as calculated from Eq. (3.15). This is nearly identical with the contribution from the simple poles at p 0 and p 2 , which is also shown in Fig. 7. In Fig. 8 we show the real part of the function u (1) 0 (x, 0). Here the simple poles do not contribute. The magnitude of the wavefunction at the origin u (1) 0 (0, 0) = −0.149 − 0.005i may be compared with that of the zeroth approximation u (0) 0 (0, 0) = 1.346. This shows that the first order correction is an order of magnitude smaller than the zeroth order approximation. Consequently we may expect that the sum u (0) 0 (x, 0) + u (1) 0 (x, 0) provides a close approximation to the exact value. In Fig. 9 we show the absolute value |F

VI. DISCUSSION
In the above we have employed an iterative scheme inspired by the exact solution of two fundamental scattering problems, reflection by a step potential in one-dimensional quantum mechanics, shown in Appendix A, and Fresnel reflection of electromagnetic radiation by a half-space, shown in Appendix B. For the planar dielectric waveguide we have implemented only the first step of the iterative scheme. In the numerical example shown in Sec. V even this first step leads to interesting results. The method is sufficiently successful that it encourages application in other situations.
In particular it will be of interest to apply the method to a circular cylindrical dielectric waveguide or optical fiber. The mathematics of the method carries over straightforwardly to this more complicated geometry, with the plane wave behavior in the transverse direction replaced by Bessel functions.
Due to symmetry the problem for both planar and cylindrical geometry can be reduced to an equation for a scalar wavefunction, so that the theory is similar to that for sound propagation. This suggests that an interesting comparison can be made with a lattice Boltzmann simulation. For a rigid circular pipe such a simulation has already been performed by da Silva and Scavone [22], with interesting results. A finite element method has been applied to a rigid open-ended duct of more general cross section [23].

Appendix A
In this Appendix we show how the iterative scheme reproduces the exact solution of the time-independent one-dimensional Schrödinger equation with a step potential. We consider the equation with potential V (z) = 0 for z < 0 and V (z) = V for z > 0. In proper units p 2 is the energy. We denote the solution for z > 0 as v(z). For a wave incident from the left the exact solution reads where p ′ = p 2 − V , and the reflection coefficient B and transmission coefficient C are given by The wavefunction and its derivative are continuous at z = 0. We apply the iterative scheme and put to zeroth order The antenna solution u A (z) solves the equation for all z. It is given by To zeroth order the source ρ is We put the first order solution equal to Note the minus sign in −K(z)ρ (0) . The value at the exit z = 0 is sufficient to calculate the coefficients B and C from the geometric series By continuation one finds for the wave function at order j for j ≥ 1  Geometry of the planar waveguide. Plot of the reduced wavenumber p n (k)/k of the lowest order guided waves for n = 0, 1, 2, as functions of kd for values of the dielectric constant given in the text. Plot of the wavefunctions ψ 0 (x) and ψ 2 (x) of the guided modes with n = 0 (no nodes) and n = 2 (two nodes) as functions of x/d. Plot of the Fourier transform φ 0 (q) and φ 2 (q) of the wavefunctions of the guided modes with n = 0 (solid curve) and n = 2 (dashed curve) as functions of qd.    Plot of the real part of the kernel K(x, 0, 0), as given by Eq. (3.8), as a function of x/d. Plot of the real part of the first order wavefunction u  Plot of the imaginary part of the first order wavefunction u (1) 0 (x, 0) at the exit plane as a function of x/d. Plot of the absolute value of the Fourier transform |F