1.3.2 计算流体力学方法

1.3.2 计算流体力学方法

计算流体力学(computational fluid dynamics,CFD)是通过计算机对流体流动等问题进行数值计算,分析的对象包括流体流动、热传导、物质输运等物理学现象。CFD进行数值分析的基本思想是通过一系列有限个离散点上的变量值的集合,来代替空间域上连续的物理场,如速度场和压力场;然后按照一定方式建立这些离散点上的场变量之间关系的代数方程组,通过求解代数方程组获得真实的变量近似值[15-16]。近年来,随着医学影像技术的进步和计算机运算能力的提高,计算流体力学方法在心脑血管系统的血液动力学模拟中也发挥着越来越重要的作用[17]

1)有限体积法

有限体积法(finite volume method,FVM)又称为控制体积法,其基本思路是将计算区域划分为网格,并使每一个网格点周围有一个互不重复的控制体积;将待求解的微分方程对每一个控制体积积分,从而得出一组离散方程。其中的未知数是网格点上的因变量[18]。为了求出控制体积的积分,必须假定因变量的值在网格点之间的变化规律。有限体积法具有以下特点:

(1)有限体积法的出发点是积分形式的控制方程,同时积分方程可以表示变量在控制体积内的守恒特性。

(2)积分方程的每一项都有明确的物理意义,从而使方程离散时,各离散项都可以得到一定的物理解释。

(3)离散后的各个节点有互不重复的控制体积,从而使整个求解域中场变量的守恒可以由各个控制体积中特征变量的守恒保证。

有限体积法常用离散格式包括以下几种:

(1)一阶迎风格式。适用于结构化网格(如四边形网格或六面体网格),它具有稳定性高、计算速度快的优点,但是在网格方向与流动方向不一致时,产生的数值误差比较大,而且当网格密度不足时,一阶迎风格式的求解精度有限。

(2)二阶迎风格式。一般适用于流动方向和网格的边缘不一致的情况。二阶格式的计算精度高于一阶格式,但是相对而言,其计算时间比较长,收敛性也相对较差。

(3)QUICK格式。在用结构化网格计算旋转流动问题时(比如对轴流血泵划分结构化网格后进行数值模拟),QUICK格式可以提供更高的计算精度。

除此之外,有限体积法的离散格式还包括中心差分格式、指数格式、乘方格式等。需要注意的是,在实际计算中,选择格式时需要兼顾精度、收敛性和计算时间等方面的要求。对于复杂流动,可以使用一阶格式进行计算,达到收敛后再设置离散格式为高阶格式。有限体积法作为求解流动和传热问题的数值计算中较成功的方法,已被多数工程流体和传热计算软件所采用。

2)格子-玻尔兹曼方法

格子-玻尔兹曼方法(lattice Boltzmann method,LBM)与建立在连续介质模型上的传统计算流体力学方法不同,LBM基于统计力学和分子运动论。LBM从介观尺度出发,建立离散的速度模型,在满足质量、动量和能量守恒的条件下,得出粒子分布函数,然后对粒子分布函数进行统计计算,得到压力、流速等宏观变量[19]。LBM作为具有显著优势的流体计算方法,已被广泛用于理论研究和处理工程问题。LBM有以下优点:

(1)算法简单。简单的线性运算加上一个松弛过程,就能模拟各种复杂的非线性宏观现象。

(2)编程容易。计算的前后处理也非常简单,具有很高的并行性。

(3)能够处理复杂的边界条件。(https://www.daowen.com)

(4)LBM中的压力可由状态方程直接求解。

(5)能直接模拟有复杂几何边界的诸如多孔介质等连通域流场,无须做计算网格的转换。

(6)LBM具备描述粒子运动的特性,使得其在处理流体与固体作用时相对直观,在解决气固和流固耦合方面具备优势。

(7)LBM不受连续介质假设的约束,能够用来解决纳/微尺度的流动和传质或稀薄气体输运等连续介质假设不适用的问题。

(8)LBM在处理多相多组分流体问题时,相比于传统计算流体力学方法,在抓取移动和变形的界面、描述组分间相互作用方面具备明显优势。

LBM具有原理简单、编程方便、可进行大规模并行计算等优点,在模拟医疗器械植入后的血液流动及血栓形成等方面具有广泛的应用。

3)光滑粒子流体动力学方法

在计算流体力学史上,欧拉观点是最常用的描述流体运动和物质传输现象的方法。基于欧拉观点的计算流体力学方法,如有限体积法、有限差分法、有限元法等,需要将物理空间离散化成网格,以求解离散偏微分方程或积分公式。使用网格方法最大的困难在于处理具有自由表面、移动界面、带有复杂几何形状变化的流体流动问题。为了保证数值模拟的精确度,在每次数值迭代时都需要根据变形几何重新划分网格件,从而大大增加了数值模拟的计算量。

光滑粒子流体动力学(smooth particle hydrodynamics,SPH)方法是一种基于拉格朗日观点的无网格方法,将具备连续性假设的流体(或固体)用相互作用的粒子来描述,每个粒子上承载各种物理量,包括质量、速度等,通过求解粒子群的动力学方程和跟踪每个粒子的运动轨迹,从而求得整个系统的力学行为。SPH方法在应用上首先将求解域离散化,物理场被离散为空间中粒子质心上的物理值,将场函数表示成一个狄拉克δ函数和任意点上函数值的积分,并用一个光滑函数(核函数)近似替代δ函数。通过这种离散方式,流体力学控制方程组中的函数积分在求解域中任意点的值转化为该点周围的所有粒子的物理量叠加求和,即每个粒子的物理量都可以由周围粒子的物理量插值得到[20-21]。SPH方法的特点如下:

(1)将求解域离散在其中任意分布的粒子,无须划分网格,粒子间也不需要连续性约束,SPH方法是一种纯粹的无网格方法。

(2)用含有核函数的积分近似方法对场函数近似,积分近似方法是一种弱形式的守恒方程,保证了SPH方法在计算上的稳定性。

(3)核函数的积分形式需要进一步转化为粒子上物理值的求和,即控制方程中场函数的积分和微分都需要通过对局部区域内离散粒子上物理值求和得到,该局部区域被称作支持域,支持域可以提高计算效率。

(4)计算过程的每个时间步上的场函数都需要当前支持域上粒子的物理值确定,这体现了SPH方法在处理大变形问题和复杂几何形状中流体运动的优势。

(5)控制方程中的每一项都由空间中位置不断变化的粒子上的物理量近似得到。

(6)SPH方法使用显式积分法求解离散后的微分方程组。

在SPH方法中,核函数非常重要,它反映了粒子间相互作用的贡献大小。随着粒子间距的增加,相互作用应该减弱,因此函数应取为距离的单调减函数。核函数不仅决定了函数近似式的形式、定义了粒子支持区域的尺寸,而且还决定了核近似和粒子近似的一致性和精度。常用的核函数有高斯型光滑函数、B样条函数等。SPH方法由于无须划分网格,适合于模拟大变形的流固耦合问题或者自由面问题,因此常被用于心脏瓣膜或者心室的流固耦合仿真中。