Matlab Physics

NUMERICAL ANALYSIS OF OPTICAL AND ELECTROMAGNETIC PHENOMENA

MODELLING THE BEHAVIOUR OF GAUSSIAN BEAMS

(PARAXIAL REGIME)

op_beams_001.m

The Script is used to model the behaviour of a Gaussian beam propagating in the Z direction with spherical wavefronts in the paraxial region. The beam is specified by its wavelength (visible part of the electromagnetic spectrum), power and the value of its waist. Also, the XY and XZ planes for the irradiance calculations and plots are also specified in the INPUT Section of the Script. The following parameters are computed in the CALCULATION Section of the Script: the Z axis grid, [m]; Rayleigh range, [m]; the ratio, ; beam divergence angle, ; beam spot size, [m]; Gouy phase angle, ; axial irradiance (), ; radial irradiance (XY plane at ), and ; irradiance displayed in XZ plane, ; phase, . The calculation outputs are displayed in a series of Figure Windows. For different simulations, most of the parameters only need to be changed in the INPUT Section of the Script. The Script calls ColorCode.m which assigns the color of a plot to the wavelength color.

INTRODUCTION

Lasers are widely used in many fields today. Having a mathematical description of laser light is essential. The light emitted from a laser is composed of a narrow band of wavelengths, such that we can assume the light is monochromatic. The light emitted from the laser is usually collimated as it propagates in a straight line as a narrow beam. As a starting point, we will assume the beam to be unpolarized and the light intensity to have a Gaussian profile in any plane that is normal to the direction of propagation. This mode is referred to as the TEM₀₀ mode and it describes the output of most lasers.

Fig. 1. Red light emitted from a HeNe laser. The light is seen due to scattering. The light propagates as a collimated beam in a straight line with slight divergence in its cross-section.

The light emitted from a laser can be modelled by solving the [3D] scalar wave equation for the electric field

(1)

We will assume that a monochromatic beam with a circular cross-section propagates in free spaces and that the electric field can be expressed as a product of its spatial part and time dependent part

(2)

Substitution of equation 2 into equation 1 gives the Helmholtz equation

(3)

We will only consider the solution for a wave propagating in the Z direction where

(4)

The term accounts for the wave oscillation along the propagation direction.

Substitution of equation 4 into equation 3 gives

(5)

We will assume that there is a slow decrease in the amplitude of the wave as the wave propagates in the Z direction and only consider the wave far from the origin but close to the

Z axis such that and . Thus, we can say that varies slowly with z, and thus we can neglect the term . Therefore, the wave equation can be expressed as

(6)

This is known as the paraxial wave equation. Without going into all the mathematical details, a solution to paraxial wave equation is

(7)

where

is the amplitude of the wave.

is the beam waist.

is the radius of the beam at position z at which the amplitude has decreased by the factor 1/e from its axial value or where the irradiance is 1/e² of its axial value. is known as the beam spot.

is the radius of curvature of the wavefront (surface that contains all points of the wave that have the same phase) at position z. As the beam propagates in the +Z direction, the width of the beam expands, and this is measured by the beam spot as it increases with z. For the expanding beam, its wavefront must have a spherical shape because a wave always propagates in a direction perpendicular to its wavefront.

is the Gouy phase.

is the Rayleigh range

(8)

The paraxial approximation is only valid when and is not accurate for strongly diverging beams.

We can better understand the solution of the scalar paraxial wave equation by running the Script op_beams_001.m with the following input parameters:

% INPUTS ==============================================================

% Max power transmitted in beam [W]

P0 = 1e-3;

% beam waist [m]

w0 = 0.5e-3;

% wavelength [m]

wL = 632.8e-9;

% Length of Z domain in multiples of zR : z = nZ * zR

nZ = 5;

% Number of grid points in calculations (must be an odd number)

N = 501;

% Z position of XY plane for radial irradiance plots: zPR = zP / zR

% [e,g, zPR = 0, 1, 2, ... ]

zPR = 4;

% Max radial distance for radial irradiance plot: Rmax = nR * w0

nR = 5;

All the figures in this document where created using the script op_beams_001.m.

Beam spot

At position z along the axis, the beam spot parameters is given by

(9)

Fig. 2. Beam spot w as a function of position z.

When z = 0, we have . This the smallest value that can have. Thus, z = 0 is a special point in the propagation of the beam. For z > 0 the beam spot increases. Also, in the negative Z direction the beam spot also increases. When the beam is focused by a lens there is always a minimum width of the beam at the focal point such that the beam spot is equal to the value of the waist .

It is important to note that the Raleigh range is an indicator of the divergence of the beam. When we have . This is a turning point in the propagation of the light as the beam spot makes the transition from being nearly constant to increasing linearly as shown in figure 2.

When the ray description of the propagation of the light breaks down as the beam spot slowly increases. When , the beam spot expands linearly as

(10)

where is the limiting value of the divergence and is called the divergence angle of the beam. The magenta line in figure 2 shows the divergence of the beam for .

(11)

Note: the smaller the value of the quicker the beam will expand in a linear manner and the smaller the value of the waist, the larger the divergence angle of the beam. Because the divergence is inversely proportional to the spot size, for a given wavelength , a Gaussian beam that is focused to a small spot diverges rapidly as it propagates away from the focus. Conversely, to minimize the divergence of a laser beam in the far field (and increase its peak intensity at large distances) it must have a large cross-section () at the waist and thus a large diameter where it is launched, since . This relationship between beam width and divergence is a fundamental characteristic of diffraction, and of the Fourier transform which describes Fraunhofer diffraction. Since the Gaussian beam model uses the paraxial approximation, it fails when wavefronts are tilted by more than about 30° from the axis of the beam. From the above expression for divergence, this means the Gaussian beam model is only accurate for beams with waists where .

The shape of a Gaussian beam of a given wavelength is governed solely by one parameter, the beam waist . This is a measure of the beam size at the point of its focus () and corresponds to the smallest value of the beam spot parameter. From waist , other parameters describing the beam geometry are determined, such as, the Rayleigh range and asymptotic beam divergence .

Radius of curvature of the wavefront

(12)

The phase term in equation (7)

gives the curvature of the wavefront which is spherical with radius . Figure 3 shows the variation in the radius of curvature as a function of z.

Fig.3. Radius of curvature of the wavefront as a function of z.

When z = 0 we have , the wavefront is a plane wave and at this point all parts of the wave are moving in the same direction.

When , then . The wavefront is a nearly a spherical surface traveling away from z = 0. This is the geometric optics limit. When we focus a beam of light we expect the light rays to go in straight lines to and from the focal point (z = 0). This is correct if . For , the wave aspect of the light applies.

When then which is the minimum value of the radius of curvature and corresponds to the turning point between ray-optics and wave-optics.

Gouy phase

(13)

The Gouy phase slightly shifts the phase of the wavefront of the wave as a whole. For a focused beam, the dependence of the Gouy phase as a function of z is displayed in figure 4.

Fig.4. Goy phase plot for a focused Gaussian beam.

As the Gouy phase asymptotes to . Thus, only a phase shift occurs from to . This results in inversion of a signal that has passed through the waist , that corresponds to the inversion obtained with the ray-tracing approximation (Geometrical Optics).

The most rapid change in phase occurs in the region from to and when , then and , then . The Gouy phase shift along the beam remains within the range (for a fundamental Gaussian beam) and is not observable in most experiments. However, it is of theoretical importance and takes on a greater range for higher-order Gaussian modes.

ENERGY: Irradiance (Intensity) S(r, z) and power P(z)

The irradiance (intensity) S is the time average flow of energy per unit time per unit area [W.m^-2]. The irradiance S can be calculated from the electric field E of the wave from equation 14

(14)

The symbol S is used for the irradiance in the Script and notes and not the more commonly used letter I. c is the speed of light and is the permittivity of free space.

For our propagating Gaussian beam, the electric field E is given by equation 7, so the irradiance S is

(15A)

(15B)

(9)

The maximum irradiance occurs at the location where

(15C)

The variation of the intensity along the Z axis is given by

(15D)

Figure 5 shows a plot of the axial irradiance S_z as a function of z.

Fig. 5. Axial irradiance as a function of position z along the axis of the beam. The maximum value of the irradiance from equation 15C is . Note, at the on-axis irradiance is one-half the maximum irradiance at the waist : .

When , then and . Therefore, the amplitude of the wave decreases with increasing z in the same manner as a spherical wave.

The radial irradiance S_r is the variation in the irradiance in the XY plane located at position z. and typical plots are shown in figure 6.

(15E)

Fig. 6. The variation in the radial irradiance S at two z positions along the axis. The profile of the plots are also a Gaussian function as is the amplitude of the wave.

Figure 7 shows the radial irradiance S_r as a [3D] plot and as image that may be viewed on a screen for the same two z positions as in figure 6.

Fig. 7. Beam profile: Radial irradiance S_r.

Figure 8 shows the beam in the XZ plane and how the beam diverges and how its intensity decreases with increasing z distance from the waist at z = 0.

Fig. 8. Profile of the beam in the XZ plane. The yellow lines show the beam spot as shown in figure 2.

The power P transmitted through a circular disk of radius r_P in an XY plane at position z_P is

(16A)

When then

(16B) total power transmitted by the beam

If we know the total power transmitted by the beam, we can then calculate the electric field amplitude and the maximum irradiance of the beam

(17A)

(17B)

Figure 9 shows the percentage power that would pass through a circular aperture of varying radius placed perpendicular to the beam at the beam waist.

Fig. 9. The percentage power that would pass through a circular aperture of varying radius placed perpendicular to the beam at the beam waist. 86% of the light passes through the aperture when the radius of the aperture is equal to the beam waist.

Phase of the wave