Matlab 學習案例 群智慧演算法篇
記錄於 2021-10-07。
差分進化演算法（DE）
離散差分進化演算法（DE）
差分進化演算法實現指數擬合（DE-ExpFit）
梯度下降法曲線擬合（GradientDescent）
灰狼演算法（GWO）
灰狼-布穀鳥演算法（GWOCS）
粒子群演算法（PSO）
離散粒子群演算法（PSO）
粒子群演算法實現指數擬合（PSO-ExpFit）
隨機遊走（RandWalk）
模擬退火演算法實現指數擬合（SA）
模擬退火演算法求解旅行商問題（SA-TSP）
樽海鞘演算法（SSA）
狀態轉移演算法（STA）
改進狀態轉移演算法（STA）
鯨魚演算法（WOA）
基礎篇：
Demo 篇：
差分進化演算法（DE）
基本流程：
初始化：
變異：
![v_{i,G+1}=x_{r1,G}+F*(x_{r2,G}-x_{r3,G}),F\in{[0,2]}](https://img.bearask.com/upload/website_attach/fiU511GvNc9_Ij1gK0cs8Hoi84MUO5VLJGN5HYMjb_RwJkhOch7ydmrSAWPgswJb1pVs6Oxb6zh4no2wISUREkYN3UfVswJbwKV_K-xaazh4no2wIS3REkYN3UfVOahvoqIz1RG-O8Nzno2Sow3g1TuT+jxi65l-oqIzvUXWOYUc-euQgk+Rx0Jlb.jpeg)
交叉：
選擇：試驗值與當前目標值比較，實驗值小則更新對應的種群值
邊界處理：超過邊界的可以用範圍內隨機值代替或者進行邊界吸收直接用邊界值代替
改進：
自適應差分進化演算法：
%%%%%%%%%%%%% 初始化 %%%%%%%%%%%%%
clear
；
clc
；
% 種群數量
NP
=
20
；
% 變數維度
D
=
2
；
% 最大進化代數
G
=
200
；
% 初始變異運算元
F0
=
0。5
；
% 交叉運算元
CR
=
0。1
；
% 上限
Xs
=
10
；
% 下限
Xx
=
-
10
；
% 閾值
yz
=
1e-8
；
%%%%%%%%%%%%% 賦初值 %%%%%%%%%%%%%
Ob
=
zeros
（
NP
，
1
）；
Ob1
=
zeros
（
NP
，
1
）；
% 變異種群
v
=
zeros
（
D
，
NP
）；
% 選擇種群
u
=
zeros
（
D
，
NP
）；
% 初始種群
x
=
rand
（
D
，
NP
）
*
（
Xs
-
Xx
）
+
Xx
；
%%%%%%%%%%%%% 計算目標函式 %%%%%%%%%%%%%
for
m
=
1
：
NP
Ob
（
m
）
=
ObjFun
（
x
（：，
m
））；
end
%%%%%%%%%%%%% 差分進化迴圈 %%%%%%%%%%%%%
trace
=
zeros
（
G
，
1
）；
trace
（
1
）
=
min
（
Ob
）；
for
gen
=
1
：
G
%%%%%%%%%%%%% 變異操作 %%%%%%%%%%%%%
%%%%%%%%%%%%% 自適應變異運算元 %%%%%%%%%
lamda
=
exp
（
1
-
G
/
（
G
+
1
-
gen
））；
F
=
F0
*
2
^
lamda
；
%%%%%%%%%% r1、r2、r3和m互不相同 %%%%%%%%%%
for
m
=
1
：
NP
r1
=
randi
（［
1
，
NP
］，
1
）；
while
r1
==
m
r1
=
randi
（［
1
，
NP
］，
1
）；
end
r2
=
randi
（［
1
，
NP
］，
1
）；
while
r2
==
m
||
r2
==
r1
r2
=
randi
（［
1
，
NP
］，
1
）；
end
r3
=
randi
（［
1
，
NP
］，
1
）；
while
r3
==
m
||
r3
==
r2
||
r3
==
r1
r3
=
randi
（［
1
，
NP
］，
1
）；
end
v
（：，
m
）
=
x
（：，
r1
）
+
F
*
（
x
（：，
r2
）
-
x
（：，
r3
））；
end
%%%%%%%%%%%%% 交叉操作 %%%%%%%%%%%%%
r
=
randi
（［
1
，
D
］，
1
）；
for
n
=
1
：
D
cr
=
rand
；
if
cr
<
=
CR
||
n
==
r
u
（
n
，
：）
=
v
（
n
，
：）；
else
u
（
n
，
：）
=
x
（
n
，
：）；
end
end
%%%%%%%%%%%%% 邊界條件的處理 %%%%%%%%%%%%%
for
n
=
1
：
D
for
m
=
1
：
NP
if
u
（
n
，
m
）
<
Xx
||
u
（
n
，
m
）
>
Xs
u
（
n
，
m
）
=
rand
*
（
Xs
-
Xx
）
+
Xx
；
end
end
end
%%%%%%%%%%%%% 邊界吸收 %%%%%%%%%%%%%
% for n = 1 ： D
% for m = 1 ： NP
% if u（n， m） < Xx
% u（n， m） = Xx；
% end
% if u（n， m） > Xs
% u（n， m） = Xs；
% end
% end
% end
%%%%%%%%%%%%% 選擇操作 %%%%%%%%%%%%%
for
m
=
1
：
NP
Ob1
（
m
）
=
ObjFun
（
u
（：，
m
））；
end
for
m
=
1
：
NP
if
Ob1
（
m
）
<
Ob
（
m
）
x
（：，
m
）
=
u
（：，
m
）；
end
end
for
m
=
1
：
NP
Ob
（
m
）
=
ObjFun
（
x
（：，
m
））；
end
trace
（
gen
+
1
）
=
min
（
Ob
）；
if
min
（
Ob
）
<
yz
break
；
end
end
［
SortOb
，
Index
］
=
sort
（
Ob
）；
x
=
x
（：，
Index
）；
% 最優變數
X
=
x
（：，
1
）；
% 最優值
Y
=
min
（
Ob
）；
disp
（［
‘最優x：’
num2str
（
X
‘
）］）；
disp
（［
’最優y：‘
num2str
（
Y
）］）；
subplot
（
1
，
2
，
1
）
plot
（
trace
，
’LineWidth‘
，
2
）；
xlabel
（
’迭代次數‘
）；
ylabel
（
’適應度值‘
）；
title
（
’適應度進化曲線‘
）；
len
=
50
；
xRange
=
linspace
（
Xx
，
Xs
，
len
）；
yRange
=
linspace
（
Xx
，
Xs
，
len
）；
［
xMap
，
yMap
］
=
meshgrid
（
xRange
，
yRange
）；
zMap
=
zeros
（
len
）；
for
i
=
1
：
len
for
j
=
1
：
len
zMap
（
i
，
j
）
=
ObjFun
（［
xMap
（
i
，
j
），
yMap
（
i
，
j
）］）；
end
end
subplot
（
1
，
2
，
2
）
surf
（
xRange
，
yRange
，
zMap
）；
colorbar
；
view
（
-
45
，
-
45
）；
shading
interp
hold
on
plot3
（
X
（
1
），
X
（
2
），
Y
，
’o‘
，
’MarkerFaceColor‘
，
’r‘
，
’MarkerSize‘
，
10
）；
hold
off
set
（
gcf
，
’Position‘
，
［
400
，
200
，
900
，
400
］）；
function
result
=
ObjFun
（
x
）
x1
=
x
（
1
）；
x2
=
x
（
2
）；
t
=
x1
^
2
+
x2
^
2
；
result
=
0。5
+
（
sin
（
sqrt
（
t
））
^
2
-
0。5
）
/
（
1
+
0。001
*
t
）
^
2
；
end
離散差分進化演算法（DE）
針對整數變數，用floor函式向下取整：
%%%%%%%%%%%%% 離散差分進化演算法求函式極值 %%%%%%%%%%%%%
%%%%%%%%%%%%% 初始化 %%%%%%%%%%%%%
clear
；
clc
；
% 種群數量
NP
=
20
；
% 變數維度
D
=
2
；
% 最大進化代數
G
=
200
；
% 初始編譯運算元
F
=
0。5
；
% 交叉運算元
CR
=
0。1
；
% 上限
Xs
=
100
；
% 下限
Xx
=
-
100
；
%%%%%%%%%%%%% 賦初值 %%%%%%%%%%%%%
Ob
=
zeros
（
NP
，
1
）；
Ob1
=
zeros
（
NP
，
1
）；
trace
=
zeros
（
G
+
1
，
1
）；
% 變異種群
v
=
zeros
（
D
，
NP
）；
% 選擇種群
u
=
zeros
（
D
，
NP
）；
% 初始種群
x
=
randi
（［
Xx
，
Xs
］，
D
，
NP
）；
%%%%%%%%%%%%% 計算目標函式 %%%%%%%%%%%%%
for
m
=
1
：
NP
Ob
（
m
）
=
ObjFunInt
（
x
（：，
m
））；
end
trace
（
1
）
=
min
（
Ob
）；
%%%%%%%%%%%%% 差分進化迴圈 %%%%%%%%%%%%%
for
gen
=
1
：
G
%%%%%%%%%%%%% 變異操作 %%%%%%%%%%%%%
%%%%%%%%%% r1、r2、r3和m互不相同 %%%%%%%%%%
for
m
=
1
：
NP
r1
=
randi
（［
1
，
NP
］，
1
）；
while
r1
==
m
r1
=
randi
（［
1
，
NP
］，
1
）；
end
r2
=
randi
（［
1
，
NP
］，
1
）；
while
r2
==
m
||
r2
==
r1
r2
=
randi
（［
1
，
NP
］，
1
）；
end
r3
=
randi
（［
1
，
NP
］，
1
）；
while
r3
==
m
||
r3
==
r2
||
r3
==
r1
r3
=
randi
（［
1
，
NP
］，
1
）；
end
v
（：，
m
）
=
floor
（
x
（：，
r1
）
+
F
*
（
x
（：，
r2
）
-
x
（：，
r3
）））；
end
%%%%%%%%%%%%% 交叉操作 %%%%%%%%%%%%%
r
=
randi
（［
1
，
D
］，
1
）；
for
n
=
1
：
D
cr
=
rand
；
if
cr
<
=
CR
||
n
==
r
u
（
n
，
：）
=
v
（
n
，
：）；
else
u
（
n
，
：）
=
x
（
n
，
：）；
end
end
%%%%%%%%%%%%% 邊界條件的處理 %%%%%%%%%%%%%
%%%%%%%%%%%%% 邊界吸收 %%%%%%%%%%%%%
for
n
=
1
：
D
for
m
=
1
：
NP
if
u
（
n
，
m
）
<
Xx
u
（
n
，
m
）
=
Xx
；
end
if
u
（
n
，
m
）
>
Xs
u
（
n
，
m
）
=
Xs
；
end
end
end
%%%%%%%%%%%%% 選擇操作 %%%%%%%%%%%%%
for
m
=
1
：
NP
Ob1
（
m
）
=
ObjFunInt
（
u
（：，
m
））；
end
for
m
=
1
：
NP
if
Ob1
（
m
）
>
Ob
（
m
）
x
（：，
m
）
=
u
（：，
m
）；
end
end
for
m
=
1
：
NP
Ob
（
m
）
=
ObjFunInt
（
x
（：，
m
））；
end
trace
（
gen
+
1
）
=
min
（
Ob
）；
end
［
SortOb
，
Index
］
=
sort
（
Ob
）；
x
=
x
（：，
Index
）；
% 最優變數
X
=
x
（：，
1
）；
% 最優值
Y
=
min
（
Ob
）；
disp
（［
’最優x：‘
num2str
（
X
’
）］）；
disp
（［
‘最優y：’
num2str
（
Y
）］）；
subplot
（
1
，
2
，
1
）
plot
（
trace
，
‘LineWidth’
，
2
）；
xlabel
（
‘迭代次數’
）；
ylabel
（
‘適應度值’
）；
title
（
‘適應度進化曲線’
）；
len
=
200
；
xRange
=
linspace
（
Xx
，
Xs
，
len
）；
yRange
=
linspace
（
Xx
，
Xs
，
len
）；
［
xMap
，
yMap
］
=
meshgrid
（
xRange
，
yRange
）；
zMap
=
zeros
（
len
）；
for
i
=
1
：
len
for
j
=
1
：
len
zMap
（
i
，
j
）
=
ObjFunInt
（［
xMap
（
i
，
j
），
yMap
（
i
，
j
）］）；
end
end
subplot
（
1
，
2
，
2
）
surf
（
xRange
，
yRange
，
zMap
）；
colorbar
；
view
（
-
45
，
45
）；
shading
interp
hold
on
plot3
（
X
（
1
），
X
（
2
），
Y
，
‘o’
，
‘MarkerFaceColor’
，
‘r’
，
‘MarkerSize’
，
10
）；
hold
off
set
（
gcf
，
‘Position’
，
［
400
，
200
，
900
，
400
］）；
function
result
=
ObjFunInt
（
x
）
result
=
-
（（
x
（
1
）
。^
2
+
x
（
2
）
-
1
）
。^
2
+
（
x
（
1
）
+
x
（
2
）
。^
2
-
7
）
。^
2
）
/
200
+
10
；
end
差分進化演算法實現指數擬合（DE-ExpFit）
clear
；
clc
；
t
=
0。2
*
（
1
：
3000
）
‘
；
data
=
400
*
exp
（
-
t
/
5
）
+
10
*
randn
（
size
（
t
））；
tic
；
p
=
DE_ExpFit
（
t
，
data
）
toc
；
fit
=
p
（
1
）
*
exp
（
-
t
/
p
（
2
））
+
p
（
3
）；
plot
（
t
，
data
，
t
，
fit
，
’LineWidth‘
，
2
）
legend
（
’model‘
，
’DE‘
）；
title
（
’差分進化演算法實現指數擬合‘
）；
function
［x_opt， y_opt］
=
DE_ExpFit
（
t， Et
）
%{
函式功能：差分進化演算法實現指數擬合：y = a*exp（-x/b） + c
輸入：
t：自變數；
Et：因變數；
輸出：
x_opt：最優解；
y_opt：適應度（目標函式值）；
%}
% 初始值
NP
=
50
；
% 種群數量
D
=
3
；
% 變數維度
G
=
200
；
% 最大進化代數
F0
=
0。5
；
% 初始變異運算元
CR
=
0。9
；
% 交叉運算元
max_sig
=
max
（
Et
）；
Xmin
=
［
0。5
*
max_sig
，
-
2
，
0
］；
% 下限
Xmax
=
［
1。5
*
max_sig
，
4
，
0。5
*
max_sig
］；
% 上限；
%%%%%%%%%%%%% 賦初值 %%%%%%%%%%%%%
Ob
=
zeros
（
NP
，
1
）；
Ob1
=
zeros
（
NP
，
1
）；
v
=
zeros
（
NP
，
D
）；
% 變異種群
u
=
zeros
（
NP
，
D
）；
% 選擇種群
x
=
rand
（
NP
，
D
）
。*
repmat
（
Xmax
-
Xmin
，
NP
，
1
）
+
repmat
（
Xmin
，
NP
，
1
）；
%%%%%%%%%%%%% 計算目標函式 %%%%%%%%%%%%%
for
m
=
1
：
NP
Ob
（
m
）
=
Obj_Fit
（
x
（
m
，
：），
t
，
Et
）；
end
y_opt
=
zeros
（
1
，
G
）；
%%%%%%%%%%%%% 差分進化迴圈 %%%%%%%%%%%%%
for
gen
=
1
：
G
%%%%%%%%%%%%% 變異操作 %%%%%%%%%%%%%
%%%%%%%%%%%%% 自適應變異運算元 %%%%%%%%%
lamda
=
exp
（
1
-
G
/
（
G
+
1
-
gen
））；
F
=
F0
*
2
^
lamda
；
%%%%%%%%%% r1、r2、r3和m互不相同 %%%%%%%%%%
for
m
=
1
：
NP
r1
=
randi
（［
1
，
NP
］，
1
）；
while
r1
==
m
r1
=
randi
（［
1
，
NP
］，
1
）；
end
r2
=
randi
（［
1
，
NP
］，
1
）；
while
r2
==
m
||
r2
==
r1
r2
=
randi
（［
1
，
NP
］，
1
）；
end
r3
=
randi
（［
1
，
NP
］，
1
）；
while
r3
==
m
||
r3
==
r2
||
r3
==
r1
r3
=
randi
（［
1
，
NP
］，
1
）；
end
v
（
m
，
：）
=
x
（
r1
，
：）
+
F
*
（
x
（
r2
，
：）
-
x
（
r3
，
：））；
end
%%%%%%%%%%%%% 交叉操作 %%%%%%%%%%%%%
r
=
randi
（［
1
，
D
］，
1
）；
for
n
=
1
：
D
cr
=
rand
；
if
cr
<
=
CR
||
n
==
r
u
（：，
n
）
=
v
（：，
n
）；
else
u
（：，
n
）
=
x
（：，
n
）；
end
end
%%%%%%%%%%%%% 邊界條件的處理 %%%%%%%%%%%%%
% ①範圍內隨機值
for
n
=
1
：
NP
for
m
=
1
：
D
if
u
（
n
，
m
）
<
Xmin
（
m
）
||
u
（
n
，
m
）
>
Xmax
（
m
）
u
（
n
，
m
）
=
rand
*
（
Xmax
（
m
）
-
Xmin
（
m
））
+
Xmin
（
m
）；
end
end
end
% ②邊界吸收
% for n = 1 ： NP
% for m = 1 ： D
% if u（n， m） < Xmin（m）
% u（n， m） = Xmin（m）；
% end
% if u（n， m） > Xmax（m）
% u（n， m） = Xmax（m）；
% end
% end
% end
%%%%%%%%%%%%% 選擇操作 %%%%%%%%%%%%%
for
m
=
1
：
NP
Ob1
（
m
）
=
Obj_Fit
（
u
（
m
，
：），
t
，
Et
）；
end
for
m
=
1
：
NP
if
Ob1
（
m
）
<
Ob
（
m
）
x
（
m
，
：）
=
u
（
m
，
：）；
end
end
for
m
=
1
：
NP
Ob
（
m
）
=
Obj_Fit
（
x
（
m
，
：），
t
，
Et
）；
end
y_opt
（
gen
）
=
min
（
Ob
）；
end
［
~
，
Index
］
=
min
（
Ob
）；
x_opt
=
x
（
Index
，
：）；
x_opt
（
2
）
=
10
^
（
x_opt
（
2
））；
end
% 目標函式
function
fitError
=
Obj_Fit
（
p0， t ，Et
）
A
=
p0
（
1
）；
B
=
p0
（
2
）；
C
=
p0
（
3
）；
f
=
A
*
exp
（
-
t
/
10
^
B
）
+
C
；
fitError
=
norm
（
Et
-
f
）；
end
梯度下降法曲線擬合（GradientDescent）
clear
；
clc
；
t
=
［
0。1350
，
0。2600
，
0。3800
，
0。5050
，
0。6450
，
0。8100
，
1。0000
，
1。2150
，
1。4600
，
1。7400
，
2。0550
，
2。4150
，
。。。
2。8150
，
3。2650
，
3。7650
，
4。3200
，
4。9350
，
5。6100
，
6。3450
，
7。1500
，
8。0200
，
8。9650
，
9。9800
，
11。0650
，
。。。
12。2200
，
13。4500
，
14。7500
，
16。1150
，
17。5450
，
19。0400
，
20。5950
，
22。1950
，
23。8450
，
25。5350
，
。。。
27。2600
，
29。2650
，
32。2700
，
36。2800
，
41。4850
，
47。8900
，
55。8950
，
65。5000
，
77。5050
，
92。3100
，
。。。
108。7150
，
128。3250
，
151。9300
，
179。9350
，
212。9400
，
252。9450
］
’
；
M
=
length
（
t
）；
% 構造模型
A
=
1800
；
B
=
15
；
qMrl
=
5
；
AA
=
50
；
C
=
10
；
rng
（
‘default’
）；
ratio
=
0
；
noise
=
ratio
*
randn
（
M
，
1
）；
Et
=
A
*
exp
（
-
t
/
B
-
0。5
*
qMrl
*
t
。^
2
）
+
AA
*
exp
（
-
t
/
B
）
+
C
+
noise
；
% 梯度下降法學習
% 初始引數
alpha
=
1
；
% 學習率大收斂快，可能有震盪
iteMax
=
1000000
；
iniA
=
1
；
iniB
=
rand
；
iniQMrl
=
rand
；
iniAA
=
rand
；
iniC
=
rand
；
initW
=
［
iniA
；
iniB
；
iniQMrl
；
iniAA
；
iniC
］；
err
=
1e-6
；
J
=
zeros
（
iteMax
，
1
）；
G
=
zeros
（
size
（
initW
））；
e
=
0。1
；
convIter
=
0
；
for
i
=
1
：
iteMax
gradA
=
2
*
iniA
*
exp
（
-
0。5
*
iniQMrl
^
2
*
t
。^
2
-
t
/
iniB
^
2
）；
gradB
=
（
2
*
iniA
^
2
*
t
。*
exp
（
-
0。5
*
iniQMrl
^
2
*
t
。^
2
-
t
/
iniB
^
2
））
/
iniB
^
3
+
（
2
*
iniAA
^
2
*
t
。*
exp
（
-
t
/
iniB
^
2
））
/
iniB
^
3
；
gradQMrl
=
-
iniA
^
2
*
iniQMrl
*
t
。^
2
。*
exp
（
-
0。5
*
iniQMrl
^
2
*
t
。^
2
-
t
/
iniB
^
2
）；
gradAA
=
2
*
iniAA
*
exp
（
-
t
/
iniB
^
2
）；
gradC
=
2
*
iniC
*
ones
（
M
，
1
）；
grad
=
1
/
M
*
［
gradA
，
gradB
，
gradQMrl
，
gradAA
，
gradC
］
‘
*
（
iniA
^
2
*
exp
（
-
t
/
iniB
^
2
-
0。5
*
iniQMrl
^
2
*
t
。^
2
）
+
iniAA
^
2
*
exp
（
-
t
/
iniB
^
2
）
+
iniC
^
2
-
Et
）；
% Adagrad 法 x = x + yita * inv（G） * grad；
G
=
G
+
grad
。^
2
；
newW
=
initW
-
alpha
*
diag
（
1
。/
sqrt
（
G
+
e
））
*
grad
；
% newW = initW - alpha * grad；
if
max
（
abs
（（
newW
-
initW
）））
<
err
convIter
=
convIter
+
1
；
if
convIter
==
10
J
（
i
+
1
：
end
）
=
［］；
disp
（［
’梯度下降法迭代次數：‘
，
num2str
（
i
）］）；
result
=
newW
’
；
result
（
1
）
=
newW
（
1
）
^
2
；
result
（
2
）
=
newW
（
2
）
^
2
；
result
（
3
）
=
newW
（
3
）
^
2
；
result
（
4
）
=
newW
（
4
）
^
2
；
result
（
5
）
=
newW
（
5
）
^
2
；
disp
（［
‘梯度下降法迭代結果：’
，
num2str
（
result
）］）；
break
；
end
else
convIter
=
0
；
initW
=
newW
；
iniA
=
newW
（
1
）；
iniB
=
newW
（
2
）；
iniQMrl
=
newW
（
3
）；
iniAA
=
newW
（
4
）；
iniC
=
newW
（
5
）；
J
（
i
）
=
1
/
（
2
*
M
）
*
norm
（
iniA
^
2
*
exp
（
-
t
/
iniB
^
2
-
0。5
*
iniQMrl
^
2
*
t
。^
2
）
+
iniAA
^
2
*
exp
（
-
t
/
iniB
^
2
）
+
iniC
^
2
-
Et
）
^
2
；
end
end
fit
=
iniA
^
2
*
exp
（
-
t
/
iniB
^
2
-
0。5
*
iniQMrl
^
2
*
t
。^
2
）
+
iniAA
^
2
*
exp
（
-
t
/
iniB
^
2
）
+
iniC
^
2
；
% 繪圖
figure
（
‘Position’
，
［
50
，
50
，
900
，
400
］）；
subplot
（
1
，
2
，
1
）
loglog
（
J
，
‘LineWidth’
，
2
）
legend
（［
‘alpha = ’
，
num2str
（
alpha
）］）；
subplot
（
1
，
2
，
2
）
semilogx
（
t
，
Et
，
‘ok’
，
‘MarkerFaceColor’
，
‘r’
，
‘MarkerSIze’
，
5
，
‘LineWidth’
，
1。5
）
hold
on
plot
（
t
，
fit
，
‘LineWidth’
，
1。5
）
hold
off
legend
（
‘Model’
，
‘Fitting’
）
title
（
‘梯度下降法曲線擬合’
）；
灰狼演算法（GWO）
灰狼-布穀鳥演算法（GWOCS）
粒子群演算法（PSO）
離散粒子群演算法（PSO）
粒子狀態更新
clear
；
clc
；
% 群體例子個數
N
=
100
；
% 粒子維度
D
=
20
；
T
=
200
；
% 最大迭代次數
c1
=
1。5
；
% 學習因子1
c2
=
1。5
；
% 學習因子2
wMax
=
0。8
；
% 慣性權重最大值
wMin
=
0。4
；
% 慣性權重最小值
xMax
=
10
；
% 位置最大值
xMin
=
-
10
；
% 位置最小值
vMax
=
10
；
% 速度最大值
vMin
=
-
10
；
% 速度最小值
%%%%%%%%%% 初始化種群個體（限定位置和速度） %%%%%%%%%
% 隨機獲得二進位制編碼的初始種群
x
=
randi
（［
0
，
1
］，
N
，
D
）；
v
=
rand
（
N
，
D
）
*
（
vMax
-
vMin
）
+
vMin
；
%%%%%%%%%% 初始化個體最優位置和最優值 %%%%%%%%%
p
=
x
；
pbest
=
ones
（
N
，
1
）；
for
i
=
1
：
N
pbest
（
i
）
=
ObjFunDiscrete
（
x
（
i
，
：），
xMax
，
xMin
）；
end
%%%%%%%%%% 初始化全域性最優位置和最優值 %%%%%%%%%
g
=
ones
（
1
，
D
）；
gbest
=
inf
；
for
i
=
1
：
N
if
pbest
（
i
）
<
gbest
g
=
p
（
i
，
：）；
gbest
=
pbest
（
i
）；
end
end
gb
=
ones
（
1
，
T
）；
vx
=
zeros
（
N
，
D
）；
%%%%%%%%%% 按照公式依次迭代直到滿足精度或者迭代次數 %%%%%%%%%
for
i
=
1
：
T
for
j
=
1
：
N
%%%%%%%%%% 更新個體最優位置和最優值 %%%%%%%%%
if
ObjFunDiscrete
（
x
（
j
，
：），
xMax
，
xMin
）
<
pbest
（
j
）
p
（
j
，
：）
=
x
（
j
，
：）；
pbest
（
j
）
=
ObjFunDiscrete
（
x
（
j
，
：），
xMax
，
xMin
）；
end
%%%%%%%%%% 更新全域性最優位置和最優值 %%%%%%%%%
if
pbest
（
j
）
<
gbest
g
=
p
（
j
，
：）；
gbest
=
pbest
（
j
）；
end
%%%%%%%%%% 動態計算慣性權重值 %%%%%%%%%
w
=
wMax
-
（
wMax
-
wMin
）
*
i
/
T
；
%%%%%%%%%% 更新位置和速度值 %%%%%%%%%
v
（
j
，
：）
=
w
*
v
（
j
，
：）
+
c1
*
rand
*
（
p
（
j
，
：）
-
x
（
j
，
：））
+
c2
*
rand
*
（
g
-
x
（
j
，
：））；
x
（
j
，
：）
=
x
（
j
，
：）
+
v
（
j
，
：）；
%%%%%%%%%% 邊界條件處理 %%%%%%%%%
for
ii
=
1
：
D
if
（
v
（
j
，
ii
）
>
vMax
）
||
（
v
（
j
，
ii
）
<
vMin
）
v
（
j
，
ii
）
=
rand
*
（
vMax
-
vMin
）
+
vMin
；
end
end
vx
（
j
，
：）
=
1
。/
（
1
+
exp
（
-
v
（
j
，
：）））；
for
jj
=
1
：
D
if
vx
（
j
，
jj
）
>
rand
x
（
j
，
jj
）
=
1
；
else
x
（
j
，
jj
）
=
0
；
end
end
end
%%%%%%%%%% 記錄歷代全域性最優值 %%%%%%%%%
gb
（
i
）
=
gbest
；
end
m
=
0
；
for
j
=
1
：
D
m
=
m
+
g
（
j
）
*
2
^
（
j
-
1
）；
end
f1
=
xMin
+
m
*
（
xMax
-
xMin
）
/
（
2
^
D
-
1
）；
disp
（［
‘最優個體：’
num2str
（
f1
）］）；
disp
（［
‘最優值：’
num2str
（
gb
（
end
））］）；
subplot
（
1
，
2
，
1
）；
plot
（
gb
，
‘LineWidth’
，
2
）；
xlabel
（
‘迭代次數’
）；
ylabel
（
‘適應度值’
）；
title
（
‘適應度進化曲線’
）；
subplot
（
1
，
2
，
2
）
xRange
=
linspace
（
xMin
，
xMax
，
100
）；
yRange
=
xRange
+
6
*
sin
（
4
*
xRange
）
+
9
*
cos
（
5
*
xRange
）；
plot
（
xRange
，
yRange
，
‘LineWidth’
，
2
）；
hold
on
；
plot
（
f1
，
gb
（
end
），
‘o’
，
‘MarkerFaceColor’
，
‘r’
，
‘MarkerSize’
，
10
）；
hold
off
；
set
（
gcf
，
‘Position’
，
［
400
，
200
，
900
，
400
］）；
function
result
=
ObjFunDiscrete
（
x， xMax， xMin
）
%%%%%%%%%% 適應度函式 %%%%%%%%%%%
m
=
0
；
D
=
length
（
x
）；
for
j
=
1
：
D
m
=
m
+
x
（
j
）
*
2
^
（
j
-
1
）；
end
% 譯碼成十進位制數
f
=
xMin
+
m
*
（
xMax
-
xMin
）
/
（
2
^
D
-
1
）；
fit
=
f
+
6
*
sin
（
4
*
f
）
+
9
*
cos
（
5
*
f
）；
result
=
fit
；
end
粒子群演算法實現指數擬合（PSO-ExpFit）
clear
；
clc
；
t
=
0。2
*
（
1
：
3000
）
‘
；
data
=
400
*
exp
（
-
t
/
5
）
+
10
*
randn
（
size
（
t
））；
tic
；
p
=
PSO_ExpFit
（
t
，
data
）
toc
；
fit
=
p
（
1
）
*
exp
（
-
t
/
p
（
2
））
+
p
（
3
）；
plot
（
t
，
data
，
t
，
fit
，
’LineWidth‘
，
2
）
legend
（
’model‘
，
’PSO‘
）；
title
（
’壓縮因子粒子群演算法實現指數擬合‘
）；
function
［x_opt， y_opt］
=
PSO_ExpFit
（
t， Et
）
%{
函式功能：壓縮因子粒子群演算法實現指數擬合：y = a*exp（-x/b） + c
輸入：
t：自變數；
Et：因變數；
輸出：
x_opt：最優解；
y_opt：適應度（目標函式值）；
%}
% 初始值
D
=
3
；
% 粒子維度，未知數個數；
IterMax
=
200
；
% 最大迭代次數；
Vmin
=
-
1
；
% 速度最小值；
Vmax
=
1
；
% 速度最大值；
c1
=
2。05
；
% 學習因子1，認知部分；
c2
=
2。05
；
% 學習因子2，社會部分；
N
=
200
；
% 種群個數；
max_sig
=
max
（
Et
）；
Xmin
=
［
0。5
*
max_sig
，
-
2
，
0
］；
% 下限
Xmax
=
［
1。5
*
max_sig
，
4
，
0。5
*
max_sig
］；
% 上限；
% 計算
phy
=
c1
+
c2
；
lamda
=
2
/
abs
（
2
-
phy
-
sqrt
（
phy
^
2
-
4
*
phy
））；
% 壓縮因子
% 初始化種群個體（限定位置和速度）
x
=
rand
（
N
，
D
）
。*
repmat
（
Xmax
-
Xmin
，
N
，
1
）
+
repmat
（
Xmin
，
N
，
1
）；
v
=
rand
（
N
，
D
）
*
（
Vmax
-
Vmin
）
+
Vmin
；
% 初始化個體最優位置和最優值
p
=
x
；
pbest
=
ones
（
1
，
N
）；
for
i
=
1
：
N
pbest
（
i
）
=
Obj_Fit
（
x
（
i
，
：），
t
，
Et
）；
end
% 初始化全域性最優位置和最優值
g
=
ones
（
1
，
D
）；
gbest
=
inf
；
for
i
=
1
：
N
if
pbest
（
i
）
<
gbest
g
=
p
（
i
，
：）；
gbest
=
pbest
（
i
）；
end
end
% 迭代直到滿足精度或者迭代次數
fx
=
zeros
（
1
，
N
）；
y_opt
=
zeros
（
1
，
IterMax
）；
for
i
=
1
：
IterMax
for
j
=
1
：
N
% 更新個體最優位置和最優值
fx
（
j
）
=
Obj_Fit
（
x
（
j
，
：），
t
，
Et
）；
if
fx
（
j
）
<
pbest
（
j
）
p
（
j
，
：）
=
x
（
j
，
：）；
pbest
（
j
）
=
fx
（
j
）；
end
% 更新全域性最優位置和最優值
if
pbest
（
j
）
<
gbest
g
=
p
（
j
，
：）；
gbest
=
pbest
（
j
）；
end
% 更新位置和速度值
v
（
j
，
：）
=
lamda
*
v
（
j
，
：）
+
c1
*
rand
*
（
p
（
j
，
：）
-
x
（
j
，
：））
+
c2
*
rand
*
（
g
-
x
（
j
，
：））；
x
（
j
，
：）
=
x
（
j
，
：）
+
v
（
j
，
：）；
% 邊界條件處理
for
ii
=
1
：
D
if
（
v
（
j
，
ii
）
>
Vmax
）
||
（
v
（
j
，
ii
）
<
Vmin
）
v
（
j
，
ii
）
=
rand
*
（
Vmax
-
Vmin
）
+
Vmin
；
end
if
（
x
（
j
，
ii
）
>
Xmax
（
ii
））
||
（
x
（
j
，
ii
）
<
Xmin
（
ii
））
x
（
j
，
ii
）
=
rand
*
（
Xmax
（
ii
）
-
Xmin
（
ii
））
+
Xmin
（
ii
）；
end
end
end
end
x_opt
=
g
；
x_opt
（
2
）
=
10
^
（
x_opt
（
2
））；
end
% 目標函式
function
fitError
=
Obj_Fit
（
p0， t ，Et
）
A
=
p0
（
1
）；
B
=
p0
（
2
）；
C
=
p0
（
3
）；
f
=
A
*
exp
（
-
t
/
10
^
B
）
+
C
；
fitError
=
norm
（
Et
-
f
）；
end
隨機遊走（RandWalk）
clear
；
clc
；
lb
=
［
-
10
，
-
10
］；
ub
=
［
10
，
10
］；
x0
=
rand
（
1
，
2
）
。*
（
ub
-
lb
）
+
lb
；
rng
（
’default‘
）；
［
mx
，
minf
］
=
OptRandWalk
（
x0
，
10
，
100
，
10
）；
disp
（［
’最優值：‘
num2str
（
mx
）］）；
disp
（［
’最優值對應目標函式：‘
num2str
（
minf
）］）；
len
=
100
；
xRange
=
linspace
（
-
10
，
10
，
len
）；
yRange
=
linspace
（
-
10
，
10
，
len
）；
［
xMap
，
yMap
］
=
meshgrid
（
xRange
，
yRange
）；
zMap
=
zeros
（
len
）；
for
i
=
1
：
len
for
j
=
1
：
len
zMap
（
i
，
j
）
=
func
（［
xMap
（
i
，
j
），
yMap
（
i
，
j
）］）；
end
end
surf
（
xRange
，
yRange
，
zMap
）；
view
（
-
45
，
-
45
）；
shading
interp
hold
on
plot3
（
mx
（
1
），
mx
（
2
），
minf
，
’o‘
，
’MarkerFaceColor‘
，
’r‘
，
’MarkerSize‘
，
10
）；
hold
off
title
（
’隨機行走法求函式的極小值‘
）；
function
［mx， minf］
=
OptRandWalk
（
x， lamda， N， n
）
%{
函式功能：隨機行走法求函式的極小值
x：初始值；
lamda：步長；
N：為了產生較好點的迭代次數；
n：單步迴圈行走次數，目的是儘可能走到全域性最優點附近
mx：最優解；
minf：最優值。
%}
F
=
zeros
（
n
，
1
）；
D
=
length
（
x
）；
X
=
zeros
（
n
，
D
）；
epsilon
=
1e-6
；
f1
=
func
（
x
）；
while
lamda
>
=
epsilon
k
=
1
；
while
（
k
<
=
N
）
u
=
2
*
rand
（
n
，
D
）
-
1
；
for
ii
=
1
：
n
X
（
ii
，
：）
=
x
+
lamda
*
（
u
（
ii
，
：）
/
norm
（
u
（
ii
，
：）））；
F
（
ii
）
=
func
（
X
（
ii
，
：））；
end
［
f11
，
kk
］
=
min
（
F
）；
if
f11
<
f1
f1
=
f11
；
x
=
X
（
kk
，
：）；
k
=
1
；
else
k
=
k
+
1
；
end
end
lamda
=
lamda
/
2
；
end
mx
=
X
（
kk
，
：）；
minf
=
f1
；
end
function
f
=
func
（
x
）
% f = -sin（sqrt（（x（1） - 100）。^2 + （x（2） - 100）。^2 ） + exp（1）） 。/ （sqrt（（x（1） - 100）。^2 + （x（2） - 100）。^2 ） + exp（1）） - 1；
% f = sum（x。^2）；
% f = 3*cos（x（1）*x（2）） + x（1） + x（2）^2；
% f = 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；
% f = 0。5 + （sin（sqrt（x（1）^2 + x（2）^2））^2 - 0。5） / （1 + 0。001*（x（1）^2 + x（2）^2））^2；
x1
=
x
（
1
）；
x2
=
x
（
2
）；
t
=
x1
^
2
+
x2
^
2
；
f
=
0。5
+
（
sin
（
sqrt
（
t
））
^
2
-
0。5
）
/
（
1
+
0。001
*
t
）
^
2
；
end
模擬退火演算法實現指數擬合（SA）
clear
；
clc
；
t
=
0。2
*
（
1
：
3000
）
’
；
data
=
400
*
exp
（
-
t
/
5
）
+
10
*
randn
（
size
（
t
））；
tic
；
p
=
SAA_ExpFit
（
t
，
data
）
toc
；
fit
=
p
（
1
）
*
exp
（
-
t
/
p
（
2
））
+
p
（
3
）；
plot
（
t
，
data
，
t
，
fit
，
‘LineWidth’
，
2
）
legend
（
‘model’
，
‘SA’
）；
title
（
‘模擬退火演算法實現指數擬合’
）；
function
p0
=
SAA_ExpFit
（
t， Et
）
%{
函式功能：模擬退火演算法實現指數擬合：y = a*exp（-x/b） + c
輸入：
t：自變數；
Et：因變數；
輸出：
x_opt：最優解；
y_opt：適應度（目標函式值）；
%}
% initial condition
T0
=
1000
；
% initial temperature
alpha
=
0。99
；
% annealing factor
N_max
=
3000
；
% max iterations
% range of varibles
max_sig
=
max
（
Et
）；
Xmin
=
［
0。5
*
max_sig
，
-
2
，
0
］；
% 下限
Xmax
=
［
1。5
*
max_sig
，
4
，
0。5
*
max_sig
］；
% 上限；
% initial value
var
=
rand
（
1
，
3
）
。*
（
Xmax
-
Xmin
）
+
Xmin
；
fun
=
Obj_Fit
（
var
，
t
，
Et
）；
% iterative
for
i
=
1
：
N_max
Tk
=
T0
*
alpha
^
（
i
-
1
）；
dvar
=
Tk
*
sign
（
rand
（
1
，
3
）
-
0。5
）
。*
（（
1
+
1
/
Tk
）
。^
abs
（
2
*
rand
（
1
，
3
）
-
1
）
-
1
）；
var1
=
var
+
dvar
。*
（
Xmax
-
Xmin
）；
% accept it or not
if
all
（
var1
>
=
Xmin
）
&&
all
（
var1
<
=
Xmax
）
new_fun
=
Obj_Fit
（
var1
，
t
，
Et
）；
if
（
new_fun
<
=
fun
）
var
=
var1
；
fun
=
new_fun
；
else
if
exp
（
-
（
new_fun
-
fun
）
/
Tk
）
>
rand
var
=
var1
；
fun
=
new_fun
；
end
end
end
end
% results
p0
=
［
var
（
1
），
10
^
var
（
2
），
var
（
3
）］
‘
；
end
function
fitError
=
Obj_Fit
（
p0， t ，Et
）
A1
=
p0
（
1
）；
B1
=
p0
（
2
）；
A0
=
p0
（
3
）；
f
=
A1
*
exp
（
-
t
/
10
^
（
B1
+
eps
））
+
A0
；
fitError
=
norm
（
Et
-
f
）；
end
模擬退火演算法求解旅行商問題（SA-TSP）
clear
；
clc
；
ratio
=
0。99
；
t_start
=
100
；
t_end
=
1
；
Markov_length
=
5000
；
coordinates
=
load
（
’CityCoordinates。txt‘
）；
amount
=
size
（
coordinates
，
1
）；
% 城市的數目
% 透過向量化的方法計算距離矩陣
coor_x_tmp1
=
coordinates
（：，
1
）
*
ones
（
1
，
amount
）；
coor_x_tmp2
=
coor_x_tmp1
’
；
coor_y_tmp1
=
coordinates
（：，
2
）
*
ones
（
1
，
amount
）；
coor_y_tmp2
=
coor_y_tmp1
‘
；
dist_matrix
=
sqrt
（（
coor_x_tmp1
-
coor_x_tmp2
）
。^
2
+
（
coor_y_tmp1
-
coor_y_tmp2
）
。^
2
）；
sol_new
=
1
：
amount
；
% 產生初始解
% sol_new 是每次產生的新解；
% sol_current 是當前解；
% sol_best 是冷卻中的最好解；
% E_current 是當前解對應的迴路距離；
% E_new 是新解的迴路距離；
% E_best 是最優解的
E_current
=
inf
；
E_best
=
inf
；
sol_current
=
sol_new
；
sol_best
=
sol_new
；
t
=
t_start
；
while
t
>
=
t_end
% 當前溫度下迴圈
for
r
=
1
：
Markov_length
% 產生隨機擾動
if
（
rand
<
0。5
）
% 隨機決定是進行兩交換還是三交換
% 兩交換
ind
=
randperm
（
amount
，
2
）；
ind1
=
ind
（
1
）；
ind2
=
ind
（
2
）；
tmp1
=
sol_new
（
ind1
）；
sol_new
（
ind1
）
=
sol_new
（
ind2
）；
sol_new
（
ind2
）
=
tmp1
；
else
% 三交換
ind
=
randperm
（
amount
，
3
）；
% 確保ind1 < ind2 < ind3
ind
=
sort
（
ind
）；
ind1
=
ind
（
1
）；
ind2
=
ind
（
2
）；
ind3
=
ind
（
3
）；
tmplist1
=
sol_new
（
ind1
+
1
：
ind2
-
1
）；
sol_new
（
ind1
+
1
：
ind1
+
ind3
-
ind2
+
1
）
=
sol_new
（
ind2
：
ind3
）；
sol_new
（
ind1
+
ind3
-
ind2
+
2
：
ind3
）
=
tmplist1
；
end
% 計算目標函式值
E_new
=
0
；
for
i
=
1
：
amount
-
1
E_new
=
E_new
+
dist_matrix
（
sol_new
（
i
），
sol_new
（
i
+
1
））；
end
% 再算上從最後一個城市到第一個城市的距離
E_new
=
E_new
+
dist_matrix
（
sol_new
（
amount
），
sol_new
（
1
））；
if
E_new
<
E_current
E_current
=
E_new
；
sol_current
=
sol_new
；
if
E_new
<
E_best
% 把冷卻過程中最好的解儲存下來
E_best
=
E_new
；
sol_best
=
sol_new
；
end
else
% 若新解的目標函式值大於當前解的，則僅以一定機率接受新解
if
rand
<
exp
（
-
（
E_new
-
E_current
）
/
t
）
E_current
=
E_new
；
sol_current
=
sol_new
；
else
sol_new
=
sol_current
；
end
end
end
t
=
t
*
ratio
；
% 降溫
end
disp
（［
’最優解為：‘
，
num2str
（
sol_best
）］）；
disp
（［
’最短距離：‘
，
num2str
（
E_best
）］）；
x_path
=
［
coordinates
（
sol_best
，
1
）；
coordinates
（
sol_best
（
1
），
1
）］；
y_path
=
［
coordinates
（
sol_best
，
2
）；
coordinates
（
sol_best
（
1
），
2
）］；
plot
（
x_path
，
y_path
，
’LineWidth‘
，
2
）
樽海鞘演算法（SSA）
演算法下載：
clear
；
clc
；
% 搜尋數量
SearchAgents_no
=
30
；
% 最大迭代次數
Max_iter
=
200
；
% 變數維度
dim
=
2
；
% 上限
ub
=
10
*
ones
（
1
，
dim
）；
% 下限
lb
=
-
10
*
ones
（
1
，
dim
）；
% 初始化位置
Positions
=
rand
（
SearchAgents_no
，
dim
）
。*
repmat
（
ub
-
lb
，
SearchAgents_no
，
1
）
+
repmat
（
lb
，
SearchAgents_no
，
1
）；
Convergence_curve
=
zeros
（
1
，
Max_iter
）；
SalpFitness
=
nan
（
1
，
Max_iter
）；
Sorted_salps
=
nan
（
SearchAgents_no
，
dim
）；
% 計算適應度值
for
i
=
1
：
SearchAgents_no
SalpFitness
（
1
，
i
）
=
ObjFun
（
Positions
（
i
，
：））；
end
［
sorted_salps_fitness
，
sorted_indexes
］
=
sort
（
SalpFitness
）；
for
i
=
1
：
SearchAgents_no
Sorted_salps
（
i
，
：）
=
Positions
（
sorted_indexes
（
i
），
：）；
end
FoodPosition
=
Sorted_salps
（
1
，
：）；
FoodFitness
=
sorted_salps_fitness
（
1
）；
Convergence_curve
（
1
）
=
FoodFitness
；
% 主迴圈
for
l
=
2
：
Max_iter
c1
=
2
*
exp
（
-
（
4
*
l
/
Max_iter
）
^
2
）；
for
i
=
1
：
SearchAgents_no
SalpPositions
=
Positions
’
；
if
i
<
=
SearchAgents_no
/
2
for
j
=
1
：
dim
c2
=
rand
（）；
c3
=
rand
（）；
if
c3
<
0。5
SalpPositions
（
j
，
i
）
=
FoodPosition
（
j
）
+
c1
*
（（
ub
（
j
）
-
lb
（
j
））
*
c2
+
lb
（
j
））；
else
SalpPositions
（
j
，
i
）
=
FoodPosition
（
j
）
-
c1
*
（（
ub
（
j
）
-
lb
（
j
））
*
c2
+
lb
（
j
））；
end
end
elseif
i
>
SearchAgents_no
/
2
point1
=
SalpPositions
（：
，
i
-
1
）；
point2
=
SalpPositions
（：
，
i
）；
SalpPositions
（：，
i
）
=
（
point2
+
point1
）
/
2
；
end
Positions
=
SalpPositions
‘
；
end
% 更新
for
i
=
1
：
SearchAgents_no
% 邊界處理：越界賦值為邊界值
Flag4ub
=
Positions
（
i
，
：）
>
ub
；
Flag4lb
=
Positions
（
i
，
：）
<
lb
；
Positions
（
i
，
：）
=
（
Positions
（
i
，
：）
。*
（
~
（
Flag4ub
+
Flag4lb
）））
+
ub
。*
Flag4ub
+
lb
。*
Flag4lb
；
SalpFitness
（
i
）
=
ObjFun
（
Positions
（
i
，
：））；
if
SalpFitness
（
i
）
<
FoodFitness
FoodPosition
=
Positions
（
i
，
：）；
FoodFitness
=
SalpFitness
（
i
）；
end
end
Convergence_curve
（
l
）
=
FoodFitness
；
end
disp
（［
’最優值：‘
num2str
（
FoodPosition
）］）；
disp
（［
’最優值對應目標函式：‘
num2str
（
FoodFitness
）］）；
subplot
（
1
，
2
，
1
）
plot
（
Convergence_curve
，
’LineWidth‘
，
2
）；
xlabel
（
’迭代次數‘
）；
ylabel
（
’適應度值‘
）；
title
（
’適應度曲線‘
）；
len
=
50
；
xRange
=
linspace
（
lb
（
1
），
ub
（
1
），
len
）；
yRange
=
linspace
（
lb
（
2
），
ub
（
2
），
len
）；
［
xMap
，
yMap
］
=
meshgrid
（
xRange
，
yRange
）；
zMap
=
zeros
（
len
）；
for
i
=
1
：
len
for
j
=
1
：
len
zMap
（
i
，
j
）
=
ObjFun
（［
xMap
（
i
，
j
），
yMap
（
i
，
j
）］）；
end
end
subplot
（
1
，
2
，
2
）
surf
（
xRange
，
yRange
，
zMap
）；
view
（
-
45
，
-
45
）；
shading
interp
hold
on
plot3
（
FoodPosition
（
1
），
FoodPosition
（
2
），
FoodFitness
，
’o‘
，
’MarkerFaceColor‘
，
’r‘
，
’MarkerSize‘
，
10
）；
hold
off
set
（
gcf
，
’Position‘
，
［
400
，
200
，
900
，
400
］）；
function
result
=
ObjFun
（
x
）
x1
=
x
（
1
）；
x2
=
x
（
2
）；
t
=
x1
^
2
+
x2
^
2
；
result
=
0。5
+
（
sin
（
sqrt
（
t
））
^
2
-
0。5
）
/
（
1
+
0。001
*
t
）
^
2
；
end
狀態轉移演算法（STA）
改進狀態轉移演算法（STA）
鯨魚演算法（WOA）
演算法下載：
The Whale Optimization Algorithm
clear
；
clc
；
% 鯨魚數量
SearchAgents_no
=
30
；
% 最大迭代次數
Max_iter
=
200
；
% 變數維度
dim
=
2
；
% 上限
ub
=
10
*
ones
（
1
，
dim
）；
% 下限
lb
=
-
10
*
ones
（
1
，
dim
）；
%%%%% 初始化
Leader_pos
=
zeros
（
1
，
dim
）；
Leader_score
=
inf
；
% 初始化位置
Positions
=
rand
（
SearchAgents_no
，
dim
）
。*
repmat
（
ub
-
lb
，
SearchAgents_no
，
1
）
+
repmat
（
lb
，
SearchAgents_no
，
1
）；
Convergence_curve
=
zeros
（
1
，
Max_iter
）；
% 主迴圈
for
t
=
1
：
Max_iter
for
i
=
1
：
SearchAgents_no
% 邊界處理：越界賦值為邊界值
Flag4ub
=
Positions
（
i
，
：）
>
ub
；
Flag4lb
=
Positions
（
i
，
：）
<
lb
；
Positions
（
i
，
：）
=
（
Positions
（
i
，
：）
。*
（
~
（
Flag4ub
+
Flag4lb
）））
+
ub
。*
Flag4ub
+
lb
。*
Flag4lb
；
% 計算目標函式
fitness
=
ObjFun
（
Positions
（
i
，
：））；
% 更新領頭位置
if
fitness
<
Leader_score
Leader_score
=
fitness
；
Leader_pos
=
Positions
（
i
，
：）；
end
end
% 線性遞減
a
=
2
-
t
*
2
/
Max_iter
；
a2
=
-
1
-
t
/
Max_iter
；
% 更新位置
for
i
=
1
：
SearchAgents_no
r1
=
rand
（）；
r2
=
rand
（）；
A
=
2
*
a
*
r1
-
a
；
C
=
2
*
r2
；
b
=
1
；
l
=
（
a2
-
1
）
*
rand
+
1
；
p
=
rand
（）；
for
j
=
1
：
dim
if
p
<
0。5
if
abs
（
A
）
>
=
1
rand_leader_index
=
floor
（
SearchAgents_no
*
rand
（）
+
1
）；
X_rand
=
Positions
（
rand_leader_index
，
：）；
D_X_rand
=
abs
（
C
*
X_rand
（
j
）
-
Positions
（
i
，
j
））；
Positions
（
i
，
j
）
=
X_rand
（
j
）
-
A
*
D_X_rand
；
elseif
abs
（
A
）
<
1
D_Leader
=
abs
（
C
*
Leader_pos
（
j
）
-
Positions
（
i
，
j
））；
Positions
（
i
，
j
）
=
Leader_pos
（
j
）
-
A
*
D_Leader
；
end
else
distance2Leader
=
abs
（
Leader_pos
（
j
）
-
Positions
（
i
，
j
））；
Positions
（
i
，
j
）
=
distance2Leader
*
exp
（
b
。*
l
）
。*
cos
（
l
。*
2
*
pi
）
+
Leader_pos
（
j
）；
end
end
end
Convergence_curve
（
t
）
=
Leader_score
；
end
disp
（［
’最優值：‘
num2str
（
Leader_pos
）］）；
disp
（［
’最優值對應目標函式：‘
num2str
（
Leader_score
）］）；
subplot
（
1
，
2
，
1
）
plot
（
Convergence_curve
，
’LineWidth‘
，
2
）；
xlabel
（
’迭代次數‘
）；
ylabel
（
’適應度值‘
）；
title
（
’適應度曲線‘
）；
len
=
50
；
xRange
=
linspace
（
lb
（
1
），
ub
（
1
），
len
）；
yRange
=
linspace
（
lb
（
2
），
ub
（
2
），
len
）；
［
xMap
，
yMap
］
=
meshgrid
（
xRange
，
yRange
）；
zMap
=
zeros
（
len
）；
for
i
=
1
：
len
for
j
=
1
：
len
zMap
（
i
，
j
）
=
ObjFun
（［
xMap
（
i
，
j
），
yMap
（
i
，
j
）］）；
end
end
subplot
（
1
，
2
，
2
）
surf
（
xRange
，
yRange
，
zMap
）；
view
（
-
45
，
-
45
）；
shading
interp
hold
on
plot3
（
Leader_pos
（
1
），
Leader_pos
（
2
），
Leader_score
，
’o‘
，
’MarkerFaceColor‘
，
’r‘
，
’MarkerSize‘
，
10
）；
hold
off
set
（
gcf
，
’Position‘
，
［
400
，
200
，
900
，
400
］）；
function
result
=
ObjFun
（
x
）
x1
=
x
（
1
）；
x2
=
x
（
2
）；
t
=
x1
^
2
+
x2
^
2
；
result
=
0。5
+
（
sin
（
sqrt
（
t
））
^
2
-
0。5
）
/
（
1
+
0。001
*
t
）
^
2
；
end
