Gaussian Beam Propagation Simulator


1.    Introduction

2.    Schrödinger equation

3.    Hankel Transform

4.    Adaptive spherical Fourier Bessel split-step method

5.    Software and user guide

6.    Results

7.    References

1.    Introduction


In this project, I will present a novel technique for solving nonlinear Schrödinger equation in three dimensions with spherical symmetry, in order to simulate the propagation of a beam in some media, linear or nonlinear with or without spatial variation of the refractive index. This will be done by using Adaptive spherical Fourier Bessel split-step method which will enable us to study many nonlinear effects including the possibility of spatio-temporal collapse.

The purpose of this project is to help us to understand the behavior of a Gaussian beam propagate through nonlinear medium whose refractive index is given as function of the Gaussian intensity. With the help of this program, we can simulate what is called ‘light bullets’ which is very important for telecommunication system due to their self confined structure.   


Back to top

2.    Schrödinger equation


The time-dependent (3 + 1)-dimensional (D = 3) paraxial wave equation which we would like to solve, in the presence of group velocity dispersion is:



Where  is the slowly varying envelope of the optical field,  is the reduced time in a moving frame of reference,  is the group velocity,  is the propagation constant,  is the group-velocity dispersion (GVD) parameter, and  is the nonlinear parameter responsible for self phase modulation (SPM). Note that  can be positive or negative depending on the nature of the dispersion in the medium, corresponding to normal dispersion (ND) and anomalous dispersion (AD), respectively. Also, for  is the intensity induced change in index where  in the case of saturable nonlinearity.

By submitting the following normalization:


We arrive at



Where  is the initial transverse spatial width, and  is the diffraction length. If  we can assume spherical symmetry of the field distribution and introduce the radial variable  to recast equation 2 as



Where D is the spatial dimension. From the linear part of equation 3 with  we have



The Hankel transform or Fourier Bessel technique can not apply directly to this operator in the case when , so we have to transform the operator from spherical coordinates to cylindrical ones by letting , where  is the order of the Fourier Bessel or Hankel transform. Eq(2.4) becomes  where


With  in the case of the cylindrical and spherical Fourier Bessel transform pair respectively, where they are related to each other by:



Where  is the  order spherical Bessel function. The above transform pair is solved by the  order finite Hankel Transform method, explained in the next section.


Back to top

3.    Hankel Transform


We use the definition of the lth order finite Hankel transform of the third kind:



in which  ’s are the roots of the transcendental equation , the inverse transform can be written as:


For b=0 and  Eq.(3.2) becomes ,where  and . Let us evaluate the radius  at  and the frequency at  and , we have , where  are the spatial and transform ranges, respectively, with . From the we can write the expansion of the function and its transform by an  order Bessel series.





Back to top

4.    Adaptive spherical Fourier Bessel split-step method


The algorithm of the adaptive spherical fourier Bessel split-step method is described in Fig1


The AFBSS algorithm is a  symmetrized version of the split-step FFT using spherical Fourier Bessel transform instead, and using adaptive longitudinal stepping and transverse gird management.


Back to top

5.    Software and user guide


The program is written in matlab code and has a friendly GUI. This enables the user to control the parameters of the problem easily and for viewing the result in such a comfort way. In order to run the program you need to download the archive file and unzipped it.




In order to run the program

·        Download the files, and then extracts them in the working directory.

·        Change the directory to your working directory (use cd) and open Matlab.

·        Type BPM and press enter


The GUI should appear like in Fig2. On the GUI you will find the following controls:

1. R1 input - the spatial vector range.

2. R2 input - the transform vector range.

3. w0 input - the beam radius

4. Amplitude of the input Gaussian.

5. Duty cycle - the period for snapshots

6. Iteration number

7. The refractive index as function of the beam intensity n(x) 

8. Checkbox to add a lens with changeable focal length

9. plot button to run the program and generate the graph.

10. defaults  some interesting cases (Defaults1-3 have long run time and Default4 has a short run time).




The generated graph simulates the propagation of a Gaussian beam through a nonlinear medium. It shows the Gaussian intensity change as it propagates the medium. The color is proportional to surface height. The run time of the simulation can take hours or even days depending on the iteration number, the duty cycle, the R1 or R2 ranges. By increasing these parameter we increase the run time of the simulation.


Back to top


6.    Results


In this section I will describe some results of special cases. The first case is a Gaussian beam go thought a lens. As we all know from basic optics, by letting the Gaussian goes through a lens, it will acquires a phase curvature which will focus it at the focal point of the lens.


Fig3 Referance[4]

And by changing the focal point we can see that the point of focusing changes in the same way.


Another case which we can check is the generation of light bullets. But I will not go to their mathematical background and will try to concentrate on the result them. As described in the paper of V.Saraka[1] and the paper of G.Nehmetallah [3], by examining a saturated nonlinear medium of the form  , a special cases emerge describing the generation of light-bullet.

 The first case describes a stable light-bullet with initial focusing followed by self trapping. We can see that by choosing the first set values when we press the defult1 button.

The second case describes a stable light-bullet with initial defocusing followed by self trapping. We can see that by choosing the second set values when we press the defult2 button.

The third case describes a stable light-bullet with initial focusing followed by self trapping. We can see that by choosing the third set values when we press the defult3 button.

The fourth case describes a Gaussian going through a lens. We can see that by choosing the fourth set values when we press the defult4 button. this set values has a short run-time.


Back to top

7.    References


[1] V. Skarka, V. Berezhiani, R. Miklaszewski, Phys. Rev. E 56 (1997) 1080.

[2] Y. Silberberg, Opt. Lett. 15 (1990) 1282

[3] G. Nehmetallah, and P. P. Banerjee,Numerical Modeling of Spatiotemporal Solitons Using an Adaptive Spherical Fourier Bessel Split Step Method,Opt. Comm. 257, 197-205 (2006).