### 1. Simulated annealing algorithm

Simulated annealing algorithm draws lessons from the idea of statistical physics. It is a simple and general heuristic optimization algorithm, and has probabilistic global optimization performance in theory. Therefore, it has been widely used in scientific research and engineering.

Annealing is a process in which the metal cools slowly from the molten state and finally reaches the equilibrium state with the lowest energy. Based on the similarity between the solution process of the optimization problem and the metal annealing process, the simulated annealing algorithm takes the optimization objective as the energy function, the solution space as the state space, and the random disturbance simulates the thermal motion of particles to solve the optimization problem ([1] KIRKPATRICK,1988).

The structure of simulated annealing algorithm is simple, which is composed of temperature update function, state generation function, state acceptance function and inner loop and outer loop termination criteria.

Temperature update function refers to the realization scheme of slow reduction of annealing temperature, also known as cooling schedule;

State generating function refers to the method of randomly generating new candidate solutions from the current solution;

State acceptance function refers to the mechanism of accepting candidate solutions, which usually adopts Metropolis criterion;

External circulation is a temperature cycle controlled by the cooling schedule;

The inner loop is the number of times that the loop iteration produces a new solution at each temperature, also known as the length of Markov chain.

The basic flow of simulated annealing algorithm is as follows:

(1) Initialization: initial temperature T, initial solution state s, number of iterations L;

(2) For each temperature state, repeat L cycles to generate and probabilistic accept a new solution:

(3) The new solution s' is generated from the current solution s by the transformation operation;

(4) Calculate the energy difference ∆ E, that is, the difference between the objective function of the new solution and the objective function of the original solution;

(5) If ∆ e < 0, accept s' as the new current solution, otherwise accept s' as the new current solution with probability exp(- ∆ E/T);

(6) After completing L internal cycles in each temperature state, reduce the temperature T until the termination temperature is reached.

### 2. Multivariable function optimization problem

The classical function optimization problem and combinatorial optimization problem are selected as test cases.

Question 1: Schwefel test function is a complex multimodal function with a large number of local extremum regions.

F(X)=418.9829×n-∑(i=1,n)〖xi* sin(√(|xi|)) 〗

In this paper, we take d = 10, x = [- 500500], and the function is the global minimum f(X)=0.0 at X=(420.9687,... 420.9687).

The basic scheme of using simulated annealing algorithm: control the temperature to decay exponentially according to T(k) = a * T(k-1), and the attenuation coefficient is a; As shown in equation (1), accept the new solution according to Metropolis criteria. For problem 1 (Schwefel function), a new solution is generated by applying a random disturbance of normal distribution to an independent variable of the current solution.

### 3. Simulated annealing algorithm Python program

# Simulated annealing algorithm program: Multivariable continuous function optimization # Program: SimulatedAnnealing_v1.py # Purpose: Simulated annealing algorithm for function optimization # Copyright 2021 YouCans, XUPT # Crated: 2021-04-30 # -*- coding: utf-8 -*- import math # Import module import random # Import module import pandas as pd # Import module import numpy as np # Import the module numpy, which is abbreviated as np import matplotlib.pyplot as plt # Import the module Matplotlib Pyplot and abbreviated as plot from datetime import datetime # Subroutine: defines the objective function of the optimization problem def cal_Energy(X, nVar): # Test function 1: Schwefel test function # -500 <= Xi <= 500 # Global extremum: (420.9687420.9687,...), f(x)=0.0 sum = 0.0 for i in range(nVar): sum += X[i] * np.sin(np.sqrt(abs(X[i]))) fx = 418.9829 * nVar - sum return fx # Subroutine: parameter setting of simulated annealing algorithm def ParameterSetting(): cName = "funcOpt" # Define question name nVar = 2 # Given the number of arguments, y = f (x1,... Xn) xMin = [-500, -500] # Lower limit of given search space, x1_min,..xn_min xMax = [500, 500] # Upper limit of given search space, x1_max,..xn_max tInitial = 100.0 # Set initial annealing temperature tFinal = 1 # Set stop temperature alfa = 0.98 # Set cooling parameters, T(k)=alfa*T(k-1) meanMarkov = 100 # Markov chain length, that is, the number of inner loop operations scale = 0.5 # Define the search step size, which can be set to a fixed value or gradually reduced return cName, nVar, xMin, xMax, tInitial, tFinal, alfa, meanMarkov, scale # Simulated annealing algorithm def OptimizationSSA(nVar,xMin,xMax,tInitial,tFinal,alfa,meanMarkov,scale): # ======Initialize random number generator====== randseed = random.randint(1, 100) random.seed(randseed) # The random number generator can set the seed or set it to a specified integer # ======Initial solution of stochastic optimization problem====== xInitial = np.zeros((nVar)) # Initialize, create array for v in range(nVar): # random.uniform(min,max) randomly generates a real number in the range of [min,max] xInitial[v] = random.uniform(xMin[v], xMax[v]) # Call sub function cal_Energy calculates the objective function value of the current solution fxInitial = cal_Energy(xInitial, nVar) # ======Simulated annealing algorithm initialization====== xNew = np.zeros((nVar)) # Initialize, create array xNow = np.zeros((nVar)) # Initialize, create array xBest = np.zeros((nVar)) # Initialize, create array xNow[:] = xInitial[:] # Initialize the current solution and set the initial solution as the current solution xBest[:] = xInitial[:] # Initialize the optimal solution and set the current solution as the optimal solution fxNow = fxInitial # Set the objective function of the initial solution to the current value fxBest = fxInitial # Set the objective function of the current solution as the optimal value print('x_Initial:{:.6f},{:.6f},\tf(x_Initial):{:.6f}'.format(xInitial[0], xInitial[1], fxInitial)) recordIter = [] # Initialization, number of external cycles recordFxNow = [] # Initialization, the objective function value of the current solution recordFxBest = [] # Initialization, the objective function value of the optimal solution recordPBad = [] # Initialization, acceptance probability of inferior solution kIter = 0 # Number of external loop iterations, number of temperature states totalMar = 0 # Total Markov chain length totalImprove = 0 # fxBest improvement times nMarkov = meanMarkov # Fixed length Markov chain # ======Start simulated annealing optimization====== # The external circulation ends when the current temperature reaches the termination temperature tNow = tInitial # Initialize current temperature while tNow >= tFinal: # The external circulation ends when the current temperature reaches the termination temperature # At the current temperature, a sufficient number of state transitions (nMarkov) are performed to achieve thermal equilibrium kBetter = 0 # Number of times to obtain high-quality solutions kBadAccept = 0 # Number of times to accept inferior solutions kBadRefuse = 0 # Number of rejections of inferior solutions # ---Inner loop, the number of cycles is the length of Markov chain for k in range(nMarkov): # Inner loop, the number of cycles is the length of Markov chain totalMar += 1 # Total Markov chain length counter # ---Generate new solutions # Generate a new solution: generate a new solution by randomly disturbing near the current solution. The new solution must be within the range of [min,max] # Scheme 1: only one of the n variables is disturbed, and the other n-1 variables remain unchanged xNew[:] = xNow[:] v = random.randint(0, nVar-1) # Generate random numbers between [0,nVar-1] xNew[v] = xNow[v] + scale * (xMax[v]-xMin[v]) * random.normalvariate(0, 1) # random. Normal variable (0, 1): generate a random real number that obeys a normal distribution with a mean of 0 and a standard deviation of 1 xNew[v] = max(min(xNew[v], xMax[v]), xMin[v]) # Ensure that the new solution is within the range of [min,max] # ---Calculate the objective function and energy difference # Call sub function cal_Energy calculates the objective function value of the new solution fxNew = cal_Energy(xNew, nVar) deltaE = fxNew - fxNow # ---Accept the new solution according to Metropolis guidelines # Acceptance judgment: decide whether to accept the new solution according to Metropolis criteria if fxNew < fxNow: # Better solution: if the objective function of the new solution is better than the current solution, the new solution is accepted accept = True kBetter += 1 else: # Tolerant solution: if the objective function of the new solution is worse than the current solution, the new solution will be accepted with a certain probability pAccept = math.exp(-deltaE / tNow) # Calculate the state transition probability of the tolerant solution if pAccept > random.random(): accept = True # Accept inferior solutions kBadAccept += 1 else: accept = False # Reject inferior solutions kBadRefuse += 1 # Save new solution if accept == True: # If the new solution is accepted, the new solution is saved as the current solution xNow[:] = xNew[:] fxNow = fxNew if fxNew < fxBest: # If the objective function of the new solution is better than the optimal solution, the new solution is saved as the optimal solution fxBest = fxNew xBest[:] = xNew[:] totalImprove += 1 scale = scale*0.99 # Variable search step size, gradually reduce the search range and improve the search accuracy # ---Data sorting after internal circulation # Complete the search of current temperature, save data and output pBadAccept = kBadAccept / (kBadAccept + kBadRefuse) # Acceptance probability of inferior solution recordIter.append(kIter) # Current external circulation times recordFxNow.append(round(fxNow, 4)) # Objective function value of current solution recordFxBest.append(round(fxBest, 4)) # Objective function value of optimal solution recordPBad.append(round(pBadAccept, 4)) # Objective function value of optimal solution if kIter%10 == 0: # Modular operation, remainder of quotient print('i:{},t(i):{:.2f}, badAccept:{:.6f}, f(x)_best:{:.6f}'.\ format(kIter, tNow, pBadAccept, fxBest)) # Slowly cool down to a new temperature. Cooling curve: T(k)=alfa*T(k-1) tNow = tNow * alfa kIter = kIter + 1 # ======End the simulated annealing process====== print('improve:{:d}'.format(totalImprove)) return kIter,xBest,fxBest,fxNow,recordIter,recordFxNow,recordFxBest,recordPBad # Result verification and output def ResultOutput(cName,nVar,xBest,fxBest,kIter,recordFxNow,recordFxBest,recordPBad,recordIter): # ======Verification and output of optimization results====== fxCheck = cal_Energy(xBest,nVar) if abs(fxBest - fxCheck)>1e-3: # Objective function test print("Error 2: Wrong total millage!") return else: print("\nOptimization by simulated annealing algorithm:") for i in range(nVar): print('\tx[{}] = {:.6f}'.format(i,xBest[i])) print('\n\tf(x):{:.6f}'.format(fxBest)) return # main program def main(): # Parameter setting, optimization problem parameter definition, simulated annealing algorithm parameter setting [cName, nVar, xMin, xMax, tInitial, tFinal, alfa, meanMarkov, scale] = ParameterSetting() # print([nVar, xMin, xMax, tInitial, tFinal, alfa, meanMarkov, scale]) # Simulated annealing algorithm [kIter,xBest,fxBest,fxNow,recordIter,recordFxNow,recordFxBest,recordPBad] \ = OptimizationSSA(nVar,xMin,xMax,tInitial,tFinal,alfa,meanMarkov,scale) # print(kIter, fxNow, fxBest, pBadAccept) # Result verification and output ResultOutput(cName, nVar,xBest,fxBest,kIter,recordFxNow,recordFxBest,recordPBad,recordIter) if __name__ == '__main__': main()

### 4. Program running results

x_Initial:-143.601793,331.160277, f(x_Initial):959.785447 i:0,t(i):100.00, badAccept:0.469136, f(x)_best:300.099320 i:10,t(i):81.71, badAccept:0.333333, f(x)_best:12.935760 i:20,t(i):66.76, badAccept:0.086022, f(x)_best:2.752498 ... i:200,t(i):1.76, badAccept:0.000000, f(x)_best:0.052055 i:210,t(i):1.44, badAccept:0.000000, f(x)_best:0.009448 i:220,t(i):1.17, badAccept:0.000000, f(x)_best:0.009448 improve:18 Optimization by simulated annealing algorithm: x[0] = 420.807471 x[1] = 420.950005 f(x):0.003352

Copyright description:

Original works

Copyright 2021 YouCans, XUPT

Crated: 2021-05-01