FDTD vs stackdipole - halfspace emission in a multilayer stack

While FDTD can be used to simulate OLED/LED designs of arbitrary geometries, these simulations tend to be very time consuming. Often, a more efficient approach is to start by using analytical methods to analyze the planar structure first. In this section we will explain how the stackdipole command (included in Lumerical's STACK optical solver) can be used to analytically calculate the farfield power density [W/(m² ×sr)] as a function of polar angle (θ).

Theory

For each dipole with given emission angle, stackdipole returns the radiance:

$$\text{stackdipole}(\theta) = \int_\lambda (j  \cdot ef \cdot st)\times \left( \frac{rd \cdot F_{rad}(\theta, \lambda)}{rd \cdot F(\lambda) + (1-rd)} \right) \cdot \left(\text{photon probability}(\lambda) \cdot E_{ph} (\lambda)\right) d\lambda$$

where:

This quantity has units of [W/(m² × sr)]. The Purcell factor (F(λ)) shows enhancement in the spontaneous emission rate of an emitter inside a microcavity and  F_rad(λ,θ) is the power density. By defining P_rad  as power radiated into the far field, P_nonrad as power lost due to absorption or captured in evanescent fields, and P₀ as the power that would be emitted in an infinite uniform medium, we have:

$$\mathrm{Purcell\ factor} = F(\lambda) = \frac{P_{rad} + P_{nonrad}}{P_0}$$

$$\mathrm{Power\ density} = F_{rad}(\theta, \lambda) = \frac{P_{rad}}{P_0}$$

$$\mathrm{Extraction\ efficiency} = \frac{F_{rad}(\theta, \lambda)}{F(\lambda)} = \frac{P_{rad}}{P_{rad} + P_{nonrad}}$$

The first factor in the integrand, j × ef × st, is the rate of singlet exciton decays, calculated for a fixed current density of 1 A/M² and user-specified exciton fraction (ef) and singlet-triple ratio (st). It is calculated as one Coulomb in atomic units (6.241e+¹⁸ electrons) multiplied by the user-specified constants. This quantity outputs the total number of generated photons.

The second is the angular density of the Quantum Yield. Calculation of this quantity in FDTD is described in the next section. The Purcell factor is calculated analytically from the stack geometry using dipole illumination. The relative decay rate (rd) denotes the proportion of singlet decays that produce a photon.

The third is derived from the user-specified intensity spectrum, normalized according to midpoint integration in wavelength. This quantity accounts for the dipole emission spectrum. The stackdipole command will automatically normalize the intensity spectrum.

FDTD calculation of quantum yield

To calculate F_rad(θ) from FDTD simulations, a Frequency-domain field and power monitor is employed to capture the near field profile. Then farfield commands are employed to project near to far field profile. The post processing step for calculating quantum yield is summarized below:

T    = transmission(monitorname);
E2   = farfield3d(monitorname);
ux   = farfieldux(monitorname);
uy   = farfielduy(monitorname);
angDistrib = E2/farfield3dintegrate(E2,ux,uy);  # normalized angular distribution
angDistrib = farfieldspherical(angDistrib,ux,uy,θ,0);  # mapped into phi=0 to eliminate azimuthal angle
Frad = T * angDistrib;
F    = dipolepower(c/λ) / sourcepower(c/λ);   # Purcell factor
quantumYieldDensity = rd*Frad/(rd*F + (1-rd));

Due to azimuthal symmetry, far field profile is mapped onto the φ=0 plane to eliminate the azimuthal angle. Note that the value φ=0 is suitable for a vertically or randomly oriented dipole.  For horizontally-oriented dipoles, either φ=0, φ=90, or average of those results can be used.

The angular distribution is normalized to one. One can check this by integrating it over a half space in spherical coordinates:

?integrate(angDistrib * sin(θ),1,θ)*2*pi;  # output should be 1

where θ is an array of values from 0 to π/2.

Dipole polarization

The stackdipole command can calculate the radiance for randomly oriented, P_vertical, or P_horizontal orientations.  Since STACK is a 1D solver we can only distinguish between horizontal and vertically oriented dipoles. Results for horizontal polarizations are normalized such that:

$$\text{random }= \frac{1}{3}P_{vertical} + \frac{2}{3} P_{horizontal} $$

If we define  F_pVertical , F_sHorizontal , F_pHorizontal as the farfield radiance calculated from FDTD simulations for vertically oriented dipole, horizontally oriented dipole with φ=90, and horizontally oriented dipole with φ=0, then the farfield radiance of randomly oriented and horizontal randomly oriented dipole are given by:

$$\text{random } = \frac{1}{3} \left( F_{pVertical} + F_{pHorizontal} + F_{sHorizontal} \right)$$

$$\text{horizontal } = \frac{1}{2} \left(F_{pHorizontal} + F_{sHorizontal} \right)$$

FDTD vs stackdipole

To compare the results of stackdipole with FDTD simulations, we consider a simple case in which the dipole is located in free space, 80nm above a medium with dielectric constant of 1.5. The FDTD simulation is performed in 3D and the power monitor is located above the dipole (upward) to capture light traveling in the +z direction. To create the plots below, download fdtd_dipole_halfspace.fsp and run fdtd_dipole_halfspace.lsf script. The script will save three different FDTD simulation files for different dipole polarizations along each Cartesian axes, and will compare it with the analytical results calculated directly from the stackdipole command. Since these simulations are for a single frequency source, the photon probability is 1 assuming that all the light is injected at the central frequency.

Important Settings

Material Fit

There are multiple script commands that can be used to get the material properties. When comparing STACK with FDTD, or another solver it is important to ensure the material values are similar if not the exact. Both getindex - Script command and getfdtdindex - Script command, will return the material properties can be used to build the stack; however, they do not entirely agree in general. The get index method will use a linear interpolation between material data points. This is common in frequency domain solvers, such as STACK, FDE or EME.

As a time domain method FDTD requires fitting the material to a superposition of functions that obey causality. For more information on how to use this FDTD's multi-coefficient material model MCM see Using the Material Explorer to view and adjust optical material models, or for a mathematical explanation see Kramers-Konig. 

When benchmarking your results against other solvers it is common to use linear interpolation. When comparing FDTD and STACK one should use the get FDTD index method, and fine tuning the material models may be required. When confirming the validity of your simulation set-up it may be easiest to stick to object defined dielectrics, and PEC for metals.

Resolving Interfaces

STACK is an analytic method, and the interface calculations will be done exactly where they are defined. FDTD on the other hand the requires a mesh on which to perform the numerical calculations. FDTD cannot resolve features smaller then the mesh and so if an interface is drawn through a mesh cells there will be some ambiguity. This sub mesh cell ambiguity can produce small discrepancies which can be overcome by using a smaller mesh and/or including mesh overrides to ensure the interfaces align with the mesh cell boundaries. 

Frequency Sampling for Farfield Analysis

When comparing FDTD with STACK you should use the farfield emission profiles. To obtain the most accurate results recommended to use Chebyshev frequency sampling. 

See also

Radiated power from incoherent isotropic point dipole sources, stackdipole, stackpurcell