Comparative Modeling of Infrared Fiber Lasers

Abstract

The modeling and design of fiber lasers facilitate the process of their practical realization. Of particular interest during the last few years is the development of lanthanide ion-doped fiber lasers that operate at wavelengths exceeding 2000 nm. There are two main host glass materials considered for this purpose, namely fluoride and chalcogenide glasses. Therefore, this study concerned comparative modeling of fiber lasers operating within the infrared wavelength region beyond 2000 nm. In particular, the convergence properties of selected algorithms, implemented within various software environments, were studied with a specific focus on the central processing unit (CPU) time and calculation residual. Two representative fiber laser cavities were considered: One was based on a chalcogenide–selenide glass step-index fiber doped with trivalent dysprosium ions, whereas the other was a fluoride step-index fiber doped with trivalent erbium ions. The practical calculation accuracy was also assessed by comparing directly the results obtained from the different models.

1. Introduction

Due to many potential applications in medicine, biology, environmental monitoring, and defense, a large research effort has been devoted to the development of fiber lasers operating at wavelengths exceeding 2000 nm. Currently available light sources for these wavelengths include gas lasers, quantum cascade lasers, interband cascade lasers, supercontinuum fiber sources, Raman fiber lasers, light emitting diodes, optical parametric oscillators, Globar©-type black-body sources, and lanthanide ion-doped fiber lasers. A particular advantage of lanthanide ion-doped fiber lasers is their high output beam quality and compact structure: So far, such fiber lasers have only been demonstrated at wavelengths <4000 nm [1,2]. Very recently, a room temperature fiber laser operation up to 3920 nm has been demonstrated [3]. For operating wavelengths of up to 2000 nm, silica glass fiber-based lasers can be used [4]. Fiber lasers operating at wavelengths from 2000 nm to 4000 nm are based on fluoride glass fibers. Lanthanide ions that have been so far applied for doping fluoride glass fibers include erbium (III), holmium (III), and dysprosium (III) [5,6,7,8,9,10,11,12]. For the development of lanthanide ion-doped fiber lasers operating at wavelengths exceeding 4000 nm, the introduction of lower phonon energy glasses is required. Particularly good candidates for this purpose are chalcogenide glasses. Chalcogenide glasses have been shown to have sufficiently good mechanical properties, chemical stability toward water and oxygen, low loss at the relevant wavelengths, good solubility for lanthanide ions, and they can be drawn into fibers [13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36]. Thus, both chalcogenide and fluoride glasses have been intensely studied, both experimentally and theoretically, for applications in fiber lasers [20,36,37,38,39,40,41].

A vital element in the development process of fiber lasers reaching long operating wavelengths is their design. Design tools are needed for the optimization of laser parameters so that a lasing action under optimal conditions can be successfully achieved. Therefore, the properties of various numerical algorithms applicable to the design and modeling of fiber lasers operating at wavelengths exceeding 2000 nm were compared in this contribution. In particular, the optical characteristics of a dysprosium ion-doped chalcogenide glass fiber laser and an erbium ion-doped fluoride glass fiber laser were numerically studied. For this purpose, several algorithms developed within various software environments were compared. The dependence of the central processing unit (CPU) time and calculation residual on the iteration number was used to assess the convergence properties of individual algorithms.

2. Materials and Methods

Figure 1 shows the configuration of the fiber laser cavity considered. The pump light was applied at one end of the fiber, whereas the signal and idler waves were collected at the other end of the fiber.

Two types of fibers were considered. The first one was a chalcogenide glass fiber core doped with trivalent dysprosium ions, whereas the other one was a fluoride fiber core doped with trivalent erbium ions. The energy level diagram for dysprosium ions, doped here into a chalcogenide–selenide glass, is shown in Figure 2. The pump laser, operating at 1710 nm, populated level 2. It was assumed simplistically that neither significant upconversion nor excited state absorption would take place, so that only the three lowest-lying energy levels needed to be included in the model. From energy level 2, a transition could take place to level 1 through either the process of spontaneous or stimulated emission, thus generating signal photons. Analogous transitions could take place between energy levels 1 and 0, accompanied by emission of idler photons. Due to the relatively long lifetime of energy level 1, measured experimentally [27], the inclusion of an idler was essential in obtaining efficient laser action. Note that the idler wave was trapped within the cavity with the help of high reflectivity mirrors, and was used to depopulate level 1. Using the energy level diagram shown in Figure 2 and the rate equation approach, one could write the following set of coupled algebraic equations that allowed for the calculation of the energy level populations of dysprosium ions:

$$\left[\begin{array}{ccc}{a}_{11}& a_{12}& a_{13}\\ a_{21}& a_{22}& a_{23}\\ 1& 1& 1\end{array}\right]\times \left[\begin{array}{c}{N}_0\\ N_1\\ N_2\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ N_{Dy}\end{array}\right],$$

(1)

where the sum of level populations N₀, N₁, and N₂ is equal to the total doping concentration N_Dy, and the coefficients a_mn are given by

$$\begin{array}{l}{a}_{11}={\sigma }_{pa}{\varphi }_p;a_{12}={\sigma }_{sa}{\varphi }_s;a_{13}=-{\sigma }_{pe}{\varphi }_p-{\sigma }_{se}{\varphi }_s-\frac{1}{{\tau }_2}\\ a_{21}={\sigma }_{ia}{\varphi }_i;a_{22}=-{\sigma }_{ie}{\varphi }_i-{\sigma }_{sa}{\varphi }_s-\frac{1}{{\tau }_1};a_{23}={\sigma }_{se}{\varphi }_s+\frac{{\beta }_{21}}{{\tau }_2}\end{array}.$$

(2)

In Equation (2), σ_xa is the absorption cross section for signal s, idler i, and pump p, whereas σ_xe gives the respective values for the emission cross section. The photon flux is ϕ_x for signal s, idler i, and pump p. The branching ratio is β₂₁ for the 2-1 transition (Figure 2), and τ₁ and τ₂ are the radiative lifetimes for levels 1 and 2, respectively (Figure 2). The rate Equation (1) is complemented by the set of six ordinary differential equations that describe the spatial evolution of the pump, idler, and signal powers for both the forward- and backward-propagating waves along the z axis:

$$\begin{array}{l}\frac{d{P}_p^{\pm }}{dz}=\mp {\mathsf{\Gamma }}_p\left[{\sigma }_{pa}{N}_0-{\sigma }_{pe}{N}_2\right]P_p^{\pm }\mp {\alpha }^{}{P}_p^{\pm }\\ \frac{d{P}_s^{\pm }}{dz}=\mp {\mathsf{\Gamma }}_s\left[{\sigma }_{sa}{N}_1-{\sigma }_{se}{N}_2\right]P_s^{\pm }\mp {\alpha }^{}{P}_s^{\pm }\\ \frac{d{P}_i^{\pm }}{dz}=\mp {\mathsf{\Gamma }}_i\left[{\sigma }_{ia}{N}_0-{\sigma }_{ie}{N}_1\right]P_i^{\pm }\mp {\alpha }^{}{P}_i^{\pm }\end{array},$$

(3)

where Γ_x is the confinement factor for signal s, idler i, and pump p; α gives the loss coefficient; and P_p, P_s, and P_i are the values of the power of the pump, signal, and idler, respectively. The numeric values of the parameters are given in

Table 1. Numerical modeling parameters used in simulations, Dy³⁺-doped chalcogenide glass fiber laser.

Quantity	Unit	Value
Dy3+ ion concentration NDy	cm−3	7 × 1019
Aeff	m2	95 × 10−12
Fiber length L	m	2.1
Fiber loss at all wavelengths α	dB/m	1
Lifetime of level 2 (Figure 2)	ms	2
Lifetime of level 1 (Figure 2)	ms	5.2
Branching ratio for 2-1 transitions		0.15
Reflectivity for idler, signal, and pump at z = 0		0.2
Reflectivity for idler, signal, and pump at z = L		0.2
Confinement factor for signal		0.8
Confinement factor for idler		0.9
Confinement factor for pump		0.034
Pump wavelength	μm	1.71
Signal wavelength (λ1)	μm	4.6
Idler wavelength (λ2)	μm	3.35
Pump emission cross section	m2	0.318 × 10−24
Pump absorption cross section	m2	0.501 × 10−24
Signal emission cross section	m2	0.912 × 10−24
Signal absorption cross section	m2	0.485 × 10−24
Idler emission cross section	m2	0.097 × 10−24
Idler absorption cross section	m2	0.016 × 10−24

. A more rigorous approach, which does not use the confinement factor approximation, involves the exact calculation of the overlapping integrals between the ion populations and the electromagnetic field by taking into account the spatial distribution of the optical propagation modes:

$$\begin{array}{l}\frac{d{P}_p^{\pm }}{dz}=\mp P_p^{\pm }{\displaystyle \underset{A_d}{\int }\left[{\sigma }_{pa}{N}_0-{\sigma }_{pe}{N}_2\right]i_{p}dA}\mp {\alpha }^{}{P}_p^{\pm }\\ \frac{d{P}_s^{\pm }}{dz}=\mp P_s^{\pm }{\displaystyle \underset{A_d}{\int }\left[{\sigma }_{sa}{N}_1-{\sigma }_{se}{N}_2\right]i_{s}dA}\mp {\alpha }^{}{P}_s^{\pm }\\ \frac{d{P}_i^{\pm }}{dz}=\mp P_i^{\pm }{\displaystyle \underset{A_d}{\int }\left[{\sigma }_{ia}{N}_0-{\sigma }_{ie}{N}_1\right]i_{i}dA}\mp {\alpha }^{}{P}_i^{\pm }\end{array},$$

(4)

where A_d is the rare earth-doped region and i_p, i_s, and i_i are the normalized intensities of the pump, signal, and idler optical modes, respectively.

Figure 3 shows the energy level diagram of erbium trivalent ions doped into fluoride glass. The pump operating at 980 nm promoted ions from the ground state to energy level 2. Also, via excited state absorption and cooperative upconversion, the ions were promoted to energy level 3. The signal operated at 2800 nm wavelength and was amplified through interaction with energy levels 1 and 2. The idler signal may have operated at approximately 1550 nm and interacted with energy levels 0 and 1. In this model, however, it was assumed that the idler signal did not build up due to the cavity loss encountered. It is noted that in the fluoride fiber laser cavity considered, the inclusion of an idler was not needed for efficient laser action. From the Er³⁺ energy level diagram in Figure 3, using the rate equations approach, one obtained consistently the following set of algebraic equations that enabled calculation of the populations of the energy levels:

$$\begin{array}{l}{W}_{22}{N}_2^2-\frac{N_4}{{\tau }_4}+R_{ESA}=0\\ \frac{{\beta }_{43}{N}_4}{{\tau }_4}-\frac{N_3}{{\tau }_3}=0\\ R_{GSA}-R_{SE}-R_{ESA}+{\displaystyle \sum _{i=3}^4\frac{{\beta }_{i2}{N}_i}{{\tau }_i}}-\frac{N_2}{{\tau }_2}-2W_{22}{N}_2^2+W_{11}{N}_1^2=0\\ R_{SE}+{\displaystyle \sum _{i=2}^4\frac{{\beta }_{i1}{N}_i}{{\tau }_i}}-\frac{N_1}{{\tau }_1}-2W_{11}{N}_1^2=0\\ -R_{GSA}+{\displaystyle \sum _{i=1}^4\frac{{\beta }_{i0}{N}_i}{{\tau }_i}}+W_{22}{N}_2^2+W_{11}{N}_1^2=0\end{array},$$

(5)

where the sum of level populations N₀, N₁, N₂, N₃, and N₄ (Figure 3) is equal to the total doping concentration N_Er. Note that τ₁, τ₂, τ₃, and τ₄ are the lifetimes of levels 1, 2, 3, and 4, respectively, whereas β_xy gives the branching ratios from level x to y. W₁₁ and W₂₂ are the cooperative upconversion coefficients for levels 1 and 2, respectively. R_GSA gives the ground state absorption rate, R_SE gives the rate of stimulated emission between levels 1 and 2, and R_ESA gives the rate of the excited state absorption from level 2 to level 4:

$$R_{GSA}=\frac{{\lambda }_p{\mathsf{\Gamma }}_p{\sigma }_{GSA}}{hc{A}_{eff}}{N}_0\left(P_p^{+}+P_p^{-}\right),$$

(6)

$$R_{SE}=\frac{{\lambda }_s{\mathsf{\Gamma }}_s{\sigma }_{se}}{hc{A}_{eff}}\left(b_2{N}_2-\frac{g_2}{g_1}{b}_1{N}_1\right)\left(P_s^{+}+P_s^{-}\right),$$

(7)

$$R_{ESA}=\frac{{\lambda }_p{\mathsf{\Gamma }}_p{\sigma }_{ESA}}{hc{A}_{eff}}{N}_2\left(P_p^{+}+P_p^{-}\right).$$

(8)

Equation (5) is complemented by a set of four ordinary differential equations that describe the evolution of the pump and signal waves. The degeneracy parameters are g₂ = g₁ = 2. The values of the relevant cross sections σ_se, σ_ESA, and σ_GSA, confinement factors Γ_x, wavelengths λ_x, effective cross section A_eff, and Boltzmann factors b_x are given in

Table 2. Numerical modeling parameters used in simulations, Er³⁺-doped fluoride fiber laser.

Quantity	Unit	Value
b1/b2		0.1/0.16
W11	m3/s	1 × 10−24
W22	m3/s	0.3 × 10−24
σGSA	m2	2.1 × 10−25
σSE	m2	4.5 × 10−25
σESA	m2	1.1 × 10−25
Γp		0.009
Γs		1.0
Er3+ ion concentration NEr	m−3	9.6 × 1026
Pump wavelength λp	Nm	976
Pump wavelength λs	Nm	2800
Fiber length L	m	2.5
Aeff	m2	314 × 10−12
αp	1/m	3 × 10−3
αs	1/m	23 × 10−3
Rp (z = 0)		0
Rp (z = L)		0.04
Rs (z = 0)		0.96
Rs (z = L)		0.04

. Aligning the fiber with the z axis of the coordinate system enabled the following four differential equations to be written in the following form:

$$\begin{array}{l}\frac{d}{dz}{P}_p^{+}=-{\mathsf{\Gamma }}_p\left({\sigma }_{GSA}{N}_0+{\sigma }_{ESA}{N}_2\right)P_p^{+}-{\alpha }_p{P}_p^{+}\\ -\frac{d}{dz}{P}_p^{-}=-{\mathsf{\Gamma }}_p\left({\sigma }_{GSA}{N}_0+{\sigma }_{ESA}{N}_2\right)P_p^{-}-{\alpha }_p{P}_p^{-}\\ \frac{d}{dz}{P}_s^{+}={\mathsf{\Gamma }}_s{\sigma }_{SE}\left(b_2{N}_2-\left(g_2/g_1\right)b_1{N}_1\right)P_s^{+}-{\alpha }_s{P}_s^{+}\\ -\frac{d}{dz}{P}_s^{-}={\mathsf{\Gamma }}_s{\sigma }_{SE}\left(b_2{N}_2-\left(g_2/g_1\right)b_1{N}_1\right)P_s^{-}-{\alpha }_s{P}_s^{-}\end{array},$$

(9)

where P_s and P_p are the powers of the signal and pump, respectively, and the superscripts + and − denote the forward- and backward-propagating waves, respectively. In Equation (9), α_x gives the value of loss.

The pump, signal, and idler powers at the terminating end-fiber faces of the laser cavity were subjected to the following boundary conditions:

$$\begin{array}{l}{P}_p^{+}\left(z=0\right)=r_p\left(z=0\right)P_p^{-}\left(z=0\right)+\left[1-r_p\left(z=0\right)\right]P_{pump}^{}\\ P_p^{-}\left(z=L\right)=r_p\left(z=L\right)P_p^{+}\left(z=L\right)\\ P_{}^{+}\left(\lambda ={\lambda }_1,z=0\right)=r_{\lambda 1}\left(z=0\right)P_{}^{-}\left(\lambda ={\lambda }_1,z=0\right)\\ P_{}^{-}\left(\lambda ={\lambda }_1,z=L\right)=r_{\lambda 1}\left(z=L\right)P_{}^{+}\left(\lambda ={\lambda }_1,z=L\right)\\ P_{}^{+}\left(\lambda ={\lambda }_2,z=0\right)=r_{\lambda 2}\left(z=0\right)P_{}^{-}\left(\lambda ={\lambda }_2,z=0\right)\\ P_{}^{-}\left(\lambda ={\lambda }_2,z=L\right)=r_{\lambda 2}\left(z=L\right)P_{}^{+}\left(\lambda ={\lambda }_2,z=L\right)\end{array}.$$

(10)

It should be noted that the boundary conditions in Equation (10) dictated the incident pump power, but not the value of the pump power after it crossed the air-fiber end interface.

Three algorithms, developed by different research groups in different computational environments, were compared here. The main characteristics considered were the CPU time and the dependence of the calculation residual on the iteration number. The calculation residual was defined as the sum of the squared differences between values from the current and previous iterations for signal, idler, and pump, calculated at z = 0. The CPU time was calculated using system functions. All three algorithms employed the relaxation method for solving the two-point boundary value problem, as follows:

• The fiber laser model developed at the Institute of Photonics and Electronics of the Czech Academy of Sciences (UFE) was implemented in C programming language (gcc 4.9.2) within the Windows 7 operating system, 64 bit Intel core i7-3930K CPU at 3.2 GHz. The UFE model is currently being developed for the study of longitudinal-mode instabilities and associated buildup of dynamic fiber Bragg gratings [40].
• The fiber laser model developed at the Politecnico di Bari (PB) was implemented in MATLAB within the Windows 10 operating system, 64-bit Intel Core i7-4790 CPU at 3.6 GHz. The numerical integration was carried out using a 4-5 Runge–Kutta algorithm, and the more rigorous overlap integrals approach was employed.
• The fiber laser model developed at the University of Nottingham and Wroclaw University of Science and Technology (NU–PWr) was implemented in MATLAB within the Windows 10 operating system, 64 bit Intel Core i5 7th Generation, CPU at 2.5 GHz. The numerical integration was carried out using a 4-5 Runge–Kutta algorithm.

3. Results

The modeling parameters for the dysprosium trivalent ion-doped chalcogenide–selenide glass fiber laser are summarized in Table 1, whereas in Table 2 the modeling parameters for the erbium trivalent ion-doped fluoride glass fiber laser are given.

It is noted that the value of 1 dB/m for chalcogenide fiber loss at 4600 nm was challenging in practical realization. However, a low value of the fiber loss is necessary for the realization of an efficient fiber laser and is thus widely used in fiber laser modeling-related literature [17,25,42,43,44].

The parameters summarized in Table 2 were extracted from experiments and verified by J.F. Li and S.D. Jackson by comparing measured and numerical results [11]. Further effort has been undertaken to verify the reliability of the numerical codes through comparing the numerical results to experimental measurements [45].

In the simulations, the value of Planck’s constant of 6.62607004 × 10⁻³⁴ J·s and the value of the speed of light in free space of 2.99792458 × 10⁸ m/s were used.

Table 3. Branching ratios and level lifetimes for erbium trivalent ions doped into a fluoride glass.

Quantity	Unit	Value
τ1	ms	9
τ2	ms	6.9
τ3	ms	0.12
τ4	ms	0.57
β21, β20		0.37, 0.63
β32, β31, β30		0.856, 0.004, 0.14
β43, β42, β41, β40		0.34, 0.04, 0.18, 0.44

shows the values of the relevant lifetimes and branching ratios for erbium trivalent ions doped into a fluoride glass.

Figure 4 shows the dependence of the residual and the CPU time in the UFE model. The values of the CPU time showed a step-wise behavior due to the quantization implemented within the C function clock that was used in the simulations. In the UFE model, the rate of residual reduction was smaller at low values of the output power. In particular, at an output power of 200 mW, one can observe that the residual decreased significantly more slowly than did the other three values of the output power.

An overall lower rate of the residual reduction in the case of the PB model, when compared to the UFE model, was observed in the results shown in Figure 5. The simulation time was also at least three orders of magnitude larger despite the application of a faster processor at PB, which shows the advantage of direct C programming. Interestingly, the PB model showed a larger rate of residual decrease for low output powers. Figure 6 shows the results obtained with the NU–PWr model. When compared to the results obtained with the PB model, one observes a much larger rate of residual reduction in the NU–PWr model. However, the CPU time in the PB model, when measured per iteration, was less. The overall calculation time for the PB model had to reach a particular value of the residual, and this took a significantly longer time than the NU–PWr model. It is noted that for the PB model, the overlap integrals between the ion populations and the optical modes of pump, signal, and idler were calculated over the rare earth-doped region according to Equation (4). These integrals were updated along the fiber length, taking into account the ion population distributions. This caused a higher calculation time, but allowed higher solution accuracy.

Finally, in

Table 4. Calculated values of signal output power values for the Dy³⁺-doped chalcogenide–selenide glass fiber laser.

Pump Power/W	Signal Power (NU–PWr)/W	Signal Power (UFE)/W	Relative Difference
0.2	4.733 × 10−3	4.731 × 10−3	0.422 × 10−3
0.4	8.744 × 10−3	8.736 × 10−3	0.915 × 10−3
1	49.13 × 10−3	49.04 × 10−3	1.833 × 10−3
5	319.1 × 10−3	318.6 × 10−3	1.568 × 10−3

,

Table 5. Calculated values of the idler output power values for the Dy³⁺-doped chalcogenide–selenide glass fiber laser.

Pump Power/W	Idler Power (NU–PWr)/W	Idler Power (UFE)/W	Relative Difference
0.2 W	0 W	4.140 × 10−6	NA
0.4 W	0 W	9.591 × 10−4	NA
1 W	55.38 × 10−3 W	55.26 × 10−3	2.169 × 10−3
5 W	426.0 × 10−3 W	425.4 × 10−3	1.409 × 10−3

and

Table 6. Calculated values of signal output power values for the Er³⁺ ion-doped fluoride glass fiber laser.

Pump Power	Signal Power (NU–PWr)/W	Signal Power (UFE)/W	Relative Difference
5 W	1.432	1.432	0
10 W	3.171	3.171	0
15 W	4.868	4.868	0
20 W	6.458	6.458	0

, the results for the output power and the idler power, calculated using the UFE and NU–PWr models, are compared for both the dysprosium trivalent ion-doped chalcogenide glass fiber laser and the erbium trivalent ion-doped fluoride glass fiber laser. For the Dy³⁺-doped chalcogenide–selenide glass fiber laser, for the results calculated using the UFE and NU–PWr models, the relative difference, defined as the ratio between the absolute value of the difference and half of the sum of the results, was then less than 0.2% for the signal and below 0.22% for the idler wave at pump powers of 1 W and 5 W, respectively. In the case of the idler wave, the small values of the idler wave power for pump powers of 0.4 W and 0.2 W made it difficult to achieve small values of the relative difference. Nonetheless, these results consistently indicated that the idler was below the lasing threshold. In the case of the Er³⁺-doped fluoride glass fiber laser, both the NU–PWr and UFE models calculated results that agreed on all four digits. It is noted that the results shown in Table 4, Table 5 and Table 6 were rounded to the nearest decimal.

4. Conclusions

In this paper, software packages developed within various environments for the modeling and design of Mid infrared MIR fiber lasers were compared. The analysis was focused on the comparison of the CPU time and the values of the computational residual. The simulation results showed an advantage to using direct encoding of the algorithm in terms of the simulation time. Also, a comparison was carried out between the results obtained by different models. Both in the case of the Dy³⁺-doped chalcogenide–selenide step-index glass fiber and in the case of the Er³⁺-doped fluoride glass fiber lasers, a very good agreement was achieved between the results calculated using the UFE and NU–PWr models.