Modeling of the
Semiconductor Laser with Optical Feedback
Introduction to the LangKabuysashi
equation and its application
The goal of this semester is to understand the behavior of the semiconductor laser system with feedback. We use LangKabuysashi 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 waveguide, 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 onedimensional case.
We want to analysis stability of the solution of LangKabuysashi. 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 steadystate 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: _{}
_{}, and _{} 

_{} 
…(3) 
_{} 

_{} 
…(3) 
Let’s consider a small perturbation from steadystate solutions
_{} 
…(3) 
_{} 

Hence
_{} 
…(3) 
_{} 

_{} 

_{} 

_{} 

Hence in terms of a matrix:
_{} 

Solving for the eigenvalue we find
_{} 

To yield
_{} 

_{} 

_{} 

It should be understood that LangKabuysashi 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 LangKobayashi equation, it is also difficult to construct the analytic solution, even in the simplest onedimensional 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 onedimensional 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 onedimensional 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 onedimensional 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
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 LK is extremely depend on the integration method we choice. We decide to use the fourthorder AdamsBashfordMoulton (ABM) PC method.
Predictor: 
_{} 

Corrector: 
_{} 

Where _{}
Phase Plot:
_{}, and_{} 

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 nonzero. (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
_{}
N/A