基于豺优化算法(Dhole Optimization Algorithm, DOA)及三次样条的机器人路径规划,50个场景任意选择,完整MATLAB代码
一、豺优化算法
豺优化算法(Dhole Optimization Algorithm, DOA)是2025年提出的一类新型种群智能元启发式优化算法,其设计灵感源于豺(Dhole,犬科豺属哺乳动物)的自然种群特征与协同狩猎行为,核心通过模拟豺群的发声驱动自适应决策机制与动态群体协作狩猎模式,实现优化过程中全局探索与局部开发能力的精准平衡,经经典基准函数及复杂优化问题验证,该算法在收敛速率、求解精度及鲁棒性方面均表现出显著优势,为复杂高维优化问题的求解提供了新的有效途径。

参考文献:
[1]Mohammed, B.O., Aghdasi, H.S. & Salehpour, P. Dhole optimization algorithm: a new metaheuristic algorithm for solving optimization problems. Cluster Comput 28, 430 (2025). https://doi.org/10.1007/s10586-024-05005-1
1、算法灵感来源
DOA的设计灵感源于豺(Dhole,犬科哺乳动物)的自然生物学特征和狩猎行为,核心参考两点:
- 豺的种群特征:豺主要分布于中亚、南亚、东亚和东南亚,与犬属动物有遗传相似性,但形态上存在独特特征(上臼齿仅1个尖头、头骨轮廓凸起、无第三个下臼齿);其种群通常由5~20只成员组成,由一对首领(alpha)领导,具备高度的协作性。
- 豺的狩猎策略:豺通过声音交流实现群体协调,狩猎活动多发生在晨昏时段,可猎取比自身大/小的猎物(鹿、野猪等),狩猎过程分为搜索猎物、包围猎物、协同攻击三个阶段,且会根据猎物大小动态调整攻击策略,这一特征成为DOA位置更新模型的核心设计依据。
2、算法核心基本定义
DOA将优化问题的求解过程映射为豺群的狩猎过程,核心概念的数学定义如下:
1. 种群与候选解
每只豺代表优化问题的一个潜在候选解,其位置对应决策变量的取值;整个豺群为一个矩阵,如公式(1):
D=[D1⋮Dt⋮DN]N⋅m=[d1,1⋯d1,j⋯d1,m⋮⋱⋮⋱⋮dt,1⋯dt,j⋯dt,m⋮⋱⋮⋱⋮dN,1⋯dN,j⋯dN,m]N⋅mD=\begin{bmatrix}D_1\\\vdots\\D_t\\\vdots\\D_N\end{bmatrix}_{N\cdot m}=\begin{bmatrix}d_{1,1}&\cdots&d_{1,j}&\cdots&d_{1,m}\\\vdots&\ddots&\vdots&\ddots&\vdots\\d_{t,1}&\cdots&d_{t,j}&\cdots&d_{t,m}\\\vdots&\ddots&\vdots&\ddots&\vdots\\d_{N,1}&\cdots&d_{N,j}&\cdots&d_{N,m}\end{bmatrix}_{N\cdot m}D=D1⋮Dt⋮DNN⋅m=d1,1⋮dt,1⋮dN,1⋯⋱⋯⋱⋯d1,j⋮dt,j⋮dN,j⋯⋱⋯⋱⋯d1,m⋮dt,m⋮dN,mN⋅m
其中:NNN为豺群总数(种群规模),mmm为决策变量数量,DiD_iDi为第iii只豺(候选解),di,jd_{i,j}di,j为第iii只豺的第jjj个决策变量值。
2. 适应度值
每只豺的适应度值通过目标函数向量评估,反映候选解的优劣,如公式(2):
[F1⋮Ft⋮FN]N⋅1=[F(D1)⋮F(Dt)⋮F(DN)]N⋅1\begin{bmatrix}F_1\\\vdots\\F_t\\\vdots\\F_N\end{bmatrix}_{N\cdot1}=\begin{bmatrix}F(D_1)\\\vdots\\F(D_t)\\\vdots\\F(D_N)\end{bmatrix}_{N\cdot1}F1⋮Ft⋮FNN⋅1=F(D1)⋮F(Dt)⋮F(DN)N⋅1
其中:FFF为目标函数向量,F(Di)F(D_i)F(Di)为第iii只豺的目标函数值,算法在迭代中持续更新最优解和最差解。
3. 群体成员数(PMN)
模拟豺群的实际规模特征,PMN随机生成于5~20之间,如公式(3):
PMN=round(rand×15+5)PMN=\text{round}(\text{rand}×15+5)PMN=round(rand×15+5)
其中:rand\text{rand}rand为[0,1]内的随机数,round(⋅)\text{round}(\cdot)round(⋅)为四舍五入取整函数。
4. 最佳狩猎时间(ps)
模拟豺晨昏狩猎的行为特征,结合种群规模、环境因子建模,如公式(4):
ps=(C11+exp(−k(PWN−μ)))2×EFps=\left( \frac{C1}{1 + \exp(-k (\text{PWN} - \mu))} \right)^2 × \text{EF}ps=(1+exp(−k(PWN−μ))C1)2×EF
其中:μ\muμ为最适群体成员数,EF\text{EF}EF为环境因子(0~1),kkk为狩猎效率系数,C1C1C1为猎物大小控制系数。
3、算法核心数学模型
DOA的核心是豺群位置更新机制,根据豺的狩猎行为分为**搜索阶段(探索)、包围阶段(开发)、攻击阶段(开发)三个阶段,通过发声信号(Vocalization,[0,1]随机数)和群体成员数(PMN)**触发不同阶段的位置更新规则,实现探索与开发的动态平衡。

阶段1:搜索阶段(全局探索)
触发条件:发声信号Vocalization<0.5\text{Vocalization}<0.5Vocalization<0.5 且 群体成员数PMN<10PMN<10PMN<10
该阶段模拟豺群分散搜索猎物的行为,核心是扩大搜索范围,增强算法的全局探索能力,避免陷入局部最优。
- 确定猎物位置:以当前种群最优位置和历史迭代最优位置的均值为目标猎物位置,如公式(5):
prey=prvNmax+prvSylmax2\text{prey} = \frac{\text{prv}_{Nmax} + \text{prv}_{Sylmax}}{2}prey=2prvNmax+prvSylmax
其中:prvNmax\text{prv}_{Nmax}prvNmax为当前种群最优位置,prvSylmax\text{prv}_{Sylmax}prvSylmax为另一迭代周期的最优位置。 - 位置更新:豺群向猎物位置靠近并随机扩大搜索范围,如公式(6):
Dijs−1=dij+C2×rand×(prvij−dij)D_{ij}^{s-1} = d_{ij} + C_2 × \text{rand} × (\text{prv}_{ij} - d_{ij})Dijs−1=dij+C2×rand×(prvij−dij)
其中:C2C_2C2为递减系数,随迭代次数增加而减小,实现探索能力的动态调整,如公式(7):
C2=1−1TC_2=1−\frac{1}{T}C2=1−T1
TTT为算法最大迭代次数,prvij\text{prv}_{ij}prvij为猎物位置的第jjj个决策变量值。
阶段2:包围阶段(局部开发)
触发条件:发声信号Vocalization<0.5\text{Vocalization}<0.5Vocalization<0.5 且 群体成员数PMN≥10PMN≥10PMN≥10
该阶段模拟豺群发现猎物后形成协作包围圈的行为,核心是收缩搜索范围,对猎物周边区域进行局部开发,逐步逼近最优解。
- 位置更新:豺的位置受其他随机个体的影响,实现群体协作包围,如公式(8):
Dijs−1=dij−dzj+prvijD_{ij}^{s-1} = d_{ij} - d_{zj} + \text{prv}_{ij}Dijs−1=dij−dzj+prvij - 随机个体索引:zzz为随机选取的豺群个体索引,且保证i≠zi≠zi=z,如公式(9):
z=1+round(rand×(N−1))z=1+\text{round}(\text{rand}×(N−1))z=1+round(rand×(N−1))
该阶段存在个体竞争,通过不同豺的位置相互影响,实现局部区域的精细搜索。
阶段3:攻击阶段(局部开发)
触发条件:发声信号Vocalization≥0.5\text{Vocalization}≥0.5Vocalization≥0.5
该阶段模拟豺群由首领发起的协同攻击行为,核心是根据猎物大小动态调整攻击策略,实现对最优解的精准逼近,是算法局部开发的核心阶段。
- 猎物大小(S)定义:结合猎物因子和适应度值建模,如公式(10):
S=C2+rand×(HwmaHwcmpop)S = C_2 + \text{rand} × \left(\frac{Hwma}{Hwcm_{pop}}\right)S=C2+rand×(HwcmpopHwma)
其中:C2C_2C2为猎物因子(常数3),HwmaHwmaHwma为当前豺的适应度值,HwcmpopHwcm_{pop}Hwcmpop为猎物的适应度值。 - 大猎物攻击策略(S≥4.5,即1.5×C2):先削弱猎物,再连续攻击,对应公式(11)和(12):
- 猎物削弱:Wprop=exp(−1S)×prvNmaxW_{prop} = \exp\left(-\frac{1}{S}\right) × \text{prv}_{Nmax}Wprop=exp(−S1)×prvNmax
- 连续攻击位置更新:Dijs−1=dij+Wprop+p×(ln(2+x+rand)−ln(2+x+rand)×Wprop+p)D_{ij}^{s-1} = d_{ij} + W_{prop} + p × \left( \ln(2 + x + \text{rand}) - \ln(2 + x + \text{rand}) × W_{prop} + p \right)Dijs−1=dij+Wprop+p×(ln(2+x+rand)−ln(2+x+rand)×Wprop+p)
- 小/虚弱猎物攻击策略(S<4.5):直接击杀猎物,位置更新如公式(13):
Dijs−1=(dij−prvSylmax)×p+p×rand×dijD_{ij}^{s-1} = (d_{ij} - \text{prv}_{Sylmax}) × p + p × \text{rand} × d_{ij}Dijs−1=(dij−prvSylmax)×p+p×rand×dij
其中:ppp为攻击强度系数,xxx为猎物状态参数,均为算法预设超参数。
4、算法流程

DOA的执行流程为参数初始化→迭代狩猎→更新最优解→收敛终止,具体步骤如下:
- 输入优化算法基础信息:包括决策变量上下界、目标函数、优化目标(最小/最大化);
- 设置核心参数:豺群规模NNN、最大迭代次数TTT;
- 初始化豺群种群:随机生成NNN个候选解,构建豺群矩阵DDD;
- 计算初始适应度值:通过目标函数评估所有豺的适应度,得到初始最优解;
- 迭代狩猎(t=1到T):
- 步骤1:通过公式(3)计算群体成员数PMNPMNPMN;
- 步骤2:通过公式(5)确定目标猎物位置;
- 步骤3:遍历每只豺(i=1到N),生成发声信号Vocalization=rand\text{Vocalization}=\text{rand}Vocalization=rand,根据触发条件执行搜索/包围/攻击阶段的位置更新;
- 步骤4:更新所有豺的适应度值,保留当前迭代的最优解;
- 步骤5:迭代次数t=t+1t=t+1t=t+1,直至达到最大迭代次数TTT;
- 输出结果:算法收敛,输出全局最优解和对应的适应度值。
二、路径规划问题的适应度-奖励函数
路径从起点 (xs,ys)(x_s, y_s)(xs,ys) 开始,到终点 (xt,yt)(x_t, y_t)(xt,yt) 结束,并包括 nnn 个航路点。路径上的每个点表示为 (xi,yi)(x_i, y_i)(xi,yi),其中 iii 从 1(起点)到 nnn(终点)。目标是求解约束优化问题,以找到起点和终点之间的最短安全路径。最短路径最小化连续点之间的总欧几里得距离,这是欧几里得空间中最直接的路线。
除了最小化距离外,路径必须满足两个约束条件。首先,所有路径点必须保持在边界内,即 xmin≤xi≤xmaxx_{\text{min}} \le x_i \le x_{\text{max}}xmin≤xi≤xmax 和 ymin≤yi≤ymaxy_{\text{min}} \le y_i \le y_{\text{max}}ymin≤yi≤ymax。其次,路径必须避开障碍物,这通过条件 g(xi,yi,b)≥0g(x_i, y_i, b) \ge 0g(xi,yi,b)≥0 强制执行,其中 g(xi,yi,b)g(x_i, y_i, b)g(xi,yi,b) 是评估与第 bbb 个障碍物碰撞条件的一般约束函数,并确保不发生碰撞。优化问题如下所示:
minimizef(x,y)=∑i=1n−1(xi+1−xi)2+(yi+1−yi)2 \begin{aligned} \text{minimize} \quad & f(x, y) = \sum_{i=1}^{n-1} \sqrt{(x_{i+1} - x_i)^2 + (y_{i+1} - y_i)^2} \end{aligned} minimizef(x,y)=i=1∑n−1(xi+1−xi)2+(yi+1−yi)2
subject toxmin≤xi≤xmax,∀i∈{1,…,n}ymin≤yi≤ymax,∀i∈{1,…,n}g(xi,yi,b)≥0,∀i∈{1,…,n},∀b∈Obstacles \begin{aligned} \text{subject to} \quad & x_{\text{min}} \le x_i \le x_{\text{max}}, \quad \forall i \in \{1, \ldots, n\} \\ & y_{\text{min}} \le y_i \le y_{\text{max}}, \quad \forall i \in \{1, \ldots, n\} \\ & g(x_i, y_i, b) \ge 0, \quad \forall i \in \{1, \ldots, n\}, \quad \forall b \in \text{Obstacles} \end{aligned} subject toxmin≤xi≤xmax,∀i∈{1,…,n}ymin≤yi≤ymax,∀i∈{1,…,n}g(xi,yi,b)≥0,∀i∈{1,…,n},∀b∈Obstacles
适应度(成本)函数通过考虑路径的总距离和与障碍物碰撞产生的任何惩罚来评估给定路径的质量。惩罚包括两个部分:点惩罚(路径点与障碍物之间的碰撞)和线段惩罚(连续航路点之间的线段与障碍物之间的碰撞)。此公式确保路径既短又安全,有效指导优化过程。
线段惩罚对于识别点惩罚无法检测到的问题至关重要。即使所有航路点都位于障碍物外部,个别路径线段仍可能与障碍物相交;因此,线段惩罚对于捕捉这些点惩罚会遗漏的违规行为至关重要。虽然线段惩罚足以检测障碍物碰撞,但也会执行点检查以评估或惩罚各个航路点。任何位于障碍物内的航路点都表明严重违规,因为该航路点无法到达,从而使路径不可行。
相比之下,如果所有航路点都位于障碍物外部,但连接它们的某些线段与障碍物相交,则该路径仍可被视为候选路径。这些路径通常只需要在后续迭代中进行 minor 优化即可消除交点。这种惩罚组合确保对不可达航路点的路径施加更严厉的惩罚,而对 minor 问题的路径仍可进行优化。此外,点惩罚允许快速初步检查,帮助跳过不必要的线段惩罚计算,当路径因航路点位于障碍物内而明显无效时,这降低了计算复杂性。
成本函数由两个主要部分组成:路径长度 LLL 和总惩罚 pTp^TpT,如下所示。路径长度 LLL 表示路径上所有点之间的欧几里得距离,定义如下:
L=∑i=1n−1(xi+1−xi)2+(yi+1−yi)2 L = \sum_{i=1}^{n-1} \sqrt{(x_{i+1} - x_i)^2 + (y_{i+1} - y_i)^2} L=i=1∑n−1(xi+1−xi)2+(yi+1−yi)2
总惩罚 pTp^TpT 考虑碰撞情况,并按权重 β\betaβ 缩放,权重设为 100 以优先考虑安全性,对不安全路径施加严厉惩罚。对于无障碍物的路径,总成本等于路径长度 LLL,代表可实现的最小成本。

总惩罚 pTp^TpT 表示所有障碍物的平均惩罚,包括基于点和线段的惩罚,如下所示,其中 BBB 是障碍物的总数。对于每个障碍物 bbb,惩罚包括两个项:∑i=1npib\sum_{i=1}^{n} p_i^b∑i=1npib,评估涉及 nnn 个路径点 (xi,yi)(x_i, y_i)(xi,yi) 的碰撞惩罚;以及 ∑j=1n−1pline,jb\sum_{j=1}^{n-1} p_{\text{line}, j}^b∑j=1n−1pline,jb,评估涉及 n−1n-1n−1 个线段 Lj=[(xj,yj),(xj+1,yj+1)]L_j = [(x_j, y_j), (x_{j+1}, y_{j+1})]Lj=[(xj,yj),(xj+1,yj+1)] 的碰撞惩罚。惩罚通过除以 BBB(障碍物总数)进行归一化,以确保惩罚反映碰撞的总体分布和严重程度。
cost=L×(1+βpT) \text{cost} = L \times (1 + \beta p^T) cost=L×(1+βpT)
对于基于网格的障碍物,两个组件 pibp_i^bpib 和 pline,jbp_{\text{line}, j}^bpline,jb 使用网格表示进行计算。搜索空间表示为基于网格的环境,网格由其原点 (xg,yg)(x_g, y_g)(xg,yg) 和分辨率 (xres,yres)(x_{\text{res}}, y_{\text{res}})(xres,yres) 定义。x 坐标从 xminx_{\text{min}}xmin 到 xmaxx_{\text{max}}xmax 以 xresx_{\text{res}}xres 为步长进行离散化,y 坐标同样从 yminy_{\text{min}}ymin 到 ymaxy_{\text{max}}ymax 以 yresy_{\text{res}}yres 为步长进行离散化。路径上的每个点 (xi,yi)(x_i, y_i)(xi,yi) 被映射到网格中的一个对应单元格,索引为 (xci,yci)(xc_i, yc_i)(xci,yci),计算方法如下:
xci=round(xi−xogxres),yci=round(yi−yogyres) xc_i = \text{round}\left( \frac{x_i - x_{og}}{x_{\text{res}}} \right), \quad yc_i = \text{round}\left( \frac{y_i - y_{og}}{y_{\text{res}}} \right) xci=round(xresxi−xog),yci=round(yresyi−yog)
网格中的每个单元格被分配一个成本 M(xci,yci)M(xc_i, yc_i)M(xci,yci),其中 0 表示空闲单元格,100 表示完全被占据的单元格,如 SLAM 等映射算法所确定。对于基于点的惩罚 pip_ipi,对应于点 (xi,yi)(x_i, y_i)(xi,yi) 的单元格的成本 M(xci,yci)M(xc_i, yc_i)M(xci,yci) 通过除以 100 进行归一化,得到范围为 [0,1][0, 1][0,1] 的惩罚,如下所示:
pi=M(xci,yci)100 p_i = \frac{M(xc_i, yc_i)}{100} pi=100M(xci,yci)
线段惩罚 pline,jp_{\text{line}, j}pline,j 通过确定由两个连续航路点连接的线段 LjL_jLj 相交的所有网格单元格来计算。使用 Bresenham 线算法确定这些单元格,其中 njn_jnj 是这些单元格的数量。


上图展示了线段与网格单元格之间的相交情况。然后通过对这些相交单元格的归一化占用成本进行平均来聚合惩罚,如下所示:
pline,j=1nj∑k=1njM(xck,yck)100 p_{\text{line}, j} = \frac{1}{n_j} \sum_{k=1}^{n_j} \frac{M(xc_k, yc_k)}{100} pline,j=nj1k=1∑nj100M(xck,yck)
方程 被修改为在基于网格的环境中计算总体路径惩罚 pTp^TpT,如下所示,其中成本直接从网格中获得。
pT=∑i=1npi+∑j=1n−1pline,j p^T = \sum_{i=1}^{n} p_i + \sum_{j=1}^{n-1} p_{\text{line}, j} pT=i=1∑npi+j=1∑n−1pline,j
成本函数支持多种障碍物类型,重点在于基于网格的障碍物,关于与圆形和矩形障碍物相关的 pibp_i^bpib 和 pline,jbp_{\text{line}, j}^bpline,jb 的详细公式以及这些类型之间的转换机制详细参考文献。
参考文献:
[1]Reda, M., Onsy, A., Haikal, A.Y. et al. A novel reinforcement learning-based multi-operator differential evolution with cubic spline for the path planning problem. Artif Intell Rev 58, 142 (2025). https://doi.org/10.1007/s10462-025-11129-6
三、部分MATLAB及结果
50个椭圆或矩阵障碍物场景任意选择
clear;
clc;
close all;
warning off all;
AlgorithmName='';%算法名称
N=40;%选择地图-可以选择1-50; 1-16 椭圆障碍物,17-50 矩形障碍物
SearchAgents_no=50; %种群大小
Max_iteration=200; %最大迭代次数
[lb,ub,dim,fobj,Fobj,model] = Get_F(N);%获取参数
[Best_score,Best_pos,cg_curve]=(SearchAgents_no,Max_iteration,lb,ub,dim,fobj);%算法求解
GlobalBest=Fobj(Best_pos);%计算最优个体的详细路径














四、完整MATLAB见下方名片
DAMO开发者矩阵,由阿里巴巴达摩院和中国互联网协会联合发起,致力于探讨最前沿的技术趋势与应用成果,搭建高质量的交流与分享平台,推动技术创新与产业应用链接,围绕“人工智能与新型计算”构建开放共享的开发者生态。
更多推荐

所有评论(0)