tsp_aco_1.m

用蚁群算法求解TSP（旅行商）问题

function [curBestDistEver, curIteration, curOutput, engRecord] = TSP_ACO(ptCorr,nPoint)
% global ptCorr;
% global 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...
    dis(i,i) = inf;
end % for i...

maxIteration = 5e2;
curIteration = 1;

% number of ants is equal to number of points...
nAnts = nPoint*5;
curBestPath = randperm(nPoint);
curBestDist = TSP_fitness(curBestPath,dist);
curBestDistPos = 0;
curBestDistEver = inf;
curDistVec = zeros(1,nAnts);
curPopulation = zeros(nAnts,nPoint);
curOutput = zeros(1,nPoint+1);
% number of nonimprovement iterations...
nBadIter = 20;
curBadIter = 0;
% initialization of pheromone...
pher = ones(nPoint,nPoint)*(1/(nPoint^2-nPoint));
for i=1:nPoint
    pher(i,i) = 0;
end % for...
% updating factor formatted as (1-\pho)...
phoStart = 0.04;
phoEnd = 0.09;
pho = 0.1;
decFactor = 0.99999;
K = 1;

while(curIteration < maxIteration && curBadIter < nBadIter)
    
    pho = phoStart - (phoStart-phoEnd)*(curIteration)/maxIteration;
    % pho = pho * decFactor;
    
    % ants go through the cities...
    for i=1:nAnts
        % since the result is a circle, let the starting point is 1
        % anyway...
        curPopulation(i,1) = 1;
        for j=2:nPoint
            curPopulation(i,j) = TSP_roulette(pher,curPopulation(i,1:j-1),i,j-1);
        end % for...
        curDistVec(i) = TSP_fitness(curPopulation(i,:),dist);
    end % for...
    [curBestDist,curBestDistPos] = min(curDistVec);
    
    % update the pheromone...
    % first, evaporation process...
    for i=1:nPoint
        for j=1:nPoint
            pher(i,j) = pher(i,j)*(1-pho);
        end % for...
    end % for...
    % second, reinforcement process...
    for i=1:nPoint-1
        pher(curPopulation(curBestDistPos,i),curPopulation(curBestDistPos,i+1)) = pher(curPopulation(curBestDistPos,i),curPopulation(curBestDistPos,i+1)) + pho/(nPoint*K);
    end % for...
    pher(curPopulation(curBestDistPos,nPoint),curPopulation(curBestDistPos,1)) = pher(curPopulation(curBestDistPos,nPoint),curPopulation(curBestDistPos,1)) + pho/(nPoint*K);
    

    if (curBestDist == curBestDistEver)
        curBadIter = curBadIter + 1;
    else
        % restart counting...
        curBestDistEver = curBestDist;
        curBadIter = 0;
    end % if...
    
    % displaying the output on the screen in every iteration...
    curOutput = [curPopulation(curBestDistPos,:) curPopulation(curBestDistPos,1)];
    
    if (curIteration == 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 = curBestDist;
        energyHandle = plot(1:size(engRecord,2),engRecord);
        title('Fitness of Travelling Salesman Problem Solution');
        xlabel('Iterations');
        ylabel('Fitness Value');
    else
        % add the most up-to-date energy state to record...
        engRecord = [engRecord curBestDist];
        set(figureHandle,'xdata',ptCorr(curOutput(:),1),'ydata',ptCorr(curOutput(:),2));
        set(energyHandle,'xdata',1:size(engRecord,2),'ydata',engRecord);
        drawnow;        
    end % if...

    
    % 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(curBestDistPos,:)));
    fprintf('[%s]Current shortest distance :%f\n',datestr(now),curBestDist);
    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...

engRecord'
save engRecord