【Matlab】智能优化算法_算数优化算法AOA
- 1.背景介绍
- 2.数学模型
- 2.1 初始化阶段
- 2.2 勘探阶段
- 2.3 开采阶段
- 3.文件结构
- 4.伪代码
- 5.详细代码及注释
- 5.1 AOA.m
- 5.2 func_plot.m
- 5.3 Get_F.m
- 5.4 initialization.m
- 5.5 main.m
- 6.运行结果
- 7.参考文献
1.背景介绍
算术是数论的基本组成部分，与几何、代数和分析一样，也是现代数学的重要组成部分。算术运算符（即乘法、除法、减法和加法）是通常用于研究数字的传统计算方法[31]。我们使用这些简单的算子作为数学优化，从一些候选备选方案（解决方案）中确定符合特定标准的最佳元素。优化问题出现在从工程、经济学、计算机科学到运筹学和工业的所有定量学科中，求解技术的改进吸引了数学界几个时代的兴趣。
所提出的AOA的主要灵感来自于在解决算术问题时使用算术运算符。在下面的小节中，将讨论算术运算符（即乘法、除法、减法和加法）的行为及其在所提出的算法中的影响。下图显示了算术运算符的层次结构及其从外到内的主导地位。然后在数学模型的基础上提出了AOA。
2.数学模型
2.1 初始化阶段
在AOA中，优化过程从矩阵（1）中所示的一组候选解（X）开始，该候选解是随机生成的，并且每次迭代中的最佳候选解被认为是迄今为止获得的最佳解或接近最佳解。
在AOA开始工作之前，它应该选择搜索阶段（即探索或开发）。，因此，数学优化器加速（MOA）函数是由式(2)计算出的系数，用于以下搜索阶段。
其中，MOA（C_Iter）表示第t次迭代时的函数值，该函数值由等式（2）计算。C_Iter表示当前迭代，其在1和最大迭代次数（M_Iter）之间。Min和Max分别表示加速函数的最小值和最大值。
2.2 勘探阶段
本节将介绍AOA的探索行为。根据算术运算符，使用除法（D）运算符或乘法（M）运算符进行的数学计算会得到高分布的值或判定（参考各种规则），这就需要探索搜索机制。然而，与其他运算符（S和A）不同，这些运算符（D和M）由于高度分散而无法轻易接近目标。图3显示了算术算子在数学计算中的影响和行为。图3显示了算术运算在数学运算中的影响和行为，使用一个基于四种数学运算的函数来显示不同运算符分布值的影响。因此，探索搜索可以检测到经过多次努力（迭代）后可能推导出的近似最优解。此外，探索算子（D和M）在这一优化阶段进行操作，以通过加强它们之间的通信来支持搜索过程中的另一阶段（开发）。
AOA的探索算子在搜索区域内随机探索多个区域，并根据两种主要搜索策略（除法搜索策略（D）和乘法搜索策略）找到更好的解决方案，其模型如式（3）所示。搜索阶段（通过执行D或M进行探索搜索，见图4）由数学优化加速（MOA）函数（见式（2））决定，条件为r 1 > MOA（r 1为随机数）。图4显示了所使用的算子如何向最优区域收敛。在此阶段，第一个算子（D）（公式(3)中的第一条规则）以r 2 < 0.5为条件，另一个算子（M）将被忽略，直到该算子完成当前任务。否则，第二个算子（M）将代替D执行当前任务（r 2为随机数）。需要注意的是，我们考虑了一个随机比例系数，以产生更多的多样化过程，并探索搜索空间的不同区域。我们采用了最简单的规则，它能够模拟算术算子的行为。本文为探索部分提出了以下位置更新方程：
其中，xi (C_Iter+1)表示下一次迭代的第 i 个解，xi, j (C_Iter)表示当前迭代第 i 个解的第 j 个位置，best(xj)是目前获得的最佳解中的第 j 个位置。ε是一个小整数，UB j和LB j分别表示第j个位置的上限值和下限值，μ是一个调整搜索过程的控制参数，根据本文的实验，将其固定为0.5。
其中，数学优化器概率（M O P）是一个系数，M O P（C Iter）表示第t次迭代时的函数值，C I ter表示当前迭代，（M I ter）表示最大迭代次数。α是一个敏感参数，定义了迭代过程中的开发精度，根据本文的实验，该精度固定为5。
2.3 开采阶段
在本节中，介绍了AOA的开发策略。根据算术算子，使用减法（S）或加法（A）的数学计算得到了高密度的结果，这参考了利用搜索机制。然而，与其他算子不同，这些算子（S和A）由于其低分散性，可以很容易地接近目标，如图3所示。因此，利用搜索检测到在几次尝试（迭代）之后可能推导出的接近最优的解决方案。此外，开发操作员（S和A）在优化的这一阶段进行了操作，以通过增强他们之间的通信来支持开发阶段。
搜索的这个阶段（通过执行S或A的利用搜索）由条件r 1不大于当前M O A（C Iter）值的MOA函数值来调节（参见等式（2））。在AOA中，AOA的利用算子（减法（S）和加法（A））在几个密集区域上深入探索搜索区域，并基于两种主要的搜索策略（即减法（S）搜索策略和加法（A）搜索策略）寻求更好的解决方案，这两种搜索策略在等式（5）中建模。
该阶段通过深度搜索来利用搜索空间，这一点在图3中非常明显。在这一阶段，第一个算子（S）（公式(5)中的第一条规则）以r 3 <0.5为条件，其他算子（A）将被忽略，直到该算子完成当前任务。否则，第二个算子（A）将代替S执行当前任务。本阶段的这些程序与上一阶段的分区类似。然而，探索搜索运算符（S和A）经常试图避免陷入局部搜索区域。这一过程有助于探索搜索策略找到最优解并保持候选解的多样性。我们精心设计了μ参数，使其在每次迭代时产生一个随机值，从而不仅在第一次迭代时保持探索性，而且在最后一次迭代时也保持探索性。这部分搜索在局部最优停滞的情况下非常有用，尤其是在最后一次迭代中。
图4解释了搜索解决方案如何在二维搜索空间中根据D、M、S和a更新其变量（位置）。可以看出，最终获得的位置可以处于由搜索范围中的D、M、S和a的位置确定的范围内的随机位置。在其他概念中，D、M、S和A估计接近最优解的位置，并且其他解在接近最优解区域周围随机地更新它们的位置。
3.文件结构
AOA.m % 算数优化算法
func_plot.m % 绘图函数
Get_F.m % 获取基准方法
initialization.m % 初始化
main.m % 主函数4.伪代码
5.详细代码及注释
5.1 AOA.m
function [Best_FF,Best_P,Conv_curve]=AOA(N,M_Iter,LB,UB,Dim,F_obj)
display('AOA Working');
%Two variables to keep the positions and the fitness value of the best-obtained solution
Best_P=zeros(1,Dim);
Best_FF=inf;
Conv_curve=zeros(1,M_Iter);
%Initialize the positions of solution
X=initialization(N,Dim,UB,LB);
Xnew=X;
Ffun=zeros(1,size(X,1));% (fitness values)
Ffun_new=zeros(1,size(Xnew,1));% (fitness values)
MOP_Max=1;
MOP_Min=0.2;
C_Iter=1;
Alpha=5;
Mu=0.499;
for i=1:size(X,1)
Ffun(1,i)=F_obj(X(i,:)); %Calculate the fitness values of solutions
if Ffun(1,i)<Best_FF
Best_FF=Ffun(1,i);
Best_P=X(i,:);
end
end
while C_Iter<M_Iter+1 %Main loop
MOP=1-((C_Iter)^(1/Alpha)/(M_Iter)^(1/Alpha)); % Probability Ratio
MOA=MOP_Min+C_Iter*((MOP_Max-MOP_Min)/M_Iter); %Accelerated function
%Update the Position of solutions
for i=1:size(X,1) % if each of the UB and LB has a just value
for j=1:size(X,2)
r1=rand();
if (size(LB,2)==1)
if r1<MOA
r2=rand();
if r2>0.5
Xnew(i,j)=Best_P(1,j)/(MOP+eps)*((UB-LB)*Mu+LB);
else
Xnew(i,j)=Best_P(1,j)*MOP*((UB-LB)*Mu+LB);
end
else
r3=rand();
if r3>0.5
Xnew(i,j)=Best_P(1,j)-MOP*((UB-LB)*Mu+LB);
else
Xnew(i,j)=Best_P(1,j)+MOP*((UB-LB)*Mu+LB);
end
end
end
if (size(LB,2)~=1) % if each of the UB and LB has more than one value
r1=rand();
if r1<MOA
r2=rand();
if r2>0.5
Xnew(i,j)=Best_P(1,j)/(MOP+eps)*((UB(j)-LB(j))*Mu+LB(j));
else
Xnew(i,j)=Best_P(1,j)*MOP*((UB(j)-LB(j))*Mu+LB(j));
end
else
r3=rand();
if r3>0.5
Xnew(i,j)=Best_P(1,j)-MOP*((UB(j)-LB(j))*Mu+LB(j));
else
Xnew(i,j)=Best_P(1,j)+MOP*((UB(j)-LB(j))*Mu+LB(j));
end
end
end
end
Flag_UB=Xnew(i,:)>UB; % check if they exceed (up) the boundaries
Flag_LB=Xnew(i,:)<LB; % check if they exceed (down) the boundaries
Xnew(i,:)=(Xnew(i,:).*(~(Flag_UB+Flag_LB)))+UB.*Flag_UB+LB.*Flag_LB;
Ffun_new(1,i)=F_obj(Xnew(i,:)); % calculate Fitness function
if Ffun_new(1,i)<Ffun(1,i)
X(i,:)=Xnew(i,:);
Ffun(1,i)=Ffun_new(1,i);
end
if Ffun(1,i)<Best_FF
Best_FF=Ffun(1,i);
Best_P=X(i,:);
end
end
%Update the convergence curve
Conv_curve(C_Iter)=Best_FF;
%Print the best solution details after every 50 iterations
if mod(C_Iter,50)==0
display(['At iteration ', num2str(C_Iter), ' the best solution fitness is ', num2str(Best_FF)]);
end
C_Iter=C_Iter+1; % incremental iteration
end5.2 func_plot.m
function func_plot(func_name)
[LB,UB,Dim,F_obj]=Get_F(func_name);
switch func_name
case 'F1'
x=-100:2:100; y=x; %[-100,100]
case 'F2'
x=-100:2:100; y=x; %[-10,10]
case 'F3'
x=-100:2:100; y=x; %[-100,100]
case 'F4'
x=-100:2:100; y=x; %[-100,100]
case 'F5'
x=-200:2:200; y=x; %[-5,5]
case 'F6'
x=-100:2:100; y=x; %[-100,100]
case 'F7'
x=-1:0.03:1; y=x %[-1,1]
case 'F8'
x=-500:10:500;y=x; %[-500,500]
case 'F9'
x=-5:0.1:5; y=x; %[-5,5]
case 'F10'
x=-20:0.5:20; y=x;%[-500,500]
case 'F11'
x=-500:10:500; y=x;%[-0.5,0.5]
case 'F12'
x=-10:0.1:10; y=x;%[-pi,pi]
case 'F13'
x=-5:0.08:5; y=x;%[-3,1]
case 'F14'
x=-100:2:100; y=x;%[-100,100]
case 'F15'
x=-5:0.1:5; y=x;%[-5,5]
case 'F16'
x=-1:0.01:1; y=x;%[-5,5]
case 'F17'
x=-5:0.1:5; y=x;%[-5,5]
case 'F18'
x=-5:0.06:5; y=x;%[-5,5]
case 'F19'
x=-5:0.1:5; y=x;%[-5,5]
case 'F20'
x=-5:0.1:5; y=x;%[-5,5]
case 'F21'
x=-5:0.1:5; y=x;%[-5,5]
case 'F22'
x=-5:0.1:5; y=x;%[-5,5]
case 'F23'
x=-5:0.1:5; y=x;%[-5,5]
end
L=length(x);
f=[];
for i=1:L
for j=1:L
if strcmp(func_name,'F15')==0 && strcmp(func_name,'F19')==0 && strcmp(func_name,'F20')==0 && strcmp(func_name,'F21')==0 && strcmp(func_name,'F22')==0 && strcmp(func_name,'F23')==0
f(i,j)=F_obj([x(i),y(j)]);
end
if strcmp(func_name,'F15')==1
f(i,j)=F_obj([x(i),y(j),0,0]);
end
if strcmp(func_name,'F19')==1
f(i,j)=F_obj([x(i),y(j),0]);
end
if strcmp(func_name,'F20')==1
f(i,j)=F_obj([x(i),y(j),0,0,0,0]);
end
if strcmp(func_name,'F21')==1 || strcmp(func_name,'F22')==1 ||strcmp(func_name,'F23')==1
f(i,j)=F_obj([x(i),y(j),0,0]);
end
end
end
surfc(x,y,f,'LineStyle','none');
end5.3 Get_F.m
function [LB,UB,Dim,F_obj] = Get_F(F)
switch F
case 'F1'
F_obj = @F1;
LB=-100;
UB=100;
Dim=10;
case 'F2'
F_obj = @F2;
LB=-10;
UB=10;
Dim=10;
case 'F3'
F_obj = @F3;
LB=-100;
UB=100;
Dim=10;
case 'F4'
F_obj = @F4;
LB=-100;
UB=100;
Dim=10;
case 'F5'
F_obj = @F5;
LB=-30;
UB=30;
Dim=10;
case 'F6'
F_obj = @F6;
LB=-100;
UB=100;
Dim=10;
case 'F7'
F_obj = @F7;
LB=-1.28;
UB=1.28;
Dim=10;
case 'F8'
F_obj = @F8;
LB=-500;
UB=500;
Dim=10;
case 'F9'
F_obj = @F9;
LB=-5.12;
UB=5.12;
Dim=10;
case 'F10'
F_obj = @F10;
LB=-32;
UB=32;
Dim=10;
case 'F11'
F_obj = @F11;
LB=-600;
UB=600;
Dim=10;
case 'F12'
F_obj = @F12;
LB=-50;
UB=50;
Dim=50;
case 'F13'
F_obj = @F13;
LB=-50;
UB=50;
Dim=10;
case 'F14'
F_obj = @F14;
LB=-65.536;
UB=65.536;
Dim=2;
case 'F15'
F_obj = @F15;
LB=-5;
UB=5;
Dim=4;
case 'F16'
F_obj = @F16;
LB=-5;
UB=5;
Dim=2;
case 'F17'
F_obj = @F17;
LB=[-5,0];
UB=[10,15];
Dim=2;
case 'F18'
F_obj = @F18;
LB=-2;
UB=2;
Dim=2;
case 'F19'
F_obj = @F19;
LB=0;
UB=1;
Dim=3;
case 'F20'
F_obj = @F20;
LB=0;
UB=1;
Dim=6;
case 'F21'
F_obj = @F21;
LB=0;
UB=10;
Dim=4;
case 'F22'
F_obj = @F22;
LB=0;
UB=10;
Dim=4;
case 'F23'
F_obj = @F23;
LB=0;
UB=10;
Dim=4;
end
end
% F1
function o = F1(x)
o=sum(x.^2);
end
% F2
function o = F2(x)
o=sum(abs(x))+prod(abs(x));
end
% F3
function o = F3(x)
dim=size(x,2);
o=0;
for i=1:dim
o=o+sum(x(1:i))^2;
end
end
% F4
function o = F4(x)
o=max(abs(x));
end
% F5
function o = F5(x)
dim=size(x,2);
o=sum(100*(x(2:dim)-(x(1:dim-1).^2)).^2+(x(1:dim-1)-1).^2);
end
% F6
function o = F6(x)
o=sum(abs((x+.5)).^2);
end
% F7
function o = F7(x)
dim=size(x,2);
o=sum([1:dim].*(x.^4))+rand;
end
% F8
function o = F8(x)
o=sum(-x.*sin(sqrt(abs(x))));
end
% F9
function o = F9(x)
dim=size(x,2);
o=sum(x.^2-10*cos(2*pi.*x))+10*dim;
end
% F10
function o = F10(x)
dim=size(x,2);
o=-20*exp(-.2*sqrt(sum(x.^2)/dim))-exp(sum(cos(2*pi.*x))/dim)+20+exp(1);
end
% F11
function o = F11(x)
dim=size(x,2);
o=sum(x.^2)/4000-prod(cos(x./sqrt([1:dim])))+1;
end
% F12
function o = F12(x)
dim=size(x,2);
o=(pi/dim)*(10*((sin(pi*(1+(x(1)+1)/4)))^2)+sum((((x(1:dim-1)+1)./4).^2).*...
(1+10.*((sin(pi.*(1+(x(2:dim)+1)./4)))).^2))+((x(dim)+1)/4)^2)+sum(Ufun(x,10,100,4));
end
% F13
function o = F13(x)
dim=size(x,2);
o=.1*((sin(3*pi*x(1)))^2+sum((x(1:dim-1)-1).^2.*(1+(sin(3.*pi.*x(2:dim))).^2))+...
((x(dim)-1)^2)*(1+(sin(2*pi*x(dim)))^2))+sum(Ufun(x,5,100,4));
end
% F14
function o = F14(x)
aS=[-32 -16 0 16 32 -32 -16 0 16 32 -32 -16 0 16 32 -32 -16 0 16 32 -32 -16 0 16 32;,...
-32 -32 -32 -32 -32 -16 -16 -16 -16 -16 0 0 0 0 0 16 16 16 16 16 32 32 32 32 32];
for j=1:25
bS(j)=sum((x'-aS(:,j)).^6);
end
o=(1/500+sum(1./([1:25]+bS))).^(-1);
end
% F15
function o = F15(x)
aK=[.1957 .1947 .1735 .16 .0844 .0627 .0456 .0342 .0323 .0235 .0246];
bK=[.25 .5 1 2 4 6 8 10 12 14 16];bK=1./bK;
o=sum((aK-((x(1).*(bK.^2+x(2).*bK))./(bK.^2+x(3).*bK+x(4)))).^2);
end
% F16
function o = F16(x)
o=4*(x(1)^2)-2.1*(x(1)^4)+(x(1)^6)/3+x(1)*x(2)-4*(x(2)^2)+4*(x(2)^4);
end
% F17
function o = F17(x)
o=(x(2)-(x(1)^2)*5.1/(4*(pi^2))+5/pi*x(1)-6)^2+10*(1-1/(8*pi))*cos(x(1))+10;
end
% F18
function o = F18(x)
o=(1+(x(1)+x(2)+1)^2*(19-14*x(1)+3*(x(1)^2)-14*x(2)+6*x(1)*x(2)+3*x(2)^2))*...
(30+(2*x(1)-3*x(2))^2*(18-32*x(1)+12*(x(1)^2)+48*x(2)-36*x(1)*x(2)+27*(x(2)^2)));
end
% F19
function o = F19(x)
aH=[3 10 30;.1 10 35;3 10 30;.1 10 35];cH=[1 1.2 3 3.2];
pH=[.3689 .117 .2673;.4699 .4387 .747;.1091 .8732 .5547;.03815 .5743 .8828];
o=0;
for i=1:4
o=o-cH(i)*exp(-(sum(aH(i,:).*((x-pH(i,:)).^2))));
end
end
% F20
function o = F20(x)
aH=[10 3 17 3.5 1.7 8;.05 10 17 .1 8 14;3 3.5 1.7 10 17 8;17 8 .05 10 .1 14];
cH=[1 1.2 3 3.2];
pH=[.1312 .1696 .5569 .0124 .8283 .5886;.2329 .4135 .8307 .3736 .1004 .9991;...
.2348 .1415 .3522 .2883 .3047 .6650;.4047 .8828 .8732 .5743 .1091 .0381];
o=0;
for i=1:4
o=o-cH(i)*exp(-(sum(aH(i,:).*((x-pH(i,:)).^2))));
end
end
% F21
function o = F21(x)
aSH=[4 4 4 4;1 1 1 1;8 8 8 8;6 6 6 6;3 7 3 7;2 9 2 9;5 5 3 3;8 1 8 1;6 2 6 2;7 3.6 7 3.6];
cSH=[.1 .2 .2 .4 .4 .6 .3 .7 .5 .5];
o=0;
for i=1:5
o=o-((x-aSH(i,:))*(x-aSH(i,:))'+cSH(i))^(-1);
end
end
% F22
function o = F22(x)
aSH=[4 4 4 4;1 1 1 1;8 8 8 8;6 6 6 6;3 7 3 7;2 9 2 9;5 5 3 3;8 1 8 1;6 2 6 2;7 3.6 7 3.6];
cSH=[.1 .2 .2 .4 .4 .6 .3 .7 .5 .5];
o=0;
for i=1:7
o=o-((x-aSH(i,:))*(x-aSH(i,:))'+cSH(i))^(-1);
end
end
% F23
function o = F23(x)
aSH=[4 4 4 4;1 1 1 1;8 8 8 8;6 6 6 6;3 7 3 7;2 9 2 9;5 5 3 3;8 1 8 1;6 2 6 2;7 3.6 7 3.6];
cSH=[.1 .2 .2 .4 .4 .6 .3 .7 .5 .5];
o=0;
for i=1:10
o=o-((x-aSH(i,:))*(x-aSH(i,:))'+cSH(i))^(-1);
end
end
function o=Ufun(x,a,k,m)
o=k.*((x-a).^m).*(x>a)+k.*((-x-a).^m).*(x<(-a));
end5.4 initialization.m
function X=initialization(N,Dim,UB,LB)
B_no= size(UB,2); % numnber of boundaries
if B_no==1
X=rand(N,Dim).*(UB-LB)+LB;
end
% If each variable has a different lb and ub
if B_no>1
for i=1:Dim
Ub_i=UB(i);
Lb_i=LB(i);
X(:,i)=rand(N,1).*(Ub_i-Lb_i)+Lb_i;
end
end5.5 main.m
clear all
clc
Solution_no=20; %Number of search solutions
F_name='F1'; %Name of the test function F1-f23
M_Iter=1000; %Maximum number of iterations
[LB,UB,Dim,F_obj]=Get_F(F_name); %Give details of the underlying benchmark function
[Best_FF,Best_P,Conv_curve]=AOA(Solution_no,M_Iter,LB,UB,Dim,F_obj); % Call the AOA
figure('Position',[454 445 694 297]);
subplot(1,2,1);
func_plot(F_name);
title('Parameter space')
xlabel('x_1');
ylabel('x_2');
zlabel([F_name,'( x_1 , x_2 )'])
subplot(1,2,2);
semilogy(Conv_curve,'Color','r','LineWidth',2)
title('Convergence curve')
xlabel('Iteration#');
ylabel('Best fitness function');
axis tight
legend('AOA')
display(['The best-obtained solution by Math Optimizer is : ', num2str(Best_P)]);
display(['The best optimal value of the objective funciton found by Math Optimizer is : ', num2str(Best_FF)]);6.运行结果
7.参考文献
[1]Laith A,Ali D,Seyedali M, et al. The Arithmetic Optimization Algorithm[J]. Computer Methods in Applied Mechanics and Engineering,2021,376.