1. Introduction

High-power narrow-linewidth laser sources with good beam quality are of considerable interest due to their utility in several applications including coherent beam combination [1] and nonlinear frequency conversion [2]. Diode-pumped rare-earth doped fiber amplifiers are highly desirable as such laser sources as they are blessed with small footprints, superior thermo-optical properties, and high electrical to optical and optical to optical efficiencies. However, fiber amplifiers are limited by the onset of optical nonlinearities including stimulated Brillouin scattering (SBS), stimulated Raman scattering (SRS) and self-phase modulation (SPM). It is well-known that the thresholds of these nonlinearities are generally proportional to the effective area of the fiber core. Consequently, a considerable amount of effort has been devoted over the past decade to increasing the fiber core size while maintaining single-mode operation. Unfortunately, several research groups have reported recently that large mode area (LMA) high-power fiber amplifiers are susceptible to the onset of modal instabilities [3,4]. These instabilities are characterized by a sudden and dramatic increase in the energy within higher-order modes (especially the LP₁₁) as well as rapid temporal oscillations.

There have been several theoretical and experimental investigations of the physics of modal instabilities [3–9]. It is now generally agreed upon that the onset of these instabilities is a result of interplay between thermal effects and modal interference. This interplay was first suggested in the literature by Jauregui et al. [5]. Subsequent work by Smith and Smith [6] emphasized the need to include the temporal dynamics to account for the energy transfer between the fundamental and higher-order modes. In order to achieve this, they introduced a frequency offset between the two modes and a steady-periodic treatment of the heat equation. It was not clear in that work how to accurately determine the frequency offset other than by trial and error. Furthermore, by allowing only one frequency offset, contributions from the other frequency components were neglected. Semi-analytical approaches by Hansen et al. [7,8] and Dong [9] that built on the frequency offset proposition provided valuable insights, but due to the difficulty of finding exact solutions did not resolve the spatial and temporal evolution of laser gain. Alternatively, the work of Ward et al. [10] showed that a full time-dependent treatment of the heat equation can capture much of the physics of the modal instabilities without assuming a frequency offset a priori. Furthermore, this work was able to capture the three regimes of temporal dynamics described in the experimental work of Otto et al. [11].

In this paper, we present extensive theoretical and computational analysis of the temporal dynamics of modal instabilities in fiber amplifiers. Our analytical treatment of laser gain, amplitude evolution equations, and the heat equation provides formal proof of the need for a dynamic approach to the modal instabilities. In this treatment, the thermal load is dictated by population inversion rather than assuming it simply follows the seed intensity profile. The time-dependent heat equation subject to a Dirichlet boundary condition is solved in polar coordinates using the alternating direction implicit (ADI) method. The implementation of this method in polar coordinates is typically avoided due to the occurrence of an artificial singularity at the center of the coordinate system and due to the added complexity of applying the appropriate boundary conditions. However, since the polar coordinate system is more suitable to describe the fiber geometry, we have proceeded with implementation of this technique. The propagation of the electromagnetic waves is described using a coupled-mode theory approach. The resultant system of coupled nonlinear equations is solved using a fourth-order Runge-Kutta solver. Low resolution of the spatial and temporal numerical grid can potentially lead to inaccurate solutions including an artificially induced transfer of energy. Here, every numerical evaluation is verified and quantum efficiency as related to energy conservation has been monitored to track any possible numerical errors for an accurate computation of the field propagation. Numerical simulations exploring dependence on pump power, seed power, and seed modal ratio are presented. In addition, two different gain-tailored fiber designs are simulated and compared. It is shown numerically that both gain-tailoring and increased seed power can mitigate the modal instabilities. Finally, the effect of modulating the seed source is investigated using the time-dependent approach. The results are in general agreement with the predictions of the steady-periodic approach and indicate that the two approaches contain similar underlying physics.

The paper is organized as follows: section 2 describes the governing equations. In section 3, we describe the analytical and numerical approaches to this phenomenon. The results of our simulations are presented in section 4 and a summary is provided in section 5.

2. Governing equations

Since we are examining the spatial and temporal evolution of the modes, we begin by considering the real part of the index of refraction to be a function of space and time, i.e. $n(\overrightarrow{r},t)$. The deviation of $n(\overrightarrow{r},t)$from the background index of refraction,$n_0$, is attributed to thermal effects. The contribution of laser action to$n(\overrightarrow{r},t)$, which is dictated by the Kramer-Krönig relations, is negligible as argued in [12]. The pertinent wave equation for an isotropic medium derived from Maxwell’s equations is then given by:

In the analysis that follows, we consider a weakly guiding, polarization-maintaining (PM) fiber with the light launched into the fiber polarized along one of its axis; chosen here to be the x-direction. Efficient coupling to the other polarization state through modal interference is not possible and thus the other polarization state is not considered here. The longitudinal component of the electric field,$E_z$, is small and is neglected in the analysis.

Since the laser gain and the variation in the index of refraction are small perturbations, we can construct the solution to the system as a superposition of the modes in the “ideal” fiber:

The variation in refractive index $\delta n=n-n_0$ is directly dependent on temperature difference, $\Delta T$, by $\delta n=\Delta T(dn/dT)$where $dn/dT=1.29\times {10}^{-5}$/°C is the thermo-optic coefficient for silica. Typically, for high-power amplifiers,$\Delta T$within the core is on the order of 100 °C. Therefore, ${(\delta n)}^2~O({10}^{-6})$while$n_0~O(1)$, and the following approximation is valid: $n^2={(n_0{}^2+\delta n)}^2=n_0{}^2+2n_0\delta n+{(\delta n)}^2~n_0{}^2+2n_0\delta n$. The polarization in the x-direction due to laser gain is given by$P_{LG,x}\approx 2i{\epsilon }_0{n}_I{n}_0{E}_x$, where $n_I$ is the imaginary index of refraction and $n_I\ll n_0$. We then proceed by relating $n_I$ to the imaginary part of the wave number,$k_I$, which corresponds to the spatial growth (or decay) of the field due to the population inversion of the two-level system:

We proceed by using Eq. (3) in the wave equation for $E_x$ and applying the differential operators. The set of ${\phi }_k$ is comprised of the eigenfunctions of the “ideal” fiber and thus:

The differential operators on the left hand side (LHS) of Eq. (2) will also yield the following terms that need further examination to determine if they are retained in the first-order system of equations: $-(2A_k{\phi }_k/n_0){\partial }_x{}^2\delta n$, $-(2A_k/n_0)({\partial }_x{\phi }_k)({\partial }_x\delta n)$, $(4i{n}_0\omega /c^2)A_k{\phi }_k{\partial }_t\delta n$, $-(4n_0/c^2){\phi }_k({\partial }_t{A}_k)({\partial }_t\delta n)$, and $-(4n_0/c^2){\phi }_k{A}_k{\partial }_t{}^2\delta n$. The last three terms are clearly negligible as ${\partial }_t\delta n$is dictated by the thermal response time which is much greater than the optical period (i.e. ${\partial }_t\delta n\ll \omega \delta n$). Of the first two terms, one would expect the second term to provide the greater contribution. Even this term is typically negligible as $({\partial }_x{\phi }_k)({\partial }_x\delta n)\ll {\beta }_k{}^2{\phi }_k\delta n$. This can be argued by noting that the length scale of the transverse profile is typically larger than the optical wavelength (i.e. ${\partial }_x{\phi }_k<{\beta }_k{\phi }_k$) and that ${\partial }_x\delta n$ which is dictated by the temperature profile in the transverse direction is even smaller as compared to ${\beta }_k\delta n$. For some of the simulations presented below, we have compared $({\partial }_x{\phi }_k)({\partial }_x\delta n)/({\beta }_k{}^2{\phi }_k\delta n)$and determined it to be <$O({10}^{-4})$.

On the right hand side (RHS) of Eq. (2), only the second derivative in time term generating${\omega }^2$is retained. Terms containing ${\partial }_t({\sigma }_s^{(e)}{N}_2-{\sigma }_s^{(a)}{N}_1)$ or ${\partial }_t^2({\sigma }_s^{(e)}{N}_2-{\sigma }_s^{(a)}{N}_1)$are neglected as the laser gain response time is much greater than the optical period. The following equation that describes the evolution of the laser field is then obtained:

Equation (7) represents a system of coupled partial differential equations (PDEs) whose size is dictated by the number of modes. The standard approach to solving this system is to use the method of characteristics. In this case, the characteristics are defined by$dt/dz=n_0/c$; leading to the characteristic lines given by $t-t_0=(z{n}_0/c)$ where $t_0$ initializes each set of characteristics. The right hand side (RHS) of Eq. (7) depends on time through$\delta n$,$N_2$and _{$N_1$}. The time scales for these quantities are much longer than the traverse time of the signal, $L{n}_0/c$ (where $L$is the length of the fiber). To be certain, for high-power amplifier experiments utilizing LMA fiber, the thermal response which dictates the variation in $\delta n$ in time is on the order of milliseconds and the population inversion occurs on a time scale of 10 µsec. In contrast, the traverse time is on the order of 10 nsec. Therefore, we can estimate $t~t_0$ which is equivalent to freezing time as each set of light rays (characteristics) traverses the fiber. Therefore the equation describing the spatial evolution of $A_j$ at a time $t_0$ is given by:

The variation of $\delta n$ in time is captured by solving the time-dependent heat equation which is solved in conjunction with the coupled amplitude equations. Neglecting longitudinal heat flow, the heat equation takes the following form in polar coordinates:

Finally, the evolution of the pump intensity is given by:

3. Analytical and numerical approach

We simulate a co-pumped 1.6 m LMA Yb-doped fiber amplifier with 50 μm core diameter, 250 μm pump cladding diameter, and 500 μm outer cladding diameter. The Yb doping concentration in the core is taken to be of 6 × 10²⁵ /m³. It is assumed that the fiber is “cold” at the beginning of the simulation and therefore evaluation starts immediately after the introduction of pump light into the fiber. Various seed powers ratios, and pump powers are examined and their effects on modal instabilities are reported in section 4. In this section we discuss the numerical and analytical details which are essential in understanding this phenomena as well as conducting accurate evaluations.

3.1 The need for time-dependent treatment

One of the requirements for efficient transfer of energy is that the phase-matching condition is met. This is realized through the oscillations created in the temperature (and consequently the index of refraction) due to modal interference [5,6]. It was pointed out based on general arguments that, in addition, a phase lag between the signal intensity oscillations and the index of refraction oscillations is required [6]. We will take a different approach by showing mathematically using the population inversion and amplitude propagation equations that a steady state approach to the heat equation will largely lead to a phase shift of the amplitudes of the fields and no transfer of energy. Therefore, in order to simulate the modal instabilities and observe the energy transfer, the dynamical behavior of this phenomenon needs to be investigated.

In Eq. (11) for the generated heat, the term containing $N_0$is a DC term that does not provide for the required phase matching for energy transfer among the modes. The generated heat is also a function of the upper state population inversion,$N_2$, which results in the dependence of the temperature and index of refraction variation on$N_2$. This term will provide for oscillations with the proper periodicity. We restrict our analysis to the fundamental LP₀₁ mode and the LP₁₁ mode. In this case, $I_s$is given by:

For brevity, we define$\xi =A_0{A}_1^{*}{\phi }_0{\phi }_1^{*}{e}^{i\Delta {\beta }_{01}z}+A_0^{*}{A}_1{\phi }_0^{*}{\phi }_1{e}^{-i\Delta {\beta }_{01}z}$. It can readily be shown that $\xi \le {\left|A_0\right|}^2{\left|{\phi }_0\right|}^2+{\left|A_1\right|}^2{\left|{\phi }_1\right|}^2$ leading to${\eta }_2\xi /{\alpha }_2<1$. By applying the binomial expansion, the following expression is obtained:

Equation (18) indicates that the most dominant term for the oscillations in $N_2$along $z$ is $\left(-\frac{{\alpha }_1{\eta }_2}{{\alpha }_2{}^2}+\frac{{\eta }_1}{{\alpha }_2}\right)\xi$ which oscillates with the period of $2\pi /\Delta \beta$. The term containing ${\xi }^2$does not oscillate at the required period with the next term providing this period being${\xi }^3$.

Note that $Q$ oscillates as $N_2$ does and is given by:

The term containing $\xi$ in the expression above accounts for the oscillation in the temperature. Consequently, the period of the oscillatory part of the average temperature is$2\pi /\Delta \beta$. The rest of the terms account for quantum defect heating in the absence of modal interference. It should be pointed out, however, that the temperature at the center of the core oscillates at ½ the period, i.e. $\pi /\Delta \beta$. This can be deduced by noting that $Q(r,\theta ,z)\approx Q(r,\theta +\pi ,z+\frac{\pi }{\Delta \beta })$and that from symmetry, the temperature at the center should be the same regardless if the higher value of $Q$at a given point in $z$ lies in the upper half ($0<\theta <\pi$) or the lower half ($0<\theta <\pi$) of the fiber.

We now proceed by examining the coupled system as described by Eq. (9). We examine, only the phase-matched terms as they are the only terms of relevance. Also, SPM due to the coupling coefficients ${\kappa }_{\delta n,0\to 0}$and ${\kappa }_{\delta n,1\to 1}$does not contribute directly to a transfer of energy and is neglected in this theoretical analysis. The function ${\phi }_0$is chosen to be real and based on symmetry arguments there is no loss of generality in also choosing in our analysis${\phi }_1$to be real. Since $\delta n$is a real quantity, ${\kappa }_{\delta n,1\to 0}$and ${\kappa }_{\delta n,0\to 1}$ are real and equal to each other. In general, one can show that the following expression holds:

Note that contribution of the term containing $A_0^{*}{A}_1{\phi }_0{\phi }_1{e}^{-i\Delta {\beta }_{01}z}$in the coupling coefficient as described by Eq. (20) is not included in the above system of equations as it does not provide phase matching. Also, note that $g_{l,0}$and $g_{l,1}$represent the laser gain for each mode in the absence of modal interference.

The terms containing $\varsigma$in Eqs. (21) are much smaller than the other terms. Therefore, the solution to system of Eqs. (21) can be approximated as:

Thus, the intensity of each mode is only affected by laser gain and there is no transfer of energy between them due to$\delta n$; i.e. only cross-phase modulation (XPM) happens. Equations (24) can be used in Eqs. (23) to obtain more explicit expressions for the amplitudes:

Consequently, a time-independent treatment of this phenomenon can only lead to a phase modulation of the amplitudes. In order for a significant energy transfer to occur, a phase lag between temperature and intensity is required which forms in time.

To illustrate how this phase lag develops, Fig. 1(a) displays the difference between core temperature and room temperature,$\Delta T$, along the fiber after 40 ms, 60 ms, and 80 ms for the fiber whose parameters are given above and at a pump power of 800 W. The numerical techniques used for the simulations are described below in subsection 3.2. The temperature inside the core is plotted at a fixed point away from the center and the close up is displayed in Fig. 1(b) to indicate the oscillations possess a period of approximately$2\pi /\Delta \beta$. The plot indicates that a phase lag occurs among the thermal gratings sampled at various times. In the inset to Fig. 1(b), we show that the oscillations at the center possess a period of $\pi /\Delta \beta$and in agreement with the conclusions given in [13].

 

Fig. 1 a) Temperature difference between fiber core and room temperature,$\Delta T$, along the fiber after 40 ms, 60 ms, and 80 ms at a point away from the center. b) Close-up of part (a) to display the phase shift and oscillations with period$2\pi /\Delta \beta$. The inset reveals the oscillations have a period$\pi /\Delta \beta$at the center.

Download Full Size | PDF

3.2 Numerical evaluation of the optical fields and temperature

We consider only the LP₀₁ and LP₁₁ in our numerical simulations. In order to solve the system of nonlinear ordinary differential Eqs. (9), we have applied the fourth-order Runge-Kutta method to compute the evolution of the optical fields along the fiber. By using a higher-order method, the solution converges with a lesser number of computational points. Note that the longitudinal length scale is dependent on the beat length; that is enough resolution is required to capture the beating between the two modes. For these computations, we used 110 computational points per beat length. The overlap integrals along the transverse direction required to solve the coupled system were evaluated to second-order using the trapezoidal rule. We have found that tracking the quantum efficiency of our system offered an effective tool to validate the numerical resolution of the grid in the longitudinal direction. Low resolution in the longitudinal direction led to a noticeable departure from the expected quantum efficiency of the system and consequently a lack of energy conservation. Unlike the analytical approach presented above, all terms in the nonlinear system of equations were retained: i.e. SPM, XPM, energy transfer as captured by${\kappa }_{\delta n,1\to 0}$, and direct coupling of the modes due to laser gain (i.e. term containing $\varsigma$). We noted that retention of the latter term is important in retaining high fidelity for conservation of energy.

The ADI method in polar coordinate was used to solve the heat equation. Due to the cylindrical shape of the fiber, it is more suitable to use polar coordinates to impose the correct boundary conditions which will result in smaller numerical error. The ADI method is typically implemented in rectangular domains to study the modal instabilities using a semi-periodic treatment of the heat equation and rectangular boundary conditions. The reason that this method is avoided in polar coordinates is the added complexity due to the occurrence of an artificial singularity at the center of the coordinate system. Nevertheless we have successfully implemented this method in polar coordinates.

In some studies of the modal instabilities, the boundary condition was artificially placed well inside the fiber. However, we believe the spatial mesh in the transverse direction must encompass a significant portion of the fiber for an accurate computation of the threshold of the modal instabilities. To illustrate this, we examine the temperature profile in a fiber of core radius of $r_{core}$ and outer radius of$r_{out}$and a constant heat deposition$Q_0$ inside the core. It is straightforward to determine the solution for a Dirichlet boundary condition to be:

We have chosen the spatial mesh in $r$ and $\theta$ so it minimized the computational cost while keeping the error in determining the temperature under 0.5%. The adaptive stepping along the radius was chosen by selecting a fine mesh inside the core and a larger spatial spacing for the pump cladding and outer cladding. For our simulations, 45 computational points in the radial direction and 90 in the azimuthal direction were used which summed up to 60 million computational points at each time step.

In order to avoid the singularity at the center of coordinate system, a grid was chosen so that the grid points are half-integrated in the radial direction and integrated in the azimuthal direction as given below:

Figure 2 displays the convergence of this method for a test problem with a known solution. This example uses$Q=(r\cos \theta +r\sin \theta -2t^2)\exp \left(rt(\sin \theta +\cos \theta )\right)$. Here $n=8,16,32$, $m=2n$, and $k_t=80n$ where $k_t$ is the time step. Table 1 and Fig. 2 confirm that the method converges with second-order accuracy as expected.

 

Fig. 2 Error in temperature versus the number of computational points in log form (Blue) for n = 8, 16 and 32 compared to a line with a slope of −2 (Red).

Download Full Size | PDF

Table 1. Error in temperature for n = 8, 16 and 32

View Table

4. Numerical results

Our numerical and analytical work has confirmed the thermal origin of modal instabilities as well as offering robust explanation on dynamical dependence of temperature effects. We believe the longitudinal dependence of heat is negligible due to the comparison of longitudinal length scale versus transverse length scale. The variation of temperature in the transverse direction is much larger than its variation in a beat length. Phase shift occurs in time without considering longitudinal heat. We have not considered any frequency offset between the two modes. Instead the time-dependent evolution has been studied whereby the Doppler shift is automatically considered. Furthermore, the polar coordinate has been used in our model where correct boundary condition has been implemented. Our numerical results are presented below. Various seed powers and ratios are examined and their effects on modal instabilities are reported. In addition, two different gain-tailored designs for modal instabilities suppression are simulated and compared. Furthermore, we investigate the effect of signal modulation on the threshold.

4.1 Power-dependent behavior

An important characteristic of the modal instabilities is power dependence and as reported by Otto et al. [11] consists of three regimes: stable, transition, and chaotic. Figure 3(a) displays the total intensity in the two modes at the output end of the fiber for various times for a seed power of 30 W of which 5% is in the LP₁₁ mode. Also shown, is the power in the two modes at the output end of the fiber as a function of time (up to 80 ms). The amplifier is co-pumped with a pump power of 720 W (below threshold). Temporal variation of intensity at the output end as well as temporal variation in modal powers along the propagation direction,$z$, is shown in Media 1. Changes occur inside the fiber below the power threshold but they are not strong enough to impact the fundamental mode after a certain time. The energy transfers from the fundamental mode to the higher-order mode back and forth, but after the temperature approaches a quasi-equilibrium phase, energy remains largely in the fundamental mode. A similar behavior is shown in Ward et al. [10]. Figure 3(b) and Media 2 show the corresponding plots for a pump power of 820 W (just past threshold). In this case, the modes keep transferring energy and no steady state is reached. Increasing power further leads to an increase in the frequency of “deep” oscillations between the two modes and this instability enters a “chaotic” regime. In this regime, the power in each mode, on average, is ~50%. Figure 3(c) and Media 3 display the corresponding plots for a pump power of 1040 W. It should be pointed out that, in all three cases, there are transient effects occurring in the first 10-20 ms related to the abrupt turn on of the pump power.

 

Fig. 3 Top picture plots the total intensity of the modes in time and the lower plot displays the fundamental mode power (Blue) and the higher order mode power (Red) in each mode at the output end of the fiber in time for pump power of 720, 820, and 1040 W. The corresponding movies (Media 1, Media 2, and Media 3) show the temporal variation of intensity at the output end as well as temporal variation in modal powers along the propagation direction.

Download Full Size | PDF

We have computed the ratio of the average power over a “long time” range in the LP₁₁ mode to the total signal power as a function of pump power. The “long time” range in our simulations was 60-80 ms. The higher-order mode content vs. pump power is shown in Fig. 4. The modal instabilities threshold represents the pump power where a sudden increase in the average power of the higher-order mode is observed. The region below the threshold is referred to in the figure as region I. Past the threshold, there exists a pump power range (region II) whereby quasi-regular temporal oscillations in the modal content occur and the average power in the higher-order mode can exceed 50% of the total signal power. In region III (chaotic regime), rapid temporal oscillations occur with almost the entire power in one mode transferring to the other mode and back and forth. These numerical simulations are in general agreement with the experimental results reported by Otto et al. [11] (refer to Fig. 4 of reference).

 

Fig. 4 Percentage of total power in higher-order mode (LP₁₁) as a function of pump power. Region I displays the “below threshold” region. Region II corresponds to a substantial content of signal power in the higher-order mode. Region III corresponds with the chaotic regime.

Download Full Size | PDF

4.2 Seed-power dependence

The work of Wirth et al. reports a 50% increase in modal instabilities threshold using a seed power of 200 W instead of 27 W [14]. To investigate this observation numerically, we varied the seed power from 25 to 175 W in increments of 25 W. The power in the LP₁₁ was set at 5% of the total seed power. We estimated the threshold by examining the modal content as a function of pump power and determining the point at which the higher-order mode power rises rapidly. Figure 5 displays the threshold as a function of seed power. We have found a linear relation over the range of seed powers investigated between the modal instabilities threshold and the seed power. The increase in threshold for a seed power changing from 25 W to 175 W is ~20%. The experimental work in [14] reports a more dramatic effect. However, the fiber used in the experiments is different than the fiber in our simulations. Also, as reported in that reference, the degree of polarization in their experiments varies as the output power changes. The effect of polarization can be one of the reasons that they observe an improvement of ~50%.

 

Fig. 5 Threshold as a function of seed power indicating a linear relation over the seed power range investigated. The power injected in the LP₁₁ mode was 5% of the total seed power.

Download Full Size | PDF

Furthermore, it has been proposed that the seed power responsible for the modal instability can be as small as quantum noise [6]. We studied the change in threshold by using different seed powers ratios in each mode in order to investigate how sensitive the threshold is to the amount of power injected into the higher-order mode. Note that here the seed power was fixed to 30 W and only the ratio of power in each mode was varied. The power injected in the higher-order mode was varied between 10^-7% of the seed power and all the way up to 10%. Threshold for all cases was determined to lie in a narrow range ~750-770 W. Therefore, the threshold is insensitive to the seeding ratio; i.e. injecting more power into the higher-order mode has little effect. Figure 6 shows the modal content at the output of the fiber vs. time for a pump power of 760 W, and a seed power of 30 W with 10^-7% injected in the LP₁₁ mode. At this pump power, we observe the beginning of an appreciable amount of energy being transferred to the higher-order mode. By pumping at slightly higher power, the amplifier will then firmly enter region II as displayed in Fig. 4.

 

Fig. 6 The fundamental mode power (Blue) and the higher order mode power (Red) at the output end of the fiber in time for seed power of 30 W whereby only 10^-7% of the seed power is injected to the higher-order mode.

Download Full Size | PDF

4.3 Gain-tailored design

The work of Robin et al. [15] reports on the observation of modal instabilities in a segmented acoustically tailored (SAT) photonic crystal fiber (PCF) amplifier. The purpose of the acoustic tailoring was to suppress SBS by creating two acoustic regions in the core: one at the center of the core and another in the six hexagons surrounding it [16]. While this design was effective in suppressing SBS, the threshold for the modal instabilities was encountered at ~530 W. A second generation SAT fiber whereby the Yb-ion concentration was tailored to have preferential overlap with the fundamental mode and reduced pump absorption was designed and fabricated [15]. The gain tailoring was achieved by ensuring that three non-adjacent hexagons lying in the outer ring were Yb-free as shown in Fig. 7. Since the fabrication of PCFs is performed using the draw and stack technique, this gain-tailored design is well-suited for PCFs.

 

Fig. 7 The core design of the gain-tailored segmented acoustic fiber. There are three acoustic regions. Three segments in the outer ring also does not contain ytterbium in order to suppress the modal instabilities.

Download Full Size | PDF

The results from the gain-tailored SAT fiber were promising and modal instabilities were mitigated appreciably. The modal instabilities threshold is increased to >990 W . The numerical model presented here is used to investigate the effect of gain-tailored designs. Two different designs have been investigated numerically and been compared to the uniformly Yb-doped core fiber. For both cases, the core and cladding dimensions are the same as the fiber described above. The first design, design A, is a similar design to the one shown in Fig. 7. The design on the numerical grid of the core is shown in Fig. 8(a). White area displays the Yb-free parts of the core. The ytterbium is removed from $\theta =\pi /3$to$\theta =2\pi /3$, $\theta =\pi$to$\theta =4\pi /3$, and $\theta =5\pi /3$to $\theta =2\pi$for $R>R_1$ where $R_1$is the radius of the inner core segment. Therefore, each doped and un-doped region in the outer ring of the core has the same area. The lengths of the gain tailored fibers were adjusted to provide for the same level of pump absorption.

 

Fig. 8 The core design of the gain tailored fibers that are numerically investigated: a) design A is similar to the gain-tailored SAT fiber whereby Yb is absent in the three white regions b) design B is such that Yb is absent in the outer ring of the core.

Download Full Size | PDF

Numerically, ~20% increase in threshold is observed for design A. Therefore, the modal instabilities were mitigated using this design but this suppression is less than the increase in the threshold observed for the gain-tailored SAT fiber used in the experiments described in [15]. We note that the simulated fiber and the fiber used in the experiments have different dimensions. In addition, the Yb doping concentration is much lower at 3.5 × 10²⁵/m³ for the latter fiber than the 6 × 10²⁵/m³ concentration used throughout our simulations. Furthermore, in the experiments, the threshold of modal instabilities is a difficult parameter to characterize due to large dependency of this threshold to system parameters which can lead to some discrepancy between experimental and theoretical results.

For design B, Yb is not present in a ring on the outer edge of the core. This design is suitable for fabrication through modified chemical vapor deposition (MCVD) process. A total of 70% of core area is doped in this design. This design is shown in Fig. 8(b) and provides a promising result of ~30% suppression of the modal instabilities threshold. Figure 9 displays the modal content of the higher-order mode averaged between 60 to 80 ms as the pump power increases for the three cases: uniformly doped core, design A, and design B. Overall, gain tailoring is an effective technique to suppress the instabilities and further investigation of various designs to optimize the suppression is warranted.

 

Fig. 9 Higher-order mode content averaged in time between 60 ms and 80 ms as a function of pump power for uniformly Yb-doped core, gain-tailored segmented acoustic design A, and gain-tailored design B whereby ytterbium is absent in annular region near edge of core.

Download Full Size | PDF

4.4 Seed modulation

We have studied the effect of seed modulation on the modal instabilities threshold and compared our time-dependent results to the work of Hansen et al. [8]. One major difference is that we do not assume in our computations a frequency offset as is assumed by the quasi-periodic models [6–9]. We assume a sinusoidal modulation such that the signal power,$P_s$, is given by:

 

Fig. 10 Higher-order mode content averaged in time as a function of pump power with modulation frequencies of 1kHz and modulation depth of 1/3000, 1/30000, and 1/300000. The solid lines correspond to best fit to the data.

Download Full Size | PDF

We also examined for the case $a=1/3000$, the modal instabilities threshold as a function of modulation frequency. Four cases were numerically investigated: $\Omega =$0.5 kHz, 1 kHz, 1.5 kHz, and 2 kHz. We found that the highest gain for the instability for the frequencies investigated occurred at $\Omega =$1 kHz, while the lowest occurred at$\Omega =$2 kHz. This would indicate that gain peak lies in the range between 0.5 kHz and 1.5 kHz. Figure 11 displays the higher-order mode content vs. pump power for the four modulation frequencies investigated.

 

Fig. 11 Higher-order mode content averaged in time as a function of pump power with modulation depth 1/3000 and modulation frequencies of 0.5 kHz, 1 kHz, 1.5 kHz, and 2 kHz. The solid lines correspond to best fit to the data.

Download Full Size | PDF

5. Summary

Using an analytical approach, we showed that a steady state solution to the heat equation could not lead to any appreciable energy transfer and a time-dependent study is required. We presented detailed numerical simulations of modal instabilities in Yb-doped fiber amplifiers using a time-dependent temperature solver in polar coordinates coupled into the optical fields and population inversion equations. Numerical convergence of the solver was presented and the importance of correct implementation of boundary conditions and computational window size was discussed. The accuracy of the spatial grid for the integration of the optical fields along the direction of propagation was confirmed trough numerical convergence and by tracking the quantum efficiency of the system as related to energy conservation. The model captures the three power-dependent regions that are characteristic of the dynamics of modal instabilities. Using this model, we found a linear relationship between threshold and seed power over a wide range. Also, we showed that the threshold is insensitive to the input modal content. Mitigation of the modal instabilities through gain tailoring was investigated by simulating two different fiber designs including a design that has been implemented experimentally in a PCF. Furthermore, we have shown that this fully time-dependent approach agrees with the quasi-periodic approach regarding the effect of signal modulation on the threshold. Finally, although this time domain approach is computationally intensive, it captures: 1) the whole range of offset frequencies that can contribute to the transfer of energy, 2) laser gain both spatially and temporally. The former has yet to be included in the quasi-periodic fully numerical approach (this approach assumes one offset frequency at a time), while the latter is not included in the semi-analytical approach that integrates over the frequency range.