In this example, we simulate a ferroelectric waveguide modulator which uses BaTiO3 - a material which exhibits a change in refractive index due to applied electric field. This device is based on the publication by Xiong et al.[1] We extract the effective index of the waveguide modulator versus voltage for a given operating frequency.
Background
The ferroelectric waveguide consists of a silicon layer and BaTiO3 (also called BTO) layer on top of a glass substrate. The BaTiO3 crystal orientation is such that the [011] direction of the crystal is parallel to the propagation direction (y-direction), and the [001] direction is along the z-direction. A strip of amorphous silicon is placed on top of the BaTiO3 layer to confine light in the lateral direction. Gold contacts are placed at a distance of 1 um from the amorphous silicon rib waveguide on either side.
In the following sections, we first use the CHARGE solver to simulation a range of applied bias voltages and export the electric field distribution across the waveguide cross section for each bias voltage. Then, we use the electric field distribution to calculate the change in refractive index over space in the BaTiO3 material and simulate the effective index of the waveguide for each bias voltage.
Electric field profile simulation with CHARGE
In the ferroelectric_modulator.ldev file's CHARGE solver object, an electric field monitor named “E_field” is set up to export the E field data over space to a file called WG_Efield.mat.
Under the Boundary Conditions tab, the cathode contact is set to 0 V, and the anode contact is set to sweep voltages over the 1 – 5 V range with a 0.5 V interval between each sweep point.
Setup comments
Effect of metal work function
If you plot the band diagram of the silicon slab, one thing that is not obvious is why the band is not flat at zero bias (fig. left). This is because in the vertical direction, we have a metal-insulator-semiconductor structure. In such a structure, the Fermi level of the metal aligns with the Fermi level of the semiconductor. If the two are not aligned to start with (i.e. the metal work function is not equal to the semiconductor work function) then there is a band bending in the semiconductor (fig. right). The amount of this band bending is equal to the difference in the metal and semiconductor work function.
Wider simulation region
In CHARGE, a semiconductor region must have a voltage boundary condition. The voltage boundary condition can come directly from a metal contact or through a neighboring semiconductor layer. In this particular device, the crystalline silicon layer is insulated from the anode and cathode contacts. Thus it has no voltage boundary condition. The solution is such case is not unique and may not be accurate. To provide a boundary condition, we use a grounding contact (gnd). However, we do not want the potential set by the grounding contact to affect the behavior of the silicon slab under the anode and cathode contacts. This is why, the ‘gnd’ contacts need to be placed far away so the potential in the silicon slab immediately under the anode/cathode contacts is controlled by the voltage applied at those contacts only. This is ensured by making sure that the band in the silicon slab is flat near the grounding contacts.
Thin BTO wire
Due to the same reason as the one mentioned above, the Si layer in the waveguide also requires a voltage boundary condition. Instead of providing the boundary condition directly from a contact (gnd) we provided the boundary condition through another semiconductor region. We have thus created a large bandgap semiconductor with bandgap and dielectric constant equal to the BTO so that in the electrostatic solution, it will act in a manner similar to the (insulator) BTO layer. However, this semiconductor version of the BTO will connect the Si layer in the waveguide to the crystalline silicon layer and thereby providing it with a voltage boundary condition.
Effective index calculation with MODE
Simulation method
We can use the FDE solver in MODE find the guided mode on a waveguide cross section and extract its effective index. To get the plot of effective index over bias voltage, first download ferroelectric_modulator.lms, and ferroelectric_neff_vs_voltage.lsf. If you have run the CHARGE simulation beforehand you will already have the WG_Efield.mat data file which contains the calculated electric field profile over space for each bias voltage. If not, you can download the file from the Associated files section above. Then load the ferroelectric_modulator.lms simulation file and the ferroelectric_neff_vs_voltage.lsf script file and run the script file. Note that this will launch one FDE simulation per bias voltage point.
Structure
The BaTiO3 crystal is a tetragonal crystal, with an ordinary refractive index, no, in the [100] and [010] directions and an extraordinary refractive index, ne, in the [001] direction. In the example files we use no=2.379, and ne=2.339, but for realistic results which accurately represent the BaTiO3 material you may need to obtain the refractive index values from measured data or a reference publication.
Tetragonal crystal directions
In the waveguide, the BaTiO3 crystal is oriented so that the [011] direction of the crystal is parallel to the direction of propagation along the waveguide, so the permittivity tensor of the crystal will need to be rotated in order to account for this.
In the coordinate system of the un-rotated crystal (the prime coordinate system where x’, y’, z’ correspond to [001], [010], [001] crystal directions), the permittivity tensor of the BaTiO3 material with no applied electric field is:
$$n=\left[\begin{array}{ccc}{n_{0}} & {0} & {0} \\ {0} & {n_{0}} & {0} \\ {0} & {0} & {n_{e}}\end{array}\right]$$
The perturbation of the inverse permittivity tensor due to applied electric field is given by the following:
$$\Delta \frac{1}{n^{2}}=\left[\begin{array}{ccc}{r_{13} E_{z^{\prime}}} & {0} & {r_{51} E_{x^{\prime}}} \\ {0} & {r_{13} E_{z^{\prime}}} & {r_{42} E_{y^{\prime}}} \\ {r_{51} E_{x^{\prime}}} & {r_{42} E_{y^{\prime}}} & {r_{33} E_{z^{\prime}}}\end{array}\right]$$
The perturbed refractive index squared is then given by
$$n^{2}=i n v\left(\left[\begin{array}{ccc}{n_{0}^{2}} & {0} & {0} \\ {0} & {n_{0}^{2}} & {0} \\ {0} & {0} & {n_{e}^{2}}\end{array}\right]^{-1}+\left[\begin{array}{ccc}{r_{13} E_{z^{\prime}}} & {0} & {r_{51} E_{x^{\prime}}} \\ {0} & {r_{13} E_{z^{\prime}}} & {r_{42} E_{y^{\prime}}} \\ {r_{51} E_{x^{\prime}}} & {r_{42} E_{y^{\prime}}} & {r_{33} E_{z^{\prime}}}\end{array}\right]\right)$$
Values of r13, r33, r42, r51 used in the example are obtained from references 1 and 2 listed in the related publication section at the top of the page.
The permittivity tensor contains off-diagonal components in the presence of electric field and we can represent this using a combination of the diagonalized permittivity and an applied matrix transform grid attribute. This method of representing a general permittivity tensor with off-diagonal components is discussed in more detail in the Matrix Transform chapter.
The analysis group named "BTO" in the ferroelectric_modulator.lms file uses a setup script to add the BTO film with permittivity perturbation due to applied electric field. The electric field profile can either be loaded from a prior CHARGE simulation if "import_E_field" is set to 1, or you can set it to 0 to specify a uniform electric field profile.
If the "import E field" setting is 1, the setup script performs the following:
1) Loads the electric field profile data from the WG_Efield.mat file for the specified bias point is loaded.
2) Interpolates the electric field data from the triangular lattice of the CHARGE simulation onto a rectangular grid of positions in the BaTiO3 thin film.
3) Calculates the electric field components in the prime coordinate system of the un-rotated crystal.
4) Loops through each position in the BaTiO3 crystal film and calculates the perturbed permittivity tensor of the material in the presence of the electric field.
5) Calculates the combination of diagonalized refractive index and matrix transform required to represent the general permittivity of the material.
6) Applies a rotation to the permittivity to account for the orientation of the crystal with respect to the waveguide direction.
7) Sets up an (n,k) import structure with spatially-varying diagonalize refractive index over space to represent the BaTiO3 layer.
8) Sets up a matrix transform grid attribute which applies the matrix transform and rotation of the permittivity tensor of BaTiO3.
Coordinates and rotation details
From the reference [1], the lattice constants of BaTiO3 are
a = 3.992 Å
c = 4.036 Å
Theta in the above diagram can be given by ϴ = atan(4.036/3.992) radians.
To rotate the permittivity of the crystal so that the [011] direction is along y, first rotate pi/2 counter clockwise around the y’-axis, then pi/2-theta clockwise around the new x’-axis as illustrated below.
To calculate the electric fields in the prime coordinates, the following relations can be used:
$$\begin{array}{l}{E_{x^{\prime}}=-E_{z}} \\ {E_{y^{\prime}}=E_{y} \sin (\theta)-E_{x} \cos (\theta)} \\ {E_{z^{\prime}}=E_{y} \cos (\theta)+E_{x} \sin (\theta)}\end{array}$$
Analysis and result
The ferroelectric_neff_vs_voltage.lsf script file loops through each bias voltage point to extract the effective index over frequency for the bias voltage. The final result is a plot of the effective index versus bias voltage for the desired operating frequency.
Note that additional convergence testing of the mesh step size over the BaTiO3 layer and the spatial resolution of the imported E field profile within the BaTiO3 layer may need to be performed to obtain higher accuracy results.
Related Publications
[1] Chi Xiong, Wolfram H. P. Pernice, Joseph H. Ngai, James W. Reiner, Divine Kumah, Frederick J. Walker, Charles H. Ahn, and Hong X. Tang, “Active Silicon Integrated Nanophotonics: Ferroelectric BaTiO3 Devices”, Nano Letters 2014 14 (3), 1419-1425
DOI: 10.1021/nl404513p
[2] M. Zgonik, P. Bernasconi, M.Duelli, R. Schlesser, P. Gunter, M.H. Garrett, D. Rytz, Y. Zhu, X. Wu, "Dielectric, elastic, piezoelectric, electro-optic, and elasto-optic tensors of BaTiO3 crystals", Phys Rev B Conens Matter. 1994 Sep 1;50(9):5941-594/.
DOI: 10.1103/PhysRevB.50.5941