Problem Description:
We have a base, longitude and latitude for (70,40). Let's say our plane is 1000 km/h.
We sent a plane from the base to scout out all enemy targets before returning to the original base. In every enemy's eye
The time spent on the aircraft (assuming that we can take a long cruise time) is not counted for the detection time of the punctuation.
This is a traveling business problem. We gave the base number 1 in turn, the enemy target was numbered 2, 3, ...,
101, finally our base repeat number is 102 (this is easy to calculate in the program). Distance matrix D = ( dij ) 102x102,
which
D ij represents the distance of I, J two points, I, j = 1,2,l, 102, where D is the real symmetric matrix. The problem is
Ask for a shortest path from point 1, which travels all the intermediate points and reaches point 102.
The algorithm is as follows:
% simulated annealing algorithm travel peddler Algorithm%ij.txt store is the latitude and longitude of each point clc;clear load c:\users\littlelala\desktop\ij.txt;% load Target data ij.txt store 25x8 Matrix X=ij (:, 1:2:8 );%x is 25x4 matrix x=x (:);%x to 10x1 Matrix Y=ij (:, 2:2:8)%y 25x4 matrix y=y (:);%y to 100x1 matrix sj=[x y];%100x2 Matrix d1=[70,40];% start point sj=[d1;sj; d1];% Add the starting point and end point sj=sj*pi/180;% angle to radians% calculate distance d=zeros (102);% build a 102x102 empty distance matrix, then fill in the corresponding distance%A (r cos X1 cos y1, r sin X1 cos y 1, R sin Y1)%B (r cos x2 cos y2, r sin x2 cos y2, r sin y2)% where R = 6370 is the radius of the Earth. %A, the actual distance of the B two point:%d = R arccos[cos (x1−x2) cos y1 cos y2 + sin y1 sin y2] for i=1:101 for j=i+1:102 D (i,j) =63 70*acos (cos (SJ (i,1)-SJ (j,1)) *cos (SJ (i,2)) *cos (SJ (j,2)) ... +sin (SJ (i,2)) *sin (SJ (j,2))); End end% The above d is the upper triangular matrix, according to symmetry, as a symmetric matrix, D matrix contains the distance information between the points d=d+d '; S0=[]; Sum=inf; Rand (' state ', sum (clock)); The% function is to define a time-varying initial value%matlab the random generating function is basically based on the rand as the basis function through the function of the relationship,% such as NORMRND,UNIDRND, Every time you restart MATLAB, run the compiled function with random number generation you will get the same result,% such as my computer restart Matlab run UNIDRND (100), each time the value is 82, because the initial values of the RAND function are the same,% So in order to avoid these problems often run before the program or add command rand (' state ', sum (clock)), so restart Matlab, run random number generationThe value is different. For i=1:1000% Loop 1000 times, as much as possible to make all the sorting situation appear once s=[1,1+randperm (), 102];%s for the Middle 2 to 101 of the permutation combination temp=0; For j=1:101 Temp=temp+d (S (j), S (j+1)); End If Temp<sum sum=temp; s0=s;% Save the sorting of the current minimum value end end e=0.1^30;% The selected terminating temperature E = 10−30 to determine if the annealing process is over at=0.999;% the selected temperature coefficient α is cooled by: T =αt l=20000; T=1; The% annealing process for k=1:l C=2+floor (100*rand); the%floor function functions as "rounding down", or "rounded down", which takes a maximum integer not greater than x%ceil function to take up the whole c=sort (c);% [A,b]=sort (X) is sorted by column from small to large, while [A,b]=sort (x,2) is by row, B is the sort case c1=c (1); C2=c (2); % calculates the cost function value Df=d (S0 (c1-1), S0 (C2)) +d (S0 (C1), S0 (c2+1))-D (S0 (c1-1), S0 (C1))-D (S0 (C2), S0 (c2+1)); If df<0% acceptance criteria S0=[s0 (1:c1-1), S0 (C2:-1:C1), S0 (c2+1:102)]; SUM=SUM+DF; ElseIf exp (-df/t) >rand (1)% with probability exp (−df/t) Accept new Path% Note When ElseIf instead of else if S0=[S0 (1:c1-1), S0 (C2:-1:C1), S0 (c2+1:102)]; SUM=SUM+DF; End t=t*at;% annealing if t<e% reach the terminating temperature break; End End % output path and path length s0,sum% draw one of the cruise path plot (SJ (s0,1)/PI*180,SJ (s0,2)/pi*180); Hold on Plot (SJ (s0,1)/PI*180,SJ (s0,2)/pi*180, ' Rx '); Axis ([ -5,75,-5,45]);
Operation Result:
Simulated annealing algorithm of MATLAB essay