1. Introduction

The Soliton self-frequency shift (SSFS)^1,2 has been demonstrated as an important mechanism for the supercontinuum (SC) generation with femtosecond pulses in optical fibers ^3,4. Specifically, many effects implied in SC generation would not be possible without the SSFS ^{5,6,7,8,9,10,11}. Its characteristic large frequency shift has been exploited for the fabrication of infrared (IR) sources ^12,13,14, even its optimization has been reported ^15,16,17.

In previous works, we have shown that a fs-pulse can simultaneously generates several pre-defined spectral peaks by means of dispersive waves in the normal GVD¹⁸ or by soliton self-frequency shift in the anomalous GVD¹⁹. It is useful for potential applications of optical coherence tomography (OCT). Now, regarding the same kind of applications, we present a computational optimization with the use of a genetic algorithm (GA) to obtain not only a tunable frequency convertor (see, e.g., Refs. ^20,21), but also achieving the maximum spectral conversion possible in the fixed channel regarding the scope of the initial set of input pulse parameters. This frequency convertor is useful for applications where high power is demanded.

We simulate pulses propagation in a commercial highly nonlinear photonic crystal fiber, NL-2.4-800 PCF (see Ref. ¹⁹ for details of cross sectional geometry), exhibiting SC generation at the Ti:Sapphire laser wavelengths. This method finds the optimal input pulse parameters, namely central wavelength, λ0, temporal width, T0, and peak power, P0, that maximize the output SC power in a spectral channel of fixed width and selected central frequency, ωc. Our typical channel width chosen here narrow enough so the optimal ωc will be around the carrier frequency of the most powerful, i.e. firstly ejected Raman soliton in the IR spectral region ²². The inverse problem, i.e. the design of PCFs to optimize the SC has indeed been previously solved satisfactorily in a wide range of situations ^{23,24,25,26,27}. Our interest in the IR region is motivated specifically by applications in OCT ^{28,29,30,31,32,33,34}.

2. Pulse propagation and genetic algorithms

Using a Fourier split-step method, we simulate the propagation of optical pulses with complex amplitude A(z, t) by integrating numerically the GNLSE ³,

where z is the axis coordinate along the fiber propagation, the dispersion coefficients βq’s (up to q = 10 in this work computed with a FEM solver) account for the linear fiber dispersion at the pump frequency ω0=2πc/λ0. Nonlinearity is included through the parameter γ and the response function RT≡1-fR Δ(T)+fRhR(T)Π(T), where fR=0.18, hR is the commonly used Raman response of silica ³⁵, and Δ(T), Π(T) are the Dirac, Heaviside functions, respectively. The dispersion coefficients βq were computed using Optiwave³⁶.

The input pulse used in this work belong to a realistic laser source, the form is A(z=0,t)=P0 sech (t/T0), where the intensity full width at half maximum TFWHM=2T0ln1+2 and P0 as the peak power. The soliton order is computed by N=T0γP0/ |β2|.

Optimization algorithms have been developed in order to solve problems involving multiple variables in which solution seems to be non trivial (see Refs. ^37,38,39 for reviews on the topic). A GA is a evolutionary computational algorithm for optimization which makes evolve an initial population of individuals in order to find global minima when a number of generations is generated ^40,41. Each individual is the result of the evaluation of a set of parameters (λ0, P0 and T0 in this work) and applies a minimization strategy to find the solutions taking the minimum value of the fitness function ϕ. The use of GAs in the fiber-design context is not new ^23,24. Previously, some works have made use of GAs for the optimization of the dispersion management ^26,27, even the design of PCF structural parameters for SC generation ²⁵. We use a GA using pre-defined functions of Matlab ⁴². In the first stage, the GA starts generating a randomly initial population of p = 50 individuals. In the second stage, the most promising individuals generated in the first stage are allowed to reproduce to determine the next generation of individuals according to the pre-established evolution rules of Mutation M^, Cross-over X^ and Random R^. The Mutation operator M^ [5%in this work] makes possible the change of one chromosome of the promising individuals. The Cross-over operator X^ [95% in this work] acts as a parent selection method, i.e. a mix of two promising individuals (parents) and their rates to obtain the offspring for the next generation. Most parent selection methods are stochastic in order to keep the diversity of the population, preventing premature to a sub-optimal convergence solution. The Random operator R^ is selected if population is smaller than the threshold value c (c < p) and by the genetic operator M^ for chromosome changes. At the end (after m = 150 in our case), the best individual is picked of the evolved population and chosen as optimal solution. The optimization process is depicted in Fig. 1.

Figure 1. Diagram of the operation of the GA. In first stage a population of possible solutions is generated randomly. In the second stage, new individuals are created by M^, X^ or R^ and each one is compared with the worst solution created in the initial population, if the new individual is better, then it is selected as a new individual in the next generation instead the worst solution, otherwise is dismissed. The process causes the mean fitness improvement in the population. At Stage 3 (after m = 150 in our case), the best individual is picked of the evolved population and chosen as optimal solution. 

3. Raman frequency conversion

It has been proved that high axial resolution in OCT systems is aimed in the spectral region of 800 nm to 1400 nm ^43,44. Aditionally, the NIR II light decrease in scattering and increase in transparency of the biological tissues over the NIR range ⁴⁵. Moreover Gaussian spectral shapes avoid spurious structures in OCT images ⁴⁶. For these reasons, IR-Raman soliton can be considered a very good option to OCT applications.

In Sec. 3.1, we search the optimal parameters, λ0, T0 and P0, by implementing the fitness function defined as

where |Ã(ω')|2 is the output pulse power amplitude, 2δ=50 nm is the chosen spectral channel width, and ωc is the central frequency, determined by inspection of the output spectrum, to minimize ϕ(ωc;δ) and, consequently to increase the efficiency of conversion given by

In order to prove the convenience of our method, it is made an exhaustive search of the best fitness value by scanning the entire ranges of parameters. It is shown in Sec. 3.2

3.1. Optimal solution using genetic algorithms

The optimization consisted in the search of parameters that originate the maximum output power in each selected spectral channel on the NIR II using the GA to vary the pulse parameters P [1, 15] kW, λ [750, 850] nm and T [30, 110] fs. The range of values used in this work are attainable in realistic Ti:Sapphire lasers. In this particular optimization problem, each individual evaluation typically required 90 s what amounted for about ∼3.75 h of CPU time to perform a single run of the GA with 150 individuals. The set of all the individuals generated by the GA in each channel optimization is shown as a example in a 3D “cloud" graphic in the space of parameters in Fig. 2(a).

Figure 2. (a) Cloud of m = 150 individuals generated by the GA for a channel of λc = 1225 nm. The color of the points shows the fitness ϕ value of each individual, lighter points have smaller (thus, better fitness values) and the optimal individual is represented as the red point sited in the coordinates of the three optimized parameters. (b) Typical fitness value evolution of executions in chronological order of generation during the optimization procedure with its mean fitness value (black line) and the minimum global value (red line), dashed magenta vertical line separate the best individual within the first random individuals (m = p _th = 50) generated in the stage 1 of the process. 

This “cloud" of individuals corresponds to all solutions generated by the GA, their fitness function being represented by the color code bar. Lighter points have smaller fitness values (thus, better) than darker ones. We observe that there exists a zone where the GA tends to accumulate points. It is precisely in this region where the best fitness value (red point) is found. It is worth mentioning that in general these regions could contain more candidates to optimal solutions than those eventually selected by the GA. Thus, keeping track of these “quasi-optimal" individuals can also be of great interest from the physical point of view since they can provide extra-local minima of the fitness function not considered in a preliminary physical analysis of the optimization scenario. Once the local minima have been detected, a more accurate search around them combining GA strategies and other optimization techniques can be performed in order to find a better minimum of the fitness function.

Figure 2(b) shows clearly the “dynamical" improvement in the fitness value as the GA evolves [referred to in Fig. 1]. The initial “optimized" value is obtained in the stage 1, when the initial population of 50 individuals is randomly generated (delimited by the vertical dashed line). After the 50th evaluation, the stage 2 of our algorithm initiates, when genetic operators start to act on the previous population. A significant improvement in the fitness of the population is apparent, the mean fitness value of the population is monotonically decreasing as new individuals are generated, as the black continuous curve. The red line shows the minimum global value until the instant of the last individual is generated in the process. Our GA has not a tendency to converge towards local optima or arbitrary points rather than the global optimum of the problem, this is caused because the operator M^ gives a lower probability of occurrence than crossover operator, this combination gives a good diversity in the generation of new individuals with the better probability to conserve the best properties of its predecessors ensuring the good convergence of the GA. The set of optimized parameters for each λc within the NIR II region and under the mentioned conditions are shown in Table I.

Table I. Optimal parameters, T0, λ0, and P0, obtained using the GA. The soliton order, N, fitness value, ϕ, output central wavelength, λ0, and efficiency of frequency conversion, ƞ, are shown as the obtained results. 

The spectral and temporal evolutions of the best individual found in the optimization process (see Fig. 1) are shown in Figs. 3(a)-(b) respectively. It clearly demonstrates the optimality of the result provided by the GA: the spectrum of the first Raman soliton (the reddest one) is accurately centered in the targeted channel (delimited by the dashed lines).

Figure 3. Spectral (a) and temporal (b) evolutions on distance z corresponding to optimized parameters T ₀ = 50.45 fs, λ0 = 829.05 nm and P ₀ = 14.54 KW for a channel centred in λc = 1225 nm (as in Table ). In the spectral window, dashed vertical lines shows the channel with maximum spectral power obtained by the GA. The black continuous line represent the zero GVD. 

The maximum spectral power found for each channel results very approximate to the Kodama and Hasegawa predictions for the soliton amplitude²².

It is known that SSFS can be made large by propagating shorter pulses with high peak powers inside highly nonlinear fibers and that the fission of higher-order solitons generates frequency-shifted pulses in form of Raman solitons ^1,35.

In order to validate of our method, it is made an exhaustive search of the best fitness value by scanning the entire ranges of parameters. It is shown in Sec. 3.2.

3.2. Optimal solution using exhaustive search

Now, a searching method is implemented to check the reliability of our results. It consists of picking the individual with the best fitness value from a systematic evaluation of all possible combinations of parameters (m = 675) within the defined ranges with defined steps. It is worth mentioning that this process requires larger capabilities of time-machine (∼16.5 h) than the shown using GAs in the previous section. The set of optimized parameters by this method is shown in Table II.

Table II. Best parameters, T0, λ0, and P0, obtained by exhaustive search with m = 675. The soliton order, N, fitness value, ϕ, output central wavelength, λc, and efficiency of frequency conversion, ƞ, are shown. 

The chart of different values for the λc=1225 nm case is shown in the Fig. 4(a). The black polygonal-line shows how the best individual of fixed power jumps randomly in the case when the use of a GA is avoided. Figure 4(b) shows the fitness value evolution in chronological order for the exhaustive search. A sawtooth behavior of the fitness value evolution is because the fitness value becomes better and worse by repeating the parameter values for each individual while they are evaluated in the search.

Figure 4. (a) Fitness value charts for different P ₀ values with m = 675 generated for λc = 1225 nm without the use of the GA (see Table ). The color of the points shows the fitness (ϕ) value of each individual, lighter points have smaller (thus, better fitness values) and the best individual is represented as the red point sited in the coordinates of the three parameters. The black polygonal-line sorts the best individuals for the specific P ₀ . (b) Fitness value evolution of executions in chronological order of generation during the optimization procedure with the minimum global value (red line). The sawtooth behavior of the fitness value evolution is shown in the inset. 

By comparing Tables and , we can see that the efficiency obtained in the optimization using the GA is better than the search made without using the GA, even using more executions in the exhaustive search. We perform a new search level by focusing on a smaller region around the neighbourhood of the previously optimized results in Table . The new conditions are: P∈P0±0.5 kW, T∈T0±5 fs and λ∈λ0±5 nm. The results of the “zoom-in” for each channel are shown in Table . The efficiency of conversion based on the exhaustive search zoom-in is still not better than the efficiency obtained using GA (see Tables and ). Although, it is possible to perform a new search level, the time required will be considerably longer (∼48 h compared with ∼3.75 h using the GA). However, this new search level still does not ensure to find the best solution like GA does.

Table III. Best parameters, T ₀ , λ0, P ₀ , obtained in the zoom-in by exhaustive search with m = 675. The soliton order, N, fitness value, ϕ, output central wavelength, λC, and efficiency of frequency conversion, ƞ, are shown. 

4. Conclusions

We have presented a well defined and efficient optimization procedure of a Ti:Saphire laser pulse parameters to obtain the maximum frequency conversion using a simple device by means of solitonic red-shift in the anomalous region. This optimization is achieved with the use of GAs. Therefore, it has been shown that efficient spectral conversion based on SSFS can be achieved using a simple PCF as a medium of generating spectral broadening pumped by a Ti:Sapphire laser just by properly controlling the input parameters of the input pulses. This scenario typically involves soliton fission and emission of dispersive waves into the normal GVD regime, situation in which precise analytical estimates are not available and therefore the use of numerical simulations in combination with GA is of great usefulness. In summary, this work results in a tool with great potential for optimization of the output of SC spectra for practical OCT applications in the NIR II region.