1. 原理
物理退火 | 模拟退火 |
---|---|
粒子状态 | 解 |
能量最低态 | 最优解 |
加温过程 | 设置初温 |
等温过程 | Metropolies采样过程 |
冷却过程 | 控制参数下降 |
能量 | 目标函数 |
2. 特点
-
以一定概率接受恶化解(Metropolies算法)
-
在每个控制参数(当前温度、当前误差..)下,由当前迭代点出发,产生邻近的随机状态
-
由控制参数确定的接受准则决定新状态的取舍
-
形成一段新的Markov链
-
缓慢降低控制参数,提高接受准则
-
-
对目标函数要求少
3. 基本参数
-
t_star :初始温度
-
t_end : 冷却温度
-
k : 衰减参数,k∈(0,1)
-
L : Markov链的长度,L∈[100,1000]
-
yz : 容差
-
s : 步长因子
4. 基本流程
-
设置初始温度和初始解,计算初始差值
-
while(终止条件)
-
for i = 1:L (Markov链)
-
产生一个范围内的新解
-
判断新解是否被接受,更新局部最优解(Metropolies算法)
-
更新全局最优解
end
-
-
温度衰减,重新计算差值
end
-
5. 改进方向
-
增加记忆功能:增加存储环节,将目前为止的最好状态储存下来、
-
增加升温、重升温的过程:避免算法陷入局部最优解
-
对于当前状态,采用多次搜索策略:在判断是否接受时,多次判断,取概率接受区间中的最优状态
-
与其他搜索机制结合
6. 经典SA算法:
clear;
close all;
clc;
% x,y的求解范围
x_max = 5;
x_min = -5;
y_max = 5;
y_min = -5;
% 参数初始值
L = 100; % 马尔可夫链长度
k = 0.99; % 衰减因子
s = 0.02; % 步长
t_star = 100; % 初始温度
t_end = 0.001; % 冷却的最终温度
yz = 1e-8; % 容差
p = 0; % 新解被接受的概率(Metropolis算法)
%% 1.设置初始值
% 初始的前一次状态
x0 = rand * (x_max - x_min) + x_min;
y0 = rand * (y_max - y_min) + y_min;
best_x0 = x0;
best_y0 = y0;
% 初始状态
x0 = rand * (x_max - x_min) + x_min;
y0 = rand * (y_max - y_min) + y_min;
best_x = x0;
best_y = y0;
deta = abs(func(best_x0,best_y0)-func(best_x,best_y));
t = t_star;
while(t > t_end && deta > yz)
for i = 1:L
%% 2.在当前位置产生一个新解
p1 = 0; % 临时变量,用来产生一个范围内的随机值
while p1==0
next_x = x0 + s*(rand*(x_max-x_min)+x_min);
next_y = y0 + s*(rand*(y_max-y_min)+y_min);
% 判断是否在范围内
if next_x >= x_min && next_x <= x_max &&...
next_y >= y_min && next_y <= y_max
p1 = 1;
end
end
%% 3.计算评价函数的增量deta
deta_e = func(next_x,next_y) - func(x0,y0);
%% 4.判断是否接受新解(Metropolis算法)
if deta_e < 0 % 接受新值
x0 = next_x;
y0 = next_y;
p = p + 1;
else
p1 = exp(-1 * deta_e / t);
if p1 > rand % 接受较差的值
x0 = next_x;
y0 = next_y;
p = p + 1;
end
end
%% 4.更新全局最优解
deta_e = func(next_x,next_y) - func(best_x,best_y);
if deta_e < 0
% 记录上一个最优解,更新当前最优解
best_x0 = best_x;
best_y0 = best_y;
best_x = next_x;
best_y = next_y;
end
trace(p+1)=func(best_x, best_y);
end
%% 5.温度衰减,重新计算差值
t = t * k;
deta = abs(func(best_x,best_y) - func(best_x0,best_y0));
end
figure
disp('最小值在点:');
best_x
best_y
disp( '最小值为:');
func(best_x, best_y)
plot(trace(2:end))
hold on;
xlabel('迭代次数')
ylabel('适应度值')
title('适应度进化曲线')
legend_str = {'SA','重升温SA','多搜索策略SA','重升温+多搜索SA'};
function results = func(x,y)
results = 5 * cos(x * y)+ x * y + y * y * y;
end
7. 重升温机制改进的SA
%% 重升温机制SA
% x,y的求解范围
x_max = 5;
x_min = -5;
y_max = 5;
y_min = -5;
% 参数初始值
L = 100; % 马尔可夫链长度
k = 0.99; % 衰减因子
s = 0.02; % 步长
t_star = 100; % 初始温度
t_end = 0.001; % 冷却的最终温度
deta_t = 0.1;
yz = 1e-8; % 容差
p = 0; % 新解被接受的概率(Metropolis算法)
%% 1.设置初始值
% 初始的前一次状态
x0 = rand * (x_max - x_min) + x_min;
y0 = rand * (y_max - y_min) + y_min;
best_x0 = x0;
best_y0 = y0;
% 初始状态
x0 = rand * (x_max - x_min) + x_min;
y0 = rand * (y_max - y_min) + y_min;
best_x = x0;
best_y = y0;
deta = abs(func(best_x0,best_y0)-func(best_x,best_y));
flag = 0;
t = t_star;
while(t > t_end && deta > yz)
%% 重升温机制
if flag == 0 && (t - t_end) < deta_t
flag = 1;
t = t_star;
end
for i = 1:L
%% 2.在当前位置产生一个新解
p1 = 0; % 临时变量,用来产生一个范围内的随机值
while p1==0
next_x = x0 + s*(rand*(x_max-x_min)+x_min);
next_y = y0 + s*(rand*(y_max-y_min)+y_min);
% 判断是否在范围内
if next_x >= x_min && next_x <= x_max &&...
next_y >= y_min && next_y <= y_max
p1 = 1;
end
end
%% 3.计算评价函数的增量deta
deta_e = func(next_x,next_y) - func(x0,y0);
%% 4.判断是否接受新解(Metropolis算法)
if deta_e < 0 % 接受新值
x0 = next_x;
y0 = next_y;
p = p + 1;
else
p1 = exp(-1 * deta_e / t);
if p1 > rand % 接受较差的值
x0 = next_x;
y0 = next_y;
p = p + 1;
end
end
%% 4.更新全局最优解
deta_e = func(next_x,next_y) - func(best_x,best_y);
if deta_e < 0
% 记录上一个最优解,更新当前最优解
best_x0 = best_x;
best_y0 = best_y;
best_x = next_x;
best_y = next_y;
end
trace(p+1)=func(best_x, best_y);
end
%% 5.温度衰减,重新计算差值
t = t * k;
deta = abs(func(best_x,best_y) - func(best_x0,best_y0));
end
disp('最小值在点:');
best_x
best_y
disp( '最小值为:');
func(best_x, best_y)
plot(trace(2:end))
hold on;
xlabel('迭代次数')
ylabel('适应度值')
title('适应度进化曲线')
function results = func(x,y)
results = 5 * cos(x * y)+ x * y + y * y * y;
end
8. 多搜索机制改进的SA
%% Metropolies算法中采取多次搜索策略SA
% x,y的求解范围
x_max = 5;
x_min = -5;
y_max = 5;
y_min = -5;
% 参数初始值
L = 100; % 马尔可夫链长度
k = 0.99; % 衰减因子
s = 0.02; % 步长
t_star = 100; % 初始温度
t_end = 0.001; % 冷却的最终温度
yz = 1e-8; % 容差
p = 0; % 新解被接受的概率(Metropolis算法)
judge_t = 10; %多次搜索策略搜索次数
%% 1.设置初始值
% 初始的前一次状态
x0 = rand * (x_max - x_min) + x_min;
y0 = rand * (y_max - y_min) + y_min;
best_x0 = x0;
best_y0 = y0;
% 初始状态
x0 = rand * (x_max - x_min) + x_min;
y0 = rand * (y_max - y_min) + y_min;
best_x = x0;
best_y = y0;
deta = abs(func(best_x0,best_y0)-func(best_x,best_y));
t = t_star;
while(t > t_end && deta > yz)
for i = 1:L
%% 2.在当前位置产生一个新解
p1 = 0; % 临时变量,用来产生一个范围内的随机值
while p1==0
next_x = x0 + s*(rand*(x_max-x_min)+x_min);
next_y = y0 + s*(rand*(y_max-y_min)+y_min);
% 判断是否在范围内
if next_x >= x_min && next_x <= x_max &&...
next_y >= y_min && next_y <= y_max
p1 = 1;
end
end
%% 3.计算评价函数的增量deta
deta_e = func(next_x,next_y) - func(x0,y0);
%% 4.判断是否接受新解(Metropolis算法)(多次判断取最值)
if deta_e < 0 % 接受新值
x0 = next_x;
y0 = next_y;
p = p + 1;
else
p1 = exp(-1 * deta_e / t);
r = 0;
for j = 1:judge_t
r = max(r,rand);
end
if p1 > r % 接受较差的值
x0 = next_x;
y0 = next_y;
p = p + 1;
end
end
%% 4.更新全局最优解
deta_e = func(next_x,next_y) - func(best_x,best_y);
if deta_e < 0
% 记录上一个最优解,更新当前最优解
best_x0 = best_x;
best_y0 = best_y;
best_x = next_x;
best_y = next_y;
end
trace(p+1)=func(best_x, best_y);
end
%% 5.温度衰减,重新计算差值
t = t * k;
deta = abs(func(best_x,best_y) - func(best_x0,best_y0));
end
disp('最小值在点:');
best_x
best_y
disp( '最小值为:');
func(best_x, best_y)
plot(trace(2:end))
hold on;
xlabel('迭代次数')
ylabel('适应度值')
title('适应度进化曲线')
function results = func(x,y)
results = 5 * cos(x * y)+ x * y + y * y * y;
end