Comparative Modeling of Infrared Fiber Lasers

: The modeling and design of ﬁber lasers facilitate the process of their practical realization. Of particular interest during the last few years is the development of lanthanide ion-doped ﬁber lasers that operate at wavelengths exceeding 2000 nm. There are two main host glass materials considered for this purpose, namely ﬂuoride and chalcogenide glasses. Therefore, this study concerned comparative modeling of ﬁber 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 speciﬁc focus on the central processing unit (CPU) time and calculation residual. Two representative ﬁber laser cavities were considered: One was based on a chalcogenide–selenide glass step-index ﬁber doped with trivalent dysprosium ions, whereas the other was a ﬂuoride step-index ﬁber doped with trivalent erbium ions. The practical calculation accuracy was also assessed by comparing directly the results obtained from the different models.


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 . 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. 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. 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 . 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. 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 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

Materials and Methods
where the sum of level populations N 0 , N 1 , and N 2 is equal to the total doping concentration N Dy , and the coefficients a mn are given by 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 β 21 for the 2-1 transition (Figure 2), and τ 1 and τ 2 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: 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. 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: 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 3+ 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: where the sum of level populations N 0 , N 1 , N 2 , N 3 , and N 4 ( Figure 3) is equal to the total doping concentration N Er . Note that τ 1 , τ 2 , τ 3 , and τ 4 are the lifetimes of levels 1, 2, 3, and 4, respectively, whereas β xy gives the branching ratios from level x to y. W 11 and W 22 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: 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 2 = g 1 = 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. Aligning the fiber with the z axis of the coordinate system enabled the following four differential equations to be written in the following form: 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: 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 The pump, signal, and idler powers at the terminating end-fiber faces of the laser cavity were subjected to the following boundary conditions: 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.

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 −34 J·s and the value of the speed of light in free space of 2.99792458 × 10 8 m/s were used.  Table 3 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 Tables 4-6, 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 3+ -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 3+ -doped fluoride glass fiber laser, both the NU-PWr and UFE 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. Tables 4-6, 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 3+ -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 3+ -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 Tables 4-6 were rounded to the nearest decimal.

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 3+ -doped chalcogenide-selenide step-index glass fiber and in the case of the Er 3+ -doped fluoride glass fiber lasers, a very good agreement was achieved between the results calculated using the UFE and NU-PWr models.