动物集群运动行为模型系列之五

1、动物集群运动模型摘要本文主要模拟了鱼群的集群运动、鱼群躲避捕食者追捕的运动情况以及鸟群觅食运动的模拟，以此研究动物个体问的信息传递机制，同时也是对群体智能的初步探索。针对问题一，需要我们给出对鱼群集群运动的模型，并编写程序将运动模拟出来，对此我们建立了Boid模型，根据模型给出的准则以及算法，我们通过matlab编程，在忽略阻力等因素下分别模拟出在平面以及空间鱼群的运动，并得出密度必须大于一定值时，鱼群才能最终达到同步。鱼群的整个集群运动从刚开始的随机产生的各个个体的不均匀无规则分布到逐渐的聚拢成群再到最后的一致方向的前进。针对问题二，我们在问题一的模型的基础上给出了鱼群躲避捕食者的模型,制定

2、了鱼个体的适度逃离区域和加速逃离区域，分析捕食者与鱼个体的关系，给出进一步的模型，通过编写程序得到模拟的结果，得到了对鱼群躲避捕食者的运动的合理的动态模拟，并且给出了模型的改进方向。针对问题三，我们更加倾向于研究鸟群的觅食行为，因此我们将问题改成鸟群的觅食模拟，将鸟群的觅食行为转化为求最优解的问题，这正好与问题中提到了有一部分个体掌握食物源位置信息相对应。针对问题，我们建立了粒子群优化模型，通过PSCM法，通过鸟群寻找食物的最短路径的最优解的问题的分析，我们利用优化算法来模拟了鸟群在山问的觅食行为，得到了鸟群可以绕过我们设定的障碍物（山峰）到达食物点。关键字：动物集群运动Boid模型PSO算法

3、鸟群觅食一、问题重述在动物界，大量集结成群进行移动或者觅食的例子并不少见，这种现象在食草动物、鸟、鱼和昆虫中都存在。这些动物群在运动过程中具有很明显的特征：群中的个体聚集性很强，运动方向、速度具有一致性。通过数学模型来模拟动物群的集群运动行为以及探索动物群中的信息传递机制一直是仿生学领域的一项重要内容。通过观察附件中给出的图片和视频资料，或者在网上搜索相关资料观察，思考动物集群运动的机理，建立数学模型刻画动物集群运动、躲避威胁等行为，例如，可以考虑以下问题的分析建模：建立数学模型模拟动物的集群运动。建立数学模型刻画鱼群躲避黑鳍礁鲨鱼的运动行为。假定动物群中有一部分个体是信息丰富者（如掌握食物源

4、位置信息，掌握迁徙路线信息），请建模分析它们对丁群运动行为的影响，解释群运动方向决策如何达成。建议与说明：在上述问题的讨论中，如果能适时分析动物群中的信息传递机制无疑是更好I如果对问题2和问题3之外的其他集群运动行为更感兴趣，也可将这两个问题替换为你所感兴趣的问题来讨论。建模过程中的数据资料可以在网上查询或者自行合理设定。若果感到在三维空间讨论问题太复杂，可以先在二维空间讨论，再推广至三维空间。最好能对你所做的机理分析模型给出计算机仿真方法以便丁实际情况对比评价。二、模型假设忽略障碍、阻力以及其它无关次要因素对丁集群运动的影响问题一鱼群中每个个体运动的速度都是包定一样的鱼群集群运动的模拟中不考

5、虑障碍物的存在忽略其它种群对本文所研究种群的影响不考虑集群中个体的体积，都按粒子处理三、符号说明N鱼群的总数&集群中每个个体的位置欠量瑚集群中每个个体的速度欠量V集群中每个个体运动的速度rr排斥区域的半径r。Bm dVminVmaxTXi一Xpdi(t 1)G、C2若、&Xi (t)=(为 i(t),Xi2(t),., Xd (t)Vi(t) =(Vi(t), V2(t).,Vid(t)P(t) = Pi(t), Pi2(t),., Pid(t)%t)= P gl(t), Pg2(t),.Pgd(t)吸引区域的半径惯性权重粒子数空间维数最小速度最大速度粒子i的位置欠M捕食者的位置欠M个体i在t

6、 +1时刻的预期方向学习因子（加速因子）均匀分布在（0,1）之间的随机数 在第t次迭代时粒子i的位置表示 在第t次迭代时粒子i的速度表示 个体极值全局极值致区域的半径四、问题分析本问题是一个动物集群运动的模型问题，动物的集群运动包括很多，其中有觅食、追尾、躲避捕食者等等运动，问题一需要我们考虑动物集群的运动模型，也就是鱼群的游弋、鸟群的飞翔等行为，是不需要考虑觅食、追尾等行为活动的，我们通过建立Boid模型进行鱼群集群运动的模拟。问题二需要我们给出鱼群躲避捕食者的运动模拟，要解决问题，那就需要我们在问题一模型的基础上给出鱼的逃逸模式，然后对逃逸运动进行模拟。问题三是需要我们模拟集群运动中存在领

7、导者时的集群运动的模拟，可以运用和问题一一样的思路，但是我们对丁鸟群的觅食行为更感兴趣，所以我们转而对鸟群觅食进行建模，我们选择PSO算法,通过模型求最优解的过程对鸟群觅食行为进行模拟，从而建立起了比较合适的模型。针对这些问题，我们主要的工作是首先建立合适的模型，通过我们建立的模型，根据模型的算法，我们可以编写程序得到对集群运动的模拟。五、模型的建立与求解鱼群集群运动的模拟模型的建立我们根据问题一的要求，通过查阅资料得出了Boid模型可以解决类似的问题，这里我们就选用Boid模型作为此文的模型。Boid群模型包括三个简单的指导规则，它们描述了一个单一的“Boid”如何基丁位置和邻近个体的速度进

8、行分离、内聚和排序活动的。模型的三个规则(Reynolds聚结规则)如下：群中心定位：试图与邻近的群个体保持接近；避免障碍：避免与邻近的群个体发生冲突；速度匹配：试图与邻近的群个体速度匹配；设系统有N个个体组成，它们的位置和速度欠量分别为C*、寸，每个个体在三维空间中按照包定的速度V运动，dj(t+1)为个体i在t+1时刻的预期方向。在每一步t,每个个体可以感知到三个不重叠的区域中其他个体的位置和角度，这些信息用丁计算di(t+1),这三个区域分别为：排斥区域，一致区域，吸引区域，也称为Three-circle模型，其模型的三个区域如图(图1)所示：图1Boid模型的三个区域示意图，zor为排

9、斥区域，zoo为一致区域，zoa为吸引区域；3600-。为视野盲区，?取值范围为200-3600个体的运动规则为：首先，每个个体尽M与排斥区域(以该个体为中心，以r为半径的球)中的其他个体保持最小距离，并记其中的个体数为nr，则个体i的预期方向按照下面的方式调整nrdi(t 1)=-、AAi = 1,., N ,t a 0j=irj (t)其中rj(t)=(Cj-G)。其次，如果n二0，贝U个体i的预期方向受“一致区域”(以个体i位中心，处丁rr和r。之间的球形区域，除去该个体后面的角度为360-ot的盲区)及“吸引区域”(以个体i位中心，处丁r。和G之间的球形区域，除去该个体后面的角度为36

10、Oo-a的盲区)中的个体的影响，记相应区域中的邻居个数分别为n。、na,则可定义do(t+1)、da(t+1)如下:do(t 1)=Vj(t)rij (t)j 主-V j (tda(t + 1)=Ej# rij (t)如果na-fi1)a=0,则djt+1)=d(t+1);同样，如果n0,则d(t+1)=d(t+1);如果两者都不为0,贝U定义di(t+1)=5,d(t+1)+da(t+1)；如果经过上面的运算后所得到的di(t+1)=0，或者在三个区域中都没有个体，贝Udi(t+1)=Vi(t)o设旋转速率为9,即每一时步个体所能转过的最大角度为9 ,如果Vi与di (t +1)之间的角度差

11、小丁6,则di(t+1)=Vi(t+1),否则，个体i向期望的方向旋转角度6，这样就得到了个体下一步的运动方向。模型的求解根据模型，我们利用matlab编程，得到二维空间(平面)和三维空间(立体)的鱼群运动的模拟。平面鱼群运动的模拟程序见附录1。我们取了个体为100，个体的旋转速率为0.5,得到了运动的模拟动画。卜面我们截取了整个动态模拟过程中的3幅图。空间鱼群运动的模拟程序见附录2。我们取了个体为100,个体的旋转速率为0.5,得到了运动的模拟动画。下面我们截取了整个动态模拟过程中的3幅图。图7个体一致向一个方向运动截图从上面的各个图中，我们可以看出，鱼群的整个集群运动从刚开始的随机产生的各

12、个个体的不均匀无规则分布到逐渐的聚拢成群再到最后的一致方向的前进。整个模拟过程是比较合理的，比较符合实际，由丁论文上无法显示动态的图,所以我们采取了截图的方式。从我们的模拟中，我们还发现，取不同的密度，得到的模拟结果是不一样的,在高密度的情景下，我们得出鱼群在经过有限的运动时间后，会最终达到同步，即运动方向达到一致。而且只有在密度大丁一定范围时系统才能最终达到一致同步，随着个体密度的减小，我们可以推测存在一个临界密度，只有在大丁临界密度的时候系统才能最终达到同步。密度越大，越容易达到同步。密度约小，鱼群只能形成多个小规模的鱼群，无法形成整个鱼群运动一致。取不同的旋转速率时，所得到的模拟动画也不

13、一样。通过输入不同的旋转速率，我们得出结论：旋转速率越大，鱼群总是在原位置附近徘徊，无法达到一致同步的趋势；旋转速率较小时，能够达到最终的一致同步，即使不能达到一致同步，也能以小规模的集群同步运动。鱼群躲避捕食者的运动的模拟模型的建立在问题一的模型的基础上，为考虑有外来捕食者的情况下群集的应急机制，设计了如下的粒子个体应急措施。如下图所示图8捕食者模型示意图粒子个体的应急区域分为两个层次：适度逃离区和加速逃离区。图阴影部分即为加速逃离区，内圆与外圆间的环形部分是适度逃离区。图下半部分显示了捕食者分别位丁两个区域时，粒子个体的受力大小和方向。为简化模型，规定加速逃离区同原模型的排斥区重合。当有捕

14、食者进入到粒子个体的加速逃离区时，“逃跑”就是个体的最高策略，个体会沿捕食者与个体所在的直线告诉逃离捕食者，速度的大小同捕食者和个体间的距离存在一定的函数关系，此时个体将不考虑吸引、排斥异己方向同步等的作用。当有捕食者位丁粒子个体的适度逃离区时，个体会受到捕食者与个体所在直线远离捕食者方向的作用,同时还将受到吸引、排斥以及方向同步等的作用。因此，个体受到的实际作用是两者作用的叠加。还应指出，当捕食者位丁粒子个体的适度逃离区时，个体会加大方向同步区域的大小，期望尽快与其他粒子的方向保持同步。粒子的逃离速度公式如下：Vi(Xi-Xp),X=norm(Xi-Xp)其中a是加速逃离区半径的大小，与排斥

15、区域半径大小一样，即a=r.b是适度逃离区的半径大小，与吸引区域半径一样，即b=rac是正常情况下同步区域半径大小，d是此半径的上界（不超过b）。当捕食者位丁粒子的加速逃离区时，粒子的实际速度就是其逃离速度，当捕食者位丁粒子的适度逃离区时，粒子的实际速度是其逃离速度与在吸引、排斥以及方向同步等的作用下的速度的合成。模型的求解通过matlab编程，我们得出了对鱼群逃逸的运动的模拟（程序见附录3）截图如下：matlab 编程得到12鱼群躲避捕食者的运动模拟截图（绿色圆圈为捕食者，红色义为鱼）结论：我们通过对鱼群躲避捕食者的模型建立，然后通过了对鱼群运动的模拟，总的可以看出鱼群的逃逸运动。对鸟群的觅

16、食运动的模拟5.3.1模型的分析与建宜通过对问题的分析，我们建立了PSO莫型对鸟群的觅食行为进行模拟5.3.1.1PSO模型的简介粒子群优化算法是基丁群体的演化算法，其思想来源丁人工生命和演化计算理论。Reynolds对鸟群飞行的研究发现，鸟仅仅是追踪它有限数量的邻居，但最终的整体结果是整个鸟群好像在一个中心的控制之下，即复杂的全局行为是由简单规则的相互作用引起的。PSO即源丁对鸟群捕食行为的研究,一群鸟在随机的搜寻食物，如果区域只有一块食物，那么找到食物的最简单有效的策略就是搜寻目前离食物最近的鸟的周围区域。人们通常是以他们自己及他人的经验来作为决策的依据，这就构成了PSO勺一个基本概念。算

17、法的原理PSO算法将群体中的每个个体视为多维搜索空间中的一个没有质量和体积的粒子，这些粒子在搜索空间中以一定的速度飞行，并根据粒子本身的飞行经验以及同伴的飞行经验对自己的飞行速度进行动态调整，即每个粒子通过统计迭代过程中自身的最优解和群体的最优值来不断修正自己的前进方向和速度大小，从而形成群体寻优的正反馈机制。PS顷法就是这样依据每个粒子对环境的适应度将个体逐步移到较优的区域，并最终搜索、寻找到问题的最优解。PSCB法具有鲜明的生物社会背景：认知过程和社会行为，即在寻求一致的认知过程中，个体往往记住它们的信念，同时考虑其它同伴的信念，当个体察觉同伴的信念较好时，将进行适应性调整。模型的建立在P

18、SCB法中，用粒子的位置表示待优化问题的解，每个粒子性能的优劣程度取决丁待优化问题目标函数确定的适应值，每个粒子由一个速度欠量决定其飞行方向和速率大小。设在一个d维的目标搜索空间中，有m个粒子组成一个群体，其中，在第t次迭代时粒子组成一个群体，其中，在第t次迭代时粒子i的位置表示为Xi(t)=(X1(t),Xi2(t),.,Xd(t),相应的飞行速度表示为Vi(t)=(V1(t),V2(t),.,Md(t)。开始执行PSO算法时，首先随机初始化m个粒子的位置和速度，然后通过迭代寻找最优解，在每一次迭代中，粒子通过跟踪两个极值来更新自己的速度和位置：一个极值是粒子本身迄今为止搜索到的最优解，称为

19、个体极值，表示为P(t)=(Pi(t),Pi2(t),.,Pid(t)；另一个极值是整个粒子群到目前为止找到的最优解，称为全局极值，表示为Pg(t)=(Pg1(t),Pg2(t),.Pgd(t)。在第t+1次迭代计算时，粒子i根据下列规则来更新自己的速度和位置:Vik(t1)=Vik(t)1C1Pk(bXk(t)22CPk(-0Xk(tXk(t1)=Xik(t)Vik(t1);(i=1,2,.,m;k=1,2,.,d)s取大值可使算法具有较强的全局搜索能力，b取小值则算法倾向丁局部搜索。一般的做法是将w初始取0.9,并使其迭代次数的增加而线性递减至0.4,部精细这样下去可以先侧重丁全局搜索，使

20、搜索空间快速收敛丁某一区域，然后采用局搜索以获得高精度的解；G、%一般取2。另外粒子的每一维的速度Vi取值范围为VminVmax。如果当前粒子的加速度导致它在某一维的速度超过该维上的最大速度Vmax，则该维的速度被限制值为最大速度。位置Xi的取值范围为XminXmax。PSO算法流程随机初始化粒子群体的位置和速度；通常是在允许的范围内随机产生的，每个粒子的pbest坐标设置为其当前位置，且计算出其相应的个体极值(即个体的适应度值)，而全局极值(即全局的适应度值)就是个体极值中最好的，记录该最好值的粒子序号，并将gbest设置为该最好粒子的当前位置；计算每个粒子的适应值；对每个粒子，将其适应值与

21、个体极值进行比较，如果较优，则更新当前的个体极值；对每个粒子，将其适应值与全局极值进行比较，如果较优，则更新当前的全局极值；根据更新位置速度的公式，更新每个粒子的位置和飞行速度；如未达到预先设定的停止准则(通常设置为最大迭代次数)，则返回步骤(2),达到则停止计算；5.3.2模型的求解种群大小取30,粒子大小2,Vmax=2,q、C2均取2。根据上述算法，我们通过matlab编程，得出了对鸟群觅食的模拟(程序见附录4)。截图如下：图15鸟群找到食物结论：通过上述的截图我们可以看出，我们作出的模拟十分合理，鸟群能够绕过山峰找到食物，我们通过PSC模型，通过鸟群寻找食物的最短路径的最优解算法得到对

22、鸟群觅食过程的模拟。六、模型评价对丁问题一，我们采用的Boid模型是有缺陷的，它给出的规则都是局部规则，每个个体仅仅根据它周围附近球域内的个体的行为调整自己的行为，只能做到局部一致性，局限性较强，并且模型是在虚拟的没有障碍物存在的空间内进行模拟的，而实际中不可能没有障碍物，基丁以上的缺点因素，我们提出了Vicsek模型，这个模型能克服上诉缺点。对丁问题二，我们在模型一的基础上建立了模型，模型得出的结果也是很合适的，不过模型对丁鱼群的加速逃离区和适度逃离区的选择不是很合理，模型的改进方向可以是：划定一个必须逃离区，在必须逃离区外危险系数与捕食者和个体的距离成正比，在必须逃离区之内就只能全力逃跑。

23、鉴丁时间关系，我们没有再深入下去。对丁问题三，我们采用PSC莫型，通过对鸟群觅食的最短路径最优化问题的分析，利用PSO模型作出对鸟群觅食的模拟。参考文献1姜启源，数学模型M.北京：高等教育出版社.1987年4月第一版；2胡守信，李柏年.基丁MATLAB勺数学实验M.北京：科学出版社.2004年6江铭炎，人工鱼群算法及其应用.科学出版社.:2012年1月1日第一版；附录附录1functionboid1(n,eta)%generatesboidsinafieldn=input(inputn);eta=input(inputeta);range=0.2;speed=0.03;TIME=300;den

24、sity=10;RUNS=1;fieldsize=sqrt(n/density);pos=fieldsize*rand(n,2);heading=2*pi*rand(n,1);relposX=zeros(n);relposY=zeros(n);neighbours=zeros(n);RelHead=zeros(n,1);B=zeros(n,1);newHeading=zeros(n,1);newPos=zeros(n,2);meanHeading=zeros(TIME,1);deviationMean=zeros(n,TIME);forrun=1:RUNSfortime=1:TIMEfori=

25、1:nforj=1:nrelposX(i,j)=abs(pos(i,1)-pos(j,1);relposY(i,j)=abs(pos(i,2)-pos(j,2);if(sqrt(relposX(i,j)人2+relposY(i,j)人2)=range)neighbours(i,j)=1;endendendfori=1:nforj=1:nif(neighbours(i,j)=1)relheading(i,j)=heading(j,1)-heading(i,1);endendendRelHead=sum(relheading,2);fori=1:nwhileRelHead(i,1)piRelHea

26、d(i,1)=RelHead(i,1)-pi;endendnoise=(rand(n,1)*eta)-eta/2;newHeading(:,1)=heading(:,1)+RelHead(:,1)./sum(neighbours,2)+noise(:,1);newPos(:,1)=pos(:,1)+cos(newHeading(:,1)*speed;newPos(:,2)=pos(:,2)+sin(newHeading(:,1)*speed;fork=1:nwhilenewPos(k,1)=fieldsizenewPos(k,1)=newPos(k,1)-fieldsize;endnewPos

27、(k,2)=newPos(k,2)+fieldsize;endwhilenewPos(k,2)=fieldsizenewPos(k,2)=newPos(k,2)-fieldsize;endendheading=newHeading;pos=newPos;meanHeading(time,1)=mean(heading,1);deviationMean(time,1)=mean(abs(meanHeading(time,1)-heading(:,1);scatter(pos(:,1),pos(:,2),xr);axis(0fieldsize0fieldsize);xlabel(time);M(t

28、ime)=getframe;endrunDev(:,run)=deviationMean(:,1);end%plot(runDev)movie2avi(M,boidtest.avi,quality,100);附录2functionboid1(n,eta)%generatesboidsinafieldn=input(inputn);eta=input(inputeta);range=0.2;speed=0.05;TIME=300;density=10;RUNS=1;fieldsize=sqrt(n/density);pos=fieldsize*rand(n,3);heading=2*pi*ran

29、d(n,1);relposX=zeros(n);relposY=zeros(n);relposZ=zeros(n);neighbours=zeros(n);RelHead=zeros(n,1);B=zeros(n,1);newHeading=zeros(n,1);newPos=zeros(n,3);meanHeading=zeros(TIME,1);fortime=1:TIMEfori=1:nforj=1:nrelposX(i,j)=abs(pos(i,1)-pos(j,1);relposY(i,j)=abs(pos(i,2)-pos(j,2);relposZ(i,j)=abs(pos(i,3

30、)-pos(j,3);if(sqrt(relposX(i,j)人2+relposY(i,j)人2+relposZ(i,j)人2)=range)neighbours(i,j)=1;endendendfori=1:nforj=1:nif(neighbours(i,j)=1)relheading(i,j)=heading(j,1)-heading(i,1);endendendRelHead=sum(relheading,2);fori=1:nwhileRelHead(i,1)piRelHead(i,1)=RelHead(i,1)-pi;endendnoise=(rand(n,1)*eta)-eta/

31、2;newHeading(:,1)=heading(:,1)+RelHead(:,1)./sum(neighbours,2)+noise(:,1);newPos(:,1)=pos(:,1)+cos(newHeading(:,1)*speed;newPos(:,2)=pos(:,2)+sin(newHeading(:,1)*speed;newPos(:,3)=pos(:,3)+tan(newHeading(:,1)*speed;fork=1:nwhilenewPos(k,1)=fieldsizenewPos(k,1)=newPos(k,1)-fieldsize;endwhilenewPos(k,

32、2)=fieldsizenewPos(k,2)=newPos(k,2)-fieldsize;endwhilenewPos(k,3)=fieldsizenewPos(k,3)=newPos(k,3)-fieldsize;endendheading=newHeading;pos=newPos;meanHeading(time,1)=mean(heading,1);deviationMean(time,1)=mean(abs(meanHeading(time,1)-heading(:,1);scatter3(pos(:,1),pos(:,2),pos(:,3),*r);axis(0fieldsize

33、0fieldsize0fieldsize);xlabel(time);ylabel(Y);zlabel(Z);M(time)=getframe;endrunDev(:,run)=deviationMean(:,1);end%plot(runDev)movie3avi(M,boidtest.avi,quality,100);附录3functionboid1(n,eta)%generatesboidsinafieldn=input(inputn);eta=input(inputeta);range=0.2;speed=0.03;TIME=300;density=10;RUNS=1;m=3;fiel

34、dsize=sqrt(n/density);pos2=fieldsize*rand(m,2);heading2=2*pi*rand(m,2);elposX2=zeros(m);relposY2=zeros(m);neighbours2=zeros(m);RelHead2=zeros(m,1);B=zeros(m,1);newHeading2=zeros(m,1);newPos2=zeros(m,2);meanHeading2=zeros(TIME,1);deviationMean2=zeros(m,TIME);forrun=1:RUNSfortime=1:TIMEfori=1:mforj=1:

35、mrelposX2(i,j)=abs(pos2(i,1)-pos2(j,1);relposY2(i,j)=abs(pos2(i,2)-pos2(j,2);if(sqrt(relposX2(i,j)人2+relposY2(i,j)人2)piRelHead2(i,1)=RelHead2(i,1)-pi;endendnoise=(rand(m,1)*eta)-eta/2;newHeading2(:,1)=heading2(:,1)+RelHead2(:,1)./sum(neighbours2,2)+noise(:,1);newPos2(:,1)=pos2(:,1)+cos(newHeading2(:

36、,1)*speed;newPos2(:,2)=pos2(:,2)+sin(newHeading2(:,1)*speed;fork=1:mwhilenewPos2(k,1)=fieldsizenewPos2(k,1)=newPos2(k,1)-fieldsize;endwhilenewPos2(k,2)=fieldsizenewPos2(k,2)=newPos2(k,2)-fieldsize;endendheading2=newHeading2;pos2=newPos2;meanHeading2(time,1)=mean(heading2,1);deviationMean2(time,1)=me

37、an(abs(meanHeading2(time,1)-heading2(:,1);fieldsize=sqrt(n/density);pos1=fieldsize*rand(n,2);heading1=2*pi*rand(n,1);relposX1=zeros(n);relposY1=zeros(n);neighbours1=zeros(n);RelHead1=zeros(n,1);B=zeros(n,1);newHeading1=zeros(n,1);newPos1=zeros(n,2);meanHeading1=zeros(TIME,1);time=1:TIME?e?2?fori=1:n

38、forj=1:nfors=1:mrelposX1(i,j)=abs(pos1(i,1)-pos1(j,1);relposY1(i,j)=abs(pos1(i,2)-pos1(j,2);relposX3(j,s)=abs(pos2(s,1)-pos1(j,1);relposY3(j,s)=abs(pos2(s,2)-pos1(j,2);if(sqrt(relposX1(i,j)人2+relposY1(i,j)人2)=range|sqrt(relposX3(j,s)人2+relposY3(j,s)人2)=2*range)neighbours1(i,j)=1;endendendendfori=1:n

39、forj=1:nif(neighbours1(i,j)=1)relheading1(i,j)=heading1(j,1)-heading1(i,1);endendendRelHead1=sum(relheading1,2);fori=1:nwhileRelHead1(i,1)piRelHead1(i,1)=RelHead1(i,1)-pi;endendnoise=(rand(n,1)*eta)-eta/2;newHeading1(:,1)=heading1(:,1)+RelHead1(:,1)./sum(neighbours1,2)+noise(:,1);newPos1(:,1)=pos1(:

40、,1)+cos(newHeading1(:,1)*speed;newPos1(:,2)=pos1(:,2)+sin(newHeading1(:,1)*speed;fork=1:nwhilenewPos1(k,1)0newPos1(k,1)=newPos1(k,1)+fieldsize;newPos1(k,1)=newPos1(k,1)-fieldsize;endwhilenewPos1(k,2)=fieldsizenewPos1(k,2)=newPos1(k,2)-fieldsize;endendheading1=newHeading1;pos1=newPos1;meanHeading1(ti

41、me,1)=mean(heading1,1);deviationMean1(time,1)=mean(abs(meanHeading1(time,1)-heading1(:,1);scatter(pos2(:,1),pos2(:,2),scatter(pos1(:,1),pos1(:,2),og );holdxr );holdonoffaxis(0 fieldsize 0 fieldsize); xlabel(time );M(time) = getframe;endrunDev(:,run) = deviationMean1(:,1);runDev(:,run) = deviationMea

42、n2(:,1);end%plot(runDev)movie2avi(M, boidtest.aviquality ,100);附录 4function pso F = pso_2D()pop_size = 30;part_size = 2;gbest = zeros(1,part_size+1);max_gen = 60;region=zeros(part_size,2);region=-3,3;-3,3;pop_size 种群大小part_size 粒子大小，% gbest 当前搜索到的最小的值% max_gen 最大迭代次数% 设定搜索空间范围rand( state ,sum(100*cl

43、ock);% 重置随机数发生器状态arr_present = ini_pos(pop_size,part_size);% present 当前位置 ， 随机初始化 ， rand() 的范围为 01v=ini_v(pop_size,part_size);% 初始化当前速度pbest = zeros(pop_size,part_size+1);% pbest 粒子以前搜索到的最优值，最后一列包括这些值的适应度% w_max 权系数最大值w_max=0.7298;w_min=0.4;v_max=2;%*最大速度，为粒子的范围宽度 TOC o 1-5 h z cl=2;%学习因子1c2=2;%学习因子

44、2best_record=zeros(1,max_gen);%best_record记录最好的粒子的适应度。arr_present(:，end)=ini_fit(arr_present，pop_size，part_size);。pbest=arr_present;值best_valuebest_index=min(arr_present(:,end);%初始化全局最优，即适应度为全局最小的值，根据需要也可以选取为最大值gbest=arr_present(best_index,:);x=-3:0.01:3;y=-3:0.01:3; TOC o 1-5 h z z=(x,y)3*(1-x).人2.*exp(-(x.人2)-(y+1).人2).-10*(x/5-x.A3-y.A5).*exp(-x.A2-y.A2).-1/3*exp(-(x+1).A2-y.A2);fori=1:max_gengridon;ezmesh(z),holdon,gridon,plot3(arr_present(:,1),arr_present(:,2),arr_present(:,3),*),text(-0.8,-1.2,6.5,食物);title(鸟群觅食);holdoff;drawnowF(i)=getframe;pause(0.01);w