群体动画一直是计算机图形学、交通仿真、虚拟现实、社会公共安全和视觉特效等领域的研究热点。其主要研究内容包括建模方法、路径规划和碰撞避免。如何模拟真实多样的人群行为成为人们关注的重点，然而人群的群组行为及其交互仿真仍是群体动画中一个亟待解决的问题。

现阶段群体动画研究方法主要分为两类：连续模型和离散模型。连续模型考虑大规模群体运动整体性，其研究内容是将群体视为连续整体，利用场论和偏微分方程描述群体流密度和平均速度的关系。离散模型视个体为拥有决策能力的自主智能体，主要研究个体与周围环境的交互作用。然而，连续模型常忽略局部运动细节，难以表现逼真的人群运动；离散模型以个体为研究对象，在模拟大规模人群运动时效率较低。

因此，本文提出一种基于群体动力学的群组行为仿真算法。该算法将全局规划和局部控制相结合，模拟组内跟随和组间避让行为。Lemercier等人^[1]定义了群组的表示方法，但未提供全局路径规划方案，并且需要利用额外的碰撞避免算法实现群组避让行为。本文根据组内和组间关系计算跟随加速度，并提出群组的凸包表示方法，通过构建局部势能场，结合全局路径规划模拟多样性群组行为。主要贡献如下：

1) 结合全局规划提出一种基于群体动力学的群组行为仿真算法。

2) 引入群体流约束条件，限制人群运动速率和人群密度的关系，避免行人陷入拥塞路径。

3) 提出群组凸包表示方法和局部势能场，实现群组的动态避让。

在宏观建模方法中，Treuille等人^[2]提出的结合碰撞避免的动态势场方法。Zhao等人^[3]提出一种基于场的人群仿真方法并使用GPU进行碰撞检测。Hughes等人^[4]将人群类比为流，通过求解随时间变化的偏微分方程组描述不同约束条件下的人群运动。Narain等人^[5]提出群体流单方面不可压缩约束(Unilateral Incompressibility Constraint，UIC)，用于模拟高密度情形下的人群运动。Xu等人^[6]通过求解N-S方程建立个体运动的速度场，利用层次结构平衡绘制的效率和质量。Min等人^[7]提出一种控制大规模群体的控制粒子方法，建立稳定状态控制场引导群体随着控制粒子运动。

在微观建模方法中，Reynolds^[8]建立了“boid”系统，该系统通过分离、列队和凝集规则模拟鸟类和鱼群的运动。Helbing^[9-10]提出的经典社会力模型通过分析个体所受吸引力和排斥力等模拟紧急状况下的人群逃生。Fu等人^[11]提出一种基于元胞自动机和自主智能体的行为模型，将每个离散元胞空间看作单独智能体，通过感知模型和决策模型反映个体的人格特质和心理特征。Kim等人^[12]考虑个体行走速度，提出一种基于速度的人群物理交互仿真算法。Bera等人^[13]提出一种数据驱动的群体动画模拟方法。

近年来，一部分学者开始研究群组行为，从介于宏观和微观的中层角度出发，通过分析人群运动特征来模拟多样的群组行为。Lemercier等人^[1]赋予行人感知环境的能力，并允许其根据环境调整自身行为，然后结合RVO2^[14]碰撞避免算法提出一种模拟跟随和群组避让的群体仿真算法。Bruneau等人^[15]将能量最低原则(PME)用于不同大小和密度的群组，模拟真实行人绕过小型密集群组而穿过大型稀疏群组的现象。He等人^[16]根据最省力原则(least effect principle)提出一种连续群组导航方法。

国内部分学者考虑路径规划的群组仿真，如Liu等人^[17]提出融入关系分组的人群运动仿真算法，考虑个体间距离和关系对其运动的影响。采用改进的K-Medoids算法将人群分组，并通过社会力模型实现路径规划和碰撞避免。该方法能模拟出更加贴近真实环境下的人群疏散行为，然而该模型没有考虑群体密度和场景出口总宽度等因素对疏散效率的影响，同时利用社会力模型解决组内个体的运动效率较低。Zhang等人^[18]提出一种小型行人群组的局部行为模型，利用环境约束构建动态群组结构。引入群体凝聚力因素和反应灵活性，通过改进ORCA (optimal reciprocal collision avoidance)碰撞避免算法^[19]来模拟行人交互的多样性。该方法虽然能实现行人间的碰撞避免，但避让决策容易陷入局部最优状态导致行人出现不自然的行为。同时该模型没有考虑人群密度对行人运动的影响，在模拟密度较高的人群运动中存在一定的局限性。因此，本文提出一种基于群体动力学的群组行为仿真算法，通过势能场方法为个体规划全局路径，结合行人运动特征对较大规模的群组运动进行仿真。

研究群组行为的算法常将固定大小的群组视为特殊个体来实现碰撞避免，适用于模拟规模较小的人群模拟，为实现大规模人群运动中的群组行为，本文提出一种基于群体动力学的群组行为仿真算法。首先，根据人群分布和环境信息构建全局动态势能场，计算个体朝目标运动的初始速度$ {\mathit{\boldsymbol{v}}_0} $。其次，在局部运动中引入群体动力学模型，模拟组内跟随和组间避让行为。部分研究群组运动的文献[1, 15]需要引入额外的碰撞避免算法控制个体局部运动，降低了算法效率。针对这一问题，本文提出群组的凸包表示方法，建立局部势能场$ {\varphi _{{\rm{local}}}} $并结合全局规划对多样性群组进行仿真。算法框架如图 1所示。

本文提出的算法能模拟多样性群组行为，引入群体流约束减少个体在高密度区域发生碰撞的次数。提出群组的凸包表示方法，将其与全局规划结合，用于模拟群组的避让行为。

本文采用与文献[2]类似的方法构建动态势能场，假设存在势函数$ \phi :{\mathit{\boldsymbol{R}}^2} \to \mathit{\boldsymbol{R}} $，使得目标位置处$ \phi $=0，其他位置满足程函方程

$ \left\| {\nabla \phi \left( \mathit{\boldsymbol{x}} \right)} \right\| = C $

(1)

式中， $C$ 表示单位成本。$ {\phi \left( \mathit{\boldsymbol{x}} \right)} $为$ \mathit{\boldsymbol{x}} $处的势能场大小，$ \nabla $为梯度算子。

引入格林希尔茨－线性模型避免个体陷入高密度区域， $t$ 时刻个体受群体流约束的速率为

$ f\left( {\rho, t} \right) = A - B\rho $

(2)

式中， $A$ =1.4 m/s为最大行走速率， $B$ =0.25 m³/s为衰减系数，$ \rho $表示群体密度。

如图 2所示，通过文献[2]的方法计算个体对周围网格g₁、g₂、g₃和g₄的密度贡献， $r$ 为个体的感知半径。在群体流约束条件下，个体通过感知周围网格的密度值调节自身行走速度，使其在高密度区域的行走速度较慢，反之行走速度较快。同时，个体根据所处位置的局部势能场信息调节自身运动，减少局部最优造成人群拥挤的现象。局部势能场的实现方法在第3.4.2节讨论。

$t$ 时刻个体受密度影响的期望速度$ {\mathit{\boldsymbol{v}}_0} $为

$ {\mathit{\boldsymbol{v}}_0} = - f\left( {\rho, t} \right)\frac{{\nabla \phi \left( \mathit{\boldsymbol{x}} \right)}}{{\left\| {\nabla \phi \left( \mathit{\boldsymbol{x}} \right)} \right\|}} $

(3)

根据式(1)－式(3)，利用快速行进法^[20]求解全局动态势能场的具体流程如下：

1) 离散化：根据个体所处位置将其转化为单位密度场并计算群体密度$ \rho $。

2) 速率场转化：基于群体密度$ \rho $，根据式(2)求解最大速率场 $f$ 。

3) 单位成本域：根据路径长度、行走时间和不适函数构建单位成本域 $C$ 。

4) 将目标网格存入Known_List，其他网格存入Unknown_List。

5) 循环以下步骤：

(1) 对Grid中所有网格Gridi( $i$ )，若其在已知列表Known_List中，则继续执行；否则跳出循环。

(2) 将Gridi( $i$ )所有非已知邻居节点存入候选列表Candidate_List。

6) 对Candidate_List中所有网格，计算其势能大小。

7) 找出Candidate_List中势能值最小的网格，将其存入Known_List。返回步骤5)，直到所有网格遍历完成。

本文对群组的定义为： $t$ 时刻拥有相似位置和速度的个体的集合。现实生活中行人拥有不同的目标位置和行走速度，对不适影响的感知能力因人而异。为模拟人群运动过程中动态形成群组的现象，本文根据个体位置和速度的差异将人群划分为不同的群组，若个体 $A$ 和个体 $B$ 属于同一群组，需要满足条件

$ \left( {\left\| {{\mathit{\boldsymbol{p}}_B} - {\mathit{\boldsymbol{p}}_A}} \right\| < d} \right) \wedge \left( {\left\| {{\mathit{\boldsymbol{v}}_B} - {\mathit{\boldsymbol{v}}_A}} \right\| < {\lambda _v}} \right) $

(4)

即组内个体的相对位置大小不超过一个关于距离的阈值 $d$ ，相对速度大小不超过一个关于速度的阈值$ {\lambda _v} $。其中，$ {\mathit{\boldsymbol{p}}_A} $和$ {\mathit{\boldsymbol{p}}_B} $分别表示个体 $A$ 和个体 $B$ 的位置，$ {\mathit{\boldsymbol{v}}_A} $和$ {\mathit{\boldsymbol{v}}_B} $分别表示个体 $A$ 和个体 $B$ 的速度。如图 3所示，行人在运动过程中动态形成群组，并根据组内跟随约束和组间避让规则实现群组行为仿真。

为了模拟现实生活中行人的跟随行为，文献[1]采用“Car-following”模型计算个体跟随加速度。“Car-following”模型是交通仿真领域用于模拟同一车道中车辆跟随行为的模型。本文采用类似方法计算跟随者的跟随加速度，即

$ {\mathit{\boldsymbol{a}}_{{\rm{follow}}}}\left( t \right) = W \times \frac{{\Delta {\mathit{\boldsymbol{v}}_x}\left( {t - \tau } \right)}}{{{{\left( {\Delta {\mathit{\boldsymbol{p}}_x}} \right)}^{\lambda \left( t \right)}}}} $

(5)

式中， $W$ > 0为比例因子，$ \lambda > 0 $是关于距离的调节因子，$ \tau $是个体反应时间。领导者与跟随者水平方向相对速度定义为$ \Delta {\mathit{\boldsymbol{v}}_x} = {\mathit{\boldsymbol{v}}_{jx}} - {\mathit{\boldsymbol{v}}_{ix}} $，$ {\mathit{\boldsymbol{v}}_{jx}} $和$ {\mathit{\boldsymbol{v}}_{ix}} $分别表示个体 $j$ 与个体 $i$ 在水平方向的行走速度。$ \Delta {\mathit{\boldsymbol{p}}_x} $表示两者的水平距离。如图 4所示，若个体 $i$ 跟随个体 $j$ 运动，需满足θ < θ₀，$ 0 < \Delta {\mathit{\boldsymbol{p}}_x} < {p_0} $和$ \Delta {\mathit{\boldsymbol{p}}_y} < {r_i} + {r_j} $。其中，θ₀、$ {p_0} $是预定义阈值，$ \Delta {\mathit{\boldsymbol{p}}_y} $表示两者在垂直水平方向的距离，$ {r_i} $和$ {r_j} $为个体半径。

结合全局规划计算的初始速度，计算个体的真实运动速度 $ \mathit{\boldsymbol{v}} $

$ \mathit{\boldsymbol{v}} = \left\{ \begin{array}{l} {\mathit{\boldsymbol{v}}_0} + {\mathit{\boldsymbol{a}}_{{\rm{follow}}}}\left( t \right) \times \tau \;\;{\rm{follow}}\\ {\mathit{\boldsymbol{v}}_0}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\rm{其他}} \end{array} \right. $

(6)

文献[1]采用RVO2碰撞避免算法实现群组避让行为，同时考虑了个体离群组较近或穿过群组的情况。该方法虽能实现群组避让行为，但额外的碰撞避免算法限制了其计算效率。因此，本文提出一种群组的凸包表示方法，并构建局部势能场，使群组能实现动态避让，且无需考虑额外碰撞避免。

对于群组$ \mathit{\boldsymbol{G}} = \left\{ {{i_1}, \cdots, {i_7}} \right\} $，横坐标最大和最小的个体一定是群组凸包上的个体，因此可采用分治法求解群组凸包，过程如图 5所示。

具体算法流程如下：

输入：群组个体$i$₁，$i$₂，$i$₃，$i$₄，$i$₅，$i$₆，$i$₇。

输出：群组凸包表示$i$₁，$i$₇，$i$₆，$i$₃，$i$₂。

算法1群组凸包算法(分治法)

1) 分别找出横坐标最小和最大的两个个体${i_m}$、${i_n}$(即图 5(a)中$i$₁、$i$₃);

2) 对上包，求距离直线${i_m}{i_n}$最远的个体$i$_max(即图 5(b)中$i$₇);

3) 连接直线$i$₁$i$₇、$i$₃$i$₇，把直线$i$₁$i$₇左侧和$i$₃$i$₇右侧视为上包;

4) 对直线$i$₁$i$₇、$i$₃$i$₇重复步骤2)3)，直到构成上包个体全部找到(如图 5(c)所示;

5) 对下包做类似操作(如图 5(d)所示)。

为有效结合全局规划和局部运动，利用群组的凸包结构构建局部势能场。本文设定群组产生的局部势能场不会对构成群组的个体产生影响，因此群组的凸包范围不会发生变化。将人群运动产生的群组存入不同链表，链表的长度即为构成某一群组的个体数。若有新的个体加入某一群组则将其加入链表，反之则将其从链表中移除。其他个体根据周围个体所属的链表标签感知群组，并通过计算该群组的局部势能场实现个体避开群组的行为。根据文献[21]中静态势能场的求解方法，将凸包边缘扩张δ距离表示局部势能场范围。该势能场产生垂直凸包边缘的排斥作用$ \mathit{\boldsymbol{V}} $，如图 6所示。

设群组凸包内部的场强大小为无穷大，势能场外部的场强大小为0，则任意一点的势能场$ {\varphi _{{\rm{local}}}} $计算公式为

$ {\varphi _{{\rm{local}}}} = \left\{ \begin{array}{l} \eta \left( {\delta - \left| \mathit{\boldsymbol{V}} \right|} \right)/\left( {\left| \mathit{\boldsymbol{V}} \right| \cdot \left| \mathit{\boldsymbol{V}} \right|} \right)\left| \mathit{\boldsymbol{V}} \right|\\ \;\;\;\;\;\;\;\;\;0 < \left| \mathit{\boldsymbol{V}} \right| \le \delta \;{\rm{且}}\;\delta {\rm{ > 0}}\\ 0\;\;\;\;\;\;\;\;\left| \mathit{\boldsymbol{V}} \right| > \delta \\ \infty \;\;\;\;\;\;\;\left| \mathit{\boldsymbol{V}} \right| = 0\;{\rm{或}}\;\delta = 0 \end{array} \right. $

(7)

式中，$ \eta $是比例因子，$ \mathit{\boldsymbol{V}} $表示从凸包边缘到该点的最短向量，$ \delta $为凸包边缘扩张的距离。

为了体现群组的运动趋势，将群组在仿真时间步 $t$ 内从$ {\mathit{\boldsymbol{P}}_1} $运动至$ {\mathit{\boldsymbol{P}}_2} $扫过的整个空间视为动态局部势能场。$ {\mathit{\boldsymbol{P}}_1} $表示群组当前位置，$ {\mathit{\boldsymbol{P}}_2} $表示群组下一时刻的期望位置。如图 7虚线部分所示。

当全局势能场值和局部势能场对个体的排斥影响夹角接近180°时，可能陷入局部最优而无法到达目标位置。式(7)仅通过基于距离的方法计算群组产生的局部势能场大小，不能完全解决局部最优问题。本文通过判断个体(或群组)与群组运动速度$ {\mathit{\boldsymbol{v}}_G} $的大小($ {\mathit{\boldsymbol{v}}_G} = \sum\limits_{i = 1}^n {{\mathit{\boldsymbol{v}}_i}/n} $，其中 $n$ 表示构成群组的个体数，$ {\mathit{\boldsymbol{v}}_i} $表示个体 $i$ 的速度)，为速度较小者添加一个与全局势能场相关的随机侧向力，避免其陷入局部最优状态。同时，通过初始条件判断个体是否到达目的地：当目标位置位于智能体感知范围内，认为个体到达目标位置，因此不再受势能场影响。

改进的群组凸包表示方法优点在于：1)受局部势能场影响，个体不会直接穿过群组；2)动态局部势能场能预测群组的运动趋势，使个体能提前对群组进行避让，降低碰撞发生的概率。

实验的硬件配置为Intel Core i7 3.4 GHz的处理器，8 GB的内存，NVIDIA GTX 860 M的显卡。在构建全局导航场时，为了验证提出算法的有效性和高效性，分别设计不同仿真场景并进行多次实验。实验过程包括3个阶段：全局势能场构建、局部势能场构建和局部运动控制阶段。

1) 全局势能场构建:

(1) 首先将环境划分为2维网格，并标记目标网格和障碍物网格。将目标网格存入已知列表，其余网格存入未知列表。

(2) 初始化人群并将其转化为单位密度场，计算其对网格的密度贡献及群体密度。

(3) 根据路径长度、行走时间和“不适影响”构建单位成本域。最后构建全局动态势能场为个体计算运动初始速度。

2) 局部势能场构建:

(1) 根据个体速度大小、运动方向和所处位置将人群划分成组。

(2) 通过分治法求解群组凸包边缘个体，将组内个体的平均位置视为群组凸包位置。

(3) 扩张凸包边缘并计算仿真步内群组扫过面积，构建局部势能场。

3) 局部运动控制:

(1) 结合全局势能场和局部势能场产生组间避让行为。

(2) 计算跟随加速度调整个体局部运动，实现组内跟随行为。

单向行人交通中常出现走停波(stop-and-go waves)，当个体穿过狭窄环境时会出现相互跟随行为。第1个实验模拟走廊环境下的单向行人流，如图 8所示。为表现走停波的传播方向。人群自左向右运动，从图 8(a)中可看出，红色部分个体处于减速状态，绿色部分个体处于加速状态，并且这种速度变化趋势沿行人运动的反方向传播，如图 8(b)－(d)所示。

人群速率大小满足关于群体密度的函数，即群体流基本图(fundamental diagram)，该图能反映人群运动过程中平均速率和平均密度的关系，如图 9所示。

1) 个体－群组避让实验。第1个场景模拟单个行人与群组相向而行，为说明个体受群组局部势能的影响，产生绕开群组的行为，本实验对比了提出算法和连续模型^[2]模拟的个体与群组的交互行为。本文算法的仿真结果如图 10(a所示，个体受局部势能场影响能动态避开群组。由于单个行人对群组影响较小，而群组势能场较大，因此行人选择直接绕开整个群组而非从群组中穿过。连续模型的仿真结果如图 10(b)所示，个体和群组在全局势能场的作用下相向而行，单个个体没有受到群组局部势能场影响而穿过群组：

2) 群组－群组避让实验。在模拟两个群组相向而行的实验中，受双方的动态局部势能场影响，两个群组预先感知到前方障碍并进行动态避让。传统碰撞避免算法常以个体为研究对象，难以模拟群组整体运动。本文提出的群组凸包表示方法能与全局势能场相结合，模拟出更加连续的群组行为。如图 11所示。

第3个实验将提出算法与传统连续模型^[2]进行对比，场景包括一组障碍物和4个群组。障碍物在场景中心构成狭窄通道，行人需无碰撞地通过障碍并到达对角线位置。仿真初始阶段设置群组目标并分别为每个群组构建全局势能场，如图 12 (a)所示，提出算法在仿真过程中能动态形成群组，个体受全局规划和局部势能场影响在通过狭窄区域时产生组内跟随和组间避让的行为。然而在相同环境中使用连续模型模拟人群运动时，因个体行走至障碍物区域使群体密度增高，个体相互影响发生拥塞现象，如图 12 (b)所示。

为了更好说明提出算法能有效避免个体陷入拥塞区域，分别记录了以上两个实验场景中行人发生碰撞的次数和仿真速率，如表 1所示。使用连续模型时，个体间发生碰撞的次数明显高于本文算法。

为了对比不同群体规模和网格精度对仿真效率的影响，分别设置了个体数为3 000、5 000、8 000及网格精度为80×80、40×40的仿真实验。实验结果表明本文算法在模拟较大规模的人群运动时，仍具有较高的仿真效率，如表 2所示。

在分辨率为40×40的网格上对8 000个个体运动仿真的运行时间为45.4 ms，约为22帧/s。但是，当网格分辨率增加到80×80时，因构建全局势能场遍历网格时间增大，仿真相同数量个体运动的时间增加至71.4 ms，帧速率下降至14帧/s。可见，提出算法在粗糙网格分辨率中仍能模拟较为平滑的人群运动。

本文算法能够逼真地模拟群组行为，在局部运动过程中，能产生组内跟随和组间避让行为。通过结合全局规划和局部势能场，提出算法能在群体密度较高时有效减少个体的碰撞次数。然而，该算法在网格分辨率为80×80，个体数为8 000时仿真效率较低，这是由于势函数的求解过程较为复杂，需要花费较多时间遍历网格。因此，提出算法还具有进一步提升的空间。

在未来的研究中，一方面，我们考虑采用并行求解的算法提高势函数求解效率。另一方面，为模拟更加真实的人群运动，考虑分析个体生理及心理特征，模拟具有社会性的群组行为，同时采用GPU加速算法绘制大规模人群运动。