[Из песочницы] Оптимизация на примере. Имитационный отжиг против муравьиного алгоритма. Часть 131.08.2016 20:33

% метод отжига ( на примере задачи коммивояжера )

%--------------------------------------------------------------------------
tic

% clearvars -except cities
clearvars

% -----------------------------ИТЕРАЦИИ------------------------------------
% кол-во итераций
m = 100000;
% -----------------------------ПАРАМЕТРЫ-----------------------------------
% начальная температура
Tstart = 10000;

% конечная температура
Tend = 0.1;

% начальная температура для вычислений
T = Tstart;

% расстояние
S = inf;

% количество городов
n = 100;

% --------------------------------ПАМЯТЬ-----------------------------------

% матрица расстояний
dist = zeros(n,n); 

% -------------------------------------------------------------------------

% генерация городов (x,y)
cities = rand(n,2)*100;

% задаем случайный маршрут
ROUTE = randperm(n);

% создаем матрицу расстояний
for i = 1:n
    
    for j = 1:n
        
        % dist ( расстояния )
        dist(i,j) = sqrt((cities(i,1) - cities(j,1))^2 + ...
           (cities(i,2) - cities(j,2))^2);       
            
    end
    
end
 
% поехали оптимизировать, время от кол-ва итераций

for k = 1:m
    
    % сбрасываем потенциальное расстояние
    Sp = 0;
    
    % здесь условие создания потенциальных маршрутов, ROUTEp -
    % потенциальный маршрут
    
    % потенциальный маршрут
    ROUTEp = ROUTE;

    % два случайных города
    transp = randperm(n,2);

    % переворот вектора, взяли два сгенерированных числа transp и
    % перевернули маршрут между ними
    if transp(1) < transp(2)
        ROUTEp(transp(1):transp(2)) = fliplr(ROUTEp(transp(1):transp(2)));
    else
        ROUTEp(transp(2):transp(1)) = fliplr(ROUTEp(transp(2):transp(1)));
    end
        
    % вычисляем энергию (расстояние) потенциального маршрута
    for i = 1:n-1
        Sp = Sp + dist(ROUTEp(i),ROUTEp(i+1));
    end
        Sp = Sp + dist(ROUTEp(1),ROUTEp(end));
        
    % если он меньше то, потенциальный маршрут становится основным ...
    % если нет, смотрим, осуществляется ли переход  
    
    if Sp < S
        S = Sp;
        ROUTE = ROUTEp;              
    else
        % случайно выбранное число от 0 до 1
        RANDONE = rand(1);
        % вычисляем вероятность перехода
        P = exp((-(Sp - S)) / T);
       
            if RANDONE <= P
                S = Sp;
                ROUTE = ROUTEp;                     
            end
        
    end
    
    	% уменьшаем температуру
        T = Tstart / k;

        % проверяем условие выхода
        if T < Tend
            break;
        end;
       
end

% рисуем графику

citiesOP(:,[1,2]) = cities(ROUTE(:),[1,2]);
plot([citiesOP(:,1);citiesOP(1,1)],[citiesOP(:,2);citiesOP(1,2)],'-r.')

% очищаем переменые
clearvars -except cities ROUTE S

% смотрим время
toc

% муравьиный алгоритм ( на примере задачи коммивояжера )

% -------------------------------------------------------------------------
tic

% clearvars -except cities
clearvars
 
% ------------------------------ИТЕРАЦИИ-----------------------------------

% кол-во итераций ( поколений )
age = 50;

% кол-во муравьев в поколении
countage = 30;

% кол-во городов
n = 30;

% ------------------------------ПАРМЕТРЫ-----------------------------------

% альфа - коэффициент запаха, при 0 будем ориентироваться только на
% кратчайший путь 
a = 1; 

% бета - коэффициент расстояния, при 0 будем
% ориентироваться только на оставляемый запах
b = 2.5;

% коэффициент испарения ( память поколений )
e = 0.5;

% количество выпускаемых феромонов 
Q = 1;

% кол-во элитных муравьев ( увеличение лучшей дистанции в elite раз )
elite = 1;

% начальный феромон
ph = 0.01;

% -------------------------------ПАМЯТЬ------------------------------------
% матрица расстояний
dist = zeros(n,n); 

% матрица обратных расстояний
returndist = zeros(n,n); 

% матрица маршрута муравьев в одном поколении
ROUTEant = zeros(countage,n);

% вектор расстояний муравьев в одном поколении
DISTant = zeros(countage,1); 

% вектор лучших дистанций на каждой итерации
bestDistVec = zeros(age,1);

% лучший начальный маршрут
bestDIST = inf; 

% оптимальные маршруты
ROUTE = zeros(1,n+1);

% перестановка городов без повторений ( для выхода муравьев )
RANDperm = randperm(n);

% матрица вероятностей
P = zeros(1,n);

% -------------------------------------------------------------------------

% генерация городов (x,y)
cities = rand(n,2)*100;

% матрица начальных феромонов
tao = ph*(ones(n,n));

% создаем матрицу расстояний и матрицу обратных расстояний
for i = 1:n
    
    for j = 1:n
        
        % dist ( расстояния )
        dist(i,j) = sqrt((cities(i,1) - cities(j,1))^2 + ...
           (cities(i,2) - cities(j,2))^2);
        
            % nn ( обратные расстояния )
            if i ~= j
            returndist(i,j) = 1/sqrt((cities(i,1) - cities(j,1))^2 + ...
                (cities(i,2) - cities(j,2))^2); 
            end
            
    end
    
end

% итерации
for iterration = 1:age
    
    
    % муравьи ( одно поколение)
    for k = 1:countage
        
    % ****************** НАЧАЛЬНОЕ РАСПОЛОЖЕНИЕ МУРАВЬЕВ ******************
    % выбирайте какой нужно
    
    % каждый муравей располагается случайно    
%       ROUTEant(k,1) = randi([1 n]);
    
    % с каждого города выходит один муравей ( без совпадений ), кол-во
    % городов и кол-во муравьев в поколении должны быть равны
      ROUTEant(k,1) = RANDperm(k);
    
    % с конкретного города выходят все муравьи в данном случа с 1-ого
%       ROUTEant(k,1) = 1;

% тут маршрут первому поколению задаем либо произвольный, либо с каждого
% города разный, либо с одного города все, а следующее поколение выходит по
% концам первых

%     if iterration == 1
%       ROUTEant(k,1) = randi([1 n]);
% %       ROUTEant(k,1) = RANDperm(k);
% %     ROUTEant(k,1) = 1;
%     else
%       ROUTEant(k,1) = lastROUTEant(k);  
%     end
    
    % *********************************************************************
    
    % путь каждого муравья, начиная со второго, так как первый выбран
    for s = 2:n  

        % полуаем индекс выбранного города
        ir = ROUTEant(k,s-1);
 
        % вероятность посещения городов ( числитель ) , в числителе у нас
        % следующее: tao^a*(1/S)^b 
        % 1/S -это returndist. 
        
        % поскольку данное значение будет повторяться (кол-во муравьев * на
        % колонию * кол-во городов) раз, то еще один цикл писать не выгодно,
        % скорость работы при таких вычислениях падает. Поэтому написал в 
        % этом моменте векторно. На обычном языке будет так: 
        
%         for c = 1:n             
%             P(1,c) = tao(ir,c).^a * returndist(ir,c).^b;          
%         end

        P = tao(ir,:).^a .* returndist(ir,:).^b;
        % получили числители (в формуле вероятности перехода к k-ому городу)
        % для n городов, однако в некоторых мы уже побывали, нужно исключить
        % их
        
        % проставляем нули в числитель туда, где уже были, чтобы
        % вероятность перехода была 0, следовательно в сумме знаменателя
        % формулы данный город учитываться не будет    
        P(ROUTEant(k,1:s-1)) = 0;

        % получаем вероятности перехода ( сумма строк должна быть = 1 )
        P = P ./ sum(P);
        
        % смотрим в какой город осуществляется переход
        RANDONE = rand;
        getcity = find(cumsum(P) >= RANDONE, 1, 'first');
        
        % если с вероятностями что-то не так, выдача ошибки
        if isempty(getcity) == 1
            disp 'Ошибка в вероятностях!'
        end

        % присваем s-ый город в путь k-ому муравью
        ROUTEant(k,s) = getcity;

    end
        
    % получаем маршрут k-ого муравья
    ROUTE = [ROUTEant(k,1:end),ROUTEant(k,1)];
    
    % сброс длины
    S = 0;
    
    % вычисляем маршрут k-ого муравья
    for i = 1:n
        S = S + dist(ROUTE(i),ROUTE(i+1));
    end
    
    % путь k-ого муравья, массив дистанций k-ых муравьев
    DISTant(k) = S;
    
    % присваевыем лучший маршрут и S     
    if DISTant(k) < bestDIST
        bestDIST = DISTant(k);
        bestROUTE = ROUTEant(k,1:end);
    end
    
    % вектор "последних" городов k-ых муравьев ( выбирается для старта
    % муравьев нового поколения с тех городов, где закончили путь
    % предыдущее поколение)
    
%     lastROUTEant = ROUTEant(1:end,end); 
    
    end
% -----------------------ИЗМЕНЯЕМ КОЛИЧЕСТВО ФЕРОМОНОВ---------------------

    % Испаряем феромоны "старого" пути е - коэффициент испарения
    tao = (1-e)*tao;
    
    for u = 1:countage
        
        % получаем маршрут u-ого муравья
        ROUTE = [ROUTEant(u,1:end), ROUTEant(u,1)];
        
        % для каждого города
        for t = 1:n
            
            x = ROUTE(t);
            y = ROUTE(t+1);
            
            % считаем новый феромон
            if DISTant(u) ~= max(DISTant)
                % тут по формуле добавки феромона в зависимости от
                % дистанции
                tao(x,y) = tao(x,y) + Q / DISTant(u);
            else
                % усиливаем ребера лучшего пути
                tao(x,y) = tao(x,y) + (elite*Q) / DISTant(u);
            end
            
        end

        
    end
               
end 

% строим графику
citiesOP(:,[1,2]) = cities(bestROUTE(:),[1,2]);
plot([citiesOP(:,1);citiesOP(1,1)],[citiesOP(:,2);citiesOP(1,2)],'.r-')

clearvars -except cities ROUTE bestDIST bestROUTE ROUTEant

toc