{'fraction_non_alphanumeric': 0.03672569773699992, 'fraction_numerical': 0.01378815448883877, 'mean_word_length': 4.275111531702041, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
A thorough study of the performance of simulated annealing with geometric cooling in correlated and long tailed spatial scenarios

Metaheuristics, as the simulated annealing used in the optimization of disordered systems, goes beyond physics, and the traveling salesman is a paradigmatic NP-complete problem that allows inferring important theoretical properties of the algorithm in different random environments. Many versions of the algorithm are explored in the literature, but so far the effects of the statistical distribution of the coordinates of the cities on the performance of the algorithm have been neglected. We propose a simple way to explore this aspect by analyzing the performance of a standard version of the simulated annealing (geometric cooling) in correlated systems with a simple and useful method based on a linear combination of independent random variables. Our results suggest that performance depends on the shape of the statistical distribution of the coordinates but not necessarily on its variance corroborated by the cases of uniform and normal distributions. On the other hand, a study with different power laws (different decay exponents) for the coordinates always produces different performances. We show that the performance of the simulated annealing, even in its best version, is not improved when the distribution of the coordinates does not have the first moment. However, surprisingly, we still observe improvements in situations where the second moment is not defined but not where the first one is not defined. Finite-size scaling fits, and universal laws support all of our results. In addition, our study show when the cost must be scaled.


Introduction
Some magnetic systems, known as spin glasses, present random couplings which lead to high levels of frustration and disordering. The metastability of many structures of such spin glasses leads to experimental physical measurements of the order of days due to slow decay of the magnetization to zero. This peculiarity is also hard to explore in computer simulations. Sherrington and Kirkpatrick [1] proposed an interesting exactly solvable spin-glass model that shows an slow dynamics of magnetization.
Based on these annealing phenomena (heating treatment of materials to enlarge its ductility) in solids, Kirkpatrick, Gelatt, and Vecchi [2] observed a deep connection between the Statistical Mechanics of systems with many degrees of freedom, which can be governed by the Metropolis prescription [3] to thermal equilibrium, and the combinatorial optimization of some problems related to obtain the minimal of functions with many parameters in complex landscapes.
Many combinatorial optimization problems belong to the complexity class known as NP-complete which is a subset of class simply known as NP class. Let us better explain this fundamental point in the complexity of algorithms. A problem belongs to the NP class if it can be decided in polynomial time by a Turing machine (a universal computer model, for more details see [4]). Here it is important to say that not necessarily it possesses an algorithm that can solve it in polynomial time. Sure, we know that there are problems which can be solved by an algorithm in polynomial time (fast solution) and not only decided! A simple example is sorting a list.
One algorithm in these conditions belongs to a complexity class known as P, and obviously P ⊂ NP. However there are some "rebel" problems in NP, the so-called NP-complete problems, that can be understood, in a non-rigorous way, as the harder problems in NP, since there is not a single known algorithm to date, that can solve one of these problems in polynomial time differently than occurs with problems in P.
Should we then believe in conjecture P = NP? Yes, and the reason is simple! If we know one NP-complete problem, any algorithm that can solve it, can be used to solve any other NP-complete problem (this is technically called polynomial reduction, and for a good reference in complexity theory, see again [4]). Therefore, if somebody discovers an algorithm in polynomial time that solves one of these NP-complete problems, by consequence, all of them will also have solution in polynomial time. However, such algorithms were never found, which makes stronger the conjecture.
But what do you do with these problems? An alternative for them is to use heuristics, and those which are physically motivated can be a good option. Thus, based on Ref. [2], which proposed a metaheuristics known as simulated annealing (SA), an algorithm which performs a random walk in the configuration space until it reaches the global equilibrium or, at worst, an interesting local minimum, the physicists, and computer scientists have explored several different combinatorial problems, including for example the study of thermodynamics of the protein folding [5], the evaluation of periodic orbits of multi-electron atomic systems [6], and many others.
Basically, the SA algorithm works in the following way: we start with an initial solution/configuration σ old , and calculate the energy (cost in a more general context) of this initial solution E(σ old ). After that, another state is randomly chosen, which is denoted by σ new , and this new solution is accepted with a probability that follows the Metropolis prescription [3]: otherwise the system remains in the same σ old state. This process is repeated until the ensemble sampled of the system reaches an equilibrium at the given temperature. Finally, the temperature T is decreased by a cooling schedule. The temperature will be decreased until a state with low enough energy is found.
Here, it is important to mention that many cooling schedules can be applied. There are cooling schedules that asymptotically converge towards the global minimum as the one that cools the system with a logarithm rule [7]. However such schedule converges very slowly and requires a long computation time.
There are good alternative schedules [8], that although without a rigorous guarantee of the convergence towards the global optimum, are computationally faster. One of them is the geometric cooling schedule [9] which considers that at time t, the temperature is given by T = α t T 0 , with 0 < α < 1. An interesting version of this heuristic works with two loops: an internal and another external which can be resumed by the Heuristic SAGCS (Simulated annealing 2 with geometric cooling schedule) resumed in the algorithm 1.

Algorithm 1: SAGCS (Standard Simulated Annealing with Geometric Cooling Schedule)
parameters: σ 0 : initial configuration T 0 : initial temperature T min : final temperature N iter : number of iterations of the internal loop α : coefficient of the geometric cooling Result: low energy state σ It is important to observe that such Heuristic works with two parameters: N iter that controls the number of internal iterations, while the external loop (iterations over different temperatures) is governed by: and sure, for appropriate N steps and N iter , the system must converge for a local minimum that is very close to the global minimum. Amongst the many interesting NP-complete problems, the traveling salesman problem (TSP) deserves much attention and it became a natural way to test the SA heuristic using different cooling schedules, variations of the method, and parallelization techniques in many different contexts (see for example [7,8,9,12,13,14,15]). The problem (TSP) is to find the optimal Hamiltonian cycle on a given graph (vertices and valued edges), i.e., the minimum cost closed path including all vertices passing through each vertex only once. Let us consider graphs with N vertices corresponding to points (x i , y i ), i = 1, .., N, distributed in a two dimensional scenario. In this case, the graph naturally has valued edges defined by the Euclidean distances: with i, j = 1, .., N, and therefore, the simplest situation is to consider the graph having all edges, since all costs d(i, j) can be defined for all pair of the points, which is denoted by a complete graph.
In this paper, we are not interested in testing cooling schedules, or even other SA heuristics, which are very well explored in the literature, but in building computer experiments to test the SAGCS in the TSP in order to understand its efficiency considering the effects of correlation and variance on the random coordinates x i and y i . It is important to mention that our method to include correlation is very simple and allows to change the correlation but keeping the same variance which is very important for making a fair comparison, since we simultaneously variate the two parameters, we do not know if the effects are only caused by the correlation.
Our results show how the performance of the algorithm transits from two-dimensional scenario (ρ = 0) to the one-dimensional situation (ρ = 1) by deforming the region where the points (cities) are distributed. Here ρ is the correlation coefficient between the coordinates of the points. We also investigate a parameter that controls the moments of the random variables x i and y i and we quantitatively show how the increase of the variance affects the performance. We particularly show that the shape of the distribution is more important than its variance in certain situations. On the other hand, when the points have coordinates that are power-law distributed, in the region where the average cannot be defined, the efficiency of the algorithm is negligible independently on scaling.
In the next section we will show some fundamental and pedagogical aspects of the SA for the TSP. In the following, in section 3 we define the different environments where we will apply the SAGCS. We will show in detail how to generate points with correlated coordinates, and how to generate points with power-law distributed coordinates. In section 4 we show our results and finally in section 5 we present some conclusions.

Pedagogical aspects of the Simulated Annealing in the context of the Travelling salesman problem
The simulated annealing with geometric schedule (SAGCS) is a very simple heuristic used to optimize combinatorial problems as the Traveling Salesman Problem. The points (x i , y i ), with i = 1, ..., N are scattered over a two dimensional region, where the coordinates follow some joint probability distribution p(x i , y i ). Denoting , whose cost is given by where v(N + 1) = v(1) and d(i, j) is given by Eq. 3. We aim to determine the best cycle σ , i.e., the cycle of minimal cost.
The new configuration can be obtained performing different mechanisms in this work. The first one is the "simple swap"(SS) , i.e., a point is randomly chosen, for example v( j), and the candidate to be the new configuration σ corresponds to σ by simply changing the positions of v( j) ≡ (x j , y j ) and v( j + 1) ≡ (x j+1 , y j+1 ), then we can perform the Metropolis prescription to optimize the problem in the context of the simulated annealing.
On the other hand, there is a more interesting way to obtain a new configuration according to [10,11] the 2-opt move. It is a popular procedure used to improve algorithms to approximate the solution of the Travelling Salesman Problem. Starting from a given cycle, it consists in exchanging two links of the cycle to construct a new one. This is performed by reversing the sequence of nodes between the selected links and then reconnecting the cycle back together. For example, if one has a route/cycle (1). It allows us to explore the state space of Hamiltonian cycles efficiently as well as it provides a simple way to calculate the energy difference between the states before and after the move, which is a very important feature for working in conjunction with the Metropolis algorithm. Now, let us explore the complexity of the problem. There are N! different cycles. Testing all of them is impossible for large N! For example for N = 2048 cities (this will be a number of cities of the most of our examples), by using the Stirling's formula N! ≈ √ 2πN N+1 e −N , one estimates This a staggering number when compared, for example, with the number of atoms of the universe O(10 80 ). Actually, among this "astronomical" number of different cycles (much more than astronomical), certainly there are a lot of them corresponding to a good local minima. It is exactly here that SA shows its usefulness. Let us understand how it works by showing the power of this algorithm. For example, let us examine the Fig. 1.
This figure shows three histograms: the red one, represented in Fig. 1 (a) shows a sample of the costs of one million of random drawn cycles of the same configuration with N = 2048 cities, denoting points whose coordinates are randomly distributed according to a uniform distribution, with  One can observe two very distinct Gaussian distributions whose means differ across 1000 units of the cost. It is important to notice that not even one among the one million costs of randomly drawn Hamiltonian cycles reached the cost of SA. Actually, it is worse than this, in the first case we have C rand = 2143.74 ± 0.02, while in the second case one has C SA = 1021.67 ± 0.74 and the overlapping probability between these two distributions is zero in the precision of the any known machine! It is important to mention that we provocatively used only N run = 400 runs for the SA against the brute force (random drawn cycles).
The final blow comes when we use the 2-opt choice to sample new configurations in the SAGCS (2opt-SAGCS), which is shown Fig. 1 (c), where we exactly used the same parameters of the Fig. 1 (b). In this case one has C SA = 77.959 ± 0.031 which is an extremely low cost when compared with brute force method.
Thus, this pedagogical explanation is only to show that exhaustive sampling is not a feasible solution for combinatorial problems. Heuristics as the SA are an important alternative to obtain good (not always the optimal, but in very disordered systems this is not be a meaningful difference) solutions to the TSP. Now, after this preparatory study, we present the details about the scenario for which we intend to analyze the performance of the SAGCS. In the next section we will present the method to generate points with correlated coordinates, and how to generate points with coordinates long-tailed distributed. In these environments, we intend to explore some effects on the performance of the SA which will be performed in section 4.

Correlated and long-tailed environments
One of the important questions addressed in this paper is how the SAGCS works considering that x i and y i , the coordinates of the cities, are random variables with a given correlation ρ. Additionally, for the particular case looking at the effect of γ on performance of the algorithm. Now, we discuss how to generate scenarios by following these prescriptions.

Generating correlated random coordinates from non-correlated random variables
In this section, we will show that we can generate correlated random variables from non-correlated random variables considering that both (correlated and non-correlated) have the same variance and average by imposing an additional constraint -considering the average equal to 0 for the uncorrelated random variables. This is exactly the same procedure used in [16] in the context of emerging of rogue waves in the superposition of electrical waves with correlated phases.
Let us consider spatial coordinates x e y of points (our cities) in a two-dimensional environment that are random variables: where z 1 and z 2 are independent and identically distributed random variables, such that: The variance of the variable x is given by: which implies that α 2 1 + α 2 2 = β 2 1 + β 2 2 = 1. It is worth noting that although z 1 and z 2 are non-correlated random variables, x and y are, and the correlation between these random variables can be calculated: Thus, again after some cancellations and combinations: ρ = (α 1 β 1 + α 2 β 2 ). And after a little algebra, the random variables can be now written as: and have the same average that are given by: Thus, if one considers x and y as ρ-correlated random variables generated from two independent random variables z 1 and z 2 , with average zero and variance σ 2 = (∆z) 2 , x and y also have average zero and the same variance σ 2 . This is an important point because if one works with different averages and dispersions between the cases ρ = 0 and ρ = 0 the results can be misleading since there will be no parameter for a fair comparison, i.e., in our problem no bias can be the source of the possible studied effects, since independently on correlation, the random variables are sampled with the same average and variance.
For example, using two uniform and identically distributed random variables z 1 = 2(ξ 1 − 1 2 ) and z 2 = 2(ξ 2 − 1 2 ) that assume values in [−1, 1] generated from two ξ 1 and ξ 2 uniform random variables assuming values in [0, 1], the Fig. 2 shows a plot of y versus x obtained from equations 8 and 9 for different values of ρ. In this paper, the idea is to study the SAGCS on these correlated scenarios by analyzing its performance.

Power-law distributions for the coordinates
Another important point of our study is to look at the effects of long tailed distributions for the coordinates of the points on the SA performance. Thus, we use a power-law probability density function to generate the coordinates of the two-dimensional points. For that, we initially propose the following distribution for the coordinates: It is important to observe that the power law distribution given by Eq. 10 has a necessary gap ∆ = 2x 0 with center in the origin, since the two branches can not touch the origin by normalization, however x 0 can be made arbitrarily small. Just for illustration, we show plots of p(x; x 0 , γ) in Fig. 3.
In order to draw a variable that follow the distribution of Eq. 10 one simply uses two uniform random (or more precisely pseudo-random variables) variables ξ 1 and ξ 2 assuming values in [0, 1]. First, one chooses if the variable is on the positive or negative branching. This is performed checking the sign of (2ξ 1 − 1). After that, one makes ξ 2 equal to area (γ−1) Therefore, in the general case, the random variable built from ξ 1 and ξ 2 is: (11) and naturally with other two random variables ξ 3 and ξ 4 , which are chosen from the uniform distribution in [0, 1], one has In Fig. 3 (b), a histogram for the points (x, y) sampled in this way is shown for γ = 5.5, only to check the agreement with the exact result. The idea is to observe the performance of the SAGCS in this scenario answering how important is γ and its effects on the final cost obtained with SA algorithm. Notice that points generated by Eqs. 11 and 12 have border effects. For instance, in Fig. 4, we present some examples of the points scattered according to different distributions by focusing the power law ones.
Using the Gaussian distribution (with average zero and variance 1), the pattern of points shown in Fig. 4 (a) is very different than the case which uses the points generated by Eqs. 11 and 12 for γ = 2.9 and x 0 = 10 −3 (Fig. 4  (b) ). Here, we will denote the points using the Eqs. 11 and 12 by power law I. We can observe some outlier points since, for this case, the distribution has no second moment defined. However, radial symmetry is absent in the power law I distribution and the border effects previously mentioned are notorious. The question is whether they are really important. Alternatively, in this paper, we also proposed a second power law type (power law II) to study the SA, which respects the radial symmetry: and in this case for γ = 2.9 and r 0 = 10 −3 one has the distribution of points represented in Fig. 4 (c). Similarly ξ 1 and ξ 2 are uniform random variables assuming values in [0, 1]. Thus relevant questions should be related to the performance of the SA in both situations: radially symmetric and asymmetric ones. By completing in Fig 4 (d), and (e) we show the distribution of points for a case where we have defined variance (γ = 5.5) and therefore few outliers are observed. The case γ = 3.4 complete our figure (Fig 4 (f) ) since this case will be particularly important for our results. It is important to mention that except by Fig. 4 (a)

Results
We performed computer experiments to analyze two effects on the optimization by the SA for the TSP: a) the correlation between the spatial coordinates, and b) the variance/distribution shape for the spatial random coordinates. We used the SAGCS which is a fast and standard way to perform optimization. It is natural to expect that such effects must be proportionally important in other variations of the SA employing other slower cooling schedules independently if the final cost obtained is better. The goal of this paper is not to compare different SA algorithms but performing a quantitative study of the SAGCS considering different spatial distributions for the coordinates of the points in the TSP.
However, before starting the main core of our results, it is important to understand some preliminary aspects: the effects of the number of external loop iterations N steps = ln(T min /T 0 ) ln α , and N iter which denotes the number of internal loop iterations on the final cost. We consider points (x i , y i ), = 1, ..., N uniformly distributed non-correlated random variables (case ρ = 0 in Fig. 2). In all of our simulations, we averaged our final cost over N run = 60 different runs of the SA.
So we performed simulations considering T 0 = 10 and T min = 10 −5 . We initially consider the simple swap scheme to generate new candidate configurations. Under these conditions, one obtains a plot of the cost in the steady state (C opt ) which can be observed in Fig. 5 (a) as function of total number of iterations N total = N steps N iter . Here we intent to show that such parameter is more important than the parameters N steps (α), and N iter separately. Our range for N total was from 5.10 5 iterations up to 10 9 iterations. The procedure is to change N iter from 2 10 up to 2 18 . Now, given the N total and the N iter we obtain the factor of the exponent α of the cooling: For example, when N total = 5.10 5 , α changes from 0.001 up to 0.9734 and when N total = 2.10 9 one obtains that α changes from 0.99831 up to 0.99999. There is actually a lot of information aggregated in this figure to grasp. First the size of the points is proportional to N iter and α → 0 corresponds to red, while α → 1 corresponds to blue. Fig. 5 shows that C opt depends on N total and not on N steps and N iter separately. For this experiment it was used N = 2048 cities. The 5 (b) shows that for T 0 = 100 the behavior is similar. The result suggests a reasonable power-law universal behavior C opt ∼ N −δ total for sufficiently large total number of points in the sample independently on the internal or external loops. For clarity, in the Fig. 6, we show a plot of C opt as function of N total for both T 0 = 10 and T 0 = 100 in the same plot, forgetting at this moment the effects of α and N iter separately. In log-log scale we obtain respectively δ = 0.079 ± 0.001 and δ = 0.083 ± 0.001, which shows a similar power law behavior.
We perform a similar analysis of Fig. 6 for the 2-opt prescription, which is shown in Fig. 7. In this case, somewhat differently from the simple swap, one observes a transition between two power laws with very different exponents: from δ ≈ 0.43 transiting to δ ≈ 0.02, which shows that for N total > O(e 16 ) ∼ O(10 6 ) the optimal cost found does not show a meaningful improvement, by changing from C opt ≈ 90 up to C opt ≈ 70. After this preparatory study, we can proceed to the main study of this work, to investigate the correlation and long range effects on the coordinates of the points.

Correlation effects
From now, all of our results were obtained by using a fixed set of parameters N iter = 8000, α = 0.99, N run = 60, T min = 10 −5 , and T 0 = 10, which results in N total = O(10 7 ). Except by the finite size scaling analysis, all of our experiments with SAGCS were applied in scenarios with N = 2048 cities.
Let us start our study of correlation effects, by observing the behavior of C as function of the iteration in the external loop (iteration of cooling schedule) n, for different values of ρ in scenarios with cities/points (x,y) built from identical uniform random variables z 1 and z 2 defined on [-1,1] as the Fig. 2. We will first study the SS-SAGCS performance.
In Fig. 8(a) we can observe a decreasing of the final cost C opt until it reaches a steady state after approximately 500 iterations. However, an interesting analysis is to consider the ratio C opt / C , given by where the average cost is calculated considering that cities have an average distance multiplied by the size of cycle, and here, O = (1/N run ) ∑ N run i=1 O i denotes an average of amount O over the different runs (different steady costs obtained by the SA). It is important to say that C calculated by Eq. 14 has no difference of the cost of any cycle which was randomly chosen. This is a typical characteristic of the TSP, since Hamiltonian cycles randomly generated are so far from the local minimum found by heuristics as the SA, that even generating billions of such cycles, the probability of reaching a good local minimum as the SA can be disregarded as we previously discussed in Sec. 2.
The lower C opt / C , the better the performance. Fig. 8 (b-I) shows the behavior of this ratio as function of ρ for several different values of N. We have an improvement of the performance as ρ enlarges. A reasonable (not so good) empirical scaling can be observed in Fig. 8(b-II). If one uses b = L max /L and adjust z ≈ 0.075, the different curves have a good agreement in the collapse. Here L max is the largest lattice size used in such study, which was L max = 2 11 .
We can observe that the performance for SS-SAGCS is even better for small number of cities for the simple swap prescription. But, an important point is to better investigate the performance of the SAGCS as function of ρ. Let us focus the extreme cases: ρ = 0 and ρ = 1. For ρ = 1, one has exactly the points scattered in a straight line (case ρ = 1, Fig. 2). In this particular case, we know in advance, the global optimal Hamiltonian cycle, which is exactly given by C (theor) glob = 4 √ 2 when N → ∞, since we have points increasingly closer to (−1, −1) and (1, 1) the higher the value of N, and therefore for finite N is only an approximation. On the other hand, we can analytically estimate the averaged cycle C N is the average displacement between two adjacent pair of points which is randomly distributed over the line. Thus, C (theor) ρ=1  Table 1, which shows a comparison of C opt for different values of ρ. One observes that C Actually, finding the optimal Hamiltonian cycle with points scattered in a straight line is equivalent to order a list and one has good algorithms in polynomial times (heapsort, quicksort...) that efficiently performs such classification.
Nevertheless, it is important to mention that we have no previous information about the topology of the points and the simulated annealing with 2-opt prescription simply works to find a similar result to the optimal one for ρ = 1 as well as in any other situation, which is very interesting indeed.
Let us analyze the ratio C glob / C ρ=0 to check if the curve C opt / C ×ρ obtained with SA is compatible with these extremes. Once C glob for ρ = 0 is not known, how do we estimate C glob / C ρ=0 ? First, let us start with a calculation that can be performed exactly (see Appendix): which amounts to 2136, approximately, for N = 2048. Just to check, for N run = 60, we numerically obtain C ρ=0 ≈ 2131. Although we do not exactly know C glob , constructive heuristics (out of scope of the Boltzmann machines) can supply an estimate which should work as a benchmark. A good suggestion is the nearest neighbor (NN) algorithm [17], a greedy algorithm that corresponds to a particular case of the tourist random walk [18], [19] with the maximal memory: N − 1 can be an candidate.
In this algorithm the walk start from an initial node and jumps to the nearest neighbor that has not yet been visited in the walk. The procedure is repeated until returning to the initial node by closing the cycle (Hamiltonian cycle). In this case, we applied this algorithm for ρ = 0 and one finds a better approximation C There exist many specially arranged city distributions which make the NN algorithm gives the worst route. Here we are observing, it is an excellent heuristic, however, the 2-opt-SAGCS is still a better benchmark (almost a technical draw between the two heuristics) and the SAGCS has a good complexity when compared with NN algorithm. The SA has a complexity that can be writte in a general form as O(c(N) · N total ) where c(N) = O(1) for SS-SAGCS and c(N) = O(N/2) at the worst case for the 2-opt-SAGCS. On the other hand the NN algorithm has a complexity O(N 2 ). Now it is interesting to consider some points: N total is not exactly related to N, since SA is a heuristic. Sure, N steps = ln(T min /T 0 ) ln α , and therefore it is a constant and does not depend on N. However, for larger values of N, it is needed to make considerations about N iter . It is natural to expect that the larger N, larger the number of iterations of the internal loop and we should consider N iter = O(N) and the complexity of SA would be, at the worst case, a complexity O(N 2 ) which is exactly the complexity of NN-algorithm considering the analysis of a single starting point. It is important to mention that if we include the dependence on the initial point, since different initial points can lead to different optimal costs, the complexity of the NN-algorithm would then be O(N 3 ).
SA has modest optimal costs only when one uses a simple swap scheme, which is a naive technique when compared with 2-opt that can applied in general scenarios. In realistic situations, a generalization of the TSP, the VRP (vehicle routing problem) is the best alternative compared to the NN algorithm, and Tabu search algorithm as suggested by the authors in [20].
But again, our proposal in this paper is to analyze possible effects of the environment on the SA and not a detailed comparison among the methods, yet we could not miss to show the 2-opt-SAGCS in comparison to SS-SAGCS and the NN algorithm.
Thus, it is interesting to similarly analyze the size effects on the 2-opt-SAGCS exactly as we performed in plot in Fig. 8 (b). Fig. 9 (a) I shows a similar decay of C opt / C as function of ρ but with values extremely lower than those of SS-SAGCS. Following what we applied to the SS-SAGCS, we performed a similar scaling for the 2-opt-SAGCS ( Fig.  9 (a) II) but not so with the same quality. In this case the value is z ≈ −0.480 which is very different than z ≈ 0.075 found to SS-SAGCS.
It is also interesting to analyze the effects on the points considering different statistical distributions for the coordinates. Thus we prepared some experiments to capture the effects on the ratio C opt / C as function of ρ considering different distributions for the coordinates and the result is also surprising as Fig. 9 (b) shows. If before, we were considering z 1 and z 2 identically distributed uniform random variables, now we also consider them as Gaussian random variables for a comparison. It is important to mention that the x and y coordinates, have, the same average (zero) and variance of the variables z 1 and z 2 independently on ρ. As we previously reported, the effects to be observed depends strictly on ρ since the average and variance do not change.
We can observe that in both cases, SS-SAGCS and 2-opt-SAGCS, respectively described by Fig. 9 (b-I) and Fig.  9 (b-II), the shape of the distribution seems to be more important than the variance once the Gaussian and the uniform random variables, with the same or different variances, lead to different curves, but gaussians with different variances produce practically the same behavior.
Finally, we also use an scaled cost: where x(y) max(min) = max(min){x(y) 1 , x(y) 2 , ..., x(y) N }, once, differently from uniform distribution, the Gaussian distribution is not supported on a finite interval. However, as the same Figs. 9 (b-I) and (b-II) show, no changes were observed in the two versions of the simulated annealing (the blue triangles correspond to Gaussian distribution with scaled cost). For SS-SAGCS the lower the N, the better the performance, on the other hand, for the 2-opt-SAGCS, the higher the N, the better the performance as also can be observed in the Figs. 8 (b-I) and 9 (b-I).

Finding good fits for C opt / C as a function of ρ
We also focus our results in finding good fits for C opt / C as function of ρ. So follow a (heuristic and qualitative) method to find suitable fits for such behavior.
First, we try a polynomial fit p(ρ) = n ∑ i=0 a i ρ i by testing n = 2, 3 and 4 (with 3, 4, and 5 parameters respectively). The idea here is only to make a preliminary investigation. After that, we test exponential decay fits with only two parameters:  Table 2: Values of coefficients found after nonlinear fitting by using Levenberg-Marquadt method. The parameter α corresponds to coefficient of determination of each fit once that e I (0) is fixed and taken as C opt / C (0) (it was visually investigated). Following, one considers a simple linear combination of exponential decays (with four parameters): e II (ρ) = e II (0) + a 1 e −ρ/ρ 1 + a 2 e −ρ/ρ 2 , also assuming empirically that e II (0) = C opt / C (0). Alternatively, we also experimented other functions with four parameters, and the one that presented a good result was the rational function: This is shown in Fig. 10, which shows the different fits. The fits obtained in the different cases are summarized in table 2.
We can observe that the larger the degree of the polynomial fitted, the better the coefficient of determination α (the closer to one, the better it is). Our analysis also shows that in both cases, a single exponential (exponential I) is not enough to nicely fit the curve but with two exponential (four parameters), an excellent fit with α ≈ 0.999 can be obtained. With the same four parameters, the rational function gives a good result despite not so good as the linear combination of exponentials. It is important to mention that the coefficient a found is exactly C opt / C (0) in both cases (SS and 2opt-SAGCS) and the fit should be performed with three parameters. After all, one can observe that a better match is obtained by considering of a linear combination of the exponentials described by the Eq. 16. It seems to be a universal fit for C opt / C ×ρ for different values of ρ for both SS-SAGCS and 2-opt-SAGCS.

Long tail effects
Finally, it is important to analyze the effects of long tailed distributions for the coordinates of the points distributed on the environment. In this case, we concentrate our analysis on ρ = 0. Firstly, one chooses the variables x and y according to equations 11 and 12 for the SS-SAGCS. Plotting C opt / C as function of γ according to Fig. 11 (a) and  (b).
As expected C opt / C decreases as function of exponent γ. Fig. 11 (a) shows this behavior for different values of x 0 considering the scaling in Eq. 15. We can observe a collapse of all curves independently on x 0 . The inset plot shows the same plot without performing the scaling. It is interesting to notice the instability in the regions without defined variance 2 ≤ γ ≤ 3 but yet with some efficiency. For γ < 2, which corresponds the region where the average cannot be defined, the SA has no efficiency which suggests that the algorithm should not be applied in these situations since C opt / C ≈ 1. Fig. 11 (b) is only a refinement of one case in Fig. 11 (a) with uncertainty bars and also showing the values obtained with Gaussian and uniform distributions for a comparison (dashed blue lines).
Performing similar simulations to the plot of the Fig. 11 (a) but for the case 2opt-SAGCS, which can be observed in Fig. 12 (a) and the same conclusions can be drawn for γ < 2 -the SA has no efficiency. In this case, SS-SAGCS and 2-opt-SAGCS are obviously equivalent since the scenario is really complicate. In Fig. 12 (b) we show the case where we used power-law coordinates with radial symmetry (Eq. 13). Surprisingly, the results do not change, showing that, independently from power law coordinates, the phenomena of the performance of the SA is universal as function of γ.

Summary and Conclusions
In this paper, we study the effects of the statistics on the coordinates of the points when we apply an standard simulated annealing algorithm to the travelling salesman problem. Our results are concerned with the long tail effects on the coordinates but also the correlation effects between these coordinates.
Our study shows that the performance of the simulated annealing increases as the correlation increases in both versions of the SA. The main reason is the dimensionality reduction, which transforms the simulated annealing at limit in an approximated sorting algorithm. The correlation effects show that the shape of distribution attributed to coordinates is more important than the variance when we compare Gaussian and uniform distributions.
Our results also suggest that the higher the exponent of the power law, the lower the simulated annealing performance, and for γ < 2, i.e., distributions without the first moment, we have no performance of the SA, i.e., C opt / C ≈ 1 even when one considers the benchmark for the SA, i.e., when the new configurations are obtained with 2-opt.
We obtained an universal behavior of C opt / C ×ρ and C opt / C ×γ never explored in the optimization problems with SA. In addition, we also show that 2-opt-SAGCS presents a better performance as larger is the system, differently from its inefficient version, the SS-SAGCS, where the reverse situation occurs.
We believe that both (correlation and long range effects) studies can bring an interesting knowledge for more technical applications in artificial intelligence, machine learning, and other areas. In special, in the search for global minimum in neural networks algorithms, maybe reviving the interests in simulating annealing as a viable alternative to stochastic gradient descent at the optimization step. If z 1 and z 2 are uniformly distributed, so the random variable u = |z 1 − z 2 | has a triangular distribution and the probability density function 2(1 − u), and thus making the change of variables u = |z 1 − z 2 |, v = |z 3 − z 4 |, w = z 2 , and t = z 3 , one obtains: Now is almost done! Making u = r cos θ and v = r sin θ . However a trick is to perform the integration for the lower triangle where 0 ≤ θ ≤ π/4 d = 16  The inset plot shows the same curves without scaling. (b) The same plot by using power law coordinates with radial symmetry (Eq. 13).