基于主动转速调制的离心旋转血泵流场分析

黄 峰 王楚晨 阮晓东

1.浙江科技学院机械与能源工程学院，杭州， 3100232.中国计量大学计量测试工程学院，杭州， 3100183.浙江大学流体动力与机电系统国家重点实验室，杭州， 310027

0 引言

根据《中国心血管病报告2017(概要)》推算，目前中国心血管病患病人数已达2.9亿例，其中心力衰竭450万例[1]。心脏移植依然是治疗终末期心力衰竭最有效的方法，但由于心脏供体严重短缺，因此需要心脏移植的病人中大约30%～50%死于合适供体的缺乏[2]。人工心脏能够对衰竭心脏的功能进行补偿，协助心脏为血液循环系统提供动力，可用于等待心脏供体的过渡治疗甚至完全代替原有心脏进行永久性辅助治疗，已经成为心衰治疗的重要手段[3-4]。

从1953年Gibbon利用体外循环进行心脏手术开始，研究者们对机械循环辅助开始了不断的探索[5]。到目前为止，心室辅助装置(ventricular assist device,VAD)的发展已经经历了三个阶段，悬浮式连续流旋转血泵逐渐成为现在的主流。第一代临床认可的脉动流(pulsating flow，PF)泵是一种产生脉动压力和流量的大排量泵，但其体积大，往复运动部件寿命短，这限制了该泵的实际使用范围。与PF-VADs相比，基于旋转驱动的连续流(continuous flow，CF)泵具有更好的耐用性和可靠性[6]，但在常规恒泵速运行时CF泵会降低血液循环的搏动性。有报道称，PF-VADs支持的患者比CF-VADs支持的患者的心肌恢复潜力大3倍[7]，验证了循环搏动性可能在心力衰竭逆转重构中发挥重要作用的推测。

CF-VADs通常以恒速运行。有学者提出对CF-VADs泵速主动调节作为一种潜在的临床操作策略，以恢复循环搏动性，同时抑制泵内血栓的形成[8-9]。ISING等[10]利用计算机仿真研究了几种候选的CF-VADs流量调制算法对提高血管脉动性的作用。SOUCY等[11]用慢性缺血性心力衰竭的牛模型研究了旋转血泵HVAD®(heartware ventricular assist device)在同步和异步泵速调制算法中对流量脉动和左心室容积卸荷的影响，结果表明转速调制能够增加脉动，改善心功能和器官灌注效果。美敦力公司[12]通过对HVAD®植入患者的比较研究，证明了旋转血泵周期性速度变化可以产生有效的临床效果，有转速调制功能的患者，其中风、败血症和右心衰的发生率明显低于无转速调制功能患者的上述病症发生率，而总体生存率和胃肠道出血率几乎没有受到影响[12]。

虽然转速调制方法有上述优点，但是也会带来一些负面影响。转速调制迫使旋转血泵速度急剧变化，这可能扰乱血泵的内部流场，使泵内剪切应力增大从而导致血液损伤。然而到目前为止，关于转速调制对血泵内部流场影响的相关研究很少，也未见不同转速调制波形在流场变化和剪应力增大方面的对比研究。因此，本研究以自主研发的离心旋转血泵为对象，通过心血管系统-旋转血泵数学模型的联合系统仿真得到匀速、正弦波、方波以及三角波转速调制下血泵与心脏并联时的进出口工作压力，并在该压力边界条件下进行旋转血泵计算流体动力学模拟，对比分析转速调制对血泵流场的速度分布以及剪切应力分布的影响，为旋转血泵的转速调制策略优化提供参考依据。

1 材料和方法

1.1 研究对象

本文的研究对象为实验室自主研发的连续流式离心旋转血泵[2]，血泵结构如图1所示。血泵的具体几何参数为：入口内径6 mm、出口内径8 mm、外壳直径70 mm、转子外径50 mm、表面积138 cm2、整体体积1.83×10-5m3。

(a)实体血泵分解图 (b)血泵实体图

(c)血泵剖面图图1 离心血泵分解图、实体图与剖面图Fig.1 Decomposition diagram, solid diagram andsection diagram of the heart pump

1.2 计算流体动力学模拟

1.2.1网格划分

首先使用SolidWorks对已有的离心血泵模型进行流域的提取，得到全部的流体流域。为了更方便地观察进出口处的血液流动，对出口管和入口管处做了延长处理，使血液流动能够充分地发展。采用商用数值模拟前处理软件ICEM对流域进行网格划分(图2)。由于血泵的结构比较复杂，将血泵流体域划分为入口管、叶轮主流道、蜗壳和出口管四个部分，每个部分分别划分网格以提高模拟结果的准确性，总体网格如图2a所示。网格划分全部采用结构化的六面体网格，网格总数为566万。进出口管以及叶轮中心部分采用经典的O形网格划分。由于叶轮结构较为复杂，将叶轮按叶片分别进行划分(图2c)，最后进行组合(图2d)。对叶片以及部分窄面进行局部加密处理以提高网格质量。将划分好的网格进行网格质量验证，结果表明，除叶轮叶片处极小部分网格外，大部分网格质量均较好，网格质量满足要求。此外，将各部分网格分别加密30%进行网格无关性验证，其余设置均保持一致，对比加密前后3000 r/min恒速条件下的扬程，结果表明，网格加密后的扬程变化在1%以内，证明该计算结果对网格个数不敏感，当前的网格数可以捕捉流场信息，进一步验证了网格划分的有效性。

(a)总体网络 (b)蜗壳处网格

(c)叶轮拆分网格 (d)叶轮总体网格图2 血泵的总体和局部网格划分Fig.2 Local and global meshing of blood pump

1.2.2数值模拟方案

数值模拟中使用的工作流体为血液，密度为1050 kg/m3，动力黏度为0.0035 kg/(m·s)。由于高剪切速率下血液的剪切稀化等非牛顿流体性质不明显，因此在本研究中将血液看作一种血细胞均匀分布的牛顿流体[13]。湍流模型采用剪切力传输模型SST，相比标准的κ-ω模型，该模型在广泛的流动领域中具有更高的精度和可信度。将叶轮定义为旋转域，使用旋转坐标系固定。同时将旋转域和静止域的动静交界面设为“Frozen rotor”，固体壁面设置为无滑移。

1.2.3边界条件

叶轮的初始转速设置为3000 r/min，速度调制的正弦波、方波、三角波三种波形如图3所示，三种调制转速均在2500～3500 r/min范围内波动。血泵速度周期与心脏脉动周期异步，即泵速以独立于自然心率(60 次/min)的速度调节，周期设置为2 s。

图3 三种转速调制的波形Fig.3 Waveforms of three speed modulations

血泵进口和出口都设置为动态压力边界条件。该边界条件由旋转血泵和心血管系统模型的联合系统仿真得到。心血管系统采用集总参数模型，由左右心房、左右心室、肺动脉和体动脉负荷、体静脉和肺静脉回流组成(图4)，旋转血泵模型则由系统辨识得到[14]。系统仿真在MATLAB/Simulink中进行，最后得到了旋转血泵在匀速以及正弦波、方波和三角波转速调制下的进出口压力。上述系统仿真的左心室都设置为衰竭情况，搏动能力为正常情况的30%。

Simulink系统仿真采用ode45求解器，采样周期为0.001 s，仿真时间为30 s。本文中LVAD的连接方式为左心室到主动脉，其进出口压力分别为左心室压力和主动脉压力。在得到心衰时不同调制转速下的左心室和主动脉压力随时间的变化数据后，通过MATLAB中曲线拟合模块中的傅里叶数据拟合方法得到压力变化的拟合函数，该拟合函数最终设置为CFD数值模拟的动态压力边界条件。血泵进出口动态压力的系统仿真结果与拟合曲线如图5所示，血泵入口和出口压力数据拟合的复相关系数均大于0.98，表明拟合结果很好。

(a)血泵出口压力

(b)血泵入口压力图5 血泵进出口动态压力的系统仿真结果与拟合曲线Fig.5 Systematic simulation results and fitting curvesof inlet and outlet dynamic pressures

1.2.4水力特性验证

图6 血泵稳态水力特性的CFD模拟与实验结果对比Fig.6 Comparison between numerical calculation andexperimental measurement results ofhydraulic characteristics

为了验证模型的有效性，通过CFD模拟得到了血泵在转速为3000 r/min时的稳态水力特性曲线，并与实验结果[2]进行了对比，对比结果如图6所示。CFD仿真结果与实验结果的误差在0.88%～6.7%之间，表明仿真结果与实验结果吻合较好，验证了数值模拟设置的湍流模型、边界条件的有效性，说明数值模拟能够较准确地反映实际的流动情况。此外，也进行了脉动流量条件下的CFD模型实验验证，结果表明CFD模型预测的结果和实验结果吻合较好，证明了本研究采用的CFD模型同样适合脉动条件下的模拟[15]。

1.3 标量剪切应力

血泵内的标量剪切应力(scalar shear stress，SSS)与血细胞破坏程度有着密切关系。标量剪切应力σ由模拟速度场计算得到：

其中，τii、τij、τjj为剪切应力张量在笛卡儿坐标系中的六个分量[16]。SSS按大小分为三大类：①σ<10 Pa,属于生理范畴；②10 Pa≤σ≤100 Pa,属于非生理范畴；③σ>100 Pa,属于极端非生理范畴。当血液暴露在非生理SSS时，会引起血小板活化和溶血；暴露在极端非生理SSS时，会引起血小板活化、受体脱落、细胞破裂等血细胞结构和功能损伤[17]。

2 结果

2.1 血泵流量

图7 CFD模拟得到的不同转速波形下的流量最大值、最小值和均值Fig.7 Flow values obtained by CFD simulation underdifferent speed profiles

在3000 r/min匀转速和三种转速调制波形条件下，CFD模拟得到的血泵流量最大值、最小值和均值如图7所示。转速调制运行时的流量均值分别为正弦波5.723 L/min、三角波5.742 L/min、方波5.849 L/min，均略高于匀速3000 r/min运行下的流量均值5.593 L/min。正如预测的那样，相比恒定转速的情况，转速调制下产生的流量脉动更大。而恒定泵速(3000 r/min)下平均流量的实验数据为5.0 L/min[2]，小于CFD模拟的结果，这一差异主要是由于实验数据是在血泵恒定的进出口压力下测得的，而本研究中即使是恒定转速条件，由于心脏残余的搏动能力，血泵的进出口压力仍具有一定的波动。

2.2 剪切应力分布

血泵上泵壳表面的壁面剪切应力(wall shear stress, WSS)分布云图见图8，该图显示了匀转速以及三种转速调制波形下，1个速度调制周期内5个离散时间点的WSS分布情况。匀转速状态下，各时间点的剪切力大小基本一致，高剪切区域集中在叶轮叶片的后缘。当转子处于减速阶段时，三种转速调制波形下的WSS均小于匀转速状态下的WSS，并在转子速度降到最小(2500 r/min)时同步减小到最小值；在转子加速阶段，WSS沿叶轮半径逐渐增大，在转速达到最大值(3500 r/min)时，叶轮尾缘处的WSS数值和区域面积都达到最大值，该情况在减速阶段逐渐缓解。总体来看，WSS大小与转速的大小成正比，四种状态的最大WSS分别为匀速256.887 Pa、正弦波292.401 Pa、三角波290.433 Pa、方波320.764 Pa。方波调制下(图8)各时间点的WSS数值大于其他波形下WSS数值，但总体趋势较为稳定。

(a)匀转速

(b)正弦波

(c)三角波

(d)方波图8 匀转速、正弦波、三角波和方波转速调制条件下血泵的WSS分布Fig.8 WSS distribution under constant, sinusoidal, triangular and square speed modulations

在匀速及三种速度调制下不同时间节点的SSS中，正弦波的最大SSS为153.9 Pa，三角波为146.3 Pa，方波为146.5 Pa，而匀速为131.125 Pa。可见四种状态均存在极端非生理范畴，但是转速调制会产生更大的SSS。使用CFD后处理软件得到大于100 Pa的极端非生理SSS区域如图9所示，极端非生理区域主要集中在叶轮叶片的前缘、后缘，以及叶轮和蜗壳的交界处，区域体积随速度的增大而增大。

(a)匀转速

(b)正弦波

(c)三角波

(d)方波图9 匀转速、正弦波、三角波和方波转速调制条件下血泵内SSS的极端非生理区域Fig.9 Extreme nonphysiological areas of the SSS inside the blood pump under constant, sinusoidal, triangular andsquare speed modulations

血泵内三种范畴的SSS的瞬时体积分布如图10所示。结合转速调制时血泵的速度变化情况可以看出，在0.5～1.5 s下，当泵速增大时，SSS为0～10 Pa的体积在减小，而SSS为10～100 Pa及100 Pa以上的体积在增大。这说明血泵内高SSS区域随着转速的增大而增大。将不同运行状态下三种范畴的SSS体积占比取均值，可以得到匀速状态下SSS在0～10 Pa之间的体积占比为70.32%，而正弦波下为70.87%，三角波下为70.82%，方波下为69.99%；匀速状态下SSS在10～100 Pa之间的体积占比为17.26%，而正弦波下为17.38%，三角波下为17.18%，方波下为17.90%；匀速状态下SSS大于100 Pa的体积占比为0.20%，而正弦波下为0.18%，三角波下为0.20%，方波下为0.19%。可见，与匀转速状态相比，三种速度调制波形对SSS的体积占比总体影响不大，因此不会显著增大血液损伤的风险，而这其中又以正弦波的极端非生理区域最小。

(a)SSS值0～10 Pa

(b)SSS值10～100 Pa

(c)SSS值大于100 Pa图10 血泵内SSS瞬时体积分布Fig.10 SSS instantaneous volume distributioninside the blood pump

2.3 速度

截取血泵内的一个截面来观察血液速度分布情况，截面的具体位置如图11所示。该截面上的血液速度分布如图12所示。通过叶轮处于最小转速阶段与匀速状态的对比可以明显看出，叶轮中心靠下的位置会在低转速时形成停滞区域，且正弦波和方波转速调制时的停滞区域更大；在加速阶段，血液速度明显增大，叶轮尾缘以及蜗壳和叶轮的边界处的漩涡变大；当叶轮转速最大时，叶轮外缘观察到了明显的流动分离区域，其中正弦波调制下的流动分离区域大于其他两种调制波形下的流动分离区域；在之后的转子减速阶段，流动分离区域面积明显减小，血液流动逐渐趋于稳定。

图11 截面位置Fig.11 Cross section position

3 讨论

CF-VAD作为一种可行的心力衰竭治疗手段，已经被广泛应用于临床中，并已证明它在可靠性等方面优于PF-VAD。然而，临床上不良事件发生的风险并未完全消除[18],潜在的影响因素包括血管搏动性降低、持续的心室容积卸荷以及CF-VAD产生的非生理剪应力等。在匀速条件下，CF-VAD可能降低血管搏动性，甚至导致主动脉瓣永久性关闭。尽管近几十年来对非生理搏动的潜在后果仍处于争论中，但目前最新的旋转血泵已经开始尝试在临床上利用转速调制等策略来增强CF-VAD旋转血泵的脉动性[19]。

脉动性指数(PI)是文献中常用来表征生理脉动的指标，通常定义为最大(峰)和最小(谷)流量之差与平均流量的比值。本研究中，对离心旋转血泵在四种运行模式中的PI值进行了分析，其中匀转速情况下的PI值为0.704，三种转速调制波形下的PI值分别为：正弦波1.466、三角波0.927、方波1.085。由此可见，正弦波调制可获得最大的PI值，且达到匀速情况PI值的2倍以上。

剪切应力是造成血液损伤的主要原因。CFD模拟结果表明，WSS和SSS随心脏周期和泵转速的变化而变化。叶轮转子附近的剪切应力水平与泵的转速大小密切相关。当泵转速增大时，WSS和SSS增大，反之则减小，即高转速会引起更大的血液损伤。对比三种转速调制波形的结果发现，方波调制时的WSS大于其他两种波形调制时的WSS，而速度调制对SSS的整体影响不大，正弦波调制具有最小的极端非生理区域。

(a)匀转速 (b)正弦波 (c)三角波 (d)方波图12 匀转速、正弦波、三角波和方波转速调制条件下血泵内的速度分布Fig.12 Velocity distribution inside the blood pump under constant, sinusoidal, triangular and square speed modulations

旋转血泵内血液的速度大小和方向也与泵转速密切相关，低转速会造成在转子中心靠下的位置形成停滞区域，高转速则会出现明显的漩涡和流动脱离区域。三种转速调制波形的对比结果表明，方波调制对泵内血液速度的不良影响最小。

4 结语

本研究的CFD模拟结果表明，旋转血泵在主动转速调制下运行是一种可行的方案，可以恢复循环搏动性，并且不会增加血液损伤的风险。三种速度调制波形中，正弦波可获得最大的流量脉动，且极端非生理SSS区域最小，是相对理想的调制波形。本研究得到的多种转速调制CFD模拟结果可为今后临床试验中旋转血泵的转速调制策略选择提供参考，甚至在未来的旋转血泵水力设计中转速调制工况可能会成为一种需要考虑的重要因素。由于受到一些客观条件的约束，本研究仍有以下问题有待进一步研究：

(1)本文采用的转速调制波形仍是一些基础波形，未来血泵生理控制条件下产生的转速波形更加复杂，可以在本研究的基础上结合血泵的生理控制算法进行系统仿真，进一步探究复杂生理控制下旋转血泵的流场变化情况。

(2)本文只研究了转速调制条件下旋转血泵内部的流场变化，从流场剪应力和定性的角度讨论了可能存在的血液破坏情况，未来可引入溶血模型，定量地研究复杂转速调制下旋转血泵的溶血情况。

(3)本研究的对象是自主研发的旋转血泵，转速调制对血泵流场影响的结果是否同样适用于其他血泵还有待进一步的验证，但是提出的方法对其他旋转血泵是同样适用的。