计算流体力学在古生物学形态功能定量分析中的应用

郭颖， 王京盈， 张福成， 雷丽， 赵艳， 王岩， 王孝理

(1. 临沂大学 地质与古生物研究所， 临沂 276000； 2. 山东大学 能源与动力工程学院， 济南 250061)

几个世纪以来，破译古生物形态与功能之间的关系，一直是重建古生物生活方式的重要内容。这对确定古生物的进化过程、分类位置及其古环境信息都具有极为重要的意义。形态功能分析不仅可以在系统分类中解释形态在进化过程中的变化，也有助于对生物生活环境的判断，特别是对于已灭绝生物对环境的适应性分析。然而古生物类群往往缺乏关于某一功能的客观数据，难以找到与研究对象相似的现生生物作对比[1-2]，成为目前制约其形态功能分析的最大障碍。幸运的是，随着古生物学数字化革命的快速推进，越来越多的计算机技术实现了三维空间中古生物标本的可视化及定量分析[2-4]，这些技术不仅可以运用到尚不明确的化石分类单元中，同时也不需要类似的现生生物作对照，为古生物形态功能定量化研究提供新的发展方向和新的契机。

CFD是流体力学问题的定量化求解工具，可以通过图表直观呈现计算区域的流场特性，并能够明确给出运动和力学量的值[2, 5]。20世纪60年代至今，CFD已经在航空航天[6-7]、能源环境[8-9]以及生物医疗[10-11]等众多领域取得了令人瞩目的成就。近十几年以来，CFD技术逐渐被应用到古生物学的相关研究之中[12-15]，开创了古生物学研究中新的研究方向，成为古生物形态功能定量分析的重要工具，然而到目前为止相关应用案例的数量还非常有限，尚有广阔的发展空间和巨大的潜力有待挖掘。

1 形态功能分析发展简史

传统的形态功能分析是建立在古生物形态和功能研究的基础之上，伴随着形态测量学研究发展而形成的一种研究方法。20世纪初，生物学开始从描述和定性研究向定量化研究过渡，出现了特征测量和对不同特征采用平均值进行类群间比较等分析方法[16]。20世纪50年代，随着统计学方法的发展，产生了更加先进、严格的形态测量定量分析手段，形态特征的数量描述综合了类群内和类群间的形态变化的统计分析[17]。到20世纪60年代，生物统计学家开始应用多变量统计学工具描述类群内和类群间的形态差异，主要应用的变量为线性距离、数量、比例和角度等[16, 18]。20世纪80年代末至今，形态结构数值化和数据分析的方法上有了重要突破，逐渐发展成为较为成熟的几何形态测量法[18-19]。

几何形态测量方法虽然结合了多变量统计和形态特征的数值化，但是对于古生物的形态功能分析而言仍然存在难以克服的问题，比如，已灭绝的生物类群经常难以找到与研究对象相似的现生生物作对比，往往缺乏关于某一特定功能的客观数据做验证[1-2]。不过，随着古生物学数字化的普及和推进，CFD等大量计算机技术实现了三维空间中古生物标本的可视化及定量分析[2-4]，这些技术既可以运用到尚不明确的化石分类单元中，也不需要类似的现生生物作对照，为古生物形态功能定量化研究提供了思路和方法。

2 CFD的基本原理

一般的流动问题在数学上均可由Navier-Stokes(N-S)方程组进行描述[5]：

(1)

其中t为时间，ρ、U和T分别表示流体密度、内能和温度，V为流体速度矢量，F为流体所受的体积力，P和S分别为流体的应力张量和变形速度张量，k为热导率，Q表示辐射或者其他物理化学过程的贡献；“”表示空间梯度算符，“·”为矢量内积，“:”表示二阶张量的二次收缩内积。方程组(1)中的3个方程依次被称为连续性方程、运动方程和能量方程，分别反映了流动必须满足的质量、动量和能量守恒定律。

CFD的基本原理是在有限的离散时刻和空间节点(网格)上，采用一定的数值方法将流动问题对应的N-S方程组转化为代数方程组并进行求解[2, 5]。CFD求解过程一般在计算机上完成，常用的离散数值方法有有限差分法(Finite difference method)、有限体积法(Finite volume methods)和有限元法(Finite element)等。因CFD具有求解适应性强、边界条件设置自由和成本低等显著优点，目前已成为流体力学研究的主要手段之一。

3 CFD在形态功能分析中的流程

得益于计算机的迅速发展，CFD已经能够轻松求解任意二维流场和绝大部分三维流场[5-11]，同时ANSYS Fluent、Autodesk CFD、COMSOL Multiphysics等CFD商业软件也被广泛用于科学研究和工业设计，其中便包括古生物学形态功能分析领域。

3.1 几何模型重建

获取化石形态数据、重建几何模型是古生物CFD工作的关键基础。化石数据一般是从独立的计算机辅助设计软件中导入，形态较简单的也可以在CFD软件中直接构建[2-3, 20-21]。三维数据的获取方式多样，归纳起来分为表面形态和内部形态两大类。获取表面形态三维图像通常有激光扫描、原力显微镜成像、摄影测量等方法，而构建包含内部形态的三维图像则选择显微CT、纳米CT、核磁共振、同步辐射断层成像、聚焦离子束成像等技术手段[3-5]。不同仪器的原理、分辨率和图像格式不尽相同，实际操作中需根据情况而具体选择。

几何模型重建是对化石数据再加工的过程，首先要进行地层去压实校正、成岩变形复原、孔洞缝的修补和拼接以及合理的推测和假设，然后去除噪音及其他人为因素，构建出所需的三维数据体[4, 22]。当模拟大部分无脊椎动物化石表面或周围的流体情况时，内部特征相对不重要，模型构建一般只需要考虑化石表面信息[12-13, 21-22]；当开展脊椎动物运动方式的模拟时，则需要在骨骼化石信息的基础上[14, 23]，参照相近类群和现生生物特征，重构皮肤、肌肉或者羽毛等难以保存的软组织结构。

3.2 计算区域确定

模拟三维化石模型的外部流动必须设置一个围绕化石的计算区域，且满足以下要求[13, 24, 35]：1)计算区域体积足够大，可以将化石模型所有的部分包含进来；2)在可以准确解析流体的前提下，尽可能选择最小的计算区域，减少计算成本；3)计算区域通常设定为一个立方体或圆柱体区域，方便运算和可视化呈现。

3.3 离散网格生成

CFD是在离散的时间和空间节点上给出流动问题的解，因此需要将计算区域划分为有限数量的离散单元，这一过程称之为网格生成。因几何边界和古生物形态的复杂程度不同，网格单元的形状、数量和大小会有很大差异。特别是网格的形状和稀疏分布，对求解精度的影响重大。计算网格按照拓扑形式可以分为结构和非结构两类。二维和三维结构网格单元可简单理解为四边形和六面体。结构网格生成速度快、生成质量好，同时数据结构简单，但是对复杂几何区域的适应性差，而非结构网格则自由的多。目前，网格生成工作可由商业软件完成，例如ICEM和Pointwise等，这些软件均具有自动生成功能，但是对于复杂的模型仍需要进行人工优化。

3.4 模型设置和数值求解

CFD需要设置流体的物理属性(密度、黏度等)，一般商业软件都提供了丰富的材料库可供选择，但远古时代大气属性则应依据相关报道的数据进行设置；同时，还需要基于流动的相似准则数选择适当的流态模型[26]，例如无黏或有黏、层流或紊流，以及可压缩或不可压缩等。根据流动特性选择不同的模型，在一定程度上也是在计算准确性和计算时间(内存)之间寻找平衡的过程[2, 5]。

CFD求解时必须在计算区域边界上给定流体变量所满足的条件，即边界条件。古生物流动问题一般属于不可压缩流动，只需在计算区域入口设定具体的流速，而在出口处设置为充分发展条件即可。计算区域入口的流速取决于古生物的特性以及对它生活习性的推测：1)对于静止的生物体的模拟，相关设置可以参考现代环境中水流或风的数据获得[27]；2)对飞行或游泳等运动情况的模拟，一般是参照现生动物理论或实验的速度设定一个真实的范围[28]。在描述流体在与古生物体接触界面上的运动情况时，一般使用无滑移边界条件，即将相对于古物体表面的流速限定为零；如果古生物体表面的流动摩擦力不是重点考虑的问题，则可采用滑移边界条件[2]。

完成模型和边界条件设置后，利用CFD软件迭代求解离散化方程的过程一般不需要人工干预，求解时间取决于模型的复杂度和网格划分量以及计算机的硬件性能等多重因素。另外，求解器的选择遵循以下原则：当流体物理量不随时间而改变(即定常流动)时，可以使用稳态求解器；当流体属性随时间显著变化(即非定常流动)时，应该使用瞬态求解器，因其要给出的每个时间节点处的流场解，因此比稳态求解器耗时长得多。

3.5 数据分析和可视化呈现

模拟的结果具有多种可视化的方式，基本上可以分成6大类：xy图、等值线图、流线图(矢量图)、散斑图、网格图以及组合图。对于非定常问题模拟，可以为不同的时间节点生成可视化流场图件，并组合起来生成动画。可视化图件可用来补充流动细节、评估流体对化石的作用力等，以便在模型之间进行比较，比如计算阻力(CD)和升力(CL)的无量纲系数。对于三维模型，阻力和升力被定义为：

(2)

式中，FD为流体施加的阻力(N)，FL为流体施加的升力(N)，ρ为流体密度(kg/m3)，U为特征速度(m/s)，A为特征区域(m2)。根据化石模型几何形态和流体条件，特征区域通常是研究对象的前部区域或者是表面区域[13, 24-25]。

4 CFD在古生物学形态功能分析中的应用

在过去的半个多世纪里CFD已经成为工程领域的重要方法手段[6-11]。最近十几年，在古生物学领域也取得了非常瞩目的成果，证实了计算流体力学在论证古生物长期存在的功能学和生态学假说方面具备极大的潜力，为化石生物形态功能的关联性分析和性能检验开创了新的研究方向。

4.1 验证结构和形态的功能学特征

古生物学形态功能分析最主要的目的即为验证古生物的某一特定结构或形态的生物学功能。例如，三叶虫的外形指示其为远洋生活还是海洋底栖[25]？腕足类石燕贝的壳体形态指示其体表寄生生物采用主动寄生方式还是被动寄生方式[28]？棘皮动物海百合、寒武纪海扁果运用主动捕食方式还是被动捕食方式[15, 29]？蛇颈龙的四肢采用同步协调型还是异步不协调型的水下划动策略[18]等。该应用方向一般用于解决是非判断问题，实验中通常需要构建一个标准的古生物几何模型，各对照实验中不修改其形态特征参数，只调整流场属性、边界条件[15, 25, 28-29]或者古生物运动模态等[14]，设计相互对立的两个或多个对照组，通过分析模拟结果及其差异，为判断古生物特定结构或形态的功能学特征提供定量化的数据支持。

4.2 对比真实与假想形态的生物学意义

从进化角度而言，新的形态特征使生物体获得了新的生态或生理表现，微小的变化通常会增强或削弱现有的生物学性能[12-13, 30]。在实际问题中，即使是古生物体形态上看似微小的变化，也有可能提供一个新的功能，而此类问题需要定量的方法进行探索。由于古生物化石的不完备性，某一形态或者结构的关键演化过程往往缺少一些过渡性的化石样本作证据。因此实验分析时，需要根据不同的目的和需求，设置一系列假想的形态结构特征，并进行数值模拟分析和仿真实验。例如，腕足类石燕贝外壳的沟槽是否因适应被动取食而进化出来的？研究者在石燕贝真实的沟槽深度的基础上，又重构出无沟槽、正常沟槽以及深沟槽的壳体形态[12]，为模拟分析石燕贝外壳的水动力特征和生态学意义建立了完整的证据链。

4.3 辅助和完善多类群分析

以多种尖头型的食肉鱼类是如何攻击猎物的研究为例[23]。一般而言，鱼类捕食时为了防止猎物逃跑，需要缩短攻击时间并且降低周围水流的干扰，以免被猎物的侧线感觉系统发现。通过CFD模拟的方法，研究者模拟比较了尖头型三叠纪辐鳍鱼纲龙鱼以及类似尖头型的几类现生食肉鱼(包括雀鳝属、魣脂鲤属、大麻哈鱼属、狗鱼属和颌针鱼属)的水动力特性。结果表明，灭绝的龙鱼属产生的水流扰动较低，与现生的颌针鱼属和雀鳝属类似，是一种有效的伏击捕食者；鳟鱼适应长距离游泳，但会产生相对较大的水流干扰；狗鱼属水流干扰程度排名第二，其极高的加速度性能(即高速攻击)和突出的嘴可以补偿扰动产生的影响，从而相对能够从更远的距离捕捉猎物。

5 总结和展望

CFD的实质是在物理直观和力学实验的基础上建立各种流体运动的有限的数值模型。当问题本身遵循的规律比较清楚、所建立的数学模型比较准确，并为现实所证明能反映问题本质时，CFD才能最大程度发挥出其优越性。同时，CFD也有一定的局限性：CFD可以完全控制实验参数，但不能完全代替物理实验或理论分析；离散化的结果不仅在数量上可能影响计算的精度，在性质上还可能会改变流动的特征，其可靠性仍需要实践的验证；CFD中关于稳定性和收敛性等分析方法对线性方程一般是有效的，但是对非线性方程来说有时并没有完整的理论支持，因此计算过程中有时还需要一定的技巧。

20世纪60年代CFD已经发展成为工程领域的主要方法之一，直到近十年才逐步应用于古生物学形态功能分析领域，这显然有相当程度的滞后性，主要原因来自3个方面：1)CFD所需软、硬件经济成本较高；2)缺少训练有素的CFD专业操作人员；3)对于CFD的模拟结果还存在一定程度的不信任。尽管CFD现在古生物领域应用相对较少，但是相对于水槽和风洞等物理实验技术，CFD具有经济、迅速和适应性强等特点，特别是方便通过大量对照模拟数据验证生物体功能相关的假设，并量化特定指标进行功能分析。

目前，CFD更多地被应用在水动力条件下的海洋无脊椎动物静态模型的功能性分析上，然而空气动力条件下古鸟类和具飞行能力的爬行动物的运动方式和运动能力等方面的研究也将是今后重要的探索方向。相信随着古生物学数字化革命的推进[2-4]，在不久的将来，CFD应会成为古生物学家求证化石生物功能性相关问题的重要方法和有力工具，并有可能改变我们对很多其他古生物类群在形态功能性方面的现有认识；同时借助于多学科的交叉渗透，使古生物学在从定性描述向定量化研究发展的道路上迈进一大步。