Modified Airy Function Method Applied to Optical Waveguides and Quantum Tunneling: A Critical Analysis
 Author: Lee KiYoung, Kim Chul Han, Park Dongwook, Kim ChangMin
 Publish: Current Optics and Photonics Volume 19, Issue2, p188~198, 25 Apr 2015

ABSTRACT
The validity of applying the Modified Airy Function (MAF) method to the problems of gradedindex optical waveguides and graded potential barrier analysis was critically examined. In the former case, the method yielded very accurate results from the derived eigenvalue equations. In the latter case, however, the same method produced results that deviated significantly from exact numerical results for barriers with a smooth peak. The causes of the discrepancies were investigated in detail.

KEYWORD
MAF method , Graded profile , Waveguide analysis , Tunneling probability , WKB method

I. INTRODUCTION
Analysis of optical devices such as waveguides and directional couplers is essential to be able to predict their characteristics, and involves solving the Helmholtz equation. Meanwhile, analysis of tunneling phenomena is necessary in various problems of particle physics, as well as in certain quantum devices, which requires solving the Schrӧdinger equation. Despite their apparent differences, the two problems are quite closely related from a theoretical standpoint, because the two equations are nearly identical in form. Whatever the structure, one must decide whether to employ a numerical method, or to rely on a mathematical approach for analysis. As to numerical methods for solving the wave equation, one may choose one of the wellknown simulation techniques such as FEM (finite element method), BPM (beam propagation method), the TM (transfer matrix) method, and so on, depending on the given application. Numerical methods provide very accurate results, but not necessarily physical insight. On the other hand, mathematical approaches provide such insight, allowing closedform solutions, but are limited in their use since exact solutions are available only for specific cases such as linear, exponential, and inverse hyperbolic cosine index profiles.
Thus, for a majority of graded structures, many have resorted to approximate mathematical methods. One of the most widely used methods for handling this class of problems is the WKB (WentzelKramersBrillouin) method. It is well known that results of this method tend to deviate from the exact solution on the whole since trial solutions diverge at the turning points [1]. Despite this deficiency the WKB method has frequently been employed by many researchers [2, 3] owing to the simple mathematical expressions involved, availability of eigenvalue equations in closed forms, and easy estimation of final results. In an effort to circumvent such inherent deficiency of the WKB method, Langer [4] proposed a method utilizing the Modified Airy functions (MAF), which was later applied to the analysis of optical waveguides and also to tunneling problems in connection with potential barriers by Ghatak et al. [5, 6]. More recently, Kim et al. applied the MAF method to a variety of optical waveguides [7, 8] and tunneling problems [9, 10] by employing MAFs in constructing the trial solution across the entire structure, demonstrating the method’s accuracy. Despite these recent advances, however, researchers have failed to clearly point out that the method does not yield satisfactory results in cases where particle energy is comparable to the barrier peak, for some classes of potential barrier profiles, and have not provided a rational explanation for such behavior. This work aims to address these shortcomings.
This paper is comprised as follows. First, the basics of the MAF method are addressed, with an optical waveguide as the platform for discussion. A trial solution is presented for each distinct region, followed by the corresponding ‘connection formula’. Then actual application of the MAF technique to the analysis of gradedindex waveguides is presented. Using the MAF trial solutions, a closedform expression of the eigenvalue equation is derived for gradedwaveguides. The eigenvalue equation is then solved for several gradedindex profiles, to yield the dispersion curves. The results are shown to be in very good agreement with exact numerical results obtained using the FEM.
The second part comprises the analysis of the tunneling problem. Based on the same MAF trial solutions postulated in the waveguide case, tunneling probability is derived in matrix form and evaluated for a number of graded potential barriers. The results are compared with exact numerical solutions obtained using the TM method, whose accuracy has already been verified in numerous studies [11, 12]. It turns out that MAF method fails for profiles with a smooth peak and that the WKB method does not work as well for profiles with an abrupt truncation or change.
Possible causes for the discrepancy in accuracy of the results for waveguide and tunneling problems upon application of the same MAF method are investigated in a systematic manner. This evaluation has not been previously addressed, to our knowledge. In the process, the corresponding results of the conventional WKB method and the FEM or TM method – the latter two providing exact numerical solutions – are also presented for comparison.
II. MAF TRIAL SOLUTIONS AND CONNECTION FORMULAE
In this section, a brief overview of the MAF method is presented. Although the discussion is conducted for optical waveguides, the formulation and basic concepts apply equally to the potentialbarrier problem.
Let us start with the onedimensional wave equation for optical waveguides:
where
In (2),
n (x ) denotes the given index profile across the waveguide structure andN is the eigenvalue to be found, the mode index.In Fig. 1 a diagram of Γ(
x ) depicting the distinct regions according to the sign of Γ(x ) is presented for clarification of the quantities involved. In theQ (x ) region where Γ(x ) is positive, the solution is oscillatory, and in theP (x ) region the solution is either monotonically growing or decaying.Let us at this point define the following auxiliary functions for subsequent discussions.
where
The MAF trial solutions to (1) are given by a combination of the modified Airy functions as follows:
where
Ai [−ξ (x )] andBi [−ξ (x )] are Airy functions of the first and second kind, respectively. Asymptotic expressions for these functions are shown in Appendix A and will be used in relating the MAF connection formula to the WKB connection formula.Note that
Ai (−ξ ) andAi (η ) are continuous at the turning pointx_{t} , as areBi (−ξ ) andBi (η ). Thus the boundary condition at each turning point is automatically satisfied under the present MAF trial solution construction.The connection feature of the trial solutions presented in (4) can be seen more clearly through the behavior of a typical wave function, shown in Fig. 2 below.
The MAF trial solutions of (4) take on the following asymptotic forms away from the turning point
x_{t} (see Appendix A).The expressions on the righthand side of (6) and (7) constitute the wellknown WKB connection formula, by which (6a) and (7a) are correspondingly transformed into (6b) and (7b), respectively.
III. APPLICATIONS OF THE MAF METHOD TO OPTICAL WAVEGUIDES
3.1. GradedIndex Waveguides
3.1.1. Derivation of Eigenvalue Equation
Consider a planar waveguide with an arbitrary gradedindex profile, as shown in Fig. 3.
The onedimensional Helmholtz equation is
where the mode index
N and the corresponding field distributionE (x ) are the quantities to be solved for. Trial solutions for the planar waveguide of Fig. 3 may be expressed by a combination of the modified Airy functions and the WKB trial solutions depending on the region, as follows [7].where
Away from both turning points, the expressions in (10) can be replaced by a linear combination of the sinusoidal functions of (6) and (7), as shown below.
One way to obtain the eigenvalue equation is simply to equate (13a) and (13b). Simple algebraic manipulations lead to the following result.
Another way is to impose boundary conditions requiring that the field and its derivative be continuous at
x =x _{0}.Both approaches lead to the same eigenvalue equation.
In the expression above,
δ _{1} andδ _{2} correspond to the phase shifts at the turning pointsx _{t1} andx _{t2} respectively, and are related to the amplitude coefficients byExpressions for
δ _{1} andδ _{2} can be derived by matching the fields and their derivatives at both sides ofx =a andx =b respectively, and can be cast in the following forms:where
For a symmetric profile,
δ _{1} =δ _{2} as expected.3.1.2. Simulation
To validate the eigenvalue equation of (15), the following parabolic index profile (shown in the inset of Fig. 4) will be assumed.
To facilitate computation, the normalized frequency
υ_{w} and propagation constantb_{w} , defined below, are utilized.(15) can then be cast in the following form, below which the corresponding WKB eigenvalue equation is also presented, for reference.
The normalized dispersion curves
υ_{w} −b_{w} for the first several modes are illustrated in Fig. 4, which also includes the exact numerical solution obtained via FEM for comparison. The figure clearly shows excellent agreement between the MAF result and the exact result, thus validating the use of the MAF method, whereas the WKB curves show some deviations near the cutoff points as expected.3.2. Truncated GradedIndex Waveguides
3.2.1. Derivation of Eigenvalue Equation
Consider a truncated gradedindex profile, as shown in the inset of Fig. 5. Notice that the left turning point now coincides with the location of the index discontinuity. MAF trial solutions appropriate for this structure can be written as follows.
where the definitions of
η (x ),ξ (x ), andP (x ) are identical to those given by the expressions forη _{2}(x ),ξ _{2}(x ), andP (x ) in (12), respectively, withx _{2} now replaced byx_{t} .By applying a procedure similar to that of the previous section, the eigenvalue equation can be obtained as follows.
where
While the
δ term in (26) represents the phase shift at the turning pointx =d , the inversetangent term represents the phase shift atx = 0, the point of index discontinuity.where
3.2.2. Simulation
The truncated parabolic index profile of Fig. 5 is described as
For this truncated profile, the eigenvalue equation of (26) can be cast in the following form. The WKB eigenvalue equation is also presented, for reference.
Calculation results based on the eigenvalue equation of (31) are displayed in Fig. 5 along with the corresponding WKB and FEM results. As was the case for the symmetric parabolic index profile (see Fig. 4), the MAF results show excellent agreement with the exact FEM results, while the WKB curves show minor deviations near the cutoff regions.
Based on this and the previous result, it may be concluded that the MAF technique provides a rather robust and highly accurate method for analysis of a wide range of gradedindex optical waveguides.
IV. APPLICATIONS OF THE MAF METHOD TO TUNNELING PROBLEMS
The MAFbased approach will now be applied to a set of barrier profiles, to test its feasibility in tunneling analysis. The onedimensional Schrӧdinger equation takes the form
Hereafter, we will work in atomic units where
m =ħ = 1, withx andE in the units of Bohr and Hartree, respectively.4.1. Graded Potential Barriers
4.1.1. Derivation of Tunneling Probability
Consider a potential barrier as shown in Fig. 6, with particle energy
E . Trial solutions similar to those employed in the optical waveguide analysis can also be used in this case [7].where
η (x ) is defined whereV (x ) is higher thanE , andξ (x ) whereV (x ) is lower thanE . Note that trial functions on both sides of each turning point satisfy the boundary condition automatically as before, via ‘the MAF connection formula.’Suppose the expressions of (35) are replaced by the following asymptotic forms (see Appendix A) much like in the earlier waveguide section but with different forms here.
Either by equating these two equations or by imposing the boundary conditions at a point between the two turning points, the following relation may be obtained.
At this point, it must be noted that the computation results based on (40) were nearly identical to those obtained from the WKB analysis. This is in clear contrast to the case of the optical waveguide, where use of asymptotic forms of the MAF trial solutions led to very accurate results that were distinctly different from the WKB results.
We shall now revert to the trial solutions of (34) to (36) as they are, and impose the boundary conditions at
x =a ,x_{p} , andb , respectively. Applying the continuity of the field and its derivative at the boundary pointx =a , we obtain the following matrix relation:where
Next, boundary condition matching at
x =x_{p} yieldswhere
Finally, application of the boundary conditions at
x =b yields the following relation:where
Equations (41), (43), and (45) can be combined to yield the following overall matrix equation relating
c _{1},c _{2} andc _{7}.It follows from (47) that the tunneling probability is given by the expression below.
Detailed descriptions of the terms ∆_{1}, ∆_{2},
A _{1}^{*},B _{1}^{*},A_{T} , andB_{T} are available in [9].Derivation of the tunneling probability using the WKB method for a graded potential barrier is presented in Appendix B1.
4.1.2. Simulation
Case 1: Symmetric exponential potential barrier
To verify the validity of the MAF formulation for graded potential barriers, a symmetric exponential profile (shown in the inset of Fig. 7) was selected as the first test case. The calculated results based on (48) are displayed in the figure, along with the results obtained by the exact transfer matrix (TM) method. The WKB results are also included, for reference.
Figure 7 clearly shows that the MAF method leads to solutions with very high accuracy over the entire region of
E /V _{0}. Meanwhile, the WKB method is seen to produce results with noticeable deviations, which is expected, given the WKB method’s inherent deficiencies.Case 2: Parabolic potential barrier
A parabolic potential barrier is taken as the second example, and is shown in the inset of Fig. 8.
Calculated results for this profile are illustrated in Fig. 8. The WKB results show behavior similar to that in Fig. 7. However, the MAF results show tunneling probabilities much lower than the true values across the entire range, contrary to expectation; part of this tendency can in fact be observed in the previous reports of [6, 11] upon careful examination. Even more disturbing, the calculation results fall to zero as
E /V _{0} approaches unity; the causes of this phenomenon are as follows. According to our investigation described in Appendix C, whenE approachesV _{0}, the denominators of the terms in (35) containingη '(x_{p} ) go to zero, causing the trial functions to become singular. As the trial functions in (34a) representing the incident and reflected components are finite, the coefficients of the trial functions in region II must go to zero, to satisfy the boundary conditions atx = −b . This in turn leads the coefficient in the transmitted region, and thus the tunneling probability, to vanish. It turns out that wheneverV '(x_{p} ) = 0,η '(x_{p} ) vanishes in the limit ofE →V _{0}, leading to the failure of the MAF method.4.2. Truncated Potential Barriers
4.2.1. Derivation of Tunneling Probability
Consider a truncated graded potential barrier, as shown in Fig. 9.
For the truncated graded barrier in Fig. 9, the trial solutions may be taken to be
where
By applying the boundary conditions at
x =a ,x_{t} , andb and manipulating the resulting equations, the following matrix equation can be obtained.It follows that the tunneling probability is then given by
Detailed descriptions of the terms ∆,
A _{1}^{*},B _{1}^{*},A _{2}^{*}, andB _{2}^{*} are available in [9]. Meanwhile, the tunneling probability for truncated potential barriers based on the WKB method has also been derived in Appendix B2.4.2.2. Simulation
Here the same two potential profiles analyzed in the previous section shall be examined, but with the left half of each profile removed.
Case 1: Truncated exponential potential barrier
Let us take the truncated exponential potential barrier of the inset of Fig. 10 as the first example.
The calculated results based on (54) are displayed in Fig. 10, which shows excellent agreement between the MAF results and the exact TM results. It is seen that the WKB results, on the other hand, fail: As
E /V _{0} approaches unity, the tunneling probability falls to zero. The causes of this behavior can be stated as follows. Atx = 0, the truncation point, the denominators of the terms in (B8) for the region II go to zero asE approachesV _{0}, thereby driving the trial functions to become divergent. It then follows that the coefficients of the trial functions in the region II must go to zero, for the same reason mentioned in the previous case of Fig. 8. This in turn leads to the coefficients of the trial function in the transmitted region, and thus the tunneling probability, to vanish. Based on this finding it may be concluded that whenever there is a truncation in the barrier, the denominator of the trial functions vanishes in the limit ofE →V _{0}, leading to the failure of the WKB method.Case 2: Truncated parabolic potential barrier
A truncated parabolic potential barrier, shown in the inset of Fig. 11, is taken as the second example.
Simulation results are illustrated in Fig. 11. The results are interesting in that both the MAF and the WKB curves deviate significantly from the exact solution, and fall to 0 as
E /V _{0} approaches unity. The reasons for the failure of each respective result are as stated previously: the presence ofV '(x_{p} ) = 0 in the case of the MAF application, and the very presence of a truncation in the case of the WKB application. Such characteristics of potential barriers can be traced to the divergence of the associated trial solutions in region II in the limit ofE →V _{0}, which finally causes the tunneling probability to fall to zero. That the MAF result fares much better (except forE /V _{0} in the range of unity) compared to the symmetric parabolic case of Fig. 8 can be attributed to the fact that the barrier width has been reduced by a factor of one half, thereby increasing the transmission probability, and that the barrier shape has been altered, leading to mitigation of the trial solution’s failure associated with a smooth peak.V. CONCLUSION
The plausibility of applying the MAF method to the analysis of graded optical waveguides and potential barriers was reviewed in this paper, with the conventional WKB method serving as a basis for comparison. Application of the MAF technique to the waveguide problem led to satisfactory results for both symmetric and truncated gradedindex distributions, in that it was possible to derive a simple closedform eigenvalue equation, and the calculated results were very accurate. When the same MAF approach was applied to the problem of tunneling probability in the presence of a graded potential barrier, however, the results were mixed, which was contrary to our initial expectation. The MAF method turns out to be inadequate in any case with
V '(x_{p} ) = 0, that is, when the potential contains a smooth peak, causing the trial functions to diverge at the peak locationx_{p} in theE →V _{0} limit. The WKB approach fails as well for profiles with truncation, because of the diverging trial functions at the truncation point in the same limit.The discrepancy between the end results from the MAF analysis for the two classes of problems, gradedindex optical waveguide and graded potential barrier, is believed to stem from the different respective natures of the two problems. In the waveguide problem only discrete eigenvalues are allowed, and any mode index
N obtained must be less than the peak value of the given index profile, thereby avoiding in the first place the problem of divergence of the trial solutions. In the tunneling problem, however, the concept of discrete eigenvalues does not apply;E the particle energyE is continuousvalued and can be chosen arbitrarily, including theE →V _{0} case, which can cause the trial solutions to become singular. The same reason remains valid with the WKB method.It may thus be concluded that the MAF method is not appropriate for analysis of tunneling problems in which graded potential barriers with the condition of
V '(x_{p} ) = 0 are involved. It also may be added that the WKB method is not suitable for tunneling problems involving abruptly truncated barriers.

[]

[]

[FIG. 1.] Index profile n(x) with turning point xt.

[]

[]

[]

[]

[]

[FIG. 2.] Connection of MAFs at the turning points xt1 and xt2.

[]

[]

[]

[]

[FIG. 3.] Gradedindex profile n2(x). x0 is an arbitrary point between xt1 and xt2.

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[FIG. 4.] υw ？ bw curves of a parabolic index profile for the lowestorder modes.

[]

[]

[]

[]

[FIG. 5.] υw？bw curves of a truncated gradedindex profile for the lowestorder modes.

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[FIG. 6.] Graded potential barrier. xp is the location of the potential peak V0.

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[]

[FIG. 7.] Tunneling probabilities for a symmetric exponential barrier with xp = 0, V0 = 1, and b = 5.

[FIG. 8.] Tunneling probabilities for a parabolic barrier with xp = 0, V0 = 1, and b = 1.

[FIG. 9.] Truncated graded potential barrier.

[]

[]

[]

[]

[]

[]

[FIG. 10.] Tunneling probabilities for a truncated exponential barrier with V0=1 and b=5.

[FIG. 11.] Tunneling probabilities for a truncated parabolic barrier with V0=1 and b=1.