1.3.3 流固耦合仿真方法
为了更好地分析血流和瓣膜或者血流与血管之间的相互作用,双向的流固耦合是必要的。在这样的仿真中,流体域和固体域之间的位移、压力等变量需要实时进行传递。流固耦合如果设置得当,可研究完整心动周期的瓣膜受力情况,但是在瓣膜仿真中应用双向流固耦合方法也是有很多挑战的,其中最大的问题是瓣叶打开时的快速运动及瓣叶关闭时的相互接触,都有可能引起流体域网格的拓扑发生剧烈变化而导致计算不收敛。因此,根据不同的问题使用不同的求解器,实现双向流固耦合的计算方法也不同。
1)任意拉格朗日-欧拉方法
对于连续介质(瓣膜、血液、血管、支架等)而言,在拉格朗日观点中,为了确定各物质点的位置,引入拉格朗日坐标系,它是和物质固结在一起并随物质一起运动的。质点在初始时刻的位置由其在拉格朗日坐标系中的位置矢量完全确定。在拉格朗日描述下,计算网格网格点与物质点在物质的变形过程中始终保持重合,这会导致在涉及特大变形问题中,物质的大幅度变形将导致计算网格的畸形而使计算失败。在欧拉描述中,为了描述当前构形,引入空间坐标系或欧拉坐标系,它是和空间固结在一起的。空间中各点的位置由其在欧拉坐标系中的位置矢量确定。在欧拉描述下,计算网格网格点固定在空间中,即计算网格在物质的变形过程中保持不变,因此可以很容易处理物质的大幅度变形,但此时对运动界面需要引入非常复杂的数学映射,这可能导致较大的误差。
任意拉格朗日-欧拉方法(arbitrary Lagrangian Eulerian method,ALE)是一种结合了拉格朗日描述和欧拉描述优势的新方法。在ALE描述中,引入了一个可以独立于初始构形和当前构形运动的参考构形。在物质的变形过程中,观察者始终跟随参考构形运动,因而对观察者而言,参考构形是固定不动的,而初始构形和当前构形则都相对于参考构形运动。为了确定参考构形中各参考点的位置,引入参考坐标系,参考构形中各点的位置由其在参考坐标系中的位置矢量确定。在ALE描述下,连续介质的计算网格可以在空间中任意运动,即可以独立于物质坐标系和空间坐标系运动。这样通过规定合适的网格运动形式,可以准确地描述物质的移动界面,并维持单元的合理形状[22-23]。纯拉格朗日和纯欧拉描述实际上是ALE描述的两个特例,即当网格的运动速度等于物质的运动速度时就退化为拉格朗日描述,而当网格固定于空间不动时就退化为欧拉描述。
ALE因为可以捕捉靠近界面部位的细致流场,获得更准确的流体动力学信息,所以在包括航空航天和车辆船舶的流固耦合仿真领域都有着非常广泛的应用[24]。
但如前所述,ALE在处理瓣膜的流固耦合问题时却很难实现稳健的求解,其最主要的困难还是来自瓣叶快速开闭过程中流体网格自适应重构不稳健带来的计算不稳定性。ALE要求流体域和固体域之间的界面完全匹配,在固体发生变形时,流体域的网格要跟随流体-固体的界面发生运动甚至发生重构。但在瓣叶快速开闭过程中,网格重构的频率和剧烈程度是非常大的,尤其当瓣膜关闭、瓣叶即将发生相互接触和交叠时,瓣叶之间的流体网格将被挤压、细化、重构到尺度非常小的水平,这样的重构过程很容易造成时间步长突然急剧减小、网格过于扭曲和计算不收敛。
2)浸没边界法
为了解决瓣膜运动这种大变形的流固耦合问题,Peskin提出了浸没边界法(immersed boundary method,IBM),他创造性地将传统流固耦合时界面传递物理学变量的过程,转化成用虚拟的弹性力按周围流固材料节点距离和质量权重重新分配的过程[25-27]。在浸没边界法中,流体在欧拉坐标系上表示,结构在拉格朗日坐标系上表示。对于不可压缩牛顿流体的运动,其满足的Navier-Stokes方程的形式为其中,f包括弹性边界传递给流场的所有体积力,这些力作为动量方程中的源项插入,其表达式为

(https://www.daowen.com)
式中 F——浸入边界产生的单位力;
X(s,t)——浸入边界的位移函数;
δ——Dirac函数。
可以将力扩展到多个维度,对应于弹性表面或三维实体建模。由于流体是黏性的且在界面附近速度连续,而且结构随着流体一起移动,所以有

以上这些方程构成了浸没边界法的控制方程。浸没边界法在每个时间步长内对这些控制方程数值离散的基本思路是:①先根据结构网格浸入流体域的位置和初始状态计算出当前位置产生的体积力;②将计算出的结构离散点的体积力分配到周围流体网格节点上;③在整个流场中求解Navier-Stocks方程获得新的速度场;④根据新的速度场求得浸入流体的结构运动速度,从而获得边界离散点的新位置。标准的浸没边界法执行步骤如图1-34所示。

图1-34 标准的浸没边界法执行步骤
浸没边界法求解过程稳健,其无须构造贴体网格、无须实时更新流体域网格的特点,非常适合处理有大变形、大位移的瞬态流固耦合问题。但经典浸没边界法在近壁面处流场刻画不够精确的缺点也限制了其在瓣膜仿真应用中精确分析瓣叶和血管壁上的壁面切应力。