Simulated annealing algorithm:

In order to solve the problem of local optimal solution, Kirkpatrick et al proposed simulated annealing algorithm (SA) in 1983, which can effectively solve the problem of local optimal solution. We know that in the world of molecules and atoms, the higher the energy, the less stable the molecules and atoms, and the lower the energy, the more stable the atoms. ‘Annealing’ is a physical term for the process of heating and cooling an object. The simulated annealing algorithm is derived from the process of crystal cooling. If the solid is not in the lowest energy state, the solid is heated and then cooled. As the temperature drops slowly, the atoms in the solid are arranged in a certain shape to form a regular crystal with high density and low energy, which corresponds to the global optimal solution in the algorithm. If the temperature drops too fast, the atoms may not have enough time to form a crystal structure, resulting in an amorphous structure with higher energy, which is the local optimal solution. Therefore, according to the annealing process, a little energy can be added to it, and then it is cooled. If the energy is increased and the local optimal solution is jumped out, the annealing is successful. Let’s talk in detail about how it jumps out from the local optimal solution to the global optimal solution.

The simulated annealing algorithm consists of two parts, namely Metropolis algorithm and annealing process. The Metropolis algorithm, which is how to get it to jump out in the case of a locally optimal solution, is the basis of annealing. In 1953, Metropolis proposed the importance sampling method, that is, to accept new states with probability rather than using completely determined rules, which was called Metropolis criterion and required relatively low computation. Here is the first image to say, and then in the mathematical formula:

Suppose start state on A, as the number of iterations updates to the local optimal solution of B, then find updates to B, the ability is lower than A, then close to the optimal solution, so hundred transfer, after state to B, find the next step up energy, if is gradient descent is not allowed to move on, and here will jump out of the pit at A certain probability, All of these probabilities depend on the state, the energy, and so on, and we’ll talk about it in detail, if B does eventually jump out and it gets to C, it will continue to jump out with a certain probability, and you might be wondering if you can jump back to B? As we’ll explain, until we get to D, it’s going to be stable. So the design of this probability is very important, so let’s explain it mathematically.

Let’s say that the former state is, the system according to a certain index (gradient descent, energy of the upper node), the state becomesCorrespondingly, the energy of the system is given byinto, defines the system byintoThe acceptance probability P is:

                                         

We can see from the type, if energy is reduced, so the transfer is accepted (probability of 1), if the energy increased, suggests further the system deviates from the value of the global optimal position, this algorithm is not immediately to abandon it, but for probability operation: first of all, in the interval [0, 1] produce a uniformly distributed random NumbersIf theP, the transfer is accepted, otherwise the transfer is rejected, and the next step is entered, and the cycle is repeated. Where P determines the probability P in terms of the change in energy and T, so this value is dynamic.

Parameter control of annealing algorithm

Metropolis algorithm is the basis of simulated annealing algorithm, but direct use of Metropolis algorithm may lead to too slow optimization speed to be used in practice. In order to ensure convergence in a limited time, it is necessary to set parameters for controlling algorithm convergence. In the formula above, the adjustable parameter is T; if T is too large, If the value is small, the calculation time will increase. In practical application, the annealing temperature table is used, and a large T value is adopted at the initial annealing stage, which will gradually decrease with the progress of annealing. The details are as follows:

(1) The initial temperature T(0) should be selected high enough to allow all transition states of. The higher the initial temperature, the higher the probability of obtaining high-quality solutions, and the longer the time.

(2) Annealing rate. The easiest way to go down is to go down exponentially:

                                                    

Among themIs a positive number less than 1, usually between 0.8 and 0.99. For each temperature, there are enough transfer attempts. The convergence rate of exponential descent is relatively slow. Other descent methods are as follows:

                                                      

                                                        

(3) Termination temperature

Annealing is complete if each new state can be updated or a user-set threshold is reached in the case of several iterations.

Simulated annealing steps:

1. Simulated annealing algorithm can be decomposed into three parts: solution space, objective function and initial solution.

2. Basic idea of simulated annealing:

(1) Initialization: initial temperature T(sufficiently large), initial solution state S(the starting point of algorithm iteration), and iteration times of each T value L

(2) for k = 1,… L Do steps (3) to 6:

(3) A new solution S ‘is produced

(4) Calculate the increment δ T=C(S ‘) -c (S), where C(S) is the cost function

(5) If δ T<0, accept S ‘as the new current solution; otherwise, accept S’ as the new current solution with the probability exp(- δ T/T).

(6) If the termination condition is met, output the current solution as the optimal solution and end the program.

The termination condition is usually taken to terminate the algorithm when several new solutions are not accepted.

(7) T gradually decreases, and T->0, then go to step 2.

The generation and acceptance of new solutions of simulated annealing algorithm can be divided into the following four steps:

The first step is to generate a new solution in the solution space from the current solution by a production function. To facilitate subsequent calculations and accept, reducing algorithm is time-consuming and often choose the current data processing after a simply way to produce the new transformation, such as to constitute all or part of the new element substitution, interchange, etc., noticed that the transformation method of generation of the new neighborhood structure determines the data processing, thus has certain influence the selection of the cooling schedule.

The second step is to calculate the difference of the objective function corresponding to the new solution. Because the objective function difference is generated only by the transformation part, it is best to calculate the objective function difference increments. It turns out that for most applications, this is the fastest way to calculate the difference of the objective function.

The third step is to judge whether the new solution is accepted based on an acceptance criterion. The most commonly used acceptance criterion is the Metropolis criterion: if δ T<0, S ‘is accepted as the new current solution S; otherwise, S’ is accepted as the new current solution S with the probability exp(- δ T/T).

The fourth step is to replace the current solution with the new solution when the new solution is determined to be accepted. This only requires the realization of the transformation part corresponding to the new solution in the current solution and the modification of the objective function value. At this point, the current solution implements an iteration. This could be the basis for the next round of testing. When the new solution is judged to be abandoned, the next round of experiment is continued on the basis of the original current solution.

The simulated annealing algorithm is independent of the initial value, and the solution obtained by the algorithm is independent of the initial solution state S(the starting point of the algorithm iteration). The simulated annealing algorithm has asymptotic convergence and has been proved to be a global optimization algorithm that converges to the global optimal solution with probability L in theory. Simulated annealing algorithm has parallelism.

Annealing algorithm program flow chart;

Job shop scheduling problem description

Job Shop Scheduling (JSP) is one of the most classic NP-hard problems. Its application fields are extremely wide, involving aircraft carrier dispatching, airport aircraft dispatching, port and wharf cargo ship dispatching, automobile processing assembly line and so on.

JSP problem description: A processing system has M machines, and N jobs are required to be processed, in which, job I contains the number of processes for Li. , then L is the total ordinal number of tasks set. Among them, the processing time of each process has been determined, and each operation must be processed in accordance with the sequence of the process. The task of scheduling is to arrange the processing scheduling of all jobs, so that the performance indicators can be optimized when the constraints are met.

Job shop scheduling needs to consider the following constraints:

Cons1: Each process is processed on the specified machine, and processing can only begin after the previous process is finished.

Cons2: At a point, one machine can only process one job.

Cons3: Each operation can only be processed 1 time on 1 machine;

Cons4: The process order and processing time of each operation are known, which does not change with the change of processing order.

Problem instances

An example of job shop scheduling problem is given below, where each process is marked with a pair of values (m, P), where M means that the current process must be processed on the MTH machine, and P means the processing time required for the MTH machine to process the current process. (Note: Machine and job numbers start from 0)

Jop0 = [(0, 3), (1, 2), (2)]

Jop1 = [(0, 2), (2, 1), (1, 4)]

Jop2 = [(1, 4), (2, 3)]

In this example, operation JOP0 has three processes: its first process is marked with (0,3), which means that the first process must be processed on the 0 machine, and it needs 3 units of processing time; Its second process is marked with (1,2), which means that the second process must be processed on the first machine, and it needs 2 units of processing time; The rest is the same. In total, there are eight processes in this example.

A feasible solution to this problem is L= a sequence of 8 process start times, which satisfies the constraint of the problem. The figure below shows an example of a viable solution (note: this solution is not optimal) :

  • CLC clear % = = = = = = = = = data entry, parameter adjustment = = = = = = = = = = = = = = = = = swarminitNum = 20; % Number of particles originally generated; MM=[1 2 3 4 5 6 6 6 6 6 6 6]; % workpiece and process number matrix, MM the first row represents the workpiece, the second row represents the number of processes of each workpiece; machineNum=6; % Number of processing machines; initT=1000; % initial simulated annealing temperature; gen=500; % Number of loop iterations; W1 = 0.35; % variation rate; changeNum=3; % variance transformation logarithm; restrictmatrixM=[3 1 2 4 6 5 2 3 5 6 1 4 3 4 6 1 2 5 2 1 3 4 5 6 3 2 5 6 1 4 2 4 6 1 5 3]; %job-shop machine constraint matrix; restrictmatrixT=[1 3 6 7 3 6 8 5 10 10 10 4 5 4 8 9 1 7 5 5 5 3 8 9 9 3 5 4 3 1 3 3 9 10 4 1]; %job-shop time constraint matrix; % = = = = = = = = = = = = = = = PSO algorithm = = = = = = = = = = = = = = = = = = = = = = = = = = swarminit = cell (1, swarminitNum); swarminitLong=sum(MM(2,:)); % The number of all processes is the particle length; for i=1:swarminitNum, swarminit{i}=randomparticle(MM) ; End % randomly generate the initial particle population [popu,s] = size(swarminit); trace = ones(1,gen); trace(1) = 10000; For I = 1:s, bestfit(I) = 10000; % Initial individual history optimum fitness set to be large enough end Bestpar = swarminit; For u=1:swarminitNum, fitList =[0]; end T=initT; for step = 1:gen, for q=1:swarminitNum, swarminit{j}=cross(bestparticle1,swarminit{j},l4,l3); % particle crossing; end end [a,b,c]=timedecode2(bestparticle,restrictmatrixM,restrictmatrixT,machineNum); Disp ([' optimization objective: [' swarminitNum ']; [' swarminitNum ']; [' swarminitNum ']; Int2str (changeNum)]) disp([' simulated annealing initial value :' int2str(initT) 'simulated annealing final value :' int2str(T)]) disp([' iterated loop value :' int2str(trace)]]) Disp ([' minimum average flow time: 'int2str (a)' maximum completion time: 'int2str' minimum clearance time: 'int2str (b) (c)]) disp ([' optimal particle' int2str (bestparticle)]) pause gant(bestparticle,swarminitLong,restrictmatrixM,restrictmatrixT,b)Copy the code