基于动网格的永磁悬浮章动血泵流场数值模拟及溶血预测

陈 刚，姚立纲

（1.三明学院机电工程学院，福建 三明 365004；2.福州大学机械工程及自动化学院，福建 福州 350108）

1 引言

心脏移植是治疗末期心血管病最有效的手段，但心脏移植手术存在供体严重不足和排异反应。针对这一问题，目前普遍认为人工心脏替代人体自然心脏是有效的解决办法[1-3]。磁悬浮血泵无机械接触、无摩擦、无发热，可降低血液损伤，已成为国内外学者广泛研究的热点[4-6]。

血泵溶血问题一直是制约人工心脏发展的难点，早期对血泵溶血的研究主要是采用实验法，通过采样化验以确定溶血值，该方法具有较高的可靠性，但是实验周期长，需要反复采样[7]。利用计算流体力学（Computational Fluid Dynamics，CFD）研究血泵流场流动情况及溶血预测，可以缩短研制周期、降低成本。文献[8]通过CFD技术研究了几何形状优化使轴对称狭窄几何中的血液损伤最小化，证明了基于CFD的设计技术的生物医学意义。文献[9]以提高叶轮泵的效率为目的，运用CFD数值分析研究了泵内部流场情况，探究高速叶轮泵流场规律。文献[10]通过计算流体动力学，以提高血泵水力性能和降低血液损伤为目标，对血泵结构进行了优化设计。上述文献中磁悬浮血泵多数为离心泵和轴流泵，离心泵转速较低，但体积较大，不利于植入性；轴流泵体积小、效率高，但转子速度普遍很高，产生较大的切应力，易造成血细胞破坏[11]。

章动传动是在天体行星运动或者陀螺仪运动原理的基础上提出的一种新型传动方式[12]。目前章动原理广泛应用于齿轮传动、机械手、电机、电动自行车、太空探测器等领域。笔者所在课题组基于章动传动原理提出了机械轴承章动式血泵[13]，如图1所示电机轴通过斜套驱动章动盘，迫使章动盘做定向的章动运动，电机轴线与斜套内孔轴线的夹角为章动角。然而，章动盘在容腔内由其中心位置的球副支承，球副与泵盖球腔是直接的机械接触，存在较大的机械磨损，容易导致血细胞的破坏和污染。

图1 章动转子驱动原理Fig.1 Driven Principle of Nutation Rotor

在机械轴承章动式血泵的研究基础上，自主研发了永磁悬浮章动血泵，分析了非接触磁力驱动和磁悬浮轴承的结构和工作原理，采用CFD流体动力学仿真分析了血泵内部流场，探究了血泵工作过程流场内的速度矢量、压力以及切应力等关键参数的大小与分布规律，为磁悬浮章动血泵的进一步优化改进和实验研究奠定了理论基础。

2 磁悬浮章动血泵工作原理及结构

章动血泵兼顾了离心式血泵转速低和轴流式血泵体积小的特点。在满足正常人体对血液流量5L∕min的要求下，章动血泵的理论转速和体积都相比其他类型血泵有明显的优势[13]。然而，由于机械轴承章动血泵在驱动和支承方式上都存在机械接触，不利于血泵的持久性工作和血液相容性。

磁悬浮章动血泵的结构模型，如图2所示。采用了非接触磁力驱动和永磁悬浮轴承的设计。在章动盘的盘缘均匀设计了若干块径向磁化的扇形永磁体，同时在旋转套筒内布置了两块成对角位置的径向磁化永磁体，旋转套筒内永磁体与章动盘盘缘永磁体相吸，在磁力作用下章动盘保持倾斜状态，倾斜角度即为章动角。随着旋转套筒的定向转动，章动盘上不同位置的永磁体一侧被吸引向上运动，另一侧被吸引向下运动，产生起伏摆动，章动盘与上、下盖始终成线接触，接触线的旋转方向与旋转套筒一致。在非接触磁力驱动下，章动盘定向起伏摆动，进口区域不断增大形成负压、出口区域逐渐减小压力增大，从而血液从进口流入、出口流出。

为了实现血泵转子的悬浮，在球副和泵盖内分别设计了永磁体，磁悬浮章动血泵转子支承结构，如图3所示。章动盘由中心位置的球副支承，球副分为上、下两部分，内部分别正交布置有永磁体，在泵盖的球腔内也正交布置了永磁体。章动盘球副中的永磁体与泵盖球腔内的永磁体位置对应，且互相排斥。在上盖与上球副、下盖与下球副内永磁体之间的磁力作用下，实现血泵转子的非接触支承，以减小机械接触磨损和血液污染。

3 章动式磁悬浮血泵流场数值模拟

3.1 血泵计算模型、网格划分及边界条件

采用三维计算机辅助设计软件创建血泵的流道模型，流道计算区域由转子部分、进出口流道、泵盖部分以及泵外壳内表面组成。整个流道中转子部分的球副和泵外壳内表面是球面，泵上、下盖内表面是锥面。计算区域分为转子区域和静止导流区域。对三维几何模型进行布尔减运算，将计算流道模型导入到CFD前处理软件，进行网格划分。

由于章动血泵转子的运动会造成流域形状随时间发生变化，因此本文血泵的流场模拟采用动网格模型。章动盘与泵盖内锥面成线接触，相切的接触线位置会形成楔形区域，导致该处网格质量较差。同时，在动网格定义运动部件时，运动部件与流道壁面发生接触，动网格在重构时会畸变出现负体积。因此，在划分网格时计算模型对章动盘的几何做了修改，章动盘上、下表面分别与泵上、下盖有0.1mm间隙，避免出现线接触。章动血泵的几何模型较为复杂，网格划分采用非结构化四面体网格。网格模型的单元数为164591个，网格节点数量为29114个。

计算区域为血泵整个流道，边界条件选取压力入口和压力出口，出口和入口的压力差为100mmHg（13333Pa），入口压力设置为标准大气压101325Pa，出口压力114658Pa。流体为血液，其密度为ρ=1.06g∕cm3，粘度为0.0035kg∕（m·s）。入口血流设置为不可压缩的非定常流，在程序接口载入UDF函数定义章动盘的运动和血液特性。

3.2 控制方程

根据流体力学理论，血泵内血液的流动必须满足三大物理守恒定律，即质量守恒定律、动量守恒定律及能量守恒定律。血液是不可压缩的黏性非牛顿液体，故血泵流体的质量守恒控制方程为：

式中：ρ—血液的密度；t—时间；u、v、w—流体速度矢量在x、y、z三个方向上的分量。

在流体动力学中动量守恒方程可描述为，流体微元体对时间的变化率与该时刻微元体所收到的外界合力相等。动量守恒的控制方程可表示为：

式中：F x、F y、F z—单位质量力在x、y、z三个方向上的分量；u x、u y、u z—x、y、z轴上的速度分量。

能量守恒定律即是热力学第一定律，由于人体体温基本是恒定值，血泵中流体内能近似认为保持不变，不考虑温度变化和热量交换。因此，能量守恒的控制方程表达式为：

其中，左边三项分别是流体微元的动能、重力势能和压力势能。

3.3 数值求解

血泵流场数值模拟选用瞬态计算，湍流模型采用标准κ-ε模型，血泵内部近壁面流动问题应用标准壁面函数法进行处理。动网格的更新使用Smoothing和Remeshing方法，流道中血泵内、外两侧的球面设置为变形面，章动盘为刚体。计算时间步长设置为10-4s，各变量的收敛残差设为10-6。

4 数值计算结果讨论

采用Fluent流体仿真软件，对磁悬浮章动血泵流场模型进行了三维流场数值模拟，计算了血泵流场内的流量、速度、压力及剪切应力的分布。

4.1 不同转速下的输出流量

稳定的流量输出是血泵良好性能的重要指标之一。正常成人的心脏血液输出量约为5L∕min，动脉压力在（80～120）mmHg。根据图5血泵在不同转速下的输出流量计算结果，当输出流量为5L∕min时，进出口压差为100mmHg下，血泵转子所需的转速约为1300r∕min。因此定义进出口压差100mmHg、流量5L∕min为血泵的标准工况。在分析血泵内部流场时，采用这一标准工况。

图5 不同转速下的输出流量Fig.5 Output Flow at Different Rotating Speeds

输出流量整体接近一条水平直线，表明血泵的流量输出是稳定的。在半周期的整数倍时，流量曲线出现下凹现象，但下凹值不大，且时间极短（在1500r∕min时，流量下凹值约为1.1L∕min，作用时间约0.0072s）。出现这一现象的原因是，由于与转子固连一体的章动盘设计有一处开口，当章动盘与泵盖的接触线旋转至开口时，无法形成接触线，导致进出口出现短暂的连通，造成这一位置输出流量的下降，形成了流量曲线的下凹现象。正常成年人平均心率约为75次∕min，完成一次射血全过程的时间约为0.8s，而磁悬浮章动血泵在转速为1300r∕min、输出流量5L∕min下的工作周期为0.046s，也就是每半个周期出现一次流量下降的时间是极短暂的，不会对供血产生影响。因此，流量曲线说明章动血泵的输出流量稳定。

4.2 血泵流场速度及压力分布

标准工况下磁悬浮章动血泵内部流场的速度矢量图，如图6所示。可以观察到血泵泵内部流体的流动趋势一致，整体流动速度较低且速度梯度较小，最大速度约为3.6m∕s，出现在入口和出口位置，这是由于入口和出口处流道横截面较小。在两段接触线处出现微小的回流现象，这是由于章动盘盘面和上下盖锥面之间的微小缝隙造成的。血泵内部流场没有出现流动死区。

磁悬浮章动血泵流场压力云图，如图7所示。可以看出，章动盘和泵盖的接触线将血泵内部分成一个高压力区和一个低压区。高压区的压力略大于设定的出口压力114658Pa，而低压区的压力则略小于设定的入口压力101325Pa，低压区与入口相连起“吸水”作用，而高压区与出口相连，起“压水”作用。在高压区，靠近流道外侧区域的压力较大，靠近内侧区域的压力较小，这是由于流体做圆周运动，受到离心力的作用形成的效果。总体来讲，高压区和低压区各自的压力分布都是非常均匀的。

图7 磁悬浮章动血泵流场压力云图Fig.7 Pressure Nephogram of Flow Field in the Maglev Nutation Blood Pump

4.3 血泵内剪切应力分布

标准工况下磁悬浮章动血泵内部流场的剪切应力云图，如图8所示。

图8 磁悬浮章动血泵流场剪切应力云图Fig.8 Shear Stress of Flow Field in the Maglev Nutation Blood Pump

血泵内较高的剪切应力分布在接触线位置，最大值约为588Pa左右，其他区域的切应力均非常小。根据文献[14]等的实验研究，血液发生溶血与切应力和红细胞暴露时间有关，在切应力小于150Pa时，即使红细胞暴露的时间无限长，红细胞也不会破坏。由图9剪切应力大小所占比例显示，剪切应力值小于150Pa的约占90%。大于150Pa的剪切应力分布在高速运动的接触线附近，该区域血细胞暴露的时间极短。因此，流场内剪切应力分布表明血泵内红细胞不易受到破坏，血泵的抗溶血性能较好。

图9 剪切应力直方图Fig.9 Shear Stress Histogram Showing Percentage Distribution of Stress

5 基于CFD的溶血预测模型及结果

5.1 溶血预测模型

导致血泵发生溶血的根本因素是机械切应力。采用DPM（Discrete Phase Model）离散相模型，将血液中红细胞作为血浆中的离散相，求出红细胞的运行轨迹、剪切应力值及暴露时间，以预测血泵对红细胞的破坏程度。

利用CFD数值求解获取每一个粒子在运动时下一时间步的位置，公式如下[15]：

其中，Δt i=t i+1-t i，τi—i时刻粒子所受的剪切应力值。所有微小时间段所对应的溶血积分则为血细胞的溶血值为：

选取N个粒子轨迹数目，将粒子所受的破坏程度进行累加并求平均，即为血泵的溶血损伤预测指标

5.2 血泵溶血预测数值计算与结果分析

根据所建立的溶血模型，采用粒子追踪法对磁悬浮章动血泵进行溶血预测。从入口开始追踪粒子至出口，获取粒子所受剪切应力流动轨迹，通过分析粒子的流动轨迹和所受剪切应力大小，以预测血液损伤程度。选择三个工况（进出口压差均为100mmHg时，n=1300r∕min；n=1800r∕min；n=2300r∕min）进行模拟计算，在进口释放100个粒子，粒子直径设置为8.5μm，体积与红细胞相当[7]。

不同工况下同一编号粒子在血泵内的运动时间与所受剪切应力情况，如图10所示。数值计算结果表明，血泵转子速度越高，粒子受到的剪切应力越大，但暴露时间越短。转速在2300r∕min时，最高剪切应力值约为600Pa，暴露时间为0.098s，即该粒子从入口至出口通过血泵的时间为0.098s。该工况的剪切力在0.03s后逐渐增大，主要是进口区域变大、出口区域减小，导致流速增大，出口区域剪切应力增大。在转速为1300r∕min的标准工况时，粒子所受最大剪切应力为450Pa，但出现的时间极短；大部分剪切应力值小于150Pa，粒子在血泵内的暴露时间约为0.12s。

图10 不同工况下编号80粒子受剪切应力情况Fig.10 Shear Stress of 80 Particles Under Different Working Conditions

不同工况下同一粒子在通过血泵时的溶血情况，如图11所示。溶血指数由式（6）求得。标准工况下，该粒子在进入血泵0.042s后开始出现溶血，通过血泵过程的累积溶血指数为4.1×10-9。随着血泵转子速度的增大，粒子的溶血指数增大。红细胞的溶血程度受剪切应力和暴露时间的综合影响，在转速较高时剪切应力明显增大，但其暴露时间短，因此血泵内红细胞溶血值相对也是较低的。

图11 不同工况下编号80粒子的溶血指数Fig.11 Hemolysis Indices of 80 Particles Under Different Working Conditions

6 结论

基于CFD动网格方法对磁悬浮章动血泵的内部流场进行了数值模拟和溶血预测，结果表明：

（1）磁悬浮章动血泵的章动盘开口大小对流量输出的稳定性有影响，接触线运动至开口位置时流量会下降，转速越高流量下降越明显。（2）血泵流场内速度较低，速递梯度小，无明显的回流现象。较大的剪切应力主要分布在章动盘与上、下盖的接触线区域，血泵内整体的剪切应力比较小。（3）血泵的溶血指数较小，溶血预测模型的计算结果与文献[16]的实测结果相符，这进一步验证了磁悬浮章动血泵良好的血液相容性。

通过采用CFD动网格对磁悬浮章动血泵的血流动力学分析，可以更好地了解血泵内部的流场情况，判断高剪切应力的位置，对血泵结构的优化改进及进一步的实验研究有重要意义。