Boundary Integral Equation Methods For Computational Photonics

Optical waveguides are structures that guide the propagation of light. They are the fundamental components in communications systems and integrated optical circuits. In recent years, many optical waveguides with complicated structures have appeared. As a special class of optical waveguides, photonic crystal fibers (PCFs) have been extensive-ly studied because of their many unique properties which are not available in traditional waveguides. The propagation of light in a PCF is strongly controlled by the geometry of its cross section. Periodic structures, such as diffraction gratings and photonic crys-tals (PhCs) are important optical components that can be used to control and manipulate light. Accurate and efficient numerical methods are essential in the analysis, design and optimization of optical waveguides and periodic structures.For a given optical waveguide or PCF, numerical methods that discretize the cross section of the structure give rise to linear matrix eigenvalue problems. The discretiza-tion can be obtained by using the finite difference method, the finite element method (FEM), the multi-domain pseudospectral method, etc. However, for PCFs with many holes and complicate geometries and for general optical waveguides with high-index contrast, sharp corners and complex micro-structures, the resulting matrices can be very large and the matrix eigenvalue problem can only be solved by iterative methods and the accuracy may be limited. A better approach is to formulate a nonlinear eigenval-ue problem of which the resulting matrix is much smaller. Numerical methods using the nonlinear approach include the film mode matching method, the multipole method and the boundary integral equation (BIE) method. The film mode matching method is quite successful, but it is only applicable to optical waveguides with vertical and hori-zontal refractive index discontinuities. The multipole method is accurate for PCFs with well-separated and circular inclusions, but it cannot be easily extended to other optical waveguides. For diffraction gratings, existing numerical methods include general-purpose meth-ods such as the finite-difference time-domain (FDTD) method and the FEM, and more special methods such as the analytic modal method, numerical modal methods, the BIE methods, etc. Although FDTD and FEM are extremely versatile, they are typically less efficient than the special methods. Analytic and numerical modal methods require that the structure consists of uniform layers. For gratings with high index-contrast and sharp corners in their profiles, all modal methods converge slowly and may even fail to converge, due to the possible field singularity at the corners. And for two dimension-al (2D) PhCs with circular cylinders, existing numerical methods such as the FDTD method, the FEM method, the multipole method, the scattering matrix method and the Dirichlet-to-Neumann map method are effective. However, if the cylinder in each unit cell contains corners, the above methods still suffer from a considerable loss of accuracy in the presence of the field singularity at corners.In this thesis, high order boundary integral equation methods are developed for an-alyzing optical waveguides including PCFs, diffraction gratings and photonic crystals of arbitrary unit cells. The methods rely on a standard Nystrom method for discretizing integral operators and they do not require analytic properties of the electromagnetic field (which are singular) at the corners. For PCFs with smooth interfaces, we develop a new high order BIE mode solver. The method solves two functions on the interfaces and is more efficient than existing BIE methods. The key step is to use the kernel-splitting technique for discretizing the hyper-singular boundary integral operators. For optical waveguides with high index-contrast and sharp corners, a new full-vectorial waveguide mode solver is developed based on a new formulation of boundary integral equations and the so-called Neumann-to-Dirichlet (NtD) maps for sub-domains of constant re-fractive index. The method uses the normal derivatives of the two transverse magnetic field components as the basic unknown functions, and it offers higher order of accu-racy where the order depends on a parameter used in a graded mesh for handling the corners. For diffraction gratings, we present a high order BIE-NtD method which is an improved-version of a BIE-NtD method in earlier works. The improvements include a revised formulation that is more stable numerically, and more accurate methods for computing tangential derivatives along material interfaces and for matching boundary conditions with the homogeneous top and bottom regions. For2D PhCs of arbitrary unit cells, a new BIE-NtD method is used to calculate the NtD map for each unit cell. We study two basic problems encountered in the analysis of2D PhCs, A projection technique is used for further reducing the size of the reduced NtD map for each unit cell, and it makes our method more effective.