Lambert W function methods in double square well and waveguide problems

Using methods related to the Lambert W function, we present solutions of two apparently different problems: (1) The one-dimensional double square well potential in quantum mechanics and (2) The transverse electric and magnetic modes for a step-index electromagnetic waveguide. The solution techniques provide insight into the bound energy states for the single and double square well problems and the allowed modes of propagation in a waveguide with varying refractive indices. The solutions can be viewed in either of two related complex plane representations. Comparison of the solution geometries suggests that interesting applications may be possible in nanostructures and devices which are designed to be sensitive to small changes in their environment.


Introduction
The one dimensional single and double square quantum potential wells of finite depth are a topic of fundamental interest. The double square well potential has been important in models of the ammonia molecule 6 . Binney and Skinner [1] offer a pedagogical account of the ammonia maser modeled with a one dimensional double square-well potential. Feynman in his celebrated lectures on physics [2] discusses the ammonia molecule and maser. Quantum well models have found many uses in nanostructures; Harrison [3] gives a good survey of this field.
The infinite depth single square well can be solved exactly, as is shown in most quantum mechanics textbooks; see, for instance, Griffiths [7]. It is not quite correctly said that finite depth single square well (FSW) models do not have exact analytic solutions and must be solved numerically or graphically. There are exact solutions which use contour integration in the complex plane, though such methods are beyond the scope of introductory quantum mechanics courses 7 .
Numerical solutions of models can be precise to whatever number of decimal places one wishes to carry out the calculations. The article by Jelic and Marsiglio [9] gives a thorough numerical solution to the double potential well problem. That article facilitates implementations through software tools such as Mathematica and Maple. However, a shortcoming of numerical solutions appears when one wishes to explore the sensitivity of the problem to parameter changes. Each contemplated parameter change can require an excess of calculation to determine derivatives empirically instead of by analysis.
Graphical solutions of models are easier to understand, but one may not be confident in the validity of the results, particularly when the problem parameters are near a transition point where a state may appear or disappear in the solution as a result of some slight change in the device properties or its environment.
Roberts and Valluri [8] outline a geometric-analytic approach that uses the conformal mapping ⟶ = w z w e w to find the solutions of the allowed energy levels of a bound particle in a quantum one dimensional finite single square well (FSW). That mapping is closely related to the Lambert W function 8 . In fact, the Lambert W function is the multi-branch inverse w=W(z) of the mapping = z w e w . In section 2 of the present article we take a similar geometric-analytic approach to analysis of the double well square potential (DSWP). Our objective is to exhibit the bound state energies of the DSWP in geometric terms, with associated smooth functions which can be used in analytic work to investigate sensitivities to parameter changes in the problem. Several of the ideas from the FSW method appear in the DSWP method: There is a map between two complex planes, the w-plane and the z-plane, which provides two visualizations of the situation. There are two equations, called the radial equation (RE) and the structural equation (SE), each of which has a very simple geometry in one of the planes, and a more complicated image in the mapping to the other plane. The solutions are exhibited as intersections of those curves, in either of the two planes. One can work in either plane to determine the number and positions of solutions. By sweeping around the radial curve in the w-plane, one can be sure that all the solutions have been identified. Tangency of solution curves suggests well strengths which may be suitable for adjusting device parameters for high sensitivity, if a detector is being designed. Analytic forms for the mappings provide one with the ability to use calculus to rapidly explore parameter choices.
However, there are also some differences from the single square well potential previously studied. The major difference is that the mapping between the two planes w and z is no longer analytic as a function of complex variables. Rather, it is a smooth function of two real variables. One way to interpret this is that the radial equation, a circle in the w-plane, has become an ellipse-that is, a stretching has occurred in one of the linked planes. Another interpretation is to simply recognize that the mapping between the two planes is not conformal, that is, not angle-preserving. If two curves become tangent in one plane, then they will also be tangent in the other plane. However, if the two curves intersect at a particular angle in one plane, then they will still intersect in the other plane, but possibly at a different angle.
The geometric-analytic two-planes method is also applicable to another type of problem, involving a stepindex planar waveguide. We consider this problem in section 3. The step-index planar waveguide is the simplest type of waveguide. Finding the propagation constants for different modes is a well known standard textbook exercise; see for instance Ghatak and Thyagarajan's book [13] on fiber optics. After preliminary definition and discussion one obtains a transcendental equation such as , along with another constraint which is analogous to the radial equation of the FSW or DSWP problems. As before one can employ various numerical methods to solve these transcendental equations using a computer.
The waveguide problem bears a striking similarity with the square well problems in quantum mechanics. One ends up with essentially the same transcendental equations. Thus, the Schrödinger equation and Maxwell's equations give the same result for the FSW problem. This leads us to think that the field distribution and refractive index for optical waveguides may perhaps be treated as the wave function and potential for a finite square well in quantum mechanics. It is very interesting to see how two different branches of physics, quantum mechanics and electrodynamics, have an analogy between them in the context of this problem. However it is also possible to explore these equations using Lambert W function methods, which will be discussed in section 4. A simple geometric-analytic approach may provide better insight into the problem. As with the earlier square well problems, an advantage of having the solution determined by smooth functions, instead of graphically or numerically, is that subsequent work can become easier to check the sensitivity of the solution to changes in the waveguide parameters.
In section 5 , we discuss the solutions in the w and z planes using the Lambert W function and indicate a few possible applications. Section 6 summarizes the conclusions.

Double square well potential
In this section we will analyze a double square well potential. We restrict to the particular case of a symmetric square well and a central potential island which is higher than the bound energy level. Thus the geometry of the potential well is as shown in figure 1.
Define the symmetric double square well potential (figure 1) as The energy of the particle is -E. Note that E, V 0 and V 1 are positive constants, and V 1 <E<V 0 . The stationary wave function for a bound particle of effective mass m satisfies the time-independent Schrödinger wave equation.
Because of the symmetry of the potential, ( ) ( ) -= V x V x for all x values, it suffices to solve only for wave functions which are either even or odd functions of x. Let ò denote the choice of even wave functions (ò=+1) or odd wave functions (ò=−1). We only need to consider the boundary conditions at x=a and x=L, where both ψ(x) and its derivative ( ) y¢ x must be continuous. Moreover, we need only consider ψ(x) for x0, since x for all x. We will not impose normalization since it has no effect on the conditions we will derive for -E to be the energy of a bound state -except, of course, we require that the wave function be continuous, differentiable, and not identically zero.
The wave function ψ 2 (x) in region 2 (-< < a x a) where the energy -E is less than 2 . Moreover, since region 2 straddles the origin and the wave function must be either even or odd, we can set = =  B B B 2 1 , say, so that equation (1) becomes In region 3 (a<x<L) the energy -E is greater than 2 . The relevant wave function ψ r (x) to the right outside the well, i.e. for x>L, is of the form In the wave function ψ r (x) we have removed the term which blows up as  +¥ x . Each of α, β, γ has been taken positive. Moreover, we can assume that γ is very large, which is equivalent to treating the boundary at x=L as an almost-infinite height potential wall, appropriate for instance in a model of the ammonia molecule as a double well. That is, the significant consideration in determining the energy level of a bound state is x Land x>L) are at zero potential, the middle barrier (region 2, -< < a x a) has a finite depth -V 1 , and the two symmetric wells are regions 1 (-< < -L x a) and region 3 (a<x<L). The particle energy -E lies between -V 1 and -V 0 , so that there are bound states. The depth of the potential in regions 1 and region 3 is -V 0 . the offset of that energy -E from the -V 0 and -V 1 well depth parameters-which are determined by α and β. We already have one of the conditions which describe the problem, namely 2 . That will show up later via the radial equation in our solution.
Next, consider the boundary conditions at x=a and x=L. The wave function and its first derivative must be continuous. That leads to four relationships among A B C , , We wish to derive a structural equation, which is one of the relationships which must be satisfied by the energy -E in terms of the well dimensions a, L, V 0 , V 1 (and m, ÿ). To do this we manipulate equations (5) to (8) in order to eliminate the A B C , , Further simplifications are in order before we obtain the structural equation in the desired form. First, we can use our knowledge that γ is very large in problems of interest, such as an ammonia molecule model. The most direct way is to assume that the boundary x=L is infinite height, and hence the right-hand side of the wave function continuity equation (7) is zero and the derivative continuity equation (8)  and substitution into equation (9) gives The other way to obtain equation (12) is to treat a g i in the power series of equation (10) as zero, in which case = - We also change to dimensionless variables. Make the following substitutions: The variables u, v are dimensionless, and the well geometry parameters f, R are also dimensionless. f represents the width of the two wells compared to the width of the central island.
R is a parameter called the 'strength' of a well in single FSW problems. The u, v variables are related by the radial equation That is, (u,v) is a point on a circle of radius R about the origin in the (u,v) plane. If we write = + w u iv as a complex number, then we can speak of the locus of the radial equation as a circle in the w-plane.
With the above substitutions, equation (12) becomes the structural equation iv u e iv u e iv u e iv u e .
There are of course a variety of equivalent structural equations; for example one can rewrite equation (14) by taking the reciprocals or complex conjugates, or by multiplying out the denominator. Now we come to the central idea of the geometric-analytic method, which might alternatively be called the two-planes method for visualization of the problem, in this case the DSWP bound-energy problem. There are two simultaneous equations to be solved, the radial equation (RE) and the structural equation (SE). Each equation defines a curve in the w-plane, the locus of points (u,v) which satisfy that equation. One of the equations, the RE, is quite simple to visualize in the w-plane as it is a circle of radius R about the origin. However the SE is a rather complicated-looking formula in the w-plane. Still, if we can identify the intersections of the RE and SE curves in the w-plane, we will have the bound state energy levels since -E is readily calculated from either of u or v.
The central idea is to find a smooth map from the w-plane to another plane, the z-plane, in which the SE has a very simple visualization. In contrast the RE in the z-plane may look complicated, for instance a closed curve wrapping several times around the origin and intersecting itself.
One has the choice of which plane to work in, depending on the objectives of one's study of the system of equations. If numerical values are to be calculated one can sweep around the circle of the RE, with a phase angle θ going from p to p + perhaps, calculating u and v accordingly and then substituting into the SE. Thus solving for intersections of the curves is a simple linear root-finding exercise. It may be that finding an RE-SE intersection in the z-plane is much better behaved than finding the corresponding intersection in the w-plane.
We proceed to give the details for the DSWP example. We require a smooth map from the (u,v) pair of real variables, where = + w u iv, to the (x,y) pair of real variables, where = + z x iy. Moreover we want the map to make the image of the locus of (u,v) points satisfying the SE, whose equation in the w-plane is equation (14), into a rather simpler locus of (x,y) points in the z-plane.
With the single FSW solution [8] the map from w-plane to z-plane is = z w e w . The inverse of that complex function is the Lambert W function, so the Lambert W function (which has multiple branches) is the inverse map from the z-plane to the w-plane. The locus of the SE in the w-plane is a set of lines, related to the quadratix of Hippias, which are called the 'Lambert lines'. They are images of the axial rays of the z-plane back in the w-plane, since the locus of the SE in the z-plane is the four axial directed rays from the origin: positive and negative real axes, and positive and negative imaginary axes.
With that example in mind, let's try to identify a mapping of (u,v) points which makes the SE conveniently simple. The SE for the DSWP is (14). We wish to make the two sides of that equation more similar, so we multiply out the denominator and also multiply each side bye ifv to obtain u iv e u iv e u iv e . 1 5 The first map to consider is ( ) ( ) ( , . In the latter form, we are describing the map Y as a function of two real variables u, v, the real and imaginary parts of the complex variable w. Y is considered to give two real variable results x, y, the real and imaginary parts of the complex variable z. The real-variables interpretation of Y is more rigorous mathematically, since z=Y(w) is not necessarily an analytic function of a complex variable w. Whatever the value of f, the function (x, y)=Y(u, v) is a smooth function of two real variables, and one can represent its domain and range within the w and z-planes respectively. If ¹ f 1 the map Y is not analytic, that is not conformal, and does not preserve angles of intersections of curves. However Y is still smooth, and it preserves intersections of curves in the two planes, and notions of closeness, relative magnitudes of angles of intersection, etc. We will use the z=Y(w) notation in the following quite freely, for convenience. Equation (15) becomes, where * denotes the complex conjugate, or in real variable notation, We can see that there are symmetries in equation (16). Swapping v and -v is complex conjugation of the variable w. Swapping u and v with -u and -v is negation of the variable w. Swapping of u and -u while leaving v unchanged is negation followed by conplex conjugation, that is w swapped with -w*. Recognizing those symmetries suggests a second map which is just a concise version of the first one. That map is , .
Using the z=G(w) map, equation (16) can be rewritten as The map z=G(w) is the desired map from the w-plane to the z-plane for the DSWP problem. The SE has become ( ) [ ( )] = G w G w *, which means that G(w) is a real number. That is, the image under the G map of the complicated locus of points (u,v) satisfying the w-plane equation (14) has become a simple locus of points (x, y) on the real line in the z-plane. Connected smooth portions of the SE locus in the w-plane become connected intervals on the real line in the z-plane.
How can we use the two-planes method to find the bound state energy levels for the DSWP? We first choose a value for ò, either ò=+1 to find energy levels for even wave functions, or ò=−1 to find energy levels for odd wave functions. For each choice of ò we proceed: The G map from w to z will transform the circle of radius R about the origin, the locus of points which satisfy the RE, into a multi-loop curve in the z-plane. For visualization we can examine the places where that curve crosses the real axis in the z-plane. Those will be simultaneous solutions of the RE and the SE, and hence give us u, v values which can be used to determine the energy levels. For calculation of numerical results we let the phase angle θ sweep around the circle of radius R in the w-plane, determine u, v via the relationship = = + q w e u iv i , calculate ( ) = z G w , and ( ) = y z Im . The solutions will be at θ values which satisfy y=0, that is where the RE locus crosses the real axis (which is the SE locus) in the z-plane. Because the computation is based upon graphs of well-behaved functions, we can be confident that our 1-dimensional search for solutions will have found all the energy levels.
The DSWP has one further important simplification opportunity. The SE (14) can be reduced to one of two equations, depending upon the choice of ò, i.e, even or odd wave functions. For even wavefunctions the reduced SE is cot and for odd wave functions the reduced SE is cot . The following example uses parameter values a=1, L=5 (hence f=4) and R=3.55. This can be compared with figure 5.7 of Binney and Skinner [1], except that we have set R=3.55 instead of R=3.50 in that book, because we wish to show how a threshold state appears. At approximately R=3.534 a new even solution for a bound energy level appears, but there is not a corresponding odd solution until R=3.603.
The structural equations for this example are shown in figure 2, which graphs left and right sides of the reduced structural equations as a function of the phase angle θ in a sweep around the semi-circle of radius R, for θ from 0 to π. The variables u and v are functions of θ, given by ( ) q = u R cos and ( ) q = v R sin . The intersections of the graph of ( ) u u tanh (purple line, dashes and dots) with the graph of ( ) -v fv cot (black solid curve) represent phase angles θ which correspond to bound energy states with even wave functions. The curve ( ) u u coth (red dashed line) is also graphed as a function of phase angle θ, and its intersections with the graph of ( ) -v fv cot represent phase angles θ which correspond to bound energy states with odd wave functions. The magnitude of the expression | ( , which equals the difference of the left-hand side and righthand sides of equation (15) is shown by a blue line with long dashes. Whenever that magnitude becomes zero, the structural equation (15) will be satisfied, and the blue line will touch the horizontal axis in the graph. Figure 2 is complicated and does not convey much to the viewer except the results of some calculations. Figure 3 shows the intersections of the structural equations in the vicinity of θ = 1.49 rad and 1.65 rad. One can see the intersections which represent the even solutions at those values of the phase angle, and that there are no corresponding odd solutions. A still further zoom-in of the vicinity of θ=1.49 rad, in figure 4, shows the intersection of the graph of ( ) u u tanh (purple, dashes and dots) with the graph of ( ) -v fv cot (black solid curve), which appears exactly above the θ value at which the magnitude | ( , (blue, long dashes) touches the θ-axis. Despite their complexity, figures 2 to 4 do serve to document a computational process. Using the reduced SE equations, together with the RE via a sweep through the phase angles, may be more numerically stable than the two-planes method.
In contrast, the w-plane representation is given in figures 5 and 6. It is much easier to see how the solutions arise as curve intersections in the w-plane. Further, the almost-tangency of the circle of the RE with the Lambert lines indicates that R values near the value R=3.55 are a region where new states may appear or disappear if the well's parameters are slightly altered by changes in the device's dimensions or environment. More precisely, at about R=3.534 a new bound energy level appears with an even wave function. At about R=3.603 (a 2 percent increase in the well strength) the corresponding bound energy level with an odd wave function will appear. Figures 7 and 8 show the z-plane curves for the even (red curve, dashes and dots) and odd (black curve, dots) solutions.

3.
Step-index waveguide 3.1. Solutions for a planar waveguide Optical waveguides are an integral aspect of optoelectronic devices such as optical fibers, reflectors, and semiconductor diodes [3]. They enable one to guide light in chosen directions.
The planar waveguide is a fundamental building block in the theory. It is primarily spatially inhomogeneous in one direction and light confinement is one dimensional. The light is confined in the core of higher refractive  index n 1 . The book [13] by Ghatak and Thyagarajan (G&T) considers the planar waveguide. In this section of our analysis we will parallel parts of their presentation. G & T obtain equations which are structurally similar to those for the single or double finite square quantum wells. Once we arrive at those equations, for the transverse electric or transverse magnetic modes respectively, we will diverge from the G & T presentations to show how the twoplanes geometric-analytic method may be used to visualize and identify the solutions.
We solve Maxwell's equations for the planar waveguide using the contour plots under the given boundary conditions. The planar waveguide under consideration has the profile shown in figure 9 with a step refractive index n(x) given by We impose the condition n 1 >n 2 to get the guided modes. Guided modes are those modes which are mainly confined to the core (| | < x a) and their field should decay in the cover (| | > x a), thus most of the energy associated with the field lies within the core.
Next we write down Maxwell's equations in the non-conducting and non-magnetic medium (μ≈μ 0 ).  where     , , , and denote the electric field, electric displacement, magnetic induction, and magnetic intensity respectively. They are related by the equations We take the curl of (19b) and (19d) to get equations (21a) and (21b): Equations (21a) and (21b) admit plane wave solutions. We choose the direction of propagation along the z-axis. Thus we may write where ω is the angular frequency = k c 0 , k 0 denotes the free space wave number, and c is the speed of light in free space. Also, β is the propagation constant of the mode under consideration. Since the refractive index n(x) depends only on the x coordinate, the field distribution  j and  j , also will depend only on the x coordinate.
For a guided node, the propagation constant β is constrained by the indices of refraction in the core and cover, as described later in this section.  We substitute equations (22a) and (22b) into equations (19b) and (19d) and separate the x, y, z components to get The above manipulations have reduced Maxwell's equations to two independent sets of equations. Equations (23a)-(23c) contain only E y , H x , and H z and equations (24a)-(24c) contain only E x , E z and H y . There are two different types of solutions for these equations: in the Transverse Electric (TE) modes, E z will be zero, making the electric field transverse. In the Transverse Magnetic (TM) modes, B z will be zero, making the magnetic field transverse.

Transverse electric modes
We carry out detailed modal analysis for TE modes in this subsection. We set E z =0 to make the electric field purely transverse and substitute H z and H x from equations (23a) and (23c) into equation (23b) to get We simplify equation (25) using the n(x) profile shown in figure 9.
The boundary condition on E y is that E y and ¢ = E dE dx y y must be continuous at =  x a. Since we want the field to decay outside the core we will choose It is worth noting that the Schrödinger equation for a particle of energy 2 . Since similar differential equations have a similar set of solutions it is more efficient to refer to the solutions of the single quantum FSW problem considered by Roberts and Valluri [8] that are discussed in the following section.
Next we write down the solution to equations (26a) and (26b). We define α and γ via a b = k n Because the distribution of the refractive index is symmetric about x=0, the solution must also be either symmetric or anti-symmetric about x=0. We may write We can simplify the equations (32a) and (32b) by using this symmetry argument. For the symmetric case E y must be of the form

Solutions for the TE and TM modes
In this section, we present the solution of the two transcendental equations (37a) and (37b). We match the solutions for the continuity equation for E y and its derivative at a. By forming a linear combination of the continuity equations we arrive at the following equations.
Next we divide equation (40a) and (40d), to see that B/A=D/C. Let us denote B/A=ò. Now we divide (40b) and (40c) to get = =  A B D C 1 . Hence ò=1/ò, which implies that ò=±1. Thus the four equations (40a)-(40d) reduce to two equations Dividing (41a) and (41b) gives We introduce variables u=γ a and v=α a to get The relation between u and v is as follows We introduce the Lambert W function. This remarkable function has created a renaissance in approaching the solutions of innumerable problems in various fields [10,12]. The Lambert W function W(z) is the analytic multi-branch solution of ( ) ( ) = = W z e we z W z w where z and w=W(z) are complex variables in general. We take the square root of the equation (45) Here, λ represents the square root of ò, as well as the ±1 factor that comes from taking the square root of e i v 2 . Thus λ can be either ±1 or i. Notice that the conjugate, reciprocal, or negative of λ is just another representation of the four alternate choices. Equation (46) can be rewritten as For λ=1, one obtains i.e.
. 50 Equation (50) is similar to equation (36) on appropriate redefinition of variables. The geometric aspect of the equation (47) is interesting. On writing w=u+iv, we observe that the factor e iv doesn't change the magnitude of w=R. It just rotates w in the counterclockwise direction by an angle v. Moreover, the right-hand side of the equation (47) can either be purely real or purely imaginary. Thus, the number cannot lie 'within' the quadrant, it can lie only on the axis.
The geometric description of (47) leads one to the Lambert W function. Consider The mapping = z we w carries thew plane to the z-plane, and the inverse map is the multi-branch Lambert W function. The two lines from the origin in the z-plane along the imaginary axis are the values of λ R when λ=±i. They correspond to the w-plane values for which The two rays along the real axis in the z-plane correspond to the values of w for which =v v u cot when λ=±1. Roberts and Valluri [8] have presented the solution in two different domains: the z-plane and the w-plane.

Solutions in the w and z-planes: two-planes method
We start with the circle in the w-plane. We plot + = u v R 2 2 2 . Then we map the real and imaginary axis of the z-plane to the w-plane using the multi-branch Lambert W function. Thus we get a family of lines in the w-plane. Figure 10 shows the intersection of the Lambert lines with the circle. The intersection of these lines with the circle are the solutions for the propagation modes. We go back to the circle in the w-plane and map it to the z-plane using the map = z we w . This will give a closed, self intersecting, multi-loop curve in the z-plane. This is shown in figure 11. Figure 12 shows a magnified region of the solution in the z-plane between 10 and 25. Figure 13 shows a magnified region of the solution in the z-plane between −1 and 2.5.
We discuss the TM modes given in (39a) and (39b) below. In this particular case, solving for the TM modes of a planar waveguide, the structural equation for the TM modes (symmetric case, or anti-symmetric case) is the same as for the TE modes. We already recognize that the TE modes can be solved using the same map between the two planes as was used for the quantum single square well. Hence we can also use the same map for the TM modes.
Let w and z be complex variables. Define the map F(w) from the w plane to the z plane via ( )  = = w z F w we w . The inverse map is the Lambert W function, ( )  = z w W z . The inverse map is multi-branch, as each z value could come from many w values. For example, the positive imaginary axis in the z-plane, consisting of the complex numbers w which have positive imaginary part and zero real part, goes to an infinite number of lines in the w-plane. Figure 14 shows the loci of points on the w-plane, which map to the the (η, ζ) coordinates of a symmetric TM mode. The propagation constant β can be calculated from the ζ value, as explained earlier in section 3.
As an illustrative example, suppose that n 1 =1.4847 and n 2 =1.4700, and m = a 25 m (micrometers). The value of n 2 for the cover is typical of updoped silica (silicon dioxide) 9 . The value of n 1 for the core film (doped silica)  . Figure 16 shows the TM modes of that example.
Correspondingly, to find the TE nodes, superimpose on the w-plane with Lambert lines the circle which is the radial equation for the TE nodes. Wherever that circle intersects an imaginary Lambert line, gives the (η, ζ) coordinates of an anti-symmetric TE mode. wherever that circle intersects a real Lambert line, gives the (η, ζ) coordinates of a symmetric TE mode. Figure 17 shows the TE modes of that example. The TE and TM modes in the w-plane are both shown in figure 18. The TE modes are shown as red dots, and the circle of the TE radial equation is shown as a dash-dot (Morse TE) line. The TM modes are shown as black dots, and the ellipse of the TM radial equation is shown as a dash-dash-dash (Morse TM) line.

Planar waveguide modes
We have used the Lambert W function to find the propagation constant for the step-index waveguide. The mapping ⟶ = w z we w is analytic and thus conformal. The Lambert W function is the multi-branched inverse of this mapping. We have presented the solution on two different domains. The real and the imaginary axis of the z-plane are mapped to the w-plane which give rise to the Lambert lines in the w-plane. Both TE and TM are shown in figure 18 in the w-plane. The intersection points of the Lambert lines with the circle in the z-plane represent the solution in the w-plane, which is conventional. To find the solution in the z-plane we mapped the circle in the w-plane to the z-plane using the ⟶ w we w mapping. In figure 18 we show both TE and TM modes in the w-plane. It is hard to see the difference between the TM and TE mode solutions, since the TM ellipse is only about 2 percent stretched out horizontally compared to the TE circle. Figure 19 shows a magnified region of figure 18, with the first few solutions in the first quadrant of the w-plane, that is solutions with small positive phase angles. Figure 20 shows the z-plane diagram corresponding to the TM and TE mode solutions depicted in figure 18.
The intersection points of the image of the circle with the axis give the solution in the z-plane. One can use either of the two approaches depending on the nature of the problem. A similar approach has been used by Roberts and Valluri [8,12]. Because the mapping is conformal, the angles between the circle and the Lambert lines in figure 10 are equal to the corresponding angles between the image of the circle and axial rays in figure 11, figure 12, and figure 13. This property can be used to design devices which are highly sensitive to a change in the  environment. In the case of TM modes the final equations are similar to equations (37a) and (37b), and the same process can be repeated to find the propagation constant.
Roberts and Valluri [8] have discussed a similar approach to solve for the energy eigenvalues of the onedimensional finite square well. Their solution also reduces to the infinite potential, when the outer walls go to infinity, as expected. They have used the transformation ⟶ z ze z and its multi-branched inverse, to obtain the solutions in the z-plane and w-plane.

Visualization vs analyticity vs computation
Three techniques are available for investigation of a double well system, such as a quantum potential double well or a planar optical waveguide. We have described the two-planes method in this article. Other techniques are the development of models which have analytic expressions for their solutions, and the detailed calculation of numerical models to high precision. Depending upon the requirements of the investigation, one or the other of these available techniques may be the most suitable.  Clarity of visualization is a strength of the two-planes representation. Looking at the curves in the w-plane and z-plane can help one understand what is going on. The role of a conceptual model of this type is perhaps somewhat similar to a concept such as valence in chemistry; that is, a useful guide to what might happen in nature, without being an exact representation. Applications include a conceptual framework for state splitting in trigonal pyramidal molecules such as ammonia, and a guide to sensitivity of a system based upon tangency of curves.
The simplification of the double square well model shown above as the transition from equation (11) to equation (12) is the removal of a small phase shift represented by the ratio ( ) ( ) a g a g + i i . The magnitude of that factor is 1, but it introduces a phase shift of about 2α / γ radians, and that can be visualized as a rotation of the plane in which the structural equation is represented. It may be interesting to characterize various trigonal pyramidal molecules by their α / γ ratios, with reference to the molecule's stability. For instance, reactive ammonia is intermediate between the stable phosphine and the highly unstable nitrogen triiodide.
Analytic descriptions of molecules which have inversional vibration modes and ground state splitting are developed in two interesting papers by Sitnitsky [14,15]. These are not square well potentials, but they are symmetric dual well potentials, and have analytic expressions for their properties. Thus they are suitable for exploration of molecular properties and assist with exploration of sensitivity to parameter changes. Use of such analytic expressions can speed up computational work as well as improve understanding of the system.
Numerical calculations are the third approach to double well situations. They are precise, to whatever number of decimal places one carries ou the calculations. However, they do not necessarily convey much understanding, if one wishes to find novel materials or to explore parameter sensitivity, other than via exhaustive computational search.
Each of these approaches has a role in double well systems.

Conclusions
We have presented the solution for the propagation constant in the step-index waveguide using the Lambert W function. This approach is a 'geometric-analytic' approach since we have used the special function and later plotted it on two different domains, the z-plane and the w-plane. Our approach facilitates understanding by providing a visualization of the equations. This is useful for parameter variation studies, especially in designing devices sensitive to small variation in parameters. We expect to find more applications of the ubiquitous Lambert W function in future. Kocabaş et al, [16] have used simple mathematical models to study metalinsulator-metal waveguides and observed the close analogy between the one dimensional Schrödinger equation and the electromagnetic wave equations in layered media. Their question whether studies of FSW solutions (and their relationship to changes in reflection spectra of quantum wells) could be useful in optics for the investigation of the effects of material interfaces can be answered in the affirmative. Equations in both problems are of the Sturm-Liouville form and dispersion equations for the TE of a dielectric slab waveguide can be mapped to those of the quantum FSW. This idea can be further extended to the quantum wells having different effective masses in alternating layers of semiconductors. Possible applications are indicated in the quantum confined semiconductors, by sandwiching a 'quantum well' layer of a semiconductor material such as InGaAs between barrier layers of another semicondcutor material Indium Phosphide (InP) [17]. Such semiconductor structures find applications in lasers for telecommunications with optical fibers. [3,17]. Minute changes in the material composition and structure can cause meaningful changes in properties of the engineered materials.
The advent of metamaterials has made the manipulation of minute changes in the environment more feasible, and has led to significant progress in the development of transformation optics (TO) [18].