clc;
clear;
f = inline('20 x.^2 y.^2 -10*(cos(2*pi*x) cos(2*pi*y))');
x = -5:0.01:5;
y = -5:0.01:5;
[X,Y] = meshgrid(x,y);
F = f(X,Y);
figure(1);
mesh(X,Y,F);
xlabel('横坐标x'); ylabel('纵坐标y'); zlabel('空间坐标z');
hold on;
lower_x = -5;
upper_x = 5;
lower_y = -5;
upper_y = 5;
ant = 1000;
times = 300;
rou = 0.9;
p0 = 0.2;
ant_x = zeros(1,ant);
ant_y = zeros(1,ant);
tau = zeros(1,ant);
Macro = zeros(1,ant);
for i=1:ant
ant_x(i) = (upper_x-lower_x)*rand() lower_x;
ant_y(i) = (upper_y-lower_y)*rand() lower_y;
tau(i) = f(ant_x(i),ant_y(i));
% plot3(ant_x(i),ant_y(i),tau(i),'k*');
hold on;
Macro = zeros(1,ant);
end
fprintf('蚁群搜索开始(找最小值,绿色):n');
T = 1;
trap = 0;
tau_best = zeros(1,times);
p = zeros(1,ant);
while T < times
lamda = 1/T;
[tau_best(T),bestindex] = min(tau);
plot3(ant_x(bestindex), ant_y(bestindex), f(ant_x(bestindex),ant_y(bestindex)), 'r*');
hold on;
for i = 1:ant
p(i) = (tau(bestindex) - tau(i))/tau(bestindex);
end
for i = 1:ant
if p(i) < p0
tempx = ant_x(i) 0.5*(2*rand-1)*lamda;
tempy = ant_y(i) 0.5*(2*rand-1)*lamda;
else
tempx = ant_x(i) (upper_x-lower_x)*(rand-2.5); % 5/2=2.5
tempy = ant_y(i) (upper_y-lower_y)*(rand-2.5);
end
if tempx < lower_x
tempx = lower_x;
end
if tempx > upper_x
tempx = upper_x;
end
if tempy < lower_y
tempy = lower_y;
end
if tempy > upper_y
tempy = upper_y;
end
if f(tempx,tempy) < tau(i)
ant_x(i) = tempx;
ant_y(i) = tempy;
Macro(i) = f(tempx,tempy);
end
end
for i = 1:ant
tau(i) = (1-rou)*tau(i) Macro(i);
end
T = T 1;
if T >= times
fprintf('蚁群搜索到的最小值点:(%.5f,%.5f,%.5f)n',...
ant_x(bestindex), ant_y(bestindex), f(ant_x(bestindex),ant_y(bestindex)));
fprintf('搜索次数:%dn',T)
if f(ant_x(bestindex),ant_y(bestindex)) > 0.8
trap = trap 1; % 进入陷阱次数 1
fprintf('发生大概率事件: 陷入局部极小值,自动再次执行程序!nn')
T = 1;
tau_best = zeros(1,times);
p = zeros(1,ant);
for i=1:ant
ant_x(i) = (upper_x-lower_x)*rand() lower_x;
ant_y(i) = (upper_y-lower_y)*rand() lower_y;
tau(i) = f(ant_x(i),ant_y(i));
Macro = zeros(1,ant);
end
end
end
end
hold on;
trap = trap 1;
syms x y;
f = 20 x^2 y^2 -10*(cos(2*pi*x) cos(2*pi*y));
fx = diff(f,x);
fy = diff(f,y);
acc = 0.0001;
study_step = 0.001;
x = ant_x(bestindex);
y = ant_y(bestindex);
k = 0;
fprintf('走出局部极值,梯度下降精确搜索开始(黄线):n');
while eval(fx)~=0 | eval(fy)~=0
ans_tmp = [x,y] - study_step*[eval(fx),eval(fy)];
acc_tmp = sqrt((ans_tmp(1)-x)^2 (ans_tmp(2)-y)^2);
if acc_tmp <= acc
fprintf('精确极值坐标为:(%.5f,%.5f,%.5f)n',ans_tmp(1),ans_tmp(2),f_tmp);
fprintf('迭代次数:%dnn',k);
plot3(ans_tmp(1),ans_tmp(2),f_tmp,'y.');
hold off
break;
end
x = ans_tmp(1);
y = ans_tmp(2);
f_tmp = eval(f);
plot3(x,y,f_tmp,'y.')
hold on;
k = k 1;
end
if trap > 0
fprintf('本模型走出陷阱概率:%.1f%%n',1/trap*100);
end
