## 3 The Solution Algorithm
For this purpose, the dependence of the maximum temperature \(T_{max}\), which is reached in the fiber after 1 ms, on the number of mesh elements in the computational domain was plotted. This dependence is shown in Figure 5. The black squares show the temperature values at different numbers of mesh elements, and the red line is the trend line constructed on the basis of the least-squares method. Figure 5 shows that the solution stabilizes around the \(T_{max}\approx\) 177,200 K with a total number of grid elements of more than 15,000. At the same time, the minimum grid cell size was 0.03 \(\upmu\)m, which allows us to perform calculations within the continuum approximation. Figure 4: View of the computational grid: (**a**) general view; (**b**) view near the heat source. Figure 6 shows a block diagram of the numerical algorithm, which shows the steps to solve the problem. Figure 5: Dependence of the maximum temperature in the calculated region in 1 ms on the total number of mesh elements. Black squares indicate temperature values at different number of mesh elements, and red line is trend line. Figure 6: Algorithm scheme representation. Zero step. Before starting the main calculation procedure, a number of physical parameters and constants were set, which were not changed in the course of a particular numerical calculation. A complete list of these parameters and constants is presented in Table 3.  
The first step is the initial conditions. As an initial approximation, a homogeneous temperature field with a value of 293 K was set. Using Formula (3), the initial power distribution of the heat source is calculated. For this purpose, the effective area of the mode spot \(A_{\mathit{off}}\) was calculated using Formulas (8) and (9), and the initial distribution of the absorption coefficient \(\alpha_{e}\) was calculated using Formula (7). The initial time \(t\) condition is set to zero, and the time step \(t\)s is 0.01 ms. The second step was to solve the evolution Equation (1). At each time step, the heat fluxes between the domains of the computational domain were calculated using Formula (2), and the heat source capacity was recalculated using Formulas (7)-(9) and (3). A constant temperature of 293 K was set for the outer boundaries of the computational domain. The numerical solution of Equation (1) was performed with a variable time step, which was selected on the basis of the stability conditions of the computational algorithm. The third step is the output and processing of calculation results. The algorithm for solving the problem was implemented in the application software package Comsol Multiphysics. The approximations of the thermal conductivity equation, interfacial boundary conditions, and boundary conditions were performed using the finite element method. The authors used standard computational modules, which are implemented in Comsol. The computation time was 22 s with the following system parameters: Intel Core i9-9900K Processor, 64 GB RAM, 1 Tb SSD. Calculations involved 8 cores, 1.5 GB of physical memory, and 1.6 GB of virtual memory. PARDISO Solver was used to solve the system of linear equations in multi-threaded mode. We used the constant Newton's method of solving linear equations. For searching and solving nonlinear equations at each time step Jacobian update was used with a maximum number of iterations equal to 40. To accelerate the convergence of nonlinear equations Anderson acceleration, the following parameters were used: dimension of iteration space = 5, mixing parameter = 0.9, iteration delay = 1.