**Robert S. Upton and R. John Koshel**

Tools and techniques for modeling the propagation of coherent fields—Fourier methods applied to Gaussian-beam decomposition and the finite-difference beam-propagation method—can help ensure accurate and efficient coupling. But all perspectives, from the macroscopic to the microscopic, must be considered.

Modeling the propagation of coherent fields, from the macroscopic (such as bulk optics) to the microscopic (such as guided optics) has typically involved two methods. For free-space propagation in bulk optics, Fourier methods have been applied to Gaussian-beam decomposition (GBD). For the guided-wave domain, methods based on finite-difference methods, either in time or space, can be employed.

For example, the finite-difference beam-propagation method (FD-BPM) is used to approximate the behavior of optical wave fields in media with microscopic variations in refractive index and in optical fibers, and it can be used to model devices that utilize evanescent coupling. Entire optical systems, from bulk optics to optical fibers and other integrated optical components, can be modeled with both GBD and FD-BPM. The major drawback of FD-BPM is that backward propagating optical fields are not calculated.

Modeling of the electric wave-field propagation in inhomogeneous optical systems using FD-BPM can be performed with Advanced Systems Analysis Program (ASAP) of Breault Research Organization. FD-BPM models the propagation of optical wave fields in media that are inhomogeneous and isotropic—that is, media with microscopic variations in refractive index where the refractive index is fully represented by a three-dimensional scalar function. The FD-BPM tool in the ASAP complements the GBD1 used to model the propagation of scalar wave fields through optical systems immersed in free space (such as bulk optics). So a complete coherent optical system can be modeled in both the free-space and guided-wave domains.

In FD-BPM, the wave field is propagated along the axial direction by dividing the three-dimensional volume into a sequence of planes each separated by a small and finite distance Dz. The wave field is defined in an initial plane and only the forward-propagating wave field is chosen. Hence, the FD-BPM technique cannot be reliably used for modeling systems that require the formation of the backward propagating wave field.

FD-BPM is more reliable than the Fourier BPM methods used,3 because the Fourier components of wave propagation require the mean refractive index in a given plane to be defined, so that the wavelength of the field in that plane can also be defined. FD-BPM uses a correction factor, similar to the Wentzel, Kramers, and Brillouin (WKB) approximation used in quantum theory, so the average refractive index is included to account for the inhomogeneity of the medium.

Other formalisms, such as the finite difference time-domain element (FDTD) method, model the propagation of electromagnetic wave fields by matching the wave-field boundary conditions at the edges of a four-dimensional (three spatial and one time) cell. The FDTD approach is the most rigorous and yields the most accurate answers of all the calculation formalisms. However, FDTD requires that the spatial volume of each cell be on the order of approximately (0.1l)2, where l is the wavelength of light. As a result, a vast number of cells are required to model the propagation of the wave field. Typically, 1000 x 1000 x 1000 cells are required, which requires a significant investment of computer time and resources.

FD-BPM is found to be a practical alternative that, under conditions of small variations in refractive index over the dimension of a wavelength, can model the propagation of vector, semivector3 (wave field has an axial component and a single transverse component only), and scalar wave fields through inhomogeneous media. But FD-BPM cannot model the propagation of fields in magnetic or birefringent (anisotropic) media.

FD-BPM is most appropriately used in modeling the general shape and phase of the propagating field in the forward propagating direction. Therefore, FD-BPM can be used to model the propagation of wave fields in fibers and integrated optical components, such as evanescent couplers and Y-branches. FD-BPM cannot calculate the reflected field components at dielectric discontinuities, so it cannot be used to model distributed feedback (DFB) lasers, fiber Bragg gratings (FBGs), distributed Bragg gratings, or other devices that utilize the reflected wave fields without complicating the model by including reflections in a perturbative method.

**FD-BPM DEVELOPMENT**

The advantages of FD-BPM are that it is relatively simple to implement, requires fewer computational resources than FDTD, and wave fields in inhomogeneous media can be modeled. The formalism of FD-BPM starts with Maxwell's equations within a nonmagnetic, isotropic, and inhomogeneous medium.

The resulting wave equation that is derived from Maxwell's equations is then simplified by assuming a scalar wave field. For media with large variations of refractive index over distances on the order of a wavelength, the implementation of BPM requires an approximation of the refractive index variation. In ASAP, where the ASAP user can define higher calculation accuracy, this approximation is more accurate compared to more classical implementations of BPM. However, increased accuracy results in slower calculation of the wave-field propagation. Practice has shown that this variation can be as great as going from air (n = 1.0) to a high-index material such as germanium (n ~ 3.5).

However, to remain accurate a larger number of terms in the FD-BPM need to be retained, thus slowing down the calculation. For example, ASAP includes three levels of accuracy (low, medium, and high) by the inclusion of additional terms in the model in order to provide a trade-off between calculation speed and accuracy of the output.

The system to be studied is the coupling of a laser diode into a single-mode fiber via a ball lens. Since the code is geometry-based, one can easily investigate the tolerancing of fiber axial placement to the position of the ball lens. The guided field is then propagated down the fiber until an X-branch is reached. At this point evanescent coupling to the neighboring fiber is calculated using FD-BPM. One can continue propagation using FD-BPM until the fiber output is reached. At which point, GBD can be utilized to model the coherent propagation in the bulk optics. In this example, a high degree of accuracy is used.

**FIBER-COUPLING CALCULATIONS**

A laser diode is modeled to emit at 1.55 µm into a Gaussian distribution with transverse emission angles of 30° by 15°. This light is focused by a ball lens with a 0.3-mm radius and index of 2.14 into a single-mode fiber. The single-mode fiber has a core radius of 2.45 µm and index of 1.456, a cladding radius of 50 µm and index of 1.449, and an outer coating radius of 62.5 µm and absorption factor of 0.8 µm^{1} (see Fig. 1). The ball lens introduces significant amounts of spherical aberration, which affects the coupling efficiency of the laser light into the fiber mode. The single-mode fiber is modeled with a gradient index (GRIN) step index medium in ASAP.

The coupling efficiency as a function of the axial position of the fiber in the spherical aberration caustic can be plotted (see Fig. 2, top). Two curves are realized. The first is the coupling efficiency calculated without the use of FD-BPM, where the light field is propagated to within 0.1 µm of the fiber face using GBD and then coupled into the approximate mode of the fiber. The second is the coupling efficiency calculated by propagating the light field 0.1 µm into the fiber using GBD and propagating the field an additional 20 µm down the fiber using FD-BPM.

The major differences between the two results are that the efficiency calculated using the FD-BPM method is higher for all axial positions, and, in both cases, the maximum coupling efficiencies do not coincide with the peak irradiance of the spherical aberration caustic. The FD-BPM method takes into account the affects of the index change while the GBD method approximates it. Therefore, the FD-BPM provides a more accurate representation of the realities of the system.

The difference in the peak efficiencies and peak irradiance in the caustic (as shown by comparing the top and bottom graphs in Fig. 2) indicates that peak coupling efficiency does not occur at best focus for minimum RMS wavefront error. Furthermore, the irradiance distribution within the fiber core and cladding can be graphically analyzed (see Fig. 3).

**EVANESCENT COUPLING**

Frustrated total internal reflection (FTIR) can be used to couple light fields between dielectric structures in close proximity by evanescent coupling; for example, ifn >η TIR occurs when θ' is greater than 90°. In other words, the light ray does not pass from medium 1 to medium 2, but stays in medium 1. Analyzing TIR with wave optics and electromagnetic theory, light is totally reflected at the interface between medium 1 and 2. However, a standing wave exists along the media boundary in medium 2, whose amplitude decays exponentially with distance from the media boundary.

This standing wave field is referred to as an evanescent wave. The evanescent wave field does not propagate. Bringing a structure made of a third dielectric medium, with refractive index n", into close proximity (within the dimension of a wavelength) of the dielectric boundary results in the coupling of the evanescent field to the third dielectric medium. This process is FTIR, which is used in many micro-optical systems and integrated optical systems, as well as near-field optical data storage.

So the field coupled into the first fiber can be "switched" over to the second fiber in an X-branch simply by placing the two fibers in close proximity. Thus, FD-BPM can be used to model the evanescent coupling of light between fibers placed in close proximity (see Fig. 4). This field can then be propagated down the fiber to investigate the effects of additional structures. The inclusion of FD-BPM means that coherent fields can be propagated from the free-space domain to the guided-wave domain. This feat can be accomplished seamlessly within the single piece of code.

**FIGURE 4. The geometry for evanescent coupling between fibers (top). The evanescent field energy coupled between the fibers calculated using the FD-BPM algorithm in the Advanced Systems Analysis Program (bottom). **

Additionally, because the code is based on geometry, investigations are not limited to ideal systems, but can be used for real systems such astolerances and potential fabrication errors.

**ACKNOWLEDGMENTS**

The authors wish to acknowledge the help and contributions of Alan Greynolds, Don Scott, and Carey Portnoy at Breault Research Organization.

**REFERENCES**

- A. Greynolds, Proc. SPIE, 679, (August 1986).
- L. Vincetti et al., J. Opt. Soc. Am., 17, 1124 (2000).
- J. Van Roey et al., J. Opt. Soc. Am., 71, 803 (1981).

**Robert Upton** is an optical engineer in charge of modeling and analyzing coherent optical systems. John Koshel is the manager of optical engineering services at Breault Research Organization, 6400 E.Grant Rd, Suite 350, Tucson, AZ 85715. They can be reached at 520-721-0500, or rupton@breault.com or jkoshel@breault.com.