16.3.5  Simulated annealing minimization

The simulated_annealing command attempts to find a point in a search space for which a real-valued cost function is minimized, by using the method of simulated annealing.

simulated_annealing requires that giac is linked to GSL.

• simulated_annealing takes five mandatory arguments and one or two optional arguments:
  • x₀, an initial configuration, which may be any object.
  • costfunc, the cost function which takes a configuration and returns a real number.
  • distfunc, the distance function which takes two configurations and returns the distance between them as a real number.
  • stepfunc, the step function which takes a configuration c and optionally a real number maxstep and returns a random configuration c′ which is a neighbor of c. If stepfunc accepts the maxstep argument, it can be used to limit the distance between c and c′.
  • csparam, a vector [k,T₀,µ,T_min] which defines the cooling schedule parameters, where k is a Boltzmann constant, T₀ is the initial temperature, µ is the damping factor for temperature and T_min is the minimal temperature.
  • Optionally, iterparam, a vector [n,N] of positive integers where n is the number of points tried before stepping and N is the number of iterations for each value of T. By default, n=10 and N=100.
  • Optionally, maxstep, a real number which limits the step size from above. If this parameter is set, then stepfunc must accept two arguments and maxstep will be passed as the second argument whenever calling stepfunc. By default, this parameter is unset.
• simulated_annealing(x₀,costfunc,distfunc,stepfunc,csparam ⟨,iterparam,maxstep ⟩) returns the best configuration found in the context of cost function minimization. The method of simulated annealing gives good results in avoiding local minima, so you usually get the global cost minimum.
• Simulated annealing is essentially a random walk in the search space. The probability of taking a step is determined by the Boltzmann distribution

  p=e

  − Ei+1−Ei k T

  if E_i+1>E_i, and p=1 when E_i+1≤ E_i. Here, E_i and E_i+1 are the configuration energies (costs). If the new energy is higher, a step is always taken, but it may also be taken if the new energy is lower then E_i (the probability of that event drops with T). The temperature T is initally set to a high value and lowered very slightly after the prescribed number of iterations N by using the formula T← T/µ. Here µ is usually slightly larger than 1. The cooling process stops when T reaches its minimum T_min.
• The functions costfunc, distfunc and stepfunc should return undef to raise an error. In that case, simulated_annealing will terminate immediately with an error message.
• Simulated annealing is frequently used to solve combinatorial optimization problems.

Examples

(The following examples are adapted from GSL documentation, see here.) In the first example we find the global minimum of a damped sine wave functionf(x)=e^−(x−1)2sin(8x). This function has many local minima but only a single global minimum in the interval [1.0,1.5]. We start from the point x₀=15.5, which is several local minima away from the solution. The following parameter values are used: k=1, T₀=0.008, µ=1.003, T_min=2× 10⁻⁶, n=20, N=100 and maxstep=1.

In the second example we solve an instance of traveling salesman problem with 12 cities.We enter the cities alongside their latitude and longitude coordinates in an Xcas table cell which we associate with the variable cities (see Figure 16.1).

---

Figure 16.1: Cities used in the TSP example

---

We attempt to find the shortest tour which visits each city exactly once. A tour is represented by a permutation of cities. In total, there are 12!=479001600 possible configurations.

First we need to create a distance matrix which stores the distances between each two cities because each distance may be called for multiple times. We do so by applying the function dist_km which computes the shortest distance (in kilometres) between two cities on a spherical model of Earth of radius R=6375 km.

Now we create the distance matrix by using the following commands.

We generate a neighbor of a given tour by exchanging two random consecutive cities. The distance between two tours (configurations) is simply the Hamming distance. The cost of a tour is equal to its length (in kilometres). The three functions we require are thus entered as follows:

The initial tour is a random permutation of cities.

A good tour can be found by the following command (we use the default values of n and N).

Finally, we visualize the tour. We swap latitudes and longitudes for a more appropriate layout.