1 简介
萤火虫（Glowworms）是一群能发光的昆虫，也被称为闪电虫，它们使用一种叫做生物发光（bioluminescence）的过程来发光。然而，已经发现了许多具有类似发光行为的生物，如水母、某些细菌、原生动物、水生动物等，事实上，大约80%到90%的海洋生物是由发光生物组成的。萤火虫流行的原因是它们很容易被发现且数量巨大。
看到萤火虫在夏夜忽明忽暗地眨眼，要比看到它们的群体行为容易得多，大量的萤火虫聚集在一起形成一个群体，同时发出闪光，在目睹了这些美丽的景象之后，人们不禁想知道萤火虫这种发光和聚集行为的“原理”和“原因”。这种发光行为的原因是为了吸引它们不知情的猎物进入陷阱，并吸引配偶进行繁殖。在繁殖过程中，既可以观察到个体的求偶也可以观察到群体的交配。 萤火虫的生命周期从卵开始，然后从卵变成蛹，蛹变成幼虫，幼虫变成成虫，萤火虫只需要几周的时间就可以成年，因此，以保留物种为目标，寻找交配对象进行繁殖是当务之急。
1.1 萤火虫群优化算法
针对萤火虫和萤火虫群的行为，Krishnanand和Ghose（2005）[1]提出了萤火虫群优化算法（Glowworm Swarm Optimization，GSO），并使用于诸多应用。最初，GSO的开发目标是提供数值优化问题的解决方案，而不是确定全局最优，但是由于GSO的分散决策和移动协议，它在机器人等领域做出了更多的贡献。
最初受萤火虫启发，GSO在随机搜索空间中随机分布一组或一群代理，代理间通过其他行为机制相互影响，而这些机制在其自然界中的对应物中是不存在的。算法的基本工作基于以下三种机制：
- 适应度广播 萤火虫有一种叫做荧光素（luciferin）的色素，可以使萤火虫发光。萤火虫体内荧光素的含量决定了它们在目标空间中的位置的适应度。
- 正趋性 萤火虫被比自己亮的邻居吸引，因此开始向它移动。当有多个这样的邻居时，它利用概率机制来选择一个。
- 自适应邻域 每个萤火虫利用一个自适应邻域来识别邻居，该邻域由一个具有可变范围rdi的局部决策域定义，该域的边界是一个硬限制的感知范围rs(0<rdi<rs)。这里可以使用一种合适的启发式方法来调节rdi。萤火虫的运动完全依赖于局部信息，每个萤火虫都会选择一个邻居，这个邻居的荧光素值大于它自己的荧光素值，然后向它靠近。这些运动建立在可用的局部信息和选择性的邻居交互的基础上，使得萤火虫群能够聚集成不相交的子群，朝着给定多峰函数的多个最优值移动并在其上相遇。
1.2 算法
虽然该算法被解释为寻找多峰函数的多个最优解，但可以通过简单修改用于最小化问题。最初，GSO将萤火虫随机地放置在搜索空间中，使它们分布得很好。初始时，每只萤火虫体内的荧光素含量为零。该算法的单位周期包括荧光素更新阶段、移动阶段和邻域范围更新阶段，如图1所示。
图1 GSO算法流程
2 部分代码
function [xtraj, ttraj, terminate_cond] = test_trajectory(start, stop, map, path, vis)
% TEST_TRAJECTORY simulates the robot from START to STOP following a PATH
% that's been planned for MAP.
% start - a 3d vector or a cell contains multiple 3d vectors
% stop - a 3d vector or a cell contains multiple 3d vectors
% map - map generated by your load_map
% path - n x 3 matrix path planned by your dijkstra algorithm
% vis - true for displaying visualization
%Controller and trajectory generator handles
controlhandle = @controller;
trajhandle = @trajectory_generator;
% Make cell
if ~iscell(start), start = {start}; end
if ~iscell(stop), stop = {stop}; end
if ~iscell(path), path = {path} ;end
% Get nquad
nquad = length(start);
% Make column vector
for qn = 1:nquad
start{qn} = start{qn}(:);
stop{qn} = stop{qn}(:);
end
% Quadrotor model
params = nanoplus();
%% **************************** FIGURES *****************************
% Environment figure
if nargin < 5
vis = true;
end
fprintf('Initializing figures...\n')
if vis
h_fig = figure('Name', 'Environment');
else
h_fig = figure('Name', 'Environment', 'Visible', 'Off');
end
if nquad == 1
plot_path(map, path{1});
else
% for qn = 1:nquad
% plot_path(map, path{qn});
% end
% you could modify your plot_path to handle cell input for multiple robots
end
h_3d = gca;
drawnow;
xlabel('x [m]'); ylabel('y [m]'); zlabel('z [m]')
quadcolors = lines(nquad);
set(gcf,'Renderer','OpenGL')
%% Trying Animation of Blocks
NoofBlocks = size(map(:,1),1);
x_0 = map(1,4);
x_1 = map(2,4);
y_0 = map(1,5);
y_1 = map(2,5);
z_0 = map(1,6);
z_1 = map(2,6);
for i=1:2:NoofBlocks
xb_0 = map(i,1);
xb_1 = map(i+1,1);
yb_0 = map(i,2);
yb_1 = map(i+1,2);
zb_0 = map(i,3);
zb_1 = map(i+1,3);
B_1 = [xb_0 yb_0 zb_0]';
B_2 = [xb_1 yb_0 zb_0]';
B_3 = [xb_0 yb_0 zb_1]';
B_4 = [xb_1 yb_0 zb_1]';
B_5 = [xb_0 yb_1 zb_0]';
B_6 = [xb_1 yb_1 zb_0]';
B_7 = [xb_0 yb_1 zb_1]';
B_8 = [xb_1 yb_1 zb_1]';
% BlockCoordinatesMatrix(j:j+7,:) = [B_1';B_2';B_3';B_4';B_5';B_6';B_7';B_8'];
% BlockCoordinatesMatrix(j:j+1,:) = [B_1';B_8'];
% BlockCoordinates(i,:) = {B_1 B_2 B_3 B_4 B_5 B_6 B_7 B_8};
% j = j+2;
S_1 = [B_1 B_2 B_4 B_3];
S_2 = [B_5 B_6 B_8 B_7];
S_3 = [B_3 B_4 B_8 B_7];
S_4 = [B_1 B_2 B_6 B_5];
S_5 = [B_1 B_3 B_7 B_5];
S_6 = [B_2 B_4 B_8 B_6];
fill3([S_1(1,:)' S_2(1,:)' S_3(1,:)' S_4(1,:)' S_5(1,:)' S_6(1,:)'],[S_1(2,:)' S_2(2,:)' S_3(2,:)' S_4(2,:)' S_5(2,:)' S_6(2,:)'],[S_1(3,:)' S_2(3,:)' S_3(3,:)' S_4(3,:)' S_5(3,:)' S_6(3,:)'],[1 0 0]);%[cell2mat(Block(i,8))/255 cell2mat(Block(i,9))/255 cell2mat(Block(i,10))/255]);
xlabel('x'); ylabel('y'); zlabel('z');
axis([min(x_0,x_1) (max(x_0,x_1)) min(y_0,y_1) (max(y_0,y_1)) min(z_0,z_1) (max(z_0,z_1))])
grid
hold on
end
%% *********************** INITIAL CONDITIONS ***********************
fprintf('Setting initial conditions...\n')
% Maximum time that the quadrotor is allowed to fly
time_tol = 50; % maximum simulation time
starttime = 0; % start of simulation in seconds
tstep = 0.01; % this determines the time step at which the solution is given
cstep = 0.05; % image capture time interval
nstep = cstep/tstep;
time = starttime; % current time
max_iter = time_tol / cstep; % max iteration
for qn = 1:nquad
% Get start and stop position
x0{qn} = init_state(start{qn}, 0);
xtraj{qn} = zeros(max_iter*nstep, length(x0{qn}));
ttraj{qn} = zeros(max_iter*nstep, 1);
end
% Maximum position error of the quadrotor at goal
pos_tol = 0.05; % m
% Maximum speed of the quadrotor at goal
vel_tol = 0.10; % m/s
x = x0; % state
flag = 0;
%% ************************* RUN SIMULATION *************************
fprintf('Simulation Running....\n')
for iter = 1:max_iter
timeint = time:tstep:time+cstep;
tic;
% Iterate over each quad
for qn = 1:nquad
% Initialize quad plot
if iter == 1
QP{qn} = QuadPlot(qn, x0{qn}, 0.1, 0.04, quadcolors(qn,:), max_iter, h_3d);
desired_state = trajhandle(time, qn);
QP{qn}.UpdateQuadPlot(x{qn}, [desired_state.pos; desired_state.vel], time);
h_title = title(sprintf('iteration: %d, time: %4.2f', iter, time));
end
% Run simulation
[tsave, xsave] = ode45(@(t,s) quadEOM(t, s, qn, controlhandle, trajhandle, params), timeint, x{qn});
x{qn} = xsave(end, :)';
% Save to traj
xtraj{qn}((iter-1)*nstep+1:iter*nstep,:) = xsave(1:end-1,:);
ttraj{qn}((iter-1)*nstep+1:iter*nstep) = tsave(1:end-1);
% Update quad plot
desired_state = trajhandle(time + cstep, qn);
if flag == 1
desired_state.vel = [0; 0 ;0];
flag = 0;
end
QP{qn}.UpdateQuadPlot(x{qn}, [desired_state.pos; desired_state.vel], time + cstep);
end
set(h_title, 'String', sprintf('iteration: %d, time: %4.2f', iter, time + cstep))
time = time + cstep; % Update simulation time
t = toc;
% Pause to make real-time
if (t < cstep)
pause(cstep - t);
end
% Check termination criteria
terminate_cond = terminate_check(x, time, stop, pos_tol, vel_tol, time_tol);
if terminate_cond == 3
flag = 1;
end
if terminate_cond ~= 0 && terminate_cond ~= 3
break
end
end
fprintf('Simulation Finished....\n')
%% ************************* POST PROCESSING *************************
% Truncate xtraj and ttraj
for qn = 1:nquad
xtraj{qn} = xtraj{qn}(1:iter*nstep,:);
ttraj{qn} = ttraj{qn}(1:iter*nstep);
end
% Plot the saved position and velocity of each robot
if vis
for qn = 1:nquad
% Truncate saved variables
QP{qn}.TruncateHist();
% Plot position for each quad
h_pos{qn} = figure('Name', ['Quad ' num2str(qn) ' : position']);
plot_state(h_pos{qn}, QP{qn}.state_hist(1:3,:), QP{qn}.time_hist, 'pos', 'vic');
plot_state(h_pos{qn}, QP{qn}.state_des_hist(1:3,:), QP{qn}.time_hist, 'pos', 'des');
% Plot velocity for each quad
h_vel{qn} = figure('Name', ['Quad ' num2str(qn) ' : velocity']);
plot_state(h_vel{qn}, QP{qn}.state_hist(4:6,:), QP{qn}.time_hist, 'vel', 'vic');
plot_state(h_vel{qn}, QP{qn}.state_des_hist(4:6,:), QP{qn}.time_hist, 'vel', 'des');
end
end
end
3 仿真结果
4 参考文献
-
Krishnanand, K.N. and D. Ghose. Detection of multiple source locations using a glowworm metaphor with applications to collective robotics. in Proceedings 2005 IEEE Swarm Intelligence Symposium, 2005. SIS 2005. 2005
-
部分理论引用网络文献，若有侵权联系博主删除。
5 MATLAB代码与数据下载地址
见博客主页
