参考文献
第5章 心室辅助泵建模仿真与应用
植入心室辅助泵是挽救晚期心衰患者生命的最后一种方案。心室辅助泵目前主要包括轴流式血泵和离心式血泵两种。轴流式血泵具有体积小、植入方便等优点，但是存在溶血和凝血严重、卒中风险高等缺点，所以更多地考虑将其用于做功要求较低的右心室辅助泵的设计中。近些年随着介入式心室辅助泵的快速发展，轴流式血泵又扩展到新的应用领域。本章首先介绍轴流式血泵的建模和计算流体力学仿真的方法，然后将其用于植入单心室心脏病患者的全腔静脉肺动脉连接（TCPC）的结构中，分析其在不同工况下的溶血和血栓形成风险，为轴流式右心辅助泵的优化设计提供参考。
5.1 背景介绍
功能性单心室是一种典型的先天性畸形，其发病率在活婴中约为1∶6 500，占先天性心脏病的1%～3%［1］。患者没有两个发育完善的心室，动、静脉血液在承担收缩功能的主导心室混合后泵出，造成全身器官缺氧，并且出现心力衰竭、心律失常等症状，严重的患者甚至出现猝死。目前常见的单心室患者救治方案为三阶段的姑息手术：Norwood手术、Glenn手术与Fontan手术。第一阶段的Norwood手术宜在患者出生2月内实施，通过手术切除房间隔、横断并扎闭远心端肺主动脉，将近心端肺动脉连接到发育不全的主动脉上，以提高入肺的血流量促进肺的早期发育。第二阶段的Glenn手术一般在6个月至2岁时开展，此时肺血管阻力下降，通过手术去除Norwood手术中构建的主动脉-肺分流，然后将上腔静脉（SVC）连接至肺上形成分流通路。第三阶段的Fontan手术则在4～15岁进行，通过手术将输送下身血回心的下腔静脉（IVC）从原有的右心处切除，用人工血管进行延长并连接到肺动脉上。三阶段姑息手术完成后，上下腔静脉被直接连接到肺主动脉，构成了全腔静脉肺动脉连接生理结构，从体循环返回的所有静脉血直接流向肺部。这种结构虽然完全分离了动静脉血，但仍旧由单心室为体、肺两个循环提供动力，长期的动力不足会导致全身血运效率低下。极大的后负荷还会致使单心室发生心室肥大、搏动能力下降等病变，最终导致单心室心衰。
为此，De Leval［2］在1988年提出了动力型Fontan的概念，腔肺循环动力辅助装置（CPAD）也随之走入科研与临床的视野。CPAD按是否与血液的接触可分为管外动力类型（搏动泵）与管内动力类型（离心泵与轴流泵）。搏动泵通过容积变化周期性地泵出血液，模拟真实心室输出，但存在不宜植入和易感染的问题。离心泵和轴流泵则利用叶片旋转做功，为血液加速增压，输出连续性血流。并且相比之下，轴流血泵的体积更小，更容易植入，已经发展较为成熟。本研究参考Jarvik 2000和Throckmorton［3］等的轴流泵设计出一种适用于Fontan患者的轴流血泵，旨在增加Fontan患者的肺部灌注，减轻单心室负荷，提高患者术后的生活质量。
为了研究轴流血泵的辅助性能及为临床提供动力辅助Fontan手术的指导，本章将轴流血泵虚拟植入不同的TCPC结构中，运用计算流体力学方法，从血泵压升、血液的能量增量、TCPC内的流场特性与流量分配、血管壁损伤风险、溶血与血栓形成风险等方面对比轴流血泵植入不同TCPC结构后的辅助性能差异。本章还探讨了使用单相流和多相流方法进行植入轴流血泵的TCPC结构数值模拟的差异。相较于单相流模拟而言，多相流模拟引入了红细胞、白细胞和血浆，可以模拟出红细胞受剪切后聚集的生理现象，并且从溶血的本质出发更精确地评估溶血风险，但计算资源消耗更大。本章使用的所有方法均可以应用于各种心室辅助装置（如HeartMatⅢ、Impella等）和体外膜肺氧合装置（ECMO）的计算流体力学建模仿真研究。
5.2 计算流体力学仿真方法
计算流体力学在心室辅助泵的设计中发挥着举足轻重的作用，通过数值仿真可以进行手术规划并缩短血泵的开发周期。血泵的模型一般通过CAD软件（如Solidworks、ProE等）进行结构设计，特定病例的TCPC解剖结构由CT、MRI等医学影像在MIMICS等影像学重建软件中进行三维重建后获得。所得的Pump-TCPC模型装配后，进行网格的划分。血泵部分一般使用结构性的六面体网格，并在叶片转角处进行网格精细化；TCPC部分采用快速生成的非结构性四面体网格，并在血管壁面附近使用三棱柱网格进行细化，可更好地刻画近壁面处血流的边界层效应。所得网格数据在Ansys CFX或Fluent中进行求解器的设置与求解，数值仿真结果获得压力-流量曲线及各种血液动力学参数，对结果进行分析后指导血泵的设计参数优化，并获得符合特定病例的血泵输出参数。
5.2.1 旋转机械的数值仿真方法
血泵包括转动部分和静止部分，转动部分高速旋转会带动血液的转动，因此需要特殊的数值方法模拟这种转动。对这类问题有两种常用的数值模拟方法：一种是多重参考系方法，另一种是滑移网格法［4］。
1）多重参考系方法
多重参考系（multiple reference frame，MRF）方法是一种定常的数值模拟方法，它假设网格单元做匀速运动，适用于区域边界上点的运动基本相同的问题。该方法简单且适用性强，因此被广泛应用于旋转机械的仿真。使用多重参考系方法时，血泵内的流域被划分为两个子域：转动域和静止域。控制方程在每个子域内各自求解，各子域的信息在交界面上通过将速度换算为绝对速度的形式进行交换，并且使得两个子域的速度保持连续。多重参考系方法不能精确计算随时间变化的物理量的数值，它所获得的定常解被认为是非定常解对时间的平均值。
2）滑移网格法
滑移网格法（sliding mesh）是一种非定常计算方法，使用滑移网格法时血泵内区域也被划分为转动域和静止域，转动域内的网格分布较密，且随着离散的时间步沿旋转轴产生转动，在不同的时刻重新生成网格；而静止域则保持不动。在转动域与静止域的交界面上，使用两层交界面将不对应的节点分别进行求解，然后形成新的交界面。通过新交界面上的通量传递，实现两区域内的流场耦合。滑移网格法可以得到随时间变化的参数值。
5.2.2 湍流模型
正常人体内血液流动时的雷诺数不超过2 300，处于层流状态。而血泵转子高速转动带动血流快速流动时的雷诺数远远超过2 300，血液流动进入湍流状态，因此对于血泵内的血液流动需要考虑湍流。湍流是不规则、多尺度、结构复杂的流动，具有很强的扩散性和耗散性。尽管湍流运动非常复杂，但是仍然可以用非稳态的连续性方程和Navier-Stokes方程描述湍流的瞬时运动，然而求解这样一个复杂的方程组会花费大量的计算成本，因此需要进行简化处理。其中最常用的简化方法是基于统计平均建立起来的雷诺时均化方法，即在时间域上对流场物理量进行雷诺平均化处理，然后求解所得到的时均化控制方程。雷诺时均化方法计算效率较高，解的精度也基本可以满足工程实际需要，是流体力学领域使用最为广泛的湍流模拟方法。利用雷诺时均化方法创建的湍流模型包括k-∈模型、k-ω模型等［5］。
1）k-∈湍流模型
标准k-∈模型是从试验现象中总结出来的，它被提出之后就变成工程流场计算湍流主要的工具。模型中的湍动能k输运方程是一个通过精确推导得出的精确方程，耗散率∈方程是一个半经验公式。它适用范围广，在计算成本和计算精度方面较为平衡。但是标准k-∈模型是针对发展非常充分的湍流建立的，是一种针对高雷诺数湍流的模型。当雷诺数较低时，不适宜直接采用该模型而需要进行特殊处理，如引入壁面函数等。
2）k-ω湍流模型
k-∈模型无法准确预测带有逆压力梯度的边界层流动。为解决这一问题，另外一些湍流模型被提出，k-ω模型就是其中之一，ω是湍动能的比耗散率。k-ω模型对于k-∈模型在对一些湍流进行模拟时不够精确的问题有显著的改善，比如强曲率的流动、分离流及射流，但它的收敛难度较高，并且对求解时选取的初始值非常敏感。
在k-∈模型和k-ω模型基础上又进一步改进出了具有更高精度的Shear Stress Transport（SST）k-ω模型。凭借着更好的适用性和更高的计算精度，SST k-ω模型越来越多地被用于叶轮机械的CFD仿真。
5.3 轴流血泵植入全腔静脉肺动脉连接理想模型的数值仿真研究
Fontan循环有多种搭建方式。传统的搭建方式是将SVC和IVC直接连接到右肺动脉（RPA），形成单侧上腔静脉TCPC结构。另一种搭建方式是上腔来流的肺动脉双侧分流术式［6］，这种术式在传统术式的基础上将SVC上游的左无名静脉（LIV）截下，连接至左侧肺动脉，形成带有双侧上腔静脉的TCPC结构。相较于传统搭建方式，上腔来流的肺动脉双侧分流术式可以更好地实现左右肺血供的平衡，但是手术更加复杂，难度更大。本节对比轴流血泵植入两种TCPC结构后的血液动力学差异，为临床的动力Fontan手术选择提供指导。
5.3.1 模型建立
参考Jarvik 2000和Throckmorton［3］等的轴流泵，设计一种适用于Fontan患者的轴流血泵。该血泵通过手术植入患者的IVC中，为患者提供循环动力。植入手术的施行时间可以选择在病患构建出完整的Fontan循环之后，也可选择与病患的Fontan手术同步进行。轴流血泵结构如图5-1所示，共分为引流段（inducer）、叶轮（impeller）、扩散段（diffuser）和整流区域（straightener）四个部分。各部分的功能如下：引流段将下腔血流顺畅地引导进入血泵；叶轮转动加快血流流动，将旋转的机械能转化为血液的动能；扩散段的叶片与叶轮叶片逆向，将血流的动能转化为压力能；整流区域进一步规整流出血流，减少泵出的血液出现二次流的可能性。
图5-1 轴流血泵结构示意图
本节所用到的TCPC结构为一理想模型，各个血管的尺寸和空间关系是根据真实病例的解剖几何得到的。首先简化病例的生理解剖结构得到单侧上腔静脉TCPC结构，而后依据此模型，将左无名静脉截下，连接至左侧肺主动脉上，构建出双侧上腔静脉TCPC结构，最后在所获得的两种TCPC结构的IVC中分别植入轴流血泵，构建出相应的植入轴流血泵的单侧上腔和双侧上腔静脉理想TCPC模型，分别标注为Si-SVC（单上腔静脉术式）模型和Bi-SVC（双侧双向上腔静脉术式）模型，如图5-2所示。两个模型的具体几何尺寸如下：IVC、LPA和RPA的直径均为9.5 mm；Si-SVC模型中RIV和LIV的直径分别为8.14 mm和6.74 mm；Bi-SVC模型中RSVC和LSVC的直径与RIV和LIV的直径保持一致，分别取为8.14 mm和6.74 mm；上、下腔静脉之间的偏移量设置为19 mm，对应IVC直径的2倍［7］。
图5-2 植入轴流血泵的单侧上腔与双侧上腔TCPC模型
5.3.2 网格划分
本研究通过Ansys ICEM生成流域网格。轴流血泵内部流域采用结构性六面体网格划分，这类网格可以准确刻画规整血泵内的血液流动［8-10］。血管流域则采用非结构性的四面体网格划分。四面体网格可以很好地刻画不规则的血管形状，减少畸形网格生成的可能性。此外，在血管壁附近生成五层边界层网格以对近壁面处的流场进行精确刻画，使壁面上的血液动力学参数的计算更为准确。经过网格独立性检验后，最终采用的轴流血泵与理想TCPC的网格数量见表5-1。
表5-1 植入轴流血泵的单侧上腔与双侧上腔静脉TCPC模型网格数量
5.3.3 边界条件、材料属性和求解设置
本节中的边界条件均为定常边界条件。轴流血泵入口以及左、右上腔静脉均设为流量入口边界条件，流量大小根据心输出量（CO）计算得到。轴流血泵入口流量和上腔静脉总流量分别为CO的60%和40%，LSVC与RSVC的流量分别为上腔静脉总流量的40%和60%［6］。考虑到不同病患的心输出量存在较大差异，本研究选取五个不同的CO流量进行后续分析，即2 L/min、2.5 L/min、3 L/min、3.5 L/min、4 L/min。为了捕捉LPA和RPA处的反流，LPA和RPA设置为开放边界条件。考虑到不同转速下轴流血泵的输出特性不同，LPA和RPA处的静压值与血泵转速一一对应，分别为8 mmHg对应于2 000 r/min、16 mmHg对应于3 000 r/min、24 mmHg对应于4 000 r/min、32 mmHg对应于5 000 r/min、40 mmHg对应于6 000 r/min。TCPC的血管壁面、轴流血泵的壳体壁面、叶片表面均设置为刚性壁，光滑，无滑移。根据叶片的螺旋方向，轴流血泵的叶轮流域被设置为逆时针旋转区域，其他流域为静止区域。转动区域和静止区域通过通用网格界面（general grid interface，GGI）进行耦合［11］。
模拟的流体介质为血液。血泵中是高剪切环境，因此使用不可压缩牛顿流体表征血液。血液的密度、动力学黏度和比热容分别取为1 050 kg/m3、0.003 5 Pa·s与3.594 kJ/（kg·K-1）［12］。
本节中包括稳态模拟和瞬态模拟，计算血栓形成风险时采用瞬态模拟，除此之外均采用稳态模拟。瞬态模拟时，时间步长取为0.000 33 s，计算150步，每10步保存一次，对最后10个保存步结果进行结果分析。两种算例选用Ansys CFX中的RANS求解器和标准k-∈湍流模型进行求解，此求解器设置在VAD的研究中运用广泛［8，13-15］，收敛判据设为速度和压力的计算残差小于1×10-4。
5.3.4 评估指标
评价轴流血泵植入单、双侧上腔静脉TCPC结构后性能差异包括能量增量、溶血风险和血栓形成风险三个方面，相应评价指标的介绍如下。
1）能量增量
采用轴流血泵最核心的目的是解决Fontan患者循环动力不足的难题，因此血液经过TCPC后的能量增量是动力性Fontan手术设计的一大关注点。本研究所用能量增量的计算公式如下［16］：
式中 Ptotal——总压；
Q——体积流量，inlet和outlet代表模型的入口和出口；
Pstatic——静压；
ρ——流体密度；
u——流速。
2）溶血风险
轴流血泵在运转过程中与血液直接接触，高速旋转的叶片导致极高的剪切应力，大大增加血细胞被破坏的风险。根据前人对多款血泵的溶血风险研究，我们选择细胞破坏因子（BDI）作为溶血风险评估指标［17］，该指标通过标量应力和暴露时间在流线上的积分评估血细胞被破坏的风险，其计算公式如下：
式中 ΔT——应力暴露时间；
σ——标量应力，由Bludszuweit应力公式计算，表述如下：
根据前人研究［8，1819］，本节中细胞破坏预警的阈值设定为：细胞破坏因子BDI=2%，血液在模型内的滞留时间PRT=0.6 s。
3）血栓形成风险
血栓的形成也是血泵运转中不可避免的一个问题。前人的研究表明相对滞留时间（relative resident time，RRT）可以定量评估扰动流的状态，其值的大小与血管内皮细胞破坏和血栓形成的风险有关［20］。因此本研究使用相对滞留时间评估血栓形成风险。计算公式如下：
式中 wssi——血管壁上的瞬时切应力矢量；
T——计算时长。
时间平均壁面剪应力（TAWSS）和振荡剪切指数（OSI）分别用以表征计算周期内壁面剪切应力的平均幅值和变化程度。
5.3.5 结果分析和比较
下面对轴流血泵植入单、双侧上腔静脉TCPC结构后的水力学性能差异以及溶血风险和血栓形成风险进行定量分析和比较。
1）压力-流量输出特性分析
表5-2为血泵植入单侧上腔和双侧上腔静脉TCPC模型后在3 000 r/min和5 000 r/min两个转速下的压力-流量输出特性。在3 000 r/min的转速下，随着血泵入口流量（即CO的60%）的增大，血泵在Si-SVC和Bi-SVC模型中压升的差距有所增加，但最大偏差未超过0.5%；而在5 000 r/min的转速下，血泵在两种模型中压升的相对差距均为0。因此，TCPC上腔术式的改良对轴流血泵的升压能力未产生明显影响。
表5-2 3 000 r/min和5 000 r/min下轴流血泵在Si-SVC和Bi-SVC模型中的压升
基于上述结论，接下来选取Bi-SVC模型进行轴流血泵压力-流量输出特性的分析，所得结果如图5-3所示。由图可知，轴流血泵在CO为2～4 L/min、转速为2 000～6 000 r/min时，可以提供1～36 mmHg的压升。图5-3a展示了轴流血泵的升压能力和转速之间的关系，随着转速的升高，其升压能力增大，符合轴流血泵设计预期。图5-3b展示了轴流血泵的升压能力和流量的关系，随着输入流量的增加，血泵的升压能力随之下降。可将轴流血泵比作电路中的电源，电源本身存在一定阻抗，因而随着电流的增大，在该阻抗上所消耗的电压逐渐升高，从而使得血泵所提供的压升略微减小。同时可以观察到，随着血泵转速的下降，随着流量升高，血泵压升的下降趋势更加明显，可估计血泵自身阻抗亦同转速相关，高转速下，阻抗减小，低转速下，阻抗增大。
图5-3 Bi-SVC模型中轴流血泵压升、入口流量和转速三者之间的关系
2）能量增量指标分析
表5-3为血泵在3 000 r/min和5 000 r/min两个转速下工作时，血液流经Si-SVC和Bi-SVC两个模型后获得的能量增量。从表5-3中可得，当血泵在两个转速下工作时，血液流经Bi-SVC后获得的能量增量比血液流经Si-SVC后获得的能量增量更大，并且相对差距随着血泵入口流量的增加而增加。当CO为4.0 L/min、转速为3 000 r/min和5 000 r/min时，血液流经两个模型后的能量增量的相对差距分别是5.25%和5.16%。由此可得，上腔分流与否对于全局能量增量的影响较为明显，且Bi-SVC模型要全面优于Si-SVC模型，并且在血泵转速较低时优势更为显著。
表5-3 3 000 r/min和5 000 r/min下血液流经Si-SVC和Bi-SVC模型后的能量增量
依照此结论，选取Bi-SVC模型讨论输入流量、血泵转速和能量增量之间的关系，所得结果如图5-4所示。随着轴流血泵转速的增加，能量增量呈显著上升趋势。随着输入流量的增大，能量增量基本呈升高趋势，但在血泵2 000 r/min和3 000 r/min两个工况中有所出入。特别是血泵以2 000 r/min低转速运转时，随着输入流量的增大，能量指标呈明显下降趋势，并且在较大CO流量4 L/min时，出现能量增量为负的情况，即血泵无法弥补TCPC结构的能量损耗。
图5-4 Bi-SVC模型能量增量与转速和入口流量的关系
3）下腔静脉流量分配分析
临床上在Fontan或Glenn手术时采用上腔分流的术式，主要为达到左右肺流量合理分配的目的。此种合理分配可将下腔引入的肝因子等有利于肺发育的微粒合理地分配进入两侧肺中，从而促进左右肺的发育。在动力辅助的Fontan循环中，此种术式是否可得到相同结论，流量分配的效果具体如何有待验证。
由上述结果可知，轴流血泵在CO流量为3 L/min、叶轮转速4 000 r/min的工况下压力-流量特性与能量指标均表现良好。因而根据此工况下的模拟结果绘制整个模型流场的流线图，并给出流量的分配情况，如图5-5所示。从Si-SVC和Bi-SVC模型流场及下腔流量分配可知，Si-SVC模型中下腔静脉来流的23%通过TCPC吻合口向右分配进入RPA，77%的下腔静脉来流进入LPA，存在54%的分配差。同时流场分离位置靠近右侧上腔，很好地缓冲了该处因流场汇聚而带来的紊流，与前人研究结论一致［9］。而Bi-SVC模型中下腔静脉流量的分配更加合理，下腔静脉来流较为均衡地进入左右两侧肺动脉，42%分配至RPA，58%流入LPA，分配差值降低至16%。同时下腔来流的分离点更靠近吻合口中心位置，来流冲击的效果能得到更为平衡的分配，此结果与现有Fontan上腔分流的临床与仿真模拟结果一致［21］。由此可得，从全局流场和下腔流量分配来看，Bi-SVC要优于Si-SVC。
图5-5 Si-SVC和Bi-SVC模型流场与下腔静脉流量分配特性
4）溶血风险分析
选取CO为3 L/min时，血泵三个工作转速（2 000 r/min、4 000 r/min、6 000 r/min）下的Si-SVC和Bi-SVC的稳态流场进行溶血风险分析。从血泵入口释放300个粒子，对这些粒子的标量应力和暴露时间在流线上进行积分，计算得到粒子BDI和PRT。所得结果见表5-4所列。
表5-4 不同转速下Si-SVC和Bi-SVC模型内的细胞破坏分析（CO=3 L/min）
表5-4显示各种情况下的BDI和PRT大部分低于细胞破坏阈值（BDI=2%，PRT=0.6 s）。唯一例外如下：在叶轮转速为6 000 r/min时，Si-SVC模型最大BDI为1.94%，低于阈值，而最大PRT为0.66 s，超出阈值10%。进一步分析可知，随着转速的增加，两个模型中的平均BDI和最大BDI会随之增大，而平均PRT和最大PRT会逐渐减小。但Si-SVC的最大PRT却呈现出截然相反的变化趋势，随着转速的增加，最大PRT反而增加。配合Si-SVC流线可以看到，在下腔来流分离区域，流场的左右分配不均，向右分离的小部分流场同上腔来流呈现明显竞争关系，导致了部分粒子通过此处时的运动速度减慢，从而增加了最大PRT这项指标。而Bi-SVC模型该区域的流速明显大于Si-SVC模型，使得粒子能以较快速度通过。对比Si-SVC和Bi-SVC模型的平均BDI和平均PRT可知，两个模型的溶血风险没有明显区别。
图5-6 不同转速Bi-SVC模型细胞破坏因子（BDI）统计
以下对Bi-SVC模型中粒子的BDI进行统计分析，结果如图5-6所示。在2 000 r/min的低转速下，有283个粒子（94.4%）的BDI取值低于0.12%；在轴流血泵4 000 r/min的转速下，276个粒子（92.0%）的BDI低于0.4%；6 000 r/min的高转速下，300个粒子中有280个粒子（93.4%）的BDI低于0.8%。粒子BDI的分布也展示出随着转速增高，溶血风险逐渐增高。
5）血栓形成风险分析
选取CO为3 L/min、血泵转速为4 000 r/min的瞬态模拟求解结果对血栓形成风险进行分析。两个模型TCPC部分的RRT分布云图如图5-7所示。高RRT区域主要出现在下腔来流向左右肺分流的区域，此外上腔静脉与肺动脉的吻合口也出现高RRT，但区域较小。相比而言，Bi-SVC模型的高RRT区域面积较Si-SVC模型要小得多，表明植入血泵的Si-SVC模型具有更低的血栓形成风险。结合图5-5分析，下腔静脉吻合口处的下腔血流分离以及上腔静脉吻合口处的上、下腔血流冲击均会导致血流流速降低。低流速区域易出现流体滞留的情况，因此导致对应壁面处RRT的升高。
图5-7 Si-SVC和Bi-SVC模型TCPC区域的血管壁面相对滞留时间（RRT）分布（血泵转速4 000 r/min，CO为3 L/min）
5.3.6 研究结论
本节从压力-流量特性、能量增量、流场特性与流量分配、溶血与血栓风险等方面对植入轴流血泵后的单、双侧上腔TCPC结构进行血液动力学分析，确定了TCPC上腔双侧双向分流术式作为Fontan动力循环辅助术式更有优势，为临床研究提供一定的参考和选择。研究主要获得以下结论：
（1）当心输出量为2～4 L/min，血泵转速为2 000～5 000 r/min时，血泵提供的压升为1～24 mmHg，可以为患者提供有效的动力辅助。
（2）Bi-SVC模型与Si-SVC模型仿真结果对比发现，血泵的压力-流量特性没有明显差别。但对于能量增量而言，Bi-SVC模型更有优势。
（3）从流场和下腔血流分配分析可得，Bi-SVC模型比Si-SVC模型有优势，Si-SVC模型和Bi-SVC模型左右肺血流的流量分配差分别为54%和16%。
（4）Bi-SVC模型与Si-SVC模型在溶血风险上没有明显差别，但Bi-SVC模型的血栓形成风险显著降低。
5.4 轴流血泵植入患者特异性全腔静脉肺动脉连接模型的单相流仿真研究
上节研究了理想TCPC模型植入轴流血泵后的血液动力学参数的变化，数值仿真采用了定常边界条件。真实TCPC几何结构与理想TCPC还是有较大差异，而且虽然Fontan患者上、下腔静脉处压力和流量的脉动性较正常人大幅度下降，但真实的压力和流量依然有一定的波动特征。本节选取两个真实Fontan病例的患者特异性TCPC几何模型，并采用实测的压力和流量波形作为边界条件进行仿真分析。
本节选取的两个真实TCPC几何模型代表两种典型的下腔静脉连接方式——心内直管术式（Fontan-straight）和心外弯管术式（Fontan-curved）。心内直管术式将人工血管穿过右心，直接连接下腔静脉和肺动脉；心外弯管术式则绕过右心，使用弯曲的人工血管连接下腔静脉与肺动脉。心内直管术式操作较为复杂，但血液流经直管的能量损失小。心外弯管术式操作相对简单，但血流的能量损失较大。本节探讨轴流血泵与两种术式搭配的差异，为动力Fontan中下腔静脉连接方式的选择提供指导意见。此外，为了模拟真实的右心输出，本节还对可变转速血泵进行了研究，以不同的转速来模拟右心室的各个时期，对比分析恒定转速泵和可变转速泵辅助Fontan循环的差异，探讨可变转速设计的可行性。
5.4.1 模型建立
本节的轴流血泵模型沿用上节，所用到的TCPC模型则来源于真实Fontan病例。直管道TCPC模型来源于接受双侧上腔静脉-肺动脉连接Fontan手术的13岁女性患者，弯管道TCPC模型来源于接受下腔外管道Fontan手术的14岁男性患者。植入轴流血泵的TCPC模型建立流程如下：首先使用西门子DSCT采集两名患者的原始DICOM数据。而后，将所得DICOM数据在MIMICS 17.0软件中进行三维重建，获得初步的TCPC几何模型。下一步在Geomagic 13.0中根据轴流血泵的植入位置将TCPC几何模型分为上腔静脉-肺动脉模型和下腔静脉-肝静脉模型。最后在Solidworks 2015中完成轴流血泵和两部分血管模型的装配，形成直、弯管Pump-TCPC模型，血泵的植入效果如图5-8所示。
图5-8 真实病例直、弯管解剖结构及直、弯管Pump-TCPC模型效果图
5.4.2 网格划分
本节中的流域网格同样通过Ansys ICEM生成，网格划分方法与上节相同。最终采用的网格数量见表5-5。
表5-5 植入轴流血泵的直、弯管道TCPC模型网格数量
5.4.3 边界条件、材料属性和求解设置
本节采用从14岁男性患者Fontan术后实测的流量和压力波形作为直、弯管Pump-TCPC中各个出入口的边界条件。SVC和IVC设置为流量入口，波形如图5-9a所示，并且直管Pump-TCPC模型中LSVC和RSVC的流量比为4∶6。为了捕捉反流情况，三根肝静脉和左、右肺动脉设置为开放边界条件，压力波形如图5-9b所示。其余边界条件以及材料属性均与上节保持一致。
图5-9 实测Fontan循环的流量-时间曲线与压力-时间曲线
本节采用瞬态模拟，病例的脉搏为75次/min，所以将心动周期定为0.8 s。经检验，初始化仅影响第一个心动周期前20%的结果，第二个心动周期的结果足够准确，因此设定总计算时长为1.6 s。时间步长为0.000 5 s，每10步保存一次，共320个保存步，选取第二个心动周期的模拟结果进行分析讨论。选用Ansys CFX中的RANS求解器和标准k-∈湍流模型进行求解，收敛判据设定为RMS＜1×10-4。
5.4.4 可变转速曲线设定及其生理意义
图5-10 一个心动周期内可变转速泵转速-时间曲线（红线为可变转速泵转速波形，蓝线为恒定转速泵转速基准线）
一个心动周期内，轴流血泵的恒定转速和可变转速输出曲线的设定由图5-10给出。根据上节的研究可得，轴流血泵在4 000 r/min转速运转下的压力-流量输出特性、能量增量指标均能符合Fontan循环的输出需求，因此将恒定转速泵的转速设为4 000 r/min。可变转速泵的转速则设置为一段梯形波，转速范围为2 000～6 000 r/min。0～0.3 s为高速旋转段，转速为6 000 r/min，模拟右心等压收缩射血；0.3～0.4 s为快速减速段，从高转速向低转速转变，模拟右心的等容舒张；0.4～0.7 s为低速旋转段，转速为2 000 r/min，模拟右心的等压舒张；0.7～0.8 s为快速加速段，从低转速向高转速转变，模拟右心的等容收缩，为下个心动周期的泵血做准备。此波形设置为Fontan循环提供较为接近生理情况的搏动性血流，并且保证了一个心动周期内恒定转速泵与可变转速泵具有相同的平均转速，使二者为Fontan循环提供的机械动力一致。
5.4.5 结果分析和比较
下面对轴流血泵植入患者特异性TCPC结构后的水力学性能差异以及溶血风险和血栓形成风险进行定量分析和比较。
1）血泵压升和整体能量增量分析
取血泵在恒定转速下工作时，直、弯管Pump-TCPC模型的模拟结果分析，得到两个模型中的血泵压升曲线和血液流经两个模型后的能量增量曲线，如图5-11所示。恒定转速泵在直管道TCPC结构中的压力输出如图5-11a蓝色虚线所示。在一个心动周期内，轴流血泵所提供的压升在4～6 mmHg。在0.3 s和0.55 s附近取得两个极大值，在0.05 s和0.45 s附近取得两个极小值。血液流经直管Pump-TCPC后的能量增量如图5-11a红色实线所示。在一个心动周期内，血液所获得的能量在20～40 mW。同时，能量增量指标取到极值的时刻与血泵压升极值的时刻保持一致：两个极大值出现在0.3 s和0.55 s附近，极小值出现在0.05 s和0.45 s附近。
图5-11 轴流血泵植入真实Fontan循环后压力-时间和能量-时间输出特性曲线（蓝色虚线为压力波形，红色实线为能量波形）
恒定转速泵在弯管道TCPC模型中提供的压升也在4～6 mmHg，如图5-11b蓝色虚线所示。极大值出现在0.3 s、0.55 s和0.65 s时刻，极小值出现在0.05 s和0.4 s附近。弯管Pump-TCPC模型内，能量增量随时间的变化规律如图5-11b红色实线所示。整个心动周期内，能量增量在±10 mW浮动，时均值为0.5 mW，极大值出现在0.3 s和0.65 s附近，极小值出现在0.05 s和0.4 s附近。
对比直、弯管Pump-TCPC模型的压升曲线和能量增量曲线可知，弯管道模型比直管道模型的能量损耗多30 mW左右。但因为本节中直、弯管道TCPC结构的上腔静脉连接方式也有所差别，因此有必要分析上腔分流与否和上、下腔吻合口处是否存在对冲，对能耗产生的影响。根据上节结果可判断上腔分流与否对能量增量造成的影响较小，约为2%；有关Fontan术式中上、下腔对接偏移量的研究表明，上、下腔直接对接造成的能量损耗会明显高于上、下腔来流存在偏移时的能量损耗，但这种差异约为3 mW［22］。因此，上腔分流与否和上、下腔吻合口是否存在对冲都不是造成直、弯管Pump-TCPC能量增量差异的主要原因。弯曲的人造血管部分能耗太大才是造成弯管Pump-TCPC模型能量增量远小于直管Pump-TCPC模型的主要原因。
2）溶血风险分析
对比血泵分别在恒定转速和可变转速下工作时，直、弯管Pump-TCPC模型内的溶血风险。选取0.2 s、0.35 s、0.6 s和0.75 s四个时间点（分别代表可变转速泵转速的四个阶段）进行分析。从每个模型的IVC入口处释放300个粒子，统计每个粒子从进入Pump-TCPC到流出所经过的时间以及在此时间内标量应力的积分，计算PRT和BDI。
直管道TCPC分别植入恒定转速和可变转速泵后不同时刻的BDI统计分布如图5-12和图5-13所示，表5-6展示了两种情况的对比。
图5-12 不同时刻恒定转速泵植入直管道TCPC结构的细胞破坏因子（BDI）统计
图5-13 不同时刻可变转速泵植入直管道TCPC结构的细胞破坏因子（BDI）统计
表5-6 不同时刻直管道模型植入恒定转速泵和可变转速泵的细胞破坏分析
分析表5-6可知，0.2 s时可变转速泵处于6 000 r/min的高转速运转状态，此时植入可变转速泵模型内血液流速快，粒子滞留时间短，因此平均PRT和最大PRT比恒定转速泵模型更小。但高转速也带来了高切变率和高应力环境，血细胞受剪切冲击增大，因而平均BDI和最大BDI较大。0.6 s时可变转速泵处于2 000 r/min的低转速运转状态，此时可变转速泵模型对应的平均PRT和最大PRT比恒定转速泵模型更大，但是平均BDI和最大BDI更小。0.35 s和0.75 s可变转速泵和恒定转速泵的转速均为4 000 r/min。0.35 s时刻可变转速泵正由高转速向低转速快速转变，可变转速泵模型内血液流速依旧较高，血细胞滞留时间短，可变转速泵模型的平均PRT和最大PRT均小于恒定转速泵模型。此外，该时刻可变转速泵模型平均BDI比恒定转速泵模型高，这也说明了此时流场正由高应力状态向低应力状态过渡。但最大BDI出现了恒定转速泵模型（0.77%）高于可变转速泵模型（0.73%）的情况。对比图5-12a和图5-13a可知，虽然恒定转速泵模型内大部分粒子处于低BDI水平，但有少量的粒子（1%）游离于较低风险区域之外，从而出现了这一反常的现象。同样对0.75 s时刻进行分析不难得到，此时可变转速泵模型内的血液正由低流速向高流速过渡，血细胞运动速度较缓，滞留时间长，因此平均PRT和最大PRT较大而平均BDI和最大BDI较小。
进一步分析表5-6可知，四个时刻恒定转速泵模型和可变转速泵模型内的平均BDI和平均PRT均小于风险阈值（2%和0.6 s），从总体上看，直管道可变转速泵模型和直管道恒定转速泵模型的溶血风险均较低。
不同时刻弯管道恒定转速泵和可变转速泵模型对应的BDI统计分布如图5-14、图5-15所示，对比表格见表5-7。同样可以发现在0.2 s和0.35 s，恒定转速泵模型的PRT指标较小而BDI指标较大。在0.6 s和0.75 s，恒定转速泵模型的PRT指标较大而BDI指标较小。仅在0.35 s的最大BDI指标和0.6 s的最大PRT指标呈现相反趋势，考虑为依然是极少数粒子导致最大指标的异常。此外，四个时刻恒定转速泵模型和可变转速泵模型内的平均BDI均小于风险阈值（2%），平均PRT大部分小于风险阈值（0.6 s）。可变转速泵模型在0.6 s和0.75 s以及恒定转速泵模型在0.75 s时对应的平均PRT大于风险阈值。从总体上看，弯管道恒定转速泵模型和弯管道可变转速泵模型的溶血风险依然较低。
图5-14 不同时刻恒定转速泵植入弯管道TCPC结构的细胞破坏因子（BDI）统计(https://www.daowen.com)
图5-15 不同时刻可变转速泵植入弯管道TCPC结构的细胞破坏因子（BDI）统计
表5-7 不同时刻弯管道模型植入恒定转速泵和可变转速泵的细胞破坏分析
3）血管损伤和血栓形成风险分析
选取TAWSS、OSI和RRT评价直、弯管道恒定转速泵和可变转速泵四个模型的血管损伤与血栓形成风险。TAWSS和OSI过大表征较高的血管损伤风险，而RRT过大则表征此处的血液滞留情况更严重，意味着更高的血栓形成风险。
直管道恒定转速泵和可变转速泵模型的TAWSS、OSI和RRT云图如图5-16～图5-18所示。TAWSS异常的区域主要分布在血管与血泵的接口处以及上腔静脉与肺动脉的吻合口区域，其中可变转速泵模型异常TAWSS区域明显更大。OSI异常的区域主要分布在人工血管与肺动脉的吻合口处、上腔静脉处以及LHV和RHV处，同样可变转速泵模型的异常区域更大。RRT分布图所显示的低流速阻滞区域主要分布于LHV和MHV分岔口处，并且恒定转速泵模型的低流速阻滞区域更大。
图5-16 血泵植入后恒定转速泵和可变转速泵两种工况下直管道血管壁上TAWSS的分布
弯管道恒定转速泵和可变转速泵模型的TAWSS、OSI和RRT云图如图5-19～图5-21所示。弯管道模型中三个指标的异常分布情况与直管道模型相似。并且相比直管道模型而言，弯管道模型中TAWSS、OSI异常的区域更大，表明弯管道模型血管壁损伤风险更高。其中直、弯管道模型上腔区域TAWSS和OSI的差别考虑也与上腔术式有关，弯管道模型采用单侧上腔，而直管道模型采用双侧上腔。结果表明双侧上腔术式很好地降低了上腔畸高的TAWSS，并且使得上腔的流动振荡性减弱。此结论再次证明了双侧上腔术式的优越性。直、弯管道模型整体RRT指标差异并不明显，表明直、弯管道模型的血栓形成风险差异不大。
图5-17 血泵植入后恒定转速泵和可变转速泵两种工况下直管道血管壁上OSI的分布
图5-18 血泵植入后恒定转速泵和可变转速泵两种工况下直管道血管壁上RRT的分布
图5-19 血泵植入后恒定转速泵和可变转速泵两种工况下弯管道血管壁上TAWSS的分布
图5-20 血泵植入后恒定转速泵和可变转速泵两种工况下弯管道血管壁上OSI的分布
图5-21 血泵植入后恒定转速泵和可变转速泵两种工况下弯管道血管壁上RRT的分布
此外，在同一TCPC结构中，与恒定转速泵模型相比，可变转速泵模型的TAWSS更大，OSI更大，流场呈现更高的振荡性，血管壁损伤风险更大。结合图5-22分析可知，在低转速时段，可变转速泵模型的肝静脉会出现明显反流，造成整个心动周期内下腔流场振荡加剧，产生流体迟滞的同时增加了下腔肝静脉处内皮细胞的破坏风险，因此需对可变转速泵的输出特性进行一定的调整与优化。
图5-22 0.65 s时刻直管道和弯管道模型中的流场（0.65 s对应可变转速泵的低转速时刻，红色流线表示上腔静脉来流，蓝色流线表示下腔静脉来流，灰色流线表示肝静脉来流）
5.4.6 研究结论
本节从血泵的压升特性、血液的能量增量、溶血风险、血管壁面损伤风险以及血栓形成风险等方面对植入恒定转速泵和可变转速泵后的直、弯管道TCPC模型进行分析，得到以下结论：
（1）轴流血泵在4 000 r/min的恒定转速下工作时，在直、弯管Pump-TCPC模型内均能为血液提供4～6 mmHg的压升。但血液流经直管道恒定转速泵模型所获得的能量比流经弯管道恒定转速泵模型所获得的能量多30 mW左右。分析表明弯曲的人工血管是造成能量消耗较大的主要原因。
（2）恒定转速泵和可变转速泵分别植入直、弯管道TCPC结构后，直管道恒定转速泵模型和可变转速泵模型内平均BDI和平均PRT均小于风险阈值（2%和0.6 s），弯管道恒定转速泵模型和可变转速泵模型内的平均BDI均小于风险阈值（2%），平均PRT大部分小于风险阈值（0.6 s）。从总体上看，恒定转速泵和可变转速泵植入直、弯管道TCPC结构后的溶血风险均较低。
（3）无论血泵处于恒定转速还是可变转速下工作，与直管道模型相比，弯管道模型的TAWSS和OSI异常区域更大，表明血管损伤风险较高，但RRT整体分布差异不明显，因此血栓形成风险差异不大。
（4）可变转速泵虽然可以模拟真实右心活动，但现有的转速调节方式有所不足。可变转速泵在低转速下工作时会导致肝静脉处出现反流，进而导致一个心动周期内的此处的流场振荡更加剧烈，血管壁损伤风险更大。可在转速的改变幅度、平缓转速变化以及匹配转速与边界条件这三个方面，进一步调整设计和优化。
5.5 轴流血泵植入患者特异性全腔静脉肺动脉连接模型的多相流仿真研究
血泵通过叶轮旋转做功，将机械能转化为动能。但叶轮的高速旋转会产生高剪切流场，进而导致红细胞聚集，表观黏度发生改变［23］。大量针对血泵的研究忽略这一现象，造成评价溶血风险时产生误差［2426］。也有少量研究关注这一现象，Bhavsar［27］等在对血泵的研究中采用不同的全血黏度进行模拟，发现高黏血流确实会对血泵的压升特性和溶血风险产生一定的影响。
然而传统的单相流仿真无法模拟红细胞的黏度改变，多相流仿真则为解决此问题提供了一个可行的思路。Jung等［23］在3D虚拟血管通路中采用多相流Euler-Euler方法，对血管中的血细胞（红细胞、白细胞和血小板）与血浆进行多相流模拟，分析血细胞的力学行为。
本节通过Euler-Euler方法进行多相流仿真，引入红细胞，模拟其在血泵内的切变率和随之改变的黏度。并且通过红细胞所受剪切效应的积累，直接分析溶血风险，对现有的未单独考虑红细胞受力情况的细胞破坏模型进行更新，提出新的溶血风险评价指标：多相流细胞破坏因子（mBDI）和多相流粒子滞留时间（mPRT）。为了证实多相流方法的可行性和优越性，本节针对5.4节所建立的Pump-TCPC模型进行多相流仿真，使用多相流细胞破坏模型评估溶血风险，将评价结果与5.4节使用单相流仿真和单相流细胞破坏模型得到的评价结果进行对比。
5.5.1 模型建立、网格划分、边界条件和求解设置
本节中的算例包括稳态模拟和瞬态模拟，所采用的几何模型及网格均沿用5.4节的设置。在进行多相流方法验证时以单独的轴流血泵为研究对象，采用稳态模拟，血泵入口设置为1.8 L/min的流量边界条件，出口设置为40 mmHg的压力边界条件，血泵转速设置为4 000 r/min。在应用多相流细胞破坏模型评估溶血风险时以直、弯管Pump-TCPC模型为研究对象，采用瞬态模拟，血泵转速同样设置为4 000 r/min，其余边界条件与5.4节保持一致。瞬态模拟时总计算时长为1.6 s，时间步长为0.000 5 s，每10步保存一次，共320个保存步，选取0.8～1.6 （s对应第二个心动周期）的模拟结果进行分析。两种算例均使用Ansys CFX中的RANS求解器，调用标准k-∈模型与多相流体积分数求解模块进行耦合求解，收敛判据设定为RMS＜1×10-4。
5.5.2 材料属性
在多相流模拟中，血液被模拟为由血浆、白细胞和红细胞组成的多相混合物。血浆为流体，白细胞和红细胞分别用直径11 μm和7.2 μm的固体球形颗粒表征。三种组分的密度、黏度、比热容以及初始体积分数见表5-8，其中红细胞黏度随切变率变化的关系使用Quemadacite［28］模型表示，计算公式如下所示：
式中 rp，rrbc——血浆和红细胞的体积分数；
——红细胞的切变率；
k0，k∞——
→0和
→∞时的固有黏度；
——标准切变率；
μp——血浆的黏度。
这些参数均由试验测量得到［29］，见表5-9。
表5-8 单相流与多相流血液的成分对比
表5-9 Quemadacite黏度模型的参数
5.5.3 多相流细胞破坏模型
单相流细胞破坏模型根据标量应力和应力暴露时间在流线上的积分得到细胞破坏因子。多相流细胞破坏模型借鉴相同的积分形式，但直接根据红细胞所受剪切效应的积累评估溶血风险。
Blackshear［30］等通过溶血测试，提出了红细胞所受切应力、切应力暴露时间和溶血风险的关系：
式中 τrbc——切应力；
t——暴露时间；
C1——比例常数。
参考式（5-3），将多相流细胞破坏因子表达为
其中，指数因子α的值与单相流细胞破坏模型保持一致，取为0.765。根据模型的试算，取放大因子C=1×10-6，使得由式（5-10）计算所得的细胞破坏因子与使用单相流细胞破坏模型所求得的细胞破坏因子一致，因此Dm的风险阈值与单相流细胞破坏因子D的风险阈值一致，仍取为2%。多相流细胞破坏因子D（mmBDI）最终表达为
式中 μrbc——用Quemadacite模型计算的红细胞黏度；
t——曝光时间；
γ·
rbc——红细胞的标量切变率，表达式为
此外，对应于单相流细胞破坏模型中的PRT，在多相流细胞破坏模型中提出多相流粒子滞留时间（mPRT）代表粒子从入口到出口所经历的时间。由于mPRT和PRT的物理含义相同，所以mPRT的阈值也取为0.6 s，与PRT的阈值保持一致。
5.5.4 结果分析和比较
取轴流血泵在4 000 r/min转速下工作时，单相流与多相流模拟得到的流场进行对比。多相流模拟用体积分数来代表多相流组分，各组分的速度完全一致，在此选择红细胞速度指标绘制多相流模型流场，并与单相流结果对比，如图5-23所示。图5-23a为以红细胞速度为标识的多相流流场，图5-23b为单相流流场，二者几乎相同，佐证了轴流血泵多相流模拟的可行性。
接下来以轴流血泵内血细胞浓度的分析验证多相流模拟的可行性。血细胞的标准浓度由该组分在血泵各处的体积分数除以该组分在血泵入口处的体积分数得到（红细胞和白细胞在血泵入口的体积分数分别是45%和0.5%），计算公式如下：
式中 rrbc，rwbc——红细胞和白细胞在各处的体积分数；
Crbc，Cwbc——红细胞和白细胞的标准浓度。
图5-23 4 000 r/min转速下多相流模型与单相流模型流场对比
图5-24 血细胞在流场壁面分布云图
红细胞与白细胞的标准浓度如图5-24所示。通过对比可知，白细胞有更高的聚集程度。红细胞在叶片表面和外壳内表面上的标准浓度范围分别是99.46%～100.3%和99.98%～100.14%，白细胞的标准浓度范围分别是86.9%～105.2%和97.3%～100.7%。由图5-24a和图5-24c可知，血细胞聚集度较高的位置出现在引流叶片前缘、扩散段和整流区域，而聚集度较低的位置主要分布在引流叶片表面及后缘，此低聚集度现象与叶片后缘的流动分离相关。同时，相比于红细胞而言，白细胞低聚集度的区域更大。红细胞与白细胞聚集程度的差异与Melka［31］等的多相流模拟结果相似。
由图5-24b和图5-24d可知，红细胞在外壳内表面上与叶片相对应的位置标准浓度较高，但白细胞在此处标准浓度较低。造成这一差异的主要原因是叶轮旋转时产生的离心力，红细胞由于密度较大，在离心力的作用下被甩到外壳上，而密度相对较小的白细胞则停留在叶轮表面附近。上述结果再次证明了多相流模拟的可行性。
继血细胞在流场内聚集程度分析之后，对关注的红细胞的切变率和黏度进行分析，图5-25展示血泵内表面上红细胞切变率和动力学黏度。由图5-25a可知，轴流血泵内的红细胞整体处于高切变率的应力状态，处于200～300 s-1的切变率（正常血管内红细胞所受的切变率）的红细胞占比较小。低切变率区域主要集中存在于引流叶片上缘、叶轮叶栅的间隔处以及整流叶片的进出口处。极高切变率区域主要分布于几何尖锐处，如引流叶片、叶轮叶片和扩散段叶片的前缘附近。图5-25b显示了血泵内表面的红细胞动力学黏度，大部分红细胞的动力学黏度在0.003 Pa·s附近，畸高的黏度主要出现在引流叶片上缘。同时对照图5-25a可知，低切变率区域处对应的红细胞黏度为0.003 5 Pa·s，与正常血管内200～300 s-1的切变率对应的红细胞黏度值一致［32］。上述结果证明多相流方法确实可以模拟出红细胞在轴流血泵各处的切变率改变以及相应的黏度改变。
图5-25 叶片表面红细胞切变率与动力学黏度分布云图
使用多相流细胞破坏模型分析恒定转速下直、弯管Pump-TCPC模型的溶血风险，并将评价结果与5.2节展示的单相流细胞破坏模型的评价结果进行对比。同样选取0.2 s、0.35 s、0.6 s和0.75 s四个时刻，在直、弯管Pump-TCPC模型的IVC入口处释放300个粒子，计算并统计多相流细胞破坏因子（mBDI）和多相流粒子滞留时间（mPRT）。
不同时刻直、弯管Pump-TCPC模型中粒子BDI的分布如图5-26所示，粒子BDI和PRT的平均值和最大值见表5-10。分析表5-10可知，各个时刻直、弯管道模型的平均BDI均小于2%，大部分平均PRT也小于0.6 s。仅在0.75 s时刻，弯管道模型的平均PRT大于0.6 s。因此可以判断，轴流血泵工作时，直、弯管Pump-TCPC模型的溶血风险都较低。
图5-26 不同时刻直、弯管Pump-TCPC模型的单相流细胞破坏因子（BDI）统计
表5-10 单相流模拟下直、弯管Pump-TCPC模型的溶血风险评价参数分析
接下来进行直、弯管Pump-TCPC模型溶血风险大小的比较，判断标准如下：BDI和PRT越小代表溶血风险越低。从平均BDI指标分析，前三个时刻直管道模型优于弯管道模型，最后一个时刻相反。从最大BDI指标分析，四个时刻弯管道模型均优于直管道模型。结合图5-26d分析可知，该现象仅由极少数粒子导致。从平均PRT指标分析，直管道模型在四个时刻都优于弯管道模型。从最大PRT指标分析，在0.2 s和0.35 s时刻，弯管道优于直管道模型，但是在另外两个时刻，直管道模型优于弯管道模型。综合四个指标分析可知，在单相流模拟下，直管Pump-TCPC模型的溶血风险更低。
不同时刻直、弯管Pump-TCPC模型中粒子mBDI的分布如图5-27所示，mBDI和mPRT的平均值和最大值见表5-11。分析表5-11可知，各个时刻直、弯管道模型的平均PRT均大于0.6 s，平均mBDI和最大mBDI均小于2%。尽管粒子的平均滞留时间均超过风险阈值，但是血细胞在血泵之外的血管内运动时所受损伤较小，并且粒子的平均细胞破坏因子和最大细胞破坏因子均小于风险阈值。因此多相流模拟下依然可以得到直、弯管道模型的溶血风险均较低。
图5-27 不同时刻直、弯管Pump-TCPC模型的多相流细胞破坏因子（mBDI）统计
表5-11 多相流模拟下直、弯管Pump-TCPC模型的溶血风险评价参数分析
同样进行多相流模拟下直、弯管Pump-TCPC模型溶血风险大小的比较。从平均mBDI和最大mBDI指标分析，大部分时刻，直管道模型更优。但有两个反例，一个是在0.75 s，直、弯管道模型的平均mBDI相同。另一个是在0.35 s，弯管道模型的最大mBDI指标更优，图5-27b显示这一现象同样是由于极少数粒子导致的。从平均mPRT指标分析，各个时刻直管道模型更优，但对于最大mPRT指标，四个时刻弯管道模型均更优。综合四个指标分析，在多相流模拟下，同样可以得出直管Pump-TCPC模型溶血风险更低的结论。
BDI与mBDI以及PRT与mPRT的对比见表5-12和表5-13。由表5-12可知，平均mBDI与最大mBDI都比相对应的平均BDI与最大BDI小，且平均值与最大值的最大相对变化率分别为17.8%和36.8%。由表5-13可知，平均mPRT和最大mPRT都大于相对应的平均PRT和最大PRT，且平均值和最大值的最大相对变化率分别为120.6%和388.8%。即多相流模拟中粒子的滞留时间更久，但细胞受到的损伤更小。因为细胞破坏因子是通过应力和暴露时间在流线上的积分得到，所以可以得出多相流模拟对细胞所受到应力的估计更小。
表5-12 单相流细胞破坏因子（BDI）与多相流细胞破坏因子（mBDI）对比分析
注：rcr表示相对变化率。
表5-13 单相流粒子滞留时间（PRT）与多相流粒子滞留时间（mPRT）对比分析
注：rcr表示相对变化率。
5.5.5 研究结论
本节在单相流细胞破坏模型的基础上提出了多相流细胞破坏模型。分别对直、弯管Pump-TCPC模型进行了多相流模拟，使用多相流细胞破坏模型进行溶血风险评估，将评价结果与单相流细胞破坏模型的评价结果对比，得出以下结论：
（1）通过多相流模拟所得流场与单相流结果一致，所得血细胞在轴流血泵内的分布与理论结果相符，证明了多相流模拟的可行性和有效性。
（2）单相流细胞破坏模型和多相流细胞破坏模型的评价结果都表明直管Pump-TCPC模型优于弯管Pump-TCPC模型，溶血风险更低。
（3）与单相流模拟相比，多相流模拟中粒子的滞留时间更久但细胞受到的损伤更小，多相流细胞破坏模型对细胞所受应力的估计更小。
5.6 总结与展望
本章运用计算机建模和计算流体力学仿真方法，通过压力-流量特性、能量指标、流场特性与流量分配、溶血与血栓风险等血液动力学参数分析，对动力辅助Fontan手术进行优化设计。对比植入动力辅助装置的全腔静脉-肺动脉连接（TCPC）几何结构的血液动力学参数，优化动力辅助Fontan手术设计；实现动力辅助装置在真实Fontan循环虚拟手术植入；对比不同TCPC外管道设计对植入轴流血泵的TCPC（Pump-TCPC）血液动力学参数的影响，优化术式选择；评估恒定转速泵和可变转速泵在Pump-TCPC中溶血、血栓和血管损伤风险，选取适合真实病例的血泵输出特性；引入多相流理论，更新细胞破坏模型，提出新的溶血风险评估标准，并检验此标准的可行性与合理性。
本章在轴流血泵的仿真中引入多相流理论，提出了新的细胞破坏评价标准。此标准依据他人对离心叶片泵溶血试验参数的拟合估计给出，其对轴流血泵的支持性存在一定的疑问，并且用于验证该标准的模型较少，缺乏数据支持。在今后的研究中，应寻求用更多其他的血泵模型设计来对该参数进行验证与修正，同时增加溶血试验的开展，更好地实现评价标准的更新与验证。本章在真实病例Fontan循环中直观地显示了轴流血泵的虚拟植入，但真实情况下，轴流血泵的植入位置和植入方式可能与之有较大差距，在后续的研究中应继续探讨合理的Fontan循环布局与植入方式，实现虚拟手术的快速可视化。
参考文献
［1］陈树宝，杨思源.小儿心脏病学［M］.北京：人民卫生出版社，2012.
［2］De Leval M R，Kilner P，Gewillig M，Bull C.Total cavopulmonary connection：a logical alternative to atriopulmonary connection for complex fontan operations，experimental studies and early clinical experienc［eJ］.Journal of Thoracic and Cardiovascular Surgery，1988，96（5）：682-695.
［3］Throckmorton A L，Lim D S，Mcculloch M A，Jiang W，Song X，Allaire P E，Wood H G，Olsen D B.Computational design and experimental performance testing of an axial-flow pediatric ventricular assist devic［eJ］.Asaio Journal，2005，51（5）：629-635.
［4］Ansys Inc.Ansys Fluent 12.0 theory guide［Z］.2011.
［5］Tu J Y，Yeoh G H，Liu C Q，计算流体力学——从实践中学习［M］.王晓冬，译.沈阳：东北大学出版社，2011.
［6］Soerensen D D，Pekkan K，Zélicourt D D，Sharma S，Kanter K，Fogel M，Yoganathan A P.Introduction of a new optimized total cavopulmonary connection［J］.Annals of Thoracic Surgery，2007，83（6）：2182-219.
［7］Ryu K，Healy T M，Ensley A E，Sharma S，Lucas C，Yoganathan A P.Importance of accurate geometry in the study of the total cavopulmonary connection：computational simulations and in vitro experiments［J］.Annals of Biomedical Engineering，2001（29）：844-853.
［8］Chang Y，Gao B.Modeling and identification of an intra-aorta pump［J］.ASAIO Journal，2010，56（6）：504-509.
［9］Chopski S G，Downs E，Haggerty C M，Yoganathan A P，Throckmorton A L.Laser flow measurements in an idealized total cavopulmonary connection with mechanical circulatory assistance［J］.Artificial Organs，2011，35（11）：1052-1064.
［10］Mitamura Y，Nakamura H，Sekine K，Kim D W，Yozu R，Kawada S，Okamoto E.Prediction of hemolysis in rotary blood pumps with computational fluid dynamics analysis［J］.Journal of Congestive Heart Failure and Circulatory Support，2000，1（4）：331-336.
［11］Riemer R K，Amir G，Reichenbach S H，Reinhartz O.Mechanical support of total cavopulmonary connection with an axial flow pump［J］.Journal of Thoracic and Cardiovascular Surgery，2005，130（2）：351-354.
［12］Ou C B，Huang W，Yuen M F，Qian Y.Hemodynamic modeling of leukocyte and erythrocyte transport and interactions in intracranial aneurysms by a multiphase approach［J］.Journal of Biomechanics，2016，49（14）：3476-3484.
［13］Throckmorton A L，Tahir S A，Lopes S P，Rangus O M，Sciolino M G.Steady and transient flow analysis of a magnetically levitated pediatric vad：time varying boundary conditions［J］.International Journal of Artificial Organs，2013，36（10）：693-699.
［14］Pauls J P，Stevens M C，Schummy E，Tansley G，Fraser J F，Timms D，Gregory S D.In vitro comparison of active and passive physiological control systems for biventricular assist device［sJ］.Annals of Biomedical Engineering，2016，44（5）：1370-1380.
［15］Wang D，Jones C，Ballard-Croft C，Zhao J，Zhao G，Topaz S，Zwischen-Berger J B.Development of a double lumen cannula for a percutaneous rvad［J］.ASAIO Journal，2015，61（4）：397-402.
［16］Throckmorton A L，Kapadia J Y，Chopski S G，Bhavsar S S，Moskowitz W B，Gullquist S D，Gangemi J J，Haggerty C M，Yoganathan A P.Numerical，hydraulic，and hemolytic evaluation of an intravascular axial flow blood pump to mechanically support fontan patient［sJ］.Annals of Biomedical Engineering，2011，39（1）：324-336.
［17］Throckmorton A L，Ballman K K，Myers C D，Frankel S H，Brown J W，Rodefeld M D.Performance of a 3-bladed propeller pump to provide cavopulmonary assist in the failing fontan circulation［J］.Annals of Thoracic Surgery，2008，86（4）：1343-1347.
［18］Gregoric I D，Bruckner B A，Jacob L，Loyalka P，Kar B，La Francesca S，Myers T，Frazier O H.Techniques and complications of tandemheart ventricular assist device insertion during cardiac procedures［J］.ASAIO Journal，2009，55（3）：251-254.
［19］Kennington J R，Frankel S H，Chen J，Koenig S C，Sobieski M A，Giridharan G A，Rodefeld M D.Design optimization and performance studies of an adult scale viscous impeller pump for powered fontan in an idealized total cavopulmonary connection［J］.Cardiovascular Engineering and Technology，2011，2（4）：237-243.
［20］Xiang J P，Natarajan S K，Tremmel M，Ma D，Mocco J，Hopkins L N，Siddiqui A H，Levy E I，Meng H.Hemodynamic-morphologic discriminants for intracranial aneurysm rupture［J］.Stroke，2011，42（1）：144-152.
［21］Kung E，Baretta A，Baker C，et al.Predictive modeling of the virtual hemi-fontan operation for second stage single ventricle palliation：two patient-specific case［sJ］.Journal of Biomechanics，2013，46（2）：423-429.
［22］Sundareswaran K S，Pekkan K，Dasi L P，Whitehead K，Sharma S，Kanter K R，Fogel M A，Yoganathan A P.The total cavopulmonary connection resistance：a significant impact on single ventricle hemodynamics at rest and exercise［J］.American Journal of Physiology：Heart and Circulatory Physiology，2008，295（6）：H2427.
［23］Jung J，Hassanein A，Lyczkowski R W.Hemodynamic computation using multiphase flow dynamics in a right coronary arter［yJ］.Annals of Biomedical Engineering，2006，34（3）：393-407.
［24］Fraser K H，Zhang T，Taskin M E，Griffith B P，Wu Z J.A quantitative comparison of mechanical blood damage parameters in rotary ventricular assist devices：shear stress，exposure time and hemolysis index［J］.Journal of Biomechanical Engineering，2012，134（8）：1-11.
［25］Yano T，Sekine K，Mitoh A，Mitamura Y，Okamoto E，Kim D，Nishimura I，Murabayashi S，Yozu R.An estimation method of hemolysis within an axial flow blood pump by computational fluid dynamics analysi［sJ］.Artificial Organs，2003，27（10）：920-925.
［26］Xia D D，Zhao C Z，Zhang X W，Bai J.Computational fluid dynamics modeling and hemolysis analysis of axial blood pumps with various impeller structure［sJ］.Progress in Natural Science，2007，16（9）：993-997.
［27］Bhavsar S S，Moskowitz W B，Throckmorton A L.Interaction of an idealized cavopulmonary circulation with mechanical circulatory assist using an intravascular rotary blood pump［J］.Artificial Organs，2010，34（10）：816-827.
［28］Quemada D.Rheology of concentrated disperse systemsⅡ：a model for non-newtonian shear viscosity in steady flow［sJ］.Rheologica Acta，1978，17（6）：632-642.
［29］Van Weert C L M.Numerical and experimental analysis of shear-induced migration in suspension flow［J］.Bollettino Della SocietàItaliana Di Biologia Sperimentale，1996，32（1-2）：1084-1089.
［30］Blackshear P L，Dorman F D，Steinbach J.Some mechanical effects that influence hemolysis［J］.Transactions American Society for Artificial Internal Organs，1965（11）：112-117.
［31］Melka B，Gracka M，Adamczyk W，Rojczyk M，Golda A，Nowak A J，Bialecki R A，Ostrowski Z.Multiphase simulation of blood flow within main thoracic arteries of 8-year-old child with coarctation of the aort［aJ］.Heat & Mass Transfer，2017，54：2405-2413.
［32］Pop G，Duncker D J，Gardien M，Vranckx P，Slager C J.The clinical significance of whole blood viscosity in（cardio）vascular medicin［eJ］.Netherlands Heart Journal Monthly，2002，10（12）：512-516.