tsp_ga1.m

用遗传算法求解TSP（旅行商）问题

function TSP_GA(ptCorr,nPoint)

% start timing...
tic;
dist = zeros(nPoint);

% calculate the distance matrix...
for i=1:nPoint
    for j=1:nPoint
        dist(i,j) = norm(ptCorr(i,:)-ptCorr(j,:));
    end % for j...
end % for i...

% number of population...
K = 4;
nPopulation = nPoint*K; % to be expected as 100...
% probability of cross...
pCross = 0.5;
% number of cross operations...
nCross = nPopulation*0.8;
% probability of single cross...
pSingleCross = 0.5;
% probability of mutation...
pMutation = 0.05;
% terminating relative error...
terminateError = 1e-4;

% bounds of fitness function...
maxFitness = 0;
minFitness = inf;
curFitness = inf;
curPopPos = 0;
% bound of iterations...
maxIteration = 2e3;
curIteration = 0;

% initialization of population using randomized permutations...
% and the curPopPos is the individual with the best fitness...
curPopulation = zeros(nPopulation,nPoint);
curFitnessVec = zeros(1,nPopulation);
for i=1:nPopulation
    curPopulation(i,:) = randperm(nPoint);
    curFitnessVec(i) = TSP_fitness(curPopulation(i,:),dist);
    if (curFitness > curFitnessVec(i))
        curFitness = curFitnessVec(i);
        curPopPos = i;
    end % if...
%     maxFitness = max(maxFitness, curStepFitness);
%     minFitness = min(minFitness, curStepFitness);
end % for...
maxFitness = max(curFitnessVec);
minFitness = min(curFitnessVec);

curPath = curPopulation(curPopPos,:);
curOutput = [curPath curPath(1)];

% initial display...
plot(ptCorr(curOutput,1),ptCorr(curOutput,2),'b*');
title('Travelling Salesman Problem Solution Using Genetic Algorithm');
hold on;
figureHandle = plot(ptCorr(curOutput,1),ptCorr(curOutput,2),'r');
hold off;
figure;
% dynamic vector to store all the energy state...
engRecord = curFitness;
energyHandle = plot(1:size(engRecord,2),engRecord);
title('Fitness of Travelling Salesman Problem Solution');
xlabel('Iterations');
ylabel('Fitness Value');

while(curIteration<maxIteration && (maxFitness-minFitness)>terminateError*maxFitness)
    % selection...
    newPopulation = zeros(nPopulation,nPoint);
    for i=1:nPopulation
        % tournament selection...
        % random select two individual and compete, the one with larger
        % fitness wins...
        competitors = round(rand(1,2)*nPopulation+0.5);
        if (TSP_fitness(curPopulation(competitors(1),:),dist) < TSP_fitness(curPopulation(competitors(2),:),dist))
            newPopulation(i,:) = curPopulation(competitors(1),:);
            curFitnessVec(i) = TSP_fitness(curPopulation(competitors(1),:),dist);
        else
            newPopulation(i,:) = curPopulation(competitors(2),:);
            curFitnessVec(i) = TSP_fitness(curPopulation(competitors(2),:),dist);
        end % if...
    end % for...
    curPopulation = newPopulation;
    
%     % cross...
%     for i=1:nCross
%         % there are two kinds of cross operations, 000|111 and 00|11|00,
%         % where 00 stands for the part to be crossed...
%         % single point cross...
%         if (rand < pSingleCross)
%             % permit to cross...
%             if (rand < pCross)
%                 crossors = round(rand(1,2)*nPopulation+0.5);
%                 crossPos = round(rand*nPoint+0.5);
%                 temp = curPopulation(crossors(1),1:crossPos);
%                 curPopulation(crossors(1),1:crossPos) = curPopulation(crossors(2),1:crossPos);;
%                 curPopulation(crossors(2),1:crossPos) = temp;
%             end % if...
%         % two points cross...
%         else
%             % permit to cross...
%             if (rand < pCross)
%                 crossors = round(rand(1,2)*nPopulation+0.5);
%                 crossPos = round(rand(1,2)*nPoint+0.5);
%                 % cross the first part...
%                 temp = curPopulation(crossors(1),1:crossPos(1));
%                 curPopulation(crossors(1),1:crossPos(1)) = curPopulation(crossors(2),1:crossPos(1));
%                 curPopulation(crossors(2),1:crossPos(1)) = temp;
%                 % leave the second part unchanged...
%                 % cross the third part...
%                 temp = curPopulation(crossors(1),crossPos(2):nPoint);
%                 curPopulation(crossors(1),crossPos(2):nPoint) = curPopulation(crossors(2),crossPos(2):nPoint);
%                 curPopulation(crossors(2),crossPos(2):nPoint) = temp;                
%             end % if...
%         end % if...
%     end % for...
    
    % cross operation, renewed version...
    % two position cross...
    for i=1:nCross
        if (rand < pCross)
            crossors = round(rand(1,2)*nPopulation+0.5);
            crossPos = round(rand(1,2)*nPoint+0.5);
            [curPopulation(crossors(1),:),curPopulation(crossors(2),:)] = TSP_cross(curPopulation(crossors(1),:),curPopulation(crossors(2),:),min(crossPos),max(crossPos));
            curFitnessVec(crossors(1)) = TSP_fitness(curPopulation(crossors(1),:),dist);
            curFitnessVec(crossors(2)) = TSP_fitness(curPopulation(crossors(2),:),dist);
        end % if...
    end % for...
    
    % mutation...
    % mutation process is simple that for every individual, try the
    % mutation process...
    
    % pMutation = pMutation+ (curIteration)/maxIteration*0.2;
    for i=1:nPopulation
        % permit and luck enough to mutation...
        if (rand < pMutation)
            mutationPos = round(rand(1,2)*nPoint+0.5);
%             temp = curPopulation(i,mutationPos(1));
%             curPopulation(i,mutationPos(1)) = curPopulation(i,mutationPos(2));
%             curPopulation(i,mutationPos(2)) = temp;       
            curPopulation(i,min(mutationPos):max(mutationPos)) = fliplr(curPopulation(i,min(mutationPos):max(mutationPos)));
            curFitnessVec(i) = TSP_fitness(curPopulation(i,:),dist);
        end % if...
    end % for...
    
    curPopPos = 0;
    curFitness = inf;
    maxFitness = 0;
    minFitness = inf;
    % curStepFitness = 0;
    % find the best individual...
%     for i=1:nPopulation
%         cutStepFitness = TSP_fitness(curPopulation(i,:),dist);
%         
%         if (curFitness > curStepFitness)
%             curFitness = curStepFitness;
%             curPopPos = i;
%         end % if...
% 
%         maxFitness = max(maxFitness, curStepFitness);
%         minFitness = min(minFitness, curStepFitness);
%     
%     end % for...
    [curFitness,curPopPos] = min(curFitnessVec);
    maxFitness = max(curFitnessVec);
    minFitness = min(curFitnessVec);
        
   
    % displaying the output on the screen in every iteration...
    curOutput = [curPopulation(curPopPos,:) curPopulation(curPopPos,1)];
    
    % add the most up-to-date energy state to record...
    engRecord = [engRecord curFitness];
    set(figureHandle,'xdata',ptCorr(curOutput(:),1),'ydata',ptCorr(curOutput(:),2));
    set(energyHandle,'xdata',1:size(engRecord,2),'ydata',engRecord);
    drawnow;
    
    % output optimization information in this iteration...
    % fprintf('[%s]Current temperature: %f\n',datestr(now),curTemp);
    % fprintf('[%s]Current path :%s\n',datestr(now),sprintf('%d ',curPopulation(curPopPos,:)));
    fprintf('[%s]Current fitness[min max] :%f[%f,%f]\n',datestr(now),curFitness,minFitness,maxFitness);
    fprintf('[%s]Time elapsed: %f s\n',datestr(now),toc);
    
    curIteration = curIteration + 1;
    
end % while...

% final marking the citis and orders...
figure(1);
for i=1:nPoint
    text(ptCorr(curOutput(i),1),ptCorr(curOutput(i),2),['\cdot',int2str(i)],'fontsize',8);
end % i...