Modeling of the
Semiconductor Laser with Optical Feedback
PoJen
Huang, Department of Mathematics and Physics, University of Arizona, Arizona
85721
Introduction to the Lang-Kabuysashi
equation and its application
The goal of this semester is to understand the behavior of the semiconductor laser system with feedback. We use Lang-Kabuysashi equation to describe the system:
|
|
…(1) |
|
|
|
Where 
is the line
width enhancement factor, 
is the laser
frequency without feedback, 
is the photon
lifetime, 
is the carrier
life time, and 
is the injected
current density. The feedback parameter 
measure the
intensity of the light reflected back in to the laser cavity, and the delay
time 
, which is the round trip time of the light in the external
cavity of length 
.
The modal gain per unit time is given by 
, where 
is the gain
constant and 
is the carrier
density at transparency. The modal gain 
includes the
linear gain and intensity reduction of the fain due to spatial and spectral
hole burning and carrier heating, with 
being the gain
saturation coefficient. Also the electric field is normalized so that 
is the total
photon number in the laser wave-guide, where 
is the volume of
the active region. The parameter 
is the lasing
threshold current density of a solitary laser and 
is the threshold
currier density. The typical values for the above parameters are given in the
table attached to the end of this report.
For
easier numerical convention, we measure the time in the units of the photon
lifetime 
, and normalize excess carrier number density 
. Then we can rewrite the Lang Kobayashi equation in the
normalized form:
|
|
…(1) |
|
|
|
Where 
, 
, 
, and 
. We start by looking at some special solutions and analysis
their stability, but even before we start, we should first review some of the tools
that we will use later in the paper.
Stability: The sign of the real part of the eigenvalue determines the stability of the dynamical system. Positive real part indicated an unstable solution, and negative real part indicated stable solution. (one can think the real part of the eigenvalue as the stretch/compress of the system) However, solving the characteristic polynomial to find the eigenvalue is often hard. Since the analytic solution is almost impossible to obtain, one must develop some numerical method to find the real part of the eigenvalue.
Numerical Method: To find the engenvluae in the complex plan is nontrivial. Perhaps the easiest way to find the root for a complex function is to calculate the residue. Recall Cauchy’s formula:
|
|
…(1) |
Where 
is the residue
of 
at
and 
is the winding
number of 
around 
. (since we always take simple loops, the winding number is
always one)
If we let 
, then at the singularity of 
is where the
root of the characteristic polynomial, ie. The value of the eigenvalue.
Since this method is highly depending on the contour we take, there is a potential danger of applying this method. If the contour contains more then one pole, the integral might be zero & we will miss that root. Hence it is very important to find the right size of the contour. (such that only one pole is allowed in one contour) We will illustrate how such analysis can be done in the simple one-dimensional case.
Analytical analysis for the steady state
We want to analysis stability of the solution of Lang-Kabuysashi. This is generally hard, so we start by working the system without feedback.
Case 1: Laser is off with no feedback: 
|
|
…(3) |
|
|
…(3) |
Stability- negative eigenvalue
Consider a small perturbation from steady-state solutions
|
|
…(3) |
|
|
…(3) |
Where 
Let
|
|
…(3) |
|
|
|
|
|
|
Hence in terms of a matrix:
|
|
|
Solve for the eigenvalue we find
|
|
|
To yield
|
|
|
|
|
|
|
|
|
For the laser to remain off, 
must be less
then zero, i.e. for the off state to be stable
|
|
|
|
|
|
Hence, the threshold value of J
|
|
|
This is 
, the laser will be off, for 
, any small perturbation will lead to the 
filed becoming
and staying nonzero, i.e. laser turns on.
CASE 2: Constant electric field without
feedback: 
|
|
|
|
|
…(3) |
|
|
|
|
|
…(3) |
Let’s consider a small perturbation from steady-state solutions
|
|
…(3) |
|
|
|
Hence
|
|
…(3) |
|
|
|
|
|
|
|
|
|
|
|
|
Hence in terms of a matrix:
|
|
|
Solving for the eigenvalue we find
|
|
|
To yield
|
|
|
|
|
|
|
|
|
It should be understood that Lang-Kabuysashi has a nontrivial behavior even with the simplest condition. To analysis the system better, we will develop our skill by studying the simplest case.
The One Dimensional Case
Consider a simple delay differential equation
|
|
…(1) |
Where 
is the change of
at time
,
is the delay time[1],
and 
is the value of 
at time 
.
It is easier to view the problem if we rescale the delay time to be one, i.e. define
|
|
…(2) |
Then from eq(1) we find
|
|
…(3) |
Or[2]
|
|
…(4) |
Like the real Lang-Kobayashi equation, it is also difficult to construct the analytic solution, even in the simplest one-dimensional case. The first few solution are below. (with the unit input step function)
|
Time |
Analytic Function |
Value at right hand endpoint |
Value at the left hand endpoint |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
A more compact way of express the solution is via Laplace transformation[3]:
N/A
The characteristic polynomial is:
|
|
…(5) |
Since 
is a complex number,
we can write it as 
, 
. Then eq(3) becomes:
|
|
|
Let’s look at some solutions for different
value of 
.
|
Case I |
|
|
|
|
|
|
||
|
|
|
||
Figures: The top left is the solution of the
one-dimensional differential equation with 
, integrate from one to twenty. Top right is the absolute log
plot of the top left figure. The slope is determined by the envelop (the
maximal value of each hump). Bottom right is the zero contour of the real and
imagery part of the characteristic function. Red indicates the real part and
blue indicated the imagery part. Notice that the imagery part (the vertical direction)
crossing of two contours are evenly spaced. (this is important later when we
try to drive the contour size for auto root finder) Bottom left shows the detail
around the origin. They cross on the left plan, i.e. negative real part of
eignevalue)
|
Case II |
|
|
|
|
|
|
|||
|
|
|
|||
Figures: The top left is the solution of the
one-dimensional differential equation with 
, integrate from one to twenty. Top right is the
absolute log plot of the top left figure. The slope is determined by the envelop
(the maximal value of each hump). Bottom right is the zero contour of the real
and imagery part of the characteristic function. Red indicates the real part and
blue indicated the imagery part. Notice that the imagery part (the vertical direction)
crossing of two contours are evenly spaced. Bottom left shows the detail around
the origin. The crossing is on the imagery axis, i.e. the real part of
eignevalue is identically zero)
|
Case III |
|
|
|
|
|
|
|||
|
|
|
|||
Figures: The top left is the solution of the
one-dimensional differential equation with 
, integrate from one to twenty. Top right is the
absolute log plot of the top left figure. The slope is determined by the envelop
(the maximal value of each hump). Bottom right is the zero contour of the real
and imagery part of the characteristic function. Red indicates the real part and
blue indicated the imagery part. Notice that the imagery part (the vertical direction)
crossing of two contours are evenly spaced. Bottom left shows the detail around
the origin. They cross on the right plan, i.e. positive real part of
eignevalue)
Notice the following two things, (base on observing the three cases):
1.) That the slope of the log plots is the same as the maximum of the real part of the eigenvalue.
2.)
The imagery part of the eigenvalue is
evenly spread, and the distance in between is about 
.
Observation one is not suppressing, since the ratio of the stretching/ compress is determined by the real part of the eigenvalue. Recall the ratio of stretching/ compressing is given by the formula:
The final solution can be written as:
the 
is the envelop
of the solution and 
is the oscillated
solution.
PFPFPFPFPF
PF
Adams-Bashford-Moulton
(ABM) PC Integration Method
We can rewrite the differential equation in
terms of discrete map of dimension 
by introducing
the step size to be
|
|
|
Hence the system becomes:
|
|
|
|
|
|
|
|
|
Since the L-K is extremely depend on the integration method we choice. We decide to use the fourth-order Adams-Bashford-Moulton (ABM) PC method.
|
Predictor: |
|
|
|
Corrector: |
|
|
Where 
Logarithm
derivative
Phase Plot:
|
|
|
Let’s rewrite the complex electric field,
|
|
|
Since the constant phase shift does not affect the dynamics of the system.
Consider the solution has the form:
|
|
|
The arctangent only give the value up to pi. (depend on the slice you are on), however the dynamics of the system changes as time varies. We need a way to describe
|
|
|
Take the ratio we yield
|
|
|
Provided that the electric field is non-zero. (if the field is zero, then the phase is undefined)
|
|
|
Integrate both sides we have
|
|
|
Then do a little shift on the argument to get the interval we are interested.
|
|
|
PHASE DELAY
Bifurcation diagram and Sisyphus effect
N/A