Journal of Modeling and Simulation of Microsystems, Vol. 1, No. 1, Pages 1-8, 1999.
MODELING SEMICONDUCTOR OPTICAL DEVICES G. E. Sartoris† Zürcher Hochschule Winterthur, 8401 Winterthur, Switzerland
J. Leuthold Bell Labs, Lucent Technologies, Holmde1 NJ 07733, U.S.A.
Abstract In this article numerical methods used to model photonic devices are presented. In particular, the governing equations describing the problem, the coupling between electrical and optical fields, and the numerical algorithms used for computing solutions are discussed. As a computational example, we present the modeling of a semi-conductor optical amplifier employed as an all-optical demultiplexer. Keywords: semiconductor, optical devices, modeling
1. Introduction
investigations.
Semiconductor modeling plays an important role in the industrial development cycle. As a result, several modeling tools are nowadays available to model semiconductor devices, however, for photonic devices this is generally not the case. Peculiar to these devices is a physical coupling between electrical and optical fields which plays an important role in the device behavior. Modeling photonic devices is in general a rather complex task. There is not a global approach, different devices may require different tools, and both electrical and optical fields must be computed together with a correct characterization of the coupling. In this work, we address this problem and a method for computing a certain class of photonic devices is presented as implemented in the program NM Seses (SEmiconductor Sensor and actuator Simulation).
As a computational example, we present the modeling of an all-optical switch. The device consists of a Semiconductor Optical Amplifier (SOA) built into the arms of Mach-Zehnder interferometer to induce the necessary change of the refractive index needed for switching. The SOA device consists of a planar dielectric rib waveguide and a three layer hetero p-i-n laser diode as depicted in Fig. 1. The upper and lower layers are made of p- and n-doped InP, and the middle layer is made of a quaternary intrinsic 1.55 µ m -InGaAsP lattice matched to InP. The figure only displays a small section of the real device which is actually much longer and is only a single component in a monolithically integrated all-optical switch. Because the lateral dimensions of the laser device are small compared with the longitudinal length, the laser device may be well characterized by studying a single two-dimensional section.
The electrical behavior of a semiconductor device may be well characterized by computing the electrical fields with the classical drift-diffusion model. Here, the electric potential and the electron and hole densities are computed by solving a system of three partial differential equations. More complex is the approximation and computation of the optical fields since they strongly depend on the class of problem to be modeled. A standard solution is generally not available and for different classes of devices, different numerical tools and algorithms must be implemented. In the present work, the optical fields will be computed as the quasi-TE or quasi-TM modes of twodimensional dielectric waveguides. These modes are obtained by solving the scalar Helmholtz problem for one of the transverse components of the electric or the magnetic field. Once the governing equations for the electrical and the optical fields are specified, it is possible to define the coupling between these fields quite generally with the introduction of recombination models, as for example when stimulated recombination or absorption of light is present. This step is generally straightforward, however, to specify the model’s parameters we had to carry out complex experimental
p-doped InP InGaAsP n-doped InP
Fig. 1: A semiconductor laser device forming a rib-waveguide and composed of three layers: (top) p-doped InP, (middle) intrinsic 1.55 µm -InGaAsP, (bottom) n-doped InP. Electrical contacts are deposited on the top and bottom of the device. For the steady-state solutions presented in the following, the laser diode is forward biased by 1.2 V.
The layout of this article is as follows. In sections 2 and 3, we review the governing equations for the electrical fields and the
†.
[email protected]
Computational Publications, ISSN 1524-2021.
JMSM, Vol. 1, No. 1, Pages 2-8, 1999. optical fields. In section 4 we describe the coupling between electrical and optical fields for the SOA device and the experimental set-up used to measure the material gain, a fundamental parameter of the coupling. In section 5 some issues related to the numerical algorithms used to compute solutions are presented. Finally in section 6 we present some modeling aspects and results of the SOA device.
2. Governing Equations For Electrical Fields We use here the classical unstationary drift-diffusion model [1, 2] as the governing equations for the electric potential ψ and the electron n and hole p densities. The model consists of the three partial differential equations:
result, the carrier densities n, p are also generally discontinuous at the interfaces between different materials. The computation of such discontinuous fields should be avoided. The implemented finite element method only computes continuous solutions and a very large number of elements would be necessary to approximate discontinuities. To circumvent this problem, the discretization of the driftdiffusion model is therefore performed with respect to the continuous variables ( ψ, n/n in , p/n ip ) . For the numerical results of this model, Fig. 2-3 show the potential field and electron density as computed for the device of Fig. 1. Psi [V] 1.5E-02 -2.5E-03 -2.0E-02
∇ ( ε ∇ ψ ) = q ( n – p – C ),
-3.8E-02 -5.6E-02
ψ /V – ψ /V T ∂n n ⁄ n in ), (1) – q ----- + ∇J n = qR ( n, p ), J n = qD n n in e T ∇ ( e ∂t – ψ /V T ψ /V ∂p q ------ + ∇J p = – q R ( n, p ), J p = – q D p n ip e ∇ ( e T p/n ip ). ∂t
-7.3E-02 -9.1E-02 -1.1E-01 -1.3E-01 -1.4E-01 -1.6E-01 -1.8E-01
The first equation of (1) is the Poisson equation which couples together the electric potential ψ with the charge density ρ = – q ( n – p – C ) composed of the electron and hole carrier densities n, p, and the bulk charge C . q represents the elementary charge and ε is the dielectric permittivity.
-2.0E-01
The second and third equations of (1) are the current continuity equations which describe the flow of the carrier densities n, p inside the semiconductor. The recombination R ( n, p ) accounts for the creation and annihilation of electronhole pairs inside the semiconductor. For a non-degenerate semiconductor, the diffusion coefficients D n, p are related to the mobility µ n, p by D n, p = V T µ n, p , where V T = k B T /q is the thermal voltage, k B the Boltzmann constant and T the temperature that is assumed to be constant throughout the device.
-3.2E-01
Because most semiconductor optical devices are made of heterogeneous materials, the modeling process starts with a correct characterization of the bandshapes at thermodynamic equilibrium when the currents J n, p are null. This is achieved by defining the intrinsic density n i to be a space-dependent entity. However, if the same value is defined for both electrons and holes, then the bandgap's center will always be aligned at the same energy level. Since this is not true for hetero structures, we allow the specification of a shift δ E i of the bandgap's center or equivalently of the intrinsic Fermi level E i . The electron and hole intrinsic densities n in, n ip in the model (1) are then given by n in = n i e
δ Ei/(kBT )
, n ip = n i e
–δ E i / ( k B T )
-2.1E-01 -2.3E-01 -2.5E-01 -2.7E-01 -2.9E-01 -3.0E-01
-3.4E-01
NM SESES
Fig. 2: The potential field distribution for one-half of the laser device of Fig. 1. .
N [1/cm**3] 4.6E+18 4.3E+18 4.1E+18 3.9E+18 3.7E+18 3.4E+18 3.2E+18 3.0E+18 2.7E+18 2.5E+18 2.3E+18 2.1E+18 1.8E+18 1.6E+18 1.4E+18 1.1E+18
.
(2)
9.2E+17 6.9E+17
Since the intrinsic densities n in, n ip undergo abrupt variations at the interfaces between different materials, large gradients for the potential ψ and carrier densities n = n in e ψ /V and p = n ip e –ψ /V are to be expected there, see [3, 4]. In our implementation of the drift-diffusion model, we further assume that in each material the intrinsic carrier densities are constant; i.e., they are piecewise constant functions. As a
4.6E+17 2.3E+17 1.8E-02
T
T
2
NM SESES
Fig. 3: The electron density for one-half of the laser device of Fig. 1.
Modeling Semiconductor Optical Devices Although for heterogeneous materials and non-degenerated semiconductors the Poisson equation yields a correct representation of the bandshapes, the current continuity equations based on the linearization of the Boltzmann transport equation may not yield correct current-voltage characteristics. At the interface between two semiconductors with different bandgaps, strong carrier-carrier interactions are present and only Monte Carlo methods are generally able to properly model the physics. A thermoionic emission model [5, 6] has been used in [7, 8] to model the current flow over hetero-junctions. Here, it is assumed that at the interface between different materials the quasi-Fermi levels are discontinuous. Since the discontinuity jump is an unknown parameter, the difficulty of this approach lies in its implementation. A numerical tool capable of defining doublevalued mesh nodes must be used; not a simple task to build into a large software package. For this reason, this approach has not been applied and therefore I-V characteristics may not be modeled accurately.
∂E y (4) yield the continuities of E y and -------∂x for TE modes and the H -y continuities of H y and n----12- ∂--------for TM-modes. These ∂x continuity requirements at the dielectric interfaces are implicitly defined when solving the following Helmholtz problem for the field components E y , H y on the domain
Ω = {x ∈
}
2 2πn 2 2 ∇ + ---------- – ( β ) E y ( x ) = 0, x ∈ Ω, λ (5) 1 2 π 2 β 2 TM modes: ∇ ----2-∇ + ------ – --- H y ( x ) = 0, x ∈ Ω. λ n n
TE modes:
For a piecewise constant refractive index n = n ( x ) , the solutions of these Helmholtz problems are exact solutions of the Maxwell equations. There are in general a finite number of bounded solutions, each one representing a single propagating mode with velocity β . TEmode []
3. Governing Equations For Optical Fields
1.8E+00
In the present work, we consider photonic devices where the optical fields can be computed as quasi-TE or quasi-TM modes of optical two-dimensional dielectric waveguides propagating along the z-axis. These modes are obtained by solving the scalar Helmholtz problem for the y-components of the electric and magnetic optical fields E y and H y .
1.7E+00 1.7E+00 1.6E+00 1.5E+00 1.4E+00 1.3E+00 1.2E+00 1.1E+00
Let us first consider TE and TM modes for a dielectric waveguide with µ = µ 0, ε = n 2 ε 0 ,n the refractive index, ∂/ ∂y = ∂/ ∂z = 0 , and the z-axis as the propagation direction, see [9]. For TE modes the non-zero components of the electric and magnetic fields E , H are E y, H x, H z and for TM modes H y, E x, E z . For an analysis with a single frequency, we use the dependency ( . ) ( z, t ) = ( . )e i ( wt – β z ) with w = 2 π c/ λ, λ the vacuum wavelength and β the (yet unknown) propagation velocity of the mode in the waveguide. If we assume the refractive index is piecewise constant, then from the Maxwell equations with zero source terms we obtain the following exact equations for the field component E y and H y inside each layer i of constant refractive index n i
TE modes:
2
E y = E y ( x )e
2
H y = H y ( x )e
8.2E-01 7.3E-01 6.4E-01 5.5E-01 4.6E-01 3.6E-01 2.7E-01 1.8E-01 8.8E-02 -4.2E-03
NM SESES
, (3)
In practice the refractive index of optical waveguides cannot be invariant along the y-axis since the optical mode must be limited laterally in some way. If the variations of the refractive index along the y-axis are small compared with the ones along the x-axis, then the optical modes will preserve the TE and TM character and the solution of the Helmholtz problem (5) with a refractive index n = n ( x, y ) will result in approximate solutions of the Maxwell equations. These solutions are called quasi-TE and quasi-TM modes.
(4)
Although for this problem it is natural to seek solutions of the equations (5) in the unbounded domain Ω = 2 , numerical solutions are generally only computed on bounded domains. This fact raises the question of suitable boundary conditions. For guided modes, the field amplitude decays exponentially to zero at infinity, therefore homogeneous Dirichlet boundary
2 π n-i – β H ( x ) = 0, ∇ + ---------y λ 2
TM modes:
i ( wt – β z )
9.2E-01
Fig. 4: The ground quasi-TE mode of the dielectric waveguide defined by the laser device of Fig. 1. Only one half of the device is shown.
π n-i – β 2 E ( x ) = 0, ∇ + 2---------y λ 2
1.0E+00
2
i ( wt – β z )
.
The other components are derived from the relations TE modes:
β H x = – -------E y , wµ
i ∂E H z = ------- --------y- , w µ ∂x
TM modes:
β E x = ------- H y , wε
i ∂H E z = ------- ---------y- . w ε ∂x
The continuity requirements for the tangential component of the fields E and H at the dielectric interfaces together with
3
JMSM, Vol. 1, No. 1, Pages 4-8, 1999. conditions are a good approximation but only if the domain is sufficiently large. Homogeneous Neumann boundary conditions are also of interest on the axes of symmetry. Fig. 4 displays the ground quasi-TE optical mode for the laser device of Fig. 1. Now that the modes are limited laterally the amount of energy transported is finite and given by 1 TE modes: Energy Flux = – --2 1 TM modes: Energy Flux = --2
*
∫Ω E y H x dµ
β 2 = ------------- ∫ E y dµ , 2w µ 0 Ω 2
H β * - -------y dµ ∫Ω H y E x dµ = ----------2wε 0 ∫Ω n 2
(6)
The reason why the waveguide modes are computed as quasiTE or quasi-TM modes is on the simplicity of the scalar formulation. Solving the Maxwell equations with zero source terms is more complex because of the intrinsic vector character. Here, spurious modes may appear if inadequate discretization methods are used.
4. Coupling Between Optical And Electrical Fields Different classes of devices have in general different types of coupling between optical and electrical fields, however, for a large class of devices this coupling can be given by defining recombination models. In particular, for the semiconductor laser in Fig. 1, the major source of coupling is determined by stimulated emission of light given by the relation
In order to perform reliable simulations, material parameters for the gain function (8) are needed. Due to the lack of available data, we have experimentally determined temperature-corrected material gains of bulk 1.55 µmInGaAsP from very low up to below threshold carrier densities. Two double-heterostructure laser diodes have been used, one with straight cleaved facets and the second one tilted by 10o versus the (100) direction. This allowed us to use the second diode as a semiconductor optical amplifier without applying an antireflection coating since tilted facets efficiently reduce facet reflectivities to below 10-4. The modal gain g averaged over the whole device was determined experimentally with the Hakki-Paoli [10] and Cassidy [11] method and the internal losses α int were determined with the transparency method [12]-[14]. From the knowledge of the mirror losses αm = -log(R1R2)/2L with L the device's length and R1, R2 the facet's reflectivities and the mode confinement factor Γ, the material gain can be obtained with g = ( g + αint + αm)/Γ. The current densities were related to the total injected currents by determination of the differential lifetime of the carriers in the active region [15]-[16]. For a detailed description of the experiment see [17]. Recomb [1/cm**3] 6.0E+27 5.7E+27 5.4E+27 5.1E+27 4.8E+27 4.5E+27 4.2E+27
R
stim
= gain × photon flux = gain × energy flux/ w (7)
3.9E+27 3.6E+27 3.3E+27
where gain = g ( n, p ) is the material gain. Different gain models are available, which are in general given as polynomials in the variables n and p with coefficients fitted from experimental data. For example, in our implementation we used the models 2
2.7E+27 2.4E+27 2.1E+27 1.8E+27 1.5E+27 1.2E+27 9.0E+26
g 1 ( n, p ) = A 0 z + A 1 z , 2
g 2 ( n, p ) = g 1 ( n, p ) – ( C 0 + C 1 z ) ( λ – ( λ 0 – B 1 z ) ) ,
(8)
with z = np – N 0 the carrier density, N 0 the transparency carrier density, and λ 0 the bandgap-wavelength at zero carrier density. The first model is used when working with a single optical wavelength, the second one, for example, when working with spontaneous emission of light. The energy flux is evaluated using the relation energy flux = const Φ2 where Φ is the computed and normed optical amplitude proportional to Ey or Hy and const is a variable defining the power of the optical beam. This stimulated recombination is then added to other recombination models as the Schottky-Read-Hall RSRH, the spontaneous emission RSpontan and Auger recombination RAuger models. The total recombination rate for the laser device in Fig. 1 forward biased by 1.2 V and for an optical beam of 10 mW is given in Fig. 5. As for the carrier densities, the recombination rate is also a discontinuous function at the interface between different materials. 4
3.0E+27
6.0E+26 3.0E+26 0.0E+00
NM SESES
Fig. 5: Recombination rate for one-half of the laser device of Fig. 1.
The material gain for the quaternary 1.55 µm layer is displayed in Fig. 6 for a layer temperature of 20oC. Since the diode's threshold currents at room temperature were 34 mA and 100 mA, to derive the material gain we used the first diode device for current biases between 10 mA and 30 mA and the semiconductor optical amplifier device for currents above 30 mA. The depicted curves are 9th order polynomial least square fits. We also display the gain curve for a current bias of 110 mA. However, as shown by the figure, current biases above 100 mA do not provide higher gains, since the material gain is clamped; i.e., the device starts lasing and the carrier population remains constant.
Modeling Semiconductor Optical Devices
Fig. 6: Material gain spectrum for a 1.55 µm quaternary layer as a function of wavelength and carrier injection densities. The active layer temperature is 20oC. The curves are derived from measurements on the laser diode device (dashed curves) and the semiconductor optical amplifier (solid lines).
5. Numerical Algorithms Once the governing equations for the optical and electrical fields are defined and the type of coupling has been specified, numerical methods are applied to solve all of the equations consistently. In a first step, one has to consider the type and magnitude of the coupling between the physical fields. For the laser example in Fig. 1 which is our modeling goal, the dynamical behavior when a beam of light traverses the laser device is of primary interest. Since the lateral dimensions of the device are small compared with the longitudinal length, we may perform the numerical computations on a two-dimensional section only. There are situations, however, where this approach may not lead to correct results and one has to perform truly 3D simulations. In this case, as an approximation one can define the optical fields in the 3D domain to be the quasi-TE and quasi-TM modes multiplied by a complex envelope function specifying the phase and the energy flux as a function of the z-coordinate. The change of amplitude and phase is then specified as a function of the averaged gain and alpha-factor, see for example [17]. If the optical signal is traveling in one single direction, as a further approximation one may compute on a 2D section the averaged values of the gain and alpha-factor, thus allowing to evaluate the envelope function a small step further along the z-direction. By repeating this process we may end up computing the whole device. Clearly, this simple approach may not always be possible or even meaningful. As noted previously, a general and simple modeling approach is not available, therefore we will not pursue this argument further and we will limit ourselves to the 2D computations on a z-section of the device. To compute a z-section of the device we proceed as follows. After specifying the power of the optical beam, the scalar Helmholtz equations for TE or TM-modes are solved to obtain the normed mode amplitude. The drift-diffusion model is then solved to obtain the electrical fields where the mode amplitude is used to specify the amount of stimulated recombination. Since the optical modes are computed for a device at rest; i.e., without an optical beam, this solution is
not really consistent. Stimulated recombination affects the density of the carriers which in turn affects the refractive index of the dielectric waveguide and thus the optical modes. It is possible to obtain a consistent solution with a Pickard iteration method, thus by restarting the whole solution process with the last computed field values. Since this type of coupling is weak, the convergence is assured in a few iterations and in most cases no iterations at all are required. This is even more true for hetero junction devices where the refractive index differences between film and cladding layers are large. With a weak coupling mechanism from the electrical to the optical fields as discussed above, the modeling of this laser device may be characterized as one directional. First one can compute the optical fields and then the electrical fields. This is true for the presented example, but if one uses the laser diode as a self-sustained laser then the coupling becomes bidirectional. Computing solutions for this case is typically more difficult and instabilities related to the existence of multiple solutions may occur, for example, when different optical modes are competing. Next the numerical algorithms used to solve the governing equations are briefly presented. The scalar Helmholtz equations are discretized with standard rectangular Q1 finite elements [18]. To find the first few eigenpairs (λ, x ) of the generalized eigenproblem A x = λM x for the positive definite matrices A and M , the Lanczos algorithm is used to compute the Krylov subspace consisting of orthogonal vectors. In this subspace, a QR-algorithm with shift is used to diagonalize the tridiagonal matrix obtained by the Lanczos algorithm yielding the approximate eigenpairs to the original problem [19]. A partial reorthogonalization algorithm is used to keep the Lanczos vectors almost orthogonal to each other [20]. The three equations of the drift-diffusion models are one of type elliptic and two of type parabolic. Following the method of lines we first discretize the space variables and deal with time as a continuous variable. Each single equation of the drift-diffusion model is discretized in space by a mixed finite element method and Q1 elements [21, 22]. For a successful implementation of this method, it is mandatory to have an exact and fast numerical evaluation of exponentially fitted integrals. These integrals arise when discretizing the electron and hole current continuity equations, see also [23] for further details on the numerical algorithms. After the discretization of the spatial variables, a system of ordinary differential equations for the time variable is obtained. Because the system can be very stiff, general Astable methods as the standard second order trapezoidal integration algorithm are not well suited and may lead to unwanted oscillations. Better suited are the so called L-stable methods, like the BDF2 method which has been used together with an adaptive time step selection algorithm [24, 25]. A composite trapezoidal-BDF2 method is also a possible choice and offers some advantages, see [26]. All linear systems arising during the solution process have been solved by LU-factorization and the minimum degree algorithm has been used to minimize the amount of fill-in during the factorization process [27]. As shown in Figs. 2-5 for steady state solutions of the laser device of Fig. 1, as well as for the unstationary solutions 5
JMSM, Vol. 1, No. 1, Pages 6-8, 1999. computed, the numerical solutions do not show spurious oscillations and this in spite of the fact the carrier densities are discontinuous and large gradients are presented in the potential field at the interface between different materials. This good numerical behavior is also reflected by a good convergence behavior of the non-linear Newton-Raphson solution algorithm used to solve the discretized and linearized drift-diffusion model.
caused by stimulated emission induces a change of the refractive index. With a π−phase shift the data-signal will appear at the output Px. The refractive index and gain coefficient of the waveguide are related by the KramersKrönig relations. For a waveguide, this dependency is expressed by the alpha-factor
6. Computational Example
with λ the vacuum wavelength, N eff = --βc- the effective refractive index, β the mode velocity, Φ the mode amplitude and g the averaged material gain
The aim of the presented example is to compute the dynamical behavior of a monolithically integrated all-optical switch. The device realized as an InGaAsP/InP-waveguide ridge structure [28] is depicted in Fig. 7. The device comprises two asymmetric multimode interferometer (MMI) couplers with 60:40 and 40:50 splitting ratios. The asymmetry is used to guarantee high extinction ratios [29]. Two two-section semiconductor optical amplifiers (SOA) and the Mach-Zehnder interferometer (MZI) arms provide the necessary non-linearities for switching the optical signal. The SOAs are divided into two sections to allow for gain and small offset phase tuning. The total length of the SOAs is 1 mm. Two MMI-converter-combiners (CC) are used to introduce the control signals Pc1,2 as a first-order mode into the SOAs and the other two are used to extract the first order mode control-signal from the signal path of the data-signal [30]. Pc1
CC
SOA1
CC
Pin 60:40 MMI Pc2
CC
Pc1 P=
40:60 MMI SOA2
CC
Px Pc2
Fig. 7: Monolithically integrated all-optical MZI-SOA switch with asymmetric MMI-splitters.
Control signal
Data signal
Fig. 8: Schematic view of control and data signals.
A data-signal Pin presented at the input of the switch is split asymmetrically in both arms of the MZI. Without a controlsignal, the split signals undergo the same phase shift in both MZI's arms, therefore reappearing at the output P=. For this device to work as an all-optical switch, a control signal Pc1 is used to change the refractive index in the MZI's upper arm. Here the reduction of carriers inside the SOA's active region 6
4 π ∂N eff α = – ------ -----------λ ∂g
(9)
2
∫Φ
2
Φ ∫ ------2- dxdy = 1. n
TE modes: g = ∫ g ( n, p )Φ dxdy and Φ TM modes: g = ∫ g ( n, p ) ------2- dxdy and n
2 2
dxdy = 1, (10)
2π∆N eff L
- we obtain the phase shift in With the relation ∆ φ = ------------------λ the arms of the MZI as ∆ φ = α∆gL/2 . For this device type, measured values [17] of the alpha-factor α are in the range (5 - 8), therefore for L = 1 mm we need a value of ∆g ≈ 10cm –1 to obtain a π−phase-shift.
By assuming a control signal much more intense than the data signal, the latter will not perturb the state of the upper SOA device too much, so that the switching behavior will uniquely depend on the refractive index's change caused by the control signal. Immediately upon arrival of the signal, carriers are depleted from the SOA's active region due to stimulated emission. The decay of the carrier population has a faster exponential decay if compared to the reverse process when the control signal Pc1 is turned-off. Here the carriers must be injected into the SOA's active region to restore the original state and thus the refractive index in the MZI's upper arm. In order to optimize the device, we have to study the slower recovery process. This can be done, by performing time dependent simulations of the laser device at the input facet. For this goal, we have considered a device similar to the one depicted in the Fig. 1. Between the active 1.55 µm-InGaAsP and the n-, p-InP layers, there are two intrinsic quaternary 1.3 µm-InGaAsP layers. The wavelength of the control signal is λ = 1.53 µm Since for switching a fixed ∆g value is needed, the device will perform better for large gain values. However, there is a maximum here, for large applied bias voltages the gain saturates and the device starts lasing, an undesirable effect. For our simulations we have therefore fixed the carrier density to a value below the laser threshold. In a first run, we have investigated the speed of the recovery process. From the steadystate solution with a CW control signal Pc1 at the input, at time t = 0 we have turned off the signal. The gain recovery to the steady-state value is exponential and Table 1 shows a faster recovery for larger signals. Let's now fix the control signal's power at 2 mW . The results from Table 1 show that to optimize our device we have to focus the optical signal to achieve a high photon density. For this goal there is an optimum in the thickness of the active layer as displayed in the Table 2. For very thin layers there is a degradation in the confinement factor of the optical mode in the active layer. It should be noted that reducing the active layer thickness
Modeling Semiconductor Optical Devices
Table 1: Averaged gain at time t = 0, t = ∞ , and time derivative at t = 0 as a function of the control signal power. The CW signal is turned-off at the time t = 0 . Power [mW]
10 – 12 × d g ⁄ dt [cm-1s-1]
g ( 0 ) [cm-1]
g ( ∞ ) [cm-1]
0.1
0.012
368.883
371.798
1
0.144
339.357
371.798
10
1.107
123.809
371.798
Table 2: Averaged gain difference for TE and TM modes and a 2 mW control signal turned on and off as function of the active layer thickness. Thickness [µm]
TE: g off – g on [cm-1]
TM: g off – g on [cm-1]
0.2
30.243
26.561
0.3
50.551
47.789
0.4
52.281
50.941
0.5
42.673
41.793
can use in the lower arm of the MZI a control signal Pc2 of the same shape as Pc1 but time delayed and eventually less powerful. This second signal will quickly restore the total phase shift in the MZI, so that a well defined switching window is created. The averaged gain difference between both arms g c1 – g c2 is displayed in Fig. 10, showing a switching window of ≈ 15 ps for a time delay of ∆t = 15 ps . The device modeled in this section has been realized [28] and measurements show at 10 GHz a switching window of 25 ps for pulses with peak values of 24 mW and FWHM of 8 ps. We point out that the numerical computations were explicitly done using rectangular pulses to show the ideal device's performance. -4
-6
-8
Delta Gain [1/cm]
increases the carrier population, an effect partially counterbalancing the reduction of the mode confinement. Another modeling goal is the polarization independence of the device. Table 2 shows that this independence is better achieved for thicker active layers. What counts here is the difference in the phase shift between TE and TM modes, thus TE TM TM the difference ( g TE off – g on ) – ( g off – g on ) . Similar computations may be done by varying the thickness of the upper and lower cladding layers and the width of the rib.
-10
-12
-14
-16
-18
-20 240
342
260
280
300
320 Time [ps]
340
360
380
400
Fig. 10: Gain response's difference for the two 20 GHz and 15 ps delayed control signals with pulses of equal 8 ps width and peak powers of 20 mW and 18 mW . The control signals were started at t = 0.
340 338 336
Gain [1/cm]
334
7. Conclusions
332 330 328 326 324 322 240
260
280
300
320 Time [ps]
340
360
380
400
Fig. 9: Gain response to a 20 GHz control signal with rectangular pulses of 8 ps width and 20 mW peak power. The control signal was started at t = 0.
In order to show the complete dynamic response of the device, we have applied at t = 0 a 20 GHz control signal Pc1 with rectangular pulses of 8 ps width and 20 mW peak power. The change in the averaged gain is shown in Fig. 9. Assuming as before that a value of ∆g ≈ 10 cm-1 is needed for switching, Fig. 9 shows the capabilities of the device for switching an input signal between the ports P= and Px at a rate of 20 GHz, however, most of the time the switch is in an undefined on-off state. As a remedy, in addition to Pc1 one
Through an example of a semiconductor laser diode, numerical methods for modeling photonic devices have been presented. If the drift-diffusion model is to be considered a general tool for computing the electrical fields, the governing equations for the optical fields and the numerical methods used for their solutions must be chosen according to the class of optical devices being considered. For semiconductor optical amplifiers, the major source of coupling between electrical and optical fields is the stimulated emission of light, and experimental measurements have been carried out to obtain the model's parameters. Numerical solutions of the governing equations without spurious oscillations were obtained within the proposed numerical algorithms implemented with the program NM Seses. A semiconductor optical amplifier has been modeled and good agreement with a prototype device has been met.
Acknowledgment Ch. Zellweger and M. Mayer are acknowledged for help in experimental characterization of the laser and SOA devices. 7
JMSM, Vol. 1, No. 1, Pages 8-8, 1999.
References [1]
P.A. Markowich, The stationary semiconductor device equations, Springer-Verlag, 1986.
[2]
S. Selberherr, Analysis and simulation of semiconductor devices, Springer-Verlag, 1984.
[3]
G.P. Agrawal, N.K. Dutta, Long-wavelength semiconductor lasers, Van Nostrand, 1986.
[4]
H.C. Jr. Casey, M.B. Panish, Heterostructure lasers, Academic Press, 1978.
[5]
S.M. Sze, Physics of semiconductor devices, John Wiley & Sons, 1981.
[6]
K. Hess, Advanced theory of semiconductor devices, Prentice Hall, 1988.
[7]
M. Grupen, K. Hess, G.H. Song, Simulation of transport over Heterojunctions, Simulation of semiconductor devices and processes, Vol. 4., Hartung-Gorre, 1991.
[8]
S. Mottet, J.E. Viallet, Thermionic emission in Heterojunctions, Simulation of semiconductor devices and processes, Vol. 3., Tecnoprint, 1991.
[9]
A. Yariv, Quantum Electronics, John Wiley & Sons, 1989.
[10]
B.W. Hakki, T.L. Paoli, Gain spectra in GaAs doubleheterostructure injection lasers, J. of Appl. Physics, Vol. 46, No. 3, pp.1299-1306, 1975.
[11]
D.T. Cassidy, Technique for measurement of the gain spectra of semiconductor diode lasers, J. of Appl. Phys., Vol. 56, No 11, pp. 3096-3099, 1984.
[12]
K.A. Andrekson, N.A. Olsson, R.A. Logan, D.L. Coblenty, H. Temkin, Novel technique for determining internal loss of individual semiconductor lasers, Electron. Lett., Vol. 28, No 2, pp. 171- 172, 1992.
[13]
St. Pajarola, Dual-polarization external cavity diode laser for the optical generation of millimeter-waves, Ph.D. Thesis No. 12633, Swiss Federal Institute of Technology, 1997.
[14]
G.E. Shtengel, D.A. Ackerman, Internal optical loss measurements in 1.3 mm InGaAsP lasers, Electron Lett., Vol. 31, No 14, pp. 1157-1159, 1995.
[15]
W. Rideout, B. Yu, J. Lacourse, P.K. York, K.J. Beernink, J.J. Coleman, Measurement of the carrier dependence of differential gain refractive index, and linewidth-enhancement factor in strained-layer quantum-well lasers, Appl. Phys. Lett., Vol. 56, No 8, pp. 706-708, 1990.
[16]
[17] 8
P. Granestrand, K. Fröjdh, O. Sahlén, B. Stoltz, J. Wallin, Gain characteristics of QW Lasers, European Conference on Optical Communication (ECOC’98), pp. 431-432, Sept. 1998. J. Leuthold, Advanced InP waveguide Mach-Zehnder interferometer all-optical switches and wavelength
converters: Theory, design guidelines & experimental results, Hartung-Gorre Verlag, Konstanz, 1999. [18]
P.G. Ciarlet, Basic error estimates for elliptic problems, Handbook of Numerical Analysis by P.G. Ciarlet and J.L Lions, Vol. II, North-Holland, 1991.
[19]
G.H. Golub, C.F. van Loan, Matrix computations, Johns Hopkins University Press, 1993.
[20]
H.D. Simon, The Lanczos algorithm with partial reorthogonalization, Math. of Computation, Vol. 42, Num. 166, pp. 115-142, 1984.
[21]
J.E. Roberts, J.M. Thomas, Mixed and hybrid methods, Handbook of Numerical Analysis by P.G. Ciarlet and J.L Lions, Vol. II, North-Holland, 1991.
[22]
G. Sartoris, A 3D rectangular mixed finite element method to solve the stationary semiconductor equations. SIAM J. Sci. Comp., Vol. 19, No. 2, 1998.
[23]
G. Sartoris, A fast algorithm to compute the moments of the exponential function for application to semiconductor modeling, COMPEL, Vol. 16, No. 1, 1997.
[24]
E. Hairer, S.P. Norsett, G. Wanner, Solving ordinary differential equations I, Nonstiff problems, SpringerVerlag, 1993.
[25]
E. Hairer, G. Wanner, Solving ordinary differential equations II, Stiff and Differential-Algebraic problems, Springer-Verlag, 1996.
[26]
R.E. Bank, W.M. Coughran, W. Fichtner, E.H. Grosse, D.J. Rose, R.K. Smith, Transient simulation of silicon device and circuits, IEEE Trans. on Computer-aided Design, Vol. CAD-4, No. 4, 1985.
[27]
A. George, J.W.H. Liu, The evolution of the minimum degree ordering algorithm, SIAM Review, Vol. 31, No. 1, pp. 1-19, 1989.
[28]
J. Leuthold, P.A. Besse, E. Gamper, M. Dülk, ST. Fischer, H. Melchior, Cascadable dual-order mode alloptical switch with integrated data- and control- signal separators, Electronics Lett., Vol. 34, No. 16, pp. 1598-1600, 1998.
[29]
J. Leuthold, P.A. Besse, J. Eckner, E. Gamper, M. Dülk, H. Melchior, All-optical space switches with gain and principally ideal extinction ratios, IEEE J. of Quantum Electronics, Vol. 34, No. 4, pp. 622-633, April 1998.
[30]
J. Leuthold, J. Eckner, E. Gamper, P.A. Besse, H. Melchior, "Multimode interference couplers for the conversion and combining of zero- and first-order modes", J. of Lightwave Tech. Vol. 16, No. 7, pp. 1228-1239, July 1998. Received in Cambridge, MA, USA, 25th January 1999 Paper 1/01669
