## Abstract

Using numerical simulations of thermally induced mode coupling we show how the instability threshold can be substantially reduced if the pump or injected signal is modulated in the kHz range. We also show how the mode coupling gain varies with the frequency offset of the parasitic mode. We model thresholds when the source of detuned light is quantum background, amplitude modulation of the pump power, and amplitude modulation of the signal seed. We suggest several key experimental and modeling tests of our model.

© 2012 OSA

## 1. Introduction

In attempts to scale fiber amplifiers to high power, maximizing the modal area has been an effective way to suppress nonlinearities such as stimulated Brillouin and Raman scattering. Due to practical limitations on refractive index control, large mode area fibers typically support several transverse modes, and this makes them susceptible to a modal instability in which some of the light propagating in the desired LP_{01} mode is transferred into mode LP_{11}. According to several reports [1–5], this modal instability has a threshold-like behavior, with the threshold falling in the range from a few hundred to a few thousand watts of signal power, depending on fiber parameters. Further, the threshold is said to be sensitive to operating conditions [6].

In an earlier paper [7] we applied to this problem a method often used to study physical instabilities. We constructed a numerical model of the amplifier and used it to first find an equilibrium solution which, in our case, is a set of signal and pump optical fields defined along the length of the fiber. Then we introduced small perturbations to the signal equilibrium state at the input to the fiber and reran the model to determine whether those permutations grew or decayed exponentially as the light propagated through the amplifier. Based on simple physical arguments we suggested one perturbation that promised to explain the observed modal instability. Numerical modeling showed a small amount of light in the LP_{11} mode, detuned to the red by Δ*ν* ≈ 1 kHz from the strong signal light in LP_{01}, would be exponentially amplified to threshold. In this paper we will indicate this red shifted mode by adding a prime, as LP′_{11}. At threshold a small change in the input signal or pump power leads to a large change in the ratio of output powers in the LP_{01} and LP′_{11} modes.

This raises the important question of the origin of the LP′_{11} seed light, and whether its strength is sensitive to operating conditions. Here we present a detailed discussion of this question, showing how our model can account for observed threshold powers as well as a sensitivity to operating conditions. Because our concern is coherent beam combining, we define threshold as 5% of the output signal in higher order modes, this being the usual limit set for beam combining systems. In this paper we do not model behavior above this threshold, although our model is capable of that in many cases. We will also briefly compare our model with similar models published recently by other authors [4, 6]. Finally, we will suggest straightforward experiments that can critically test our proposed mode coupling mechanism and our model.

Our previous model [7] demonstrated exponential gain of LP′_{11} and allowed us to estimate thresholds, but it was too slow to permit the full length integrations that are needed to explore in detail the performance near, and especially slightly above the instability threshold. We have since greatly increased its speed so we can now perform the desired full length fiber computations. Here we present details of our improved model along with some key modeling results showing dependence of the threshold on operating conditions - namely on modulation of the input signal seed, modulation of the pump, and accidental launch of signal seed light into mode LP_{11}. These factors can act to initially populate LP′_{11}, and our model can make quantitative predictions of the resulting threshold reductions.

## 2. Mode coupling mechanism

We briefly review the gain mechanism built into our model. The LP_{01} and LP′_{11} modes interfere to produce an irradiance grating that moves along the fiber with a velocity proportional to Δ*ν*, the frequency detuning between them. This irradiance grating creates a moving temperature grating which in turn produces a moving refractive index grating via the thermo optic effect. The temperature and index gratings have a phase lag relative to the irradiance grating because of the time lag introduced by thermal diffusion. The moving, phase-shifted index grating couples light from the strong LP_{01} mode into the weak LP′_{11} mode.

This mode coupling mechanism is a form of stimulated (near) forward thermal Rayleigh scattering. It is categorized as such because it involves a diffusive temperature wave rather than a propagating wave such as the acoustic wave in stimulated Brillouin scattering (SBS), and it does not involve a material resonance like stimulated Raman scattering (SRS) [8]. However, like SBS or SRS it can produce the large exponential gain usually associated with a sharp threshold. Like other stimulated Rayleigh processes, the one we propose has a dispersion-like gain profile centered at zero offset frequency [8, 9]. The time scale, and thus the inverse of the frequency offset at maximum gain, is approximated by the thermal diffusion time across the core radius, (*r*^{2}*ρC*/*K*) where *r* is the radius, *ρ* the density, *C* the specific heat, and *K* the thermal conductivity. For typical fiber amplifiers this time falls in the range 0.3–3 ms. An approximate semi-analytic version of this model was recently published by Hansen *et al* [10]. They illustrated the variation of gain with frequency offset for several fiber designs, and also estimated the instability thresholds.

One perhaps subtle point is that the signal spectrum can be broad compared with the 1 kHz scale frequency shifts incorporated into our model. In fact most signal seeds will be much broader than 1 kHz, for example when the signal consists of a train of short pulses or when it is cw light that is phase modulated to suppress SBS. The important point is that any signal spectrum, broad or narrow, will be shifted in its entirety by approximately 1 kHz when it scatters from the moving index grating. Further, interference between a broad signal spectrum and its frequency shifted replica will produce the strong mode beats responsible for creating the moving temperature gratings. Therefore, spectrally broad signals behave exactly like single frequency signals in those respects that matter to our model. If the two replicas acquire a time offset comparable to their inverse line width due to differing group velocities for the two modes, the interference will weaken. Mode dispersion delays are on the order of 1 ps/m so the line width limit is in the range of several hundred GHz, much broader than the line width limit imposed by coherent beam combining. Hence our model simulates equally well both narrow and broad signal spectra, provided they are frequency shifted and have limited line widths.

## 3. Model details

The methods we use for this paper are conceptually identical to those described earlier [7]. We use a fast Fourier transform (FFT) beam propagation method (BPM) [11, 12] in order to make the model as general as possible. This allows us to include a wide variety of physical effects. For example, it automatically includes all transverse modes, and it accounts for all thermal lensing. The split-step propagation procedure first updates the time varying optical field, which includes all frequency components and all modes, by adding the time dependent propagation phases due to the guiding index step and due to the thermo optic phase computed from the transverse, time varying temperature profile. The laser gain or loss is also added to the fields in this part of the split step. In the second part of the split step each time slice of the wave is propagated by Fourier transforming the field, adding the appropriate propagation phases, and performing the reverse Fourier transformation. At each *z* location we solve the time-dependent temperature equation over one beat cycle of duration 1/Δ*ν*. The heat *Q*(*x*, *y*,*t*) at each location is computed by finding the absorbed pump power from the steady state solutions to the population/gain equations. The resulting time-dependent temperature profile is converted to a time-dependent refractive index change Δ*n*(*x*, *y*,*t*) by multiplying by the thermo optic coefficient. This index change is used in the beam propagation integration. When we analyze the signal field into modes for each time slice we use the low power modes as the basis set. The time varying amplitude of each mode is further analyzed into frequencies (0, ±Δ*ν*, ±2Δ*ν*, ...) by Fourier transformation. This model can accommodate mode-specific losses due to leakage from the core or due to scattering. It can also include heating due to photo darkening or other types of absorption, but they are not included here.

In computing the laser gain and the deposited heat we use the steady state solution for the upper state population fraction

*σ*s are absorption and emission cross sections for the pump and signal,

*ν*and

_{p}*ν*are the pump and signal optical frequencies,

_{s}*I*(

_{s}*x*,

*y*,

*t*) and

*I*(

_{p}*x*,

*y*,

*t*) are irradiances, and

*τ*is the upper state lifetime. The two dimensional heat equation

*T*(

*x*,

*y*,

*t*).

We make certain approximations for the sake of computational speed. As mentioned, we use a steady-state population/heating model in which the population instantaneously follows the local signal and pump irradiances. We justify this by noting that at the power levels of interest the effective lifetime of the upper state Yb^{3+} ions is in the range of 10 *μ*s, which is much shorter than the 1 ms beat times of interest. We also compute the temperature in two dimensions rather than three. This means we cannot account for longitudinal heat flow, but we believe this is a small effect. Additionally, we place the thermal boundary approximately two core diameters from the core center rather than at the actual fiber edge. We justify this placement by noting that moving the thermal boundary out by a factor of three does not change our model results significantly but slows the model by a large factor. We also set the temperature rise at the boundary to zero. In reality, thermal transport from this boundary to the physical boundary results in a temperature rise at the model boundary, and this adds a nearly uniform temperature rise in the vicinity of the core. Such a uniform rise would affect neither thermal lensing nor mode coupling directly, but it could still have a small influence on amplifier performance because the Yb^{3+} absorption and emission cross sections as well as the upper state lifetime change slightly with temperature over the range of interest [13, 14]. The absorption and emission cross sections at 1064 nm are relatively insensitive to temperature. However, the pump cross sections at 976 nm fall by 10–20% over 0–400°C, and the lifetime falls by 10% over 0–400°C, and much faster above 500°C. Temperature dependent cross sections and lifetimes could be included in our model, but presently they are not.

In our previous model we integrated the temperature equation in time using an explicit integration method. The time step required for convergence was short, making the model slow. We have since implemented an alternating direction implicit (ADI) method [15], and also a steady-periodic Green’s function (GF) method [10, 16], both of which allow larger time steps and provide integration times of a few hours per meter of fiber. The ADI method has the advantage that any signal and pump input modulations can be used. The GF is faster but requires steady-periodic behavior. Some mathematical details of the two new methods are presented in the Appendices. Mode coupling gains computed using these two methods agree within 1%.

#### 3.1. Using the ADI method

In the results presented here based on computations using the ADI method we use a single cycle of the time-dependent heating. We string together copies of the heat cycle to make a periodic time sequence containing many cycles. This is used as the heat term in the ADI computation. We then run the ADI integration until there is negligible change in the temperature response from one cycle to the next, and we keep the results for the last cycle. This guarantees a steady-periodic response is being modeled. This allows periodic modulation of the signal and pump inputs with frequency components (0, Δ*ν*, 2Δ*ν*, ...). If desired, this method can also be used to study transient effects that follow sudden changes in the injected signal or pump light. In such cases we would not enforce the periodicity or stability criteria.

#### 3.2. Using the Green’s function method

Unlike the ADI method, the Green’s function method is restricted to steady-periodic heating. As described in Appendix B, Green’s functions are calculated for sinusoidal heating at frequencies *ω _{m}* = 2

*π*(0, Δ

*ν*, 2Δ

*ν*, ...). In applying this method, the heat at each grid point is analyzed into frequency, amplitude, and phase, and the appropriate Green’s function is used to find the 2D temperature response. This is repeated at all desired frequencies for all heated (all Yb

^{3+}-doped) grid points. Because mode coupling is directly driven predominantly by the component of temperature oscillations at frequency 1Δ

*ν*, usually the threshold computations can be performed using only the zero frequency Green’s function to account for thermal lensing, plus the 1Δ

*ν*Green’s function to handle mode coupling. For computations of above threshold behavior it is sometimes necessary to include additional, higher frequency Green’s Functions.

## 4. Gain spectrum

We have mentioned that the gain of a stimulated Rayleigh process such as we propose is expected to have a dispersion-like gain profile. Figure 1 illustrates this with gains we compute using the Green’s function method for a co-pumped amplifier with 30 *μ*m core diameter. Further details of this amplifier are given in [7] except here we use 1 kW of pump and compute the gains over 15 mode beat lengths at the input end of the amplifier. These curves show the total gain of LP′_{11} including laser gain and mode coupling gain. The laser gain, which accounts for the gain at zero frequency shift, is relatively small compared with the peak mode coupling gain. The shape of the gain curve changes with *z* position due to evolving gain saturation effects and due to thermal lensing. However, the exact details of the fiber are not important here as our point is to illustrate a typical dispersion-like gain profile, demonstrating gain on the Stokes side and loss on the anti-Stokes side. Note that in this example the mode coupling gain of LP′_{11} is approximately twice that of LP′_{02}. In fact, we find the LP′_{11} gain is generally greater than that of any other mode, explaining why, to our knowledge, LP′_{11} is the intruding mode that is always observed in the laboratory.

## 5. Model parameters

The parameters we used in the modeling presented below are chosen to be nearly equal to those used by Ward *et al* [4], and are listed in Table 1. The spatial grid is 64 × 64, the time cycle is gridded into 64 points, and the propagation step is 6 *μ*m. Unlike Ward *et al* [4] we make the doping diameter equal to the core diameter, and we shorten the fiber to 0.8 m because there is little gain in the first half of their 1.6 m version.

## 6. Sources of LP′_{11} seed light

The gain of LP′_{11} that is required to reach threshold depends on the power in that mode early in the fiber, and this can vary with operating conditions. In our previous paper we did not discuss in detail the source of the initial light. Here we consider four possible sources.

#### 6.1. Quantum noise

One possibility is that LP′_{11} is populated solely by quantum noise. If this noise is added at the input to the final amplifier stage, this would give the lowest possible starting level and thus the highest possible instability threshold. Quantum noise provides approximately one photon per beat cycle, corresponding to a power of (*hν*Δ*ν*) which, for a 1 kHz offset and linewidth, is approximately 2 × 10^{−16} W. We use 10^{−16} W in simulations aimed at finding this highest possible instability threshold.

In practice the laser and preamplifiers supplying the signal seed will almost certainly have a higher noise level than this, in which case the threshold will be somewhat lower. Some of the noisy signal input will find its way into LP_{11} either through imperfect mode matching or through scattering at fiber imperfections. The amount of light injected into LP′_{11}, and thus the threshold, will be sensitive to the level of seed noise and to the launch conditions.

#### 6.2. Spontaneous thermal Rayleigh scattering

A second possible source of LP′_{11} seed light is thermally induced fluctuations in the refractive index near the input end of the amplifier. Such fluctuations could scatter light from LP_{01} or LP_{11} into LP′_{11}. This is a type of spontaneous thermal Rayleigh scattering. However, based on measured attenuations of less than 1 dB/km for temperature independent Rayleigh scattering, combined with the fact that a temperature dependent contribution is unobservable in the background of the temperature independent part [17], we estimate that this process provides at most only slightly more seed light than does quantum noise. A signature of this contribution would be a reduction in the threshold with an increase in the temperature of the signal input end of the amplifier.

#### 6.3. Pump modulation

A third possible source of LP′_{11} seed light, and one we think may dominate, is pump modulation. If the pump light is slightly amplitude or spectrum modulated in the 1 kHz band, and if the signal is launched with a fraction of its power in LP_{11} at the fiber input, gain modulation of the signal due to the pump modulation produces frequency side bands on the signal, including the portion in LP_{11}. Two symmetric frequency sidebands are created. The blue shifted one creates an optical interference and temperature wave that moves upstream; the red shifted one creates an optical interference and temperature wave that moves downstream. Each sideband interacts strongly only with the temperature wave that co-travels with its interference wave. The phase shift between the temperature wave and the optical interference wave is opposite for the two sidebands, leading to amplification for the red shifted sideband (LP′_{11}) and deamplification for the blue shifted sideband as illustrated in Fig. 1. The threshold would thus be reduced in proportion to the degree of pump modulation and in proportion to the power inadvertently launched into or scattered into LP_{11}.

#### 6.4. Signal modulation

Similarly, if the injected signal is slightly amplitude or phase modulated with a side band in the −1 kHz range at input to the fiber, the side band on the fraction that is inadvertently injected into LP_{11} will populate LP′_{11}. The strong carrier light in LP01 will act as the pump for the weak LP′_{11} light, as usual, again leading to a reduced threshold.

## 7. Threshold calculations

As stated above, we use 5% in higher order modes as the threshold criterion of modal instability. We first varied the frequency offset to find that the maximum mode coupling gain occurs near Δ*ν* = −525 Hz, and we use this frequency offset in all subsequent threshold computations. We perform multiple simulations to accurately determine the threshold powers. The signal light is launched with the profile of the low power modes rather than the slightly thermally lensed modes. The propagating and output light is analyzed into low power modes as well. The following calculations illustrate how the threshold is strongly influenced by the source of the LP′_{11} light, whether it is from quantum noise, pump amplitude modulation, or signal amplitude modulation.

#### 7.1. Quantum noise

Figure 2 shows model results when quantum noise is used as the seed to the amplifier specified in Table 1. The LP′_{11} mode is seeded with 1 × 10^{−16} W and mode LP_{11} is unseeded. The counter pump input power is 353 W of which 7 W is unabsorbed. The input LP_{01} power is 50 W; the total signal out including light in all modes is 356 W with 323 W in LP_{01} and 8.4 W in LP′_{11}. This case yields the highest possible threshold at 356 W of signal. At this power level there is significant thermal lensing. The effective area of LP_{01} changes from 2770 *μ*m^{2} at the signal input end to approximately 1700 *μ*m^{2} at the pumped end. Our model includes this effect.

#### 7.2. Pump modulation

Figure 3 shows modal powers at threshold when the pump is sinusoidally amplitude modulated at 300 parts per million peak-to-peak at 525 Hz. The total signal input is again 50 W, but 250 mW of that is in LP_{11} with 49.75 W in LP_{01}. The pump input power at threshold is reduced from 353 W to 288 W with 2.3 W unabsorbed pump. The total signal output is reduced from 356 W to 302 W with 12 W in LP ′_{11}. As Fig. 3 shows, LP′_{11} is quickly populated to a power of approximately 10^{−12} W as a side band of LP_{11} due to gain modulation. This is much larger than the 10^{−16} W used to simulate quantum noise, so the threshold is reduced. Obviously, larger pump modulations would lead to yet lower thresholds. The blue shifted LP_{11} light is also continuously created by gain modulation, but it is simultaneously deamplified by mode coupling. Consequently, this mode remains at least 10^{6} times weaker than LP′_{11}. We think the kink in the LP′_{11} curve is due to two interfering contributions to the LP′_{11} power increments; one contribution is from mode coupling gain that transfers power from LP_{01}, the other is from pump modulation induced transfer from LP_{11}. The rapid oscillations evident in some curves of Fig. 3 are related to thermal lensing and are not an essential feature of mode coupling gain.

#### 7.3. Signal modulation

The case of signal amplitude modulation without pump input modulation is difficult to model precisely since modulation of the signal induces modulation of the pump. Because we integrate from the pump output end to the pump input end we must precompensate by specifying the pump output modulation that will lead to an unmodulated input. We have not succeeded in perfect precompensation but we have come close. The results are shown in Fig. 4. The input signal is sinusoidally amplitude modulated at 30 parts per million peak-to-peak. In this simulation we apply this amplitude modulation to the input light of both the 49.75 W in LP_{01} and the 0.25 W in LP_{11}. The AM sideband of LP_{11} forms the 10^{−12} W input LP′_{11} which is amplified to threshold by the carrier LP_{01} light. We find the threshold signal output power is 275 W total, with 5.1 W in LP′_{11}, corresponding to a pump input power of 257 W and output power of 1.4 W. This is a large reduction from the 353 W pump threshold for quantum noise alone.

## 8. Comparison with measured thresholds

We have computed the highest threshold for one particular amplifier, as well as demonstrating that substantial threshold reductions can be caused by the introduction of light into LP′_{11} by modulating the pump or the signal seed, in combination with a slight population of LP_{11}. It is interesting to compare these results with some reported measured thresholds.

Several experimental studies have reported instability threshold powers. They seem to fall in the range 100–300 W for amplifiers ≈ 1 m long with core diameters of 40 *μ*m or greater [1, 18, 19], and in the range 500–1500 W for ≈ 10 m long amplifiers with core diameters less than 40 *μ*m [1, 20]. The thresholds computed above for one example of a short, large diameter fiber ranged from 275 to 350 W of signal, in reasonable agreement with the measured values. Furthermore, the model qualitatively matches the observation of sharp power thresholds [1, 4]. An increase of 10% in pump power is usually sufficient to fully switch the output signal from LP_{01} to LP′_{11}.

Experiments have also measured the modulation frequency of the modal powers in amplifiers operated slightly above threshold [18, 21]. The oscillation frequencies were approximately 500 Hz for a fiber with mode field diameter of 75 *μ*m and approximately 2 kHz for a fiber with mode field diameter of 29 *μ*m. Ward *et al* [4] report an oscillation frequency of 2 kHz for a fiber with a mode field diameter of 30 *μ*m. This trend is in accord with our model, and the values are quite close to the values predicted by the model. The oscillations also have measured spectral widths in reasonable agreement with the gain bandwidths given by our model. The measured trend of reduced threshold in the presence of a small amount of photo darkening [22] is also matched qualitatively by our model, although we do not present those results here.

Our model results account for LP′_{11} having the lowest threshold, in agreement with the universal observation of only LP_{01} and LP_{11} above threshold. Finally, we have explained how instability thresholds can be sensitive to the operating conditions, specifically pump or signal modulations combined with slight initial population of LP_{11}.

## 9. Suggested experiments

Here we propose some critical experimental tests of the gain mechanism and of our model.

#### 9.1. Look for frequency shifts

The most direct way to look for the frequency shifts that are essential to our model is to study in detail the behavior of the output beam near the modulation instability threshold. One could look at the output beam simultaneously at two symmetric locations along the line of symmetry. They could be positioned near the peak irradiance points of the LP_{11} mode. The pump power should be turned up gradually until mode distortion is barely noticeable. Out-of-phase amplitude modulation would indicate beating between frequency offset modes with even and odd symmetry in the plane of the detectors - presumably modes LP_{01} and LP′_{11}. This measurement should determine whether or not LP′_{11} is present. Similar experiments have been performed using a single detector, and they showed sinusoidal modulation [4]. Better documentation, including a thorough catalog of fiber parameters and detector placement and aperture, plus detection with dual detectors would be an important improvement.

#### 9.2. Measure/modify the pump modulation

As we illustrated above, a small amount of pump modulation could significantly reduce the instability threshold. One could measure the degree of pump power modulation in the 1 kHz band. To see whether this is an important contributor, the pump power could be deliberately modulated. One could also alter the position of the seed beam to vary the fraction injected into LP_{11} as a further test of this idea. If power modulation is found to affect the threshold, the modulation frequency could be varied to map the instability threshold and the mode coupling gain.

It may be that pump spectral modulation rather than power modulation is important. In that case one could look for fluorescence modulation to determine whether pump modulation near 1 kHz is significant.

#### 9.3. Measure/modify the signal modulation

We also described how a small sideband on the injected signal, combined with inadvertent seeding of LP_{11} will lead to the population and amplification of LP′_{11}. One could look for amplitude or phase modulation of the signal before injection, or one could deliberately add such modulation, to test this possibility. Tuning the modulation frequency would map the gain line shape.

#### 9.4. Numerical experiments

The light transit time through the amplifier is much shorter than the mode beat time so the order of integration in our model (*t* then *z*) can be reversed. Rather than integrating over the full *t* range followed by stepping in *z*, one could integrate over the full *z* range followed by stepping in *t*. Use of this reversed order has the advantage that the heat diffusion equation can be integrated in three dimensions rather than two dimensions as we have done. This would make it possible to test whether longitudinal heat flow has any impact on mode instability. Two recent papers [4, 6] described such a model. The physics included in those models is nearly identical to that in our model. However, they were applied only to studies of transient mode coupling following abrupt turn on of the pump [6] or abrupt changes to the seed light [4]. Applying them to a case with a steady-periodic perturbation would allow a test of the role of longitudinal heat flow in the instability mechanism we are touting.

## 10. Conclusions

We have presented a detailed model of mode instability that is capable of predicting modulation frequencies and threshold powers for step index fibers. Replacing the FFT-based BPM with a more general BPM would make the model suitable for accurate modeling of photonic crystal fibers as well. The model predictions in every case that we are aware of match qualitatively the laboratory results. Unfortunately, detailed quantitative comparisons are hampered by incomplete documentation of the experiments. We hope that in the future a more thorough catalog of fiber parameters will be included, and that this standard will be applied to numerical models as well.

In this study we used narrow linewidth perturbations to compute the gain spectrum. The chosen frequency offsets were those giving the highest gain. In the laboratory, the width of the gain profile will depend on the nature of the input signal and pump light. If the source of LP′_{11} is quantum noise, the Stokes spectrum should be broad. If it is pump or signal modulation the Stokes spectrum will reflect the relevant modulation spectrum. By showing how small amounts of modulation of the pump or the signal seed can substantially reduce the instability threshold we hope to motivate measurements of such modulations and studies of their influence on thresholds. If such modulations are present, the path to maximizing thresholds is clearly their minimization.

Although our model has been successful in qualitative comparison with experiments, stronger confirmation is desired. We suggested key experiments to test the model. If it is validated, we can proceed with parametric studies to learn how to best suppress the mode instability while still combating other nonlinearities such as SBS and SRS. This will require further speed up of the model. We anticipate that the Green’s function approach in particular can made much faster using a computer with a large number of processors in place of our single processor model.

## 11. Appendix A: ADI method

ADI is short for alternating direction implicit. This method is outlined in Numerical Recipes [15], and we present our implementation here. In this method the temperature is updated over Δ*t* in two steps of Δ*t*/2 each. The temperature on the boundary is maintained at zero, so we operate on the (*N _{x}* − 2) × (

*N*− 2) interior points of the

_{y}*T*grid. In the first half-step we use explicit integration in the

*y*–direction and implicit integration in the

*x*–direction. This yields the first half-step difference equation

*A*is the (

_{x}*N*− 2)×(

_{x}*N*− 2) tridiagonal matrix with diagonal elements (1 +

_{x}*λ*) and off-diagonal elements −

_{x}*λ*/2, and

_{x}*B*is the (

_{y}*N*− 2) × (

_{y}*N*− 2) tridiagonal matrix with diagonal elements (1 −

_{y}*λ*) and off-diagonal elements

_{y}*λ*/2, where

_{y}The second half-step in time reverses the implicit/explicit order to yield a similar difference equation,

*A*and

_{y}*B*are the same as above except redimensioned as appropriate, and with

_{x}*λ*and

_{x}*λ*exchanged.

_{y}Combining these two half-step equations gives the expression used to update the temperature over the full time step

*Q*

^{t+Δt/4}and

*Q*

^{t+3Δt/4}from the heats at

*t*and

*t*+ Δ

*t*by linear interpolation.

## 12. Appendix B: GF method

The idea of the Green’s function method is that the two-dimensional (*x*, *y*) temperature profile due to one unit of heat at a chosen grid point (*x*′, *y*′) can be computed once for (*x*′, *y*′) point and used repeatedly with the proper weighting and phase throughout the integration. The temperature profile due to an actual heating profile is constructed by summing the properly weighted and phased Green’s function values from all heated grid points. This procedure can be implemented in two dimensions for steady state heating, but also for steady-periodic heating [16] where the heat at each grid point oscillates with a frequency *ω*, phase *ϕ*, and amplitude *Q*.

The steady-periodic Green’s function for a clamped temperature (*T* = 0) on the (0 = *y* = *W*) and (0 = *x* = *L*) boundaries is [16]

In our model we use a fixed time interval of (1/Δ*ν*) so steady-periodic heating can have frequencies of *ω _{m}* = 2

*πm*Δ

*ν*where

*m*= 0, 1, 2.... The functions

*G*(

*x*,

*y*,

*ω*|

_{m}*x*′,

*y*′) are computed for each (

*x*′,

*y*′) location that will be heated, which in our case means all that are Yb doped, and for each

*ω*of interest. Calculation of

_{m}*G*(

*x*,

*y*,

*ω*|

_{m}*x*′,

*y*′) for each

*ω*of interest is performed once before starting the model run.

_{m}In applying the GF we Fourier analyze the heat at each grid location to find complex Fourier coefficients *q*(*x*′, *y*′, *ω _{m}*) and use them as inputs to the appropriate GF. We compute

*T*(

*x*,

*y*,

*t*) using

*q*(

*x*′,

*y*′,

*ω*) are the complex Fourier coefficients of the heat at point (

_{m}*x*′,

*y*′), and

*max*is the highest frequency component

*ω*to be included.

_{max}## Acknowledgments

We thank Shadi Naderi of Air Force Research Laboratory for pointing out an important improvement in the ADI code. This effort was supported under funding from the Air Force Research Laboratory Directed Energy Directorate.

## References and links

**1. **T. Eidam, C. Wirth, C. Jauregui, F. Stutzki, F. Jansen, H.-J. Otto, O. Schmidt, T. Schreiber, J. Limpert, and A. Tünnermann, “Experimental observations of the threshold-like onset of mode instabilities in high power fiber amplifiers,” Opt. Express **19**, 13218–13224 (2011). [CrossRef]

**2. **T. Eidam, S. Hanf, E. Seise, T. V. Andersen, T. Gabler, C. Wirth, T. Schreiber, J. Limpert, and A. Tünnermann, “Femtosecond fiber CPA system emitting 830 W average output power,” Opt. Lett. **35**, 94–96 (2010). [CrossRef]

**3. **C. Jauregui, T. Eidam, J. Limpert, and A. Tünnermann, “Impact of modal interference on the beam quality of high-power fiber amplifiers,” Opt. Express **19**, 3258–3271 (2011). [CrossRef]

**4. **B. Ward, C. Robin, and I. Dajani, “Origin of thermal modal instabilities in large mode area fiber amplifiers,” Opt. Express **20**, 11407–11422 (2012). [CrossRef]

**5. **C. Jauregui, T. Eidam, H.-J. Otto, F. Stutzki, F. Jansen, J. Limpert, and A. Tünnermann, “Temperature-induced index gratings and their impact on mode instabilities in high-power fiber laser systems,” Opt. Express **20**, 440–451 (2012). [CrossRef]

**6. **C. Jauregui, T. Eidam, H.-J. Otto, F. Stutzki, F. Jansen, J. Limpert, and A. Tünnermann, “Physical origin of mode instabilities in high-power fiber laser systems,” Opt. Express **20**, 12912–12925 (2012). [CrossRef]

**7. **A. V. Smith and J. J. Smith, “Mode instability in high power fiber amplifiers,” Opt. Express **19**, 10180–10192 (2011). [CrossRef]

**8. **R.W. Boyd, *Nonlinear Optics*, 2nd ed. (Academic Press, 2003).

**9. **O. Okusaga, J. Cahill, A. Docherty, W. Zhou, and C. R. Menyuk, “Guided entropy mode Rayleigh scattering in optical fibers,” Opt. Lett. **37**, 683–685 (2012). [CrossRef]

**10. **K. R. Hansen, T. T. Alkeskjold, J. Broeng, and J. Laegsgaard, “Thermally induced mode coupling in rare-earth doped fiber amplifiers,” Opt. Lett. **37**, 2382–2384 (2012). [CrossRef]

**11. **M. D. Feit and J. A. Fleck, “Computation of mode eigenfunctions in graded-index optical fibers by the propagating beam method,” Appl. Opt. **19**, 2240–2246 (1980). [CrossRef]

**12. **A. V. Smith and M. S. Bowers, “Phase distortions in sum- and difference-frequency mixing in crystals,” J. Opt. Soc. Am. B **12**, 49–57 (1995). [CrossRef]

**13. **S. W. Moore, T. Barnett, T. A. Reichardt, and R. L. Farrow, “Optical properties of Yb^{3+}-doped fibers and fiber lasers at high temperature,” Opt. Commun. **284**, 5774–5780 (2011). [CrossRef]

**14. **T. C. Newell, P. Peterson, A. Gavrielides, and M. P. Sharma, “Temperature effects on the emission properties of Yb-doped optical fibers,” Opt. Commun. **273**, 256–259 (2007). [CrossRef]

**15. **W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery, *Numerical Recipes in C*, 2nd ed. (Cambridge Univ. Press, 1992).

**16. **K. D. Cole, “Steady-periodic Green’s functions and thermal-measurement applications in rectangular coordinates,” Mech. Eng. Faculty Pubs., Univ. Nebraska, Lincoln, Paper 49, http://digitalcommons.unl.edu/mechengfacpub/49.

**17. **P. C. Wait and T. P. Newson, “Landau-Placzek ratio applied to distributed fiber sensing,” Opt. Commun. **122**, 141–146 (1996). [CrossRef]

**18. **F. Stutzki, F. Jansen, T. Eidam, A. Steinmetz, C. Jauregui, J. Limpert, and A. Tünnermann, “High average power large-pitch fiber amplifier with robust single-mode operation,” Opt. Lett. **36**, 689–691 (2011). [CrossRef]

**19. **F. Stutzki, H.-J. Otto, F. Jansen, C. Gaida, C. Jauregui, J. Limpert, and A. Tünnermann, “High-speed modal decomposition of mode instabilities in high-power fiber lasers,” Opt. Lett. **36**, 4572–4574 (2011). [CrossRef]

**20. **N. Haarlammert, O. de Vries, A. Liem, A. Kliner, T. Peschel, T. Schreiber, R. Eberhardt, and A. Tüunnermann, “Build up and decay of mode instability in a high power fiber amplifier,” Opt. Express **20**, 13274–13283 (2012). [CrossRef]

**21. **H.-J. Otto, F. Stutzki, F. Jansen, T. Eidam, C. Jauregui, J. Limpert, and A. Tünnermann, “Temporal dynamics of mode instabilities in high-power fiber lasers and amplifiers,” Opt. Express **20**, 15710–15722 (2012). [CrossRef]

**22. **M. Laurila, T. T. Alkeskjold, M. M. Jørgensen, S. R. Petersen, J. Broeng, and J. Laegsgaard, “Highly efficient high power single-mode fiber amplifier utilizing the distributed mode filtering bandgap rod fiber,” in *Fiber Lasers IX: Technology, Systems, and Applications*, E.C. Honea and S.T. Hendow eds., Proc. SPIE 8237, 8237–8265 (2012).