OptiFDTD is a powerful, highly integrated, user-friendly software that allows computer

aided design and simulation of advanced passive photonic components.

The OptiFDTD software package is based on the finite-difference time-domain

(FDTD) method. The FDTD method has been established as a powerful engineering

tool for integrated and diffractive optics device simulations. This is due to its unique

combination of features, such as the ability to model light propagation, scattering and

diffraction, and reflection and polarization effects. It can also model material

anisotropy and dispersion without any pre-assumption of field behavior such as the

slowly varying amplitude approximation. The method allows for the effective and

powerful simulation and analysis of sub-micron devices with very fine structural

details. A sub-micron scale implies a high degree of light confinement and

correspondingly, the large refractive index difference of the materials (mostly

semiconductors) to be used in a typical device design.

## 2D FDTD Equations

The FDTD approach is based on a direct numerical solution of the time-dependent

Maxwell’s curl equations. The first version of OptiFDTD is in 2D. The photonic device

is laid out in the X-Z plane. The propagation is along Z. The Y-direction is assumed to

be infinite. This assumption removes all the ∂/ ∂*y* derivatives from Maxwell’s

equations and splits them into two (TE and TM) independent sets of equations.

The 2D computational domain is shown in Figure 1. The space steps in the X and Z

directions are Δx and Δy, respectively. Each mesh point is associated with a specific type of material and contains information about its properties such as refractive index, and dispersion parameters.

**Figure 1: Numerical representation of the 2D computational domain**

## TE waves

In the 2D TE case (H_{x}, E_{y}, H_{z} – nonzero components, propagation along Z, transverse

field variations along X) in lossless media, Maxwell’s equations take the following

form:

where ε = ε_{0}ε_{r} is the dielectric permittivity and is the magnetic permeability of

the vacuum. The refractive index is n = √ε_{r}.

Each field is represented by a 2D array — E_{y}(i,k), H_{x}(i,k) and H_{z}(i,k) — corresponding

to the 2D mesh grid given in Figure 1. The indices i and k account for the number of

space steps in the X and Z direction, respectively. In the case of TE, the location of

the fields in the mesh is shown in Figure 2.

**Figure 2: Location of the TE fields in the computational domain**

The TE fields stencil can be explained as follows. The E_{y} field locations coincide with

the mesh nodes given in Figure 1. In Figure 2, the solid lines represent the mesh given

in Figure 1. The E_{y} field is considered to be the center of the FDTD space cell. The

dashed lines form the FDTD cells. The magnetic fields H_{x} and H_{z} are associated with

cell edges. The locations of the electric fields are associated with integer values of the

indices i and k. The H_{x} field is associated with integer i and (k + 0.5) indices. The H_{z}

field is associated with (i + 0.5) and integer k indices. The numerical analog in

Equation 1 can be derived from the following relation:

The total set of numerical Equation 1 takes the form:

The superscript n labels the time steps while the indices i and k label the space steps and Δx and Δz along the x and z directions, respectively. This is the so-called Yee’s numerical scheme applied to the 2D TE case. It uses central difference approximations for the numerical derivatives in space and time, both having second order accuracy. The sampling in space is on a sub-wavelength scale. Typically, 10 to 20 steps per wavelength are needed. The sampling in time is selected to ensure numerical stability of the algorithm. The time step is determined by the Courant limit:

## TM waves

In the 2D TM case (E_{x}, H_{y}, E_{z} — nonzero components, propagation along Z,

transverse field variations along X) in lossless media, Maxwell’s equations take the

following form:

The location of the TM fields in the computational domain follows the same

philosophy and is shown in Figure 3.

**Figure 3: Location of the TM fields in the computational domain**

Now, the electric field components Ex and Ez are associated with the cell edges, while

the magnetic field Hy is located at the cell center. The TM algorithm can be presented

in a way similar to Equation 3.

## 3D FDTD Equations

In 3D simulations, the simulation domain is a cubic box, the space steps are Dx, Dy, and Dz in x, y, and z directions respectively.

Each field components is presented by a 3D array —**Ex(i,j,k), Ey(i,j,k), Ez(i,j,k),**

**Hx(i,j,k), Hy(i,j,k), Hz(i,j,k)**. The field components position in Yee’s Cell are shown

in Figure 4. These placements and the notation show that the E and H components are interleaved at intervals of 1 ⁄ 2Dh in space and 1 ⁄ 2Dt for the purpose of implementing a leapfrog algorithm.

**Figure 4: Displacement of the electric and magnetic field vector components about a cubic unit cell of the Yee space lattice**

In general, the time domain Maxwell’s equations are given in differential form by

When the second-order finite difference method, field notations, and field

displacements are applied to the above Maxwell’s equation for time domain and

space domain derivatives, the 3D-FDTD formulas can be written as:

## Space Step and Time Step

The fundamental constraint of FDTD method is the step size both for the time and space. Space and time steps relate to the accuracy, numerical dispersion, and the stability of the FDTD method. Many references and books have discussed these problems. In general, to keep the results as accurate as possible, with a low numerical dispersion, the mesh size often quoted is “10 cells per wavelength”, meaning that the side of each cell should be 1/10λ or less at the highest frequency (shortest wavelength).

Please note that FDTD is a volumetric computational method, so that if some portion

of the computational space is filled with penetrable material, you must use the

wavelength in the material to determine the maximum cell size.

The following equation is for the suggested mesh size:

where n_{max} is the maximum refractive index value in the computational domain.

Once the cell size is determined, the maximum size for the time step Δt immediately

follows the Courant-Friedrichs-Levy (CFL) condition.

For 3D FDTD simulation, the CFL condition is:

where v is the speed of the light in medium.

## OptiFDTD Simulation Procedures

The following is the flow chart for the FDTD simulation in OptiFDTD. It also details the

work flow in OptiFDTD.

**Figure 5: FDTD Simulation Flow Chart in OptiFDTD**

## Output data

The fields propagated by the FDTD algorithm are the time domain fields. At each

location of the computational domain they have a form similar to that given in Equation

15:

where B is the amplitude of the field at that particular location, G is the wave profile,

and φ_{i} is the corresponding phase. However, the values of B and φi are not

accessible from the time domain field values.

In order to get the full amplitude/phase wave information, we need the stationary

complex fields that correspond to the waveform Equation 15. The complex fields are

the source of all useful information, such as output and reflected powers, overlap

integrals with modal fields, etc. Those complex fields are calculated by a run time

Fourier transform performed in the last time period of the simulation. The final

complex fields can be visualized at specific Output Planes located properly in the

computational domain.

**Figure 6: Output Planes**

OptiFDTD uses TF/SF (total field/scattering field) technique for the incident plane.

User can specify the incident wave direction. Behind the incident plane, it is the pure

reflection field region, when the observation detectors are placed in this region, the

reflection function can be calculated. When the Observation detectors are placed in

the field transmission region, the transmission function can be calculated.