一种基于孔隙介质理论的冠状动脉血流储备分数计算方法与流程

本发明涉及生物力学、血液动力学和计算流体力学领域，特别是涉及一种基于孔隙介质理论的冠状动脉血流储备分数(ffr)计算方法。

背景技术：

介入性冠状动脉造影(ica)是一种已经应用于临床诊断冠状动脉粥样硬化性心脏病(简称冠心病)的有效方法，现阶段该方法被认为是诊断冠心病的“金标准”。但是，最新研究表明介入性冠状动脉造影(ica)只能够从图像中识别患者冠脉血管的狭窄程度，不能提供直接的临床生理学诊断依据，即无法精确诊断某一处冠脉狭窄能否导致心肌缺血性损伤，具有一定的误诊概率。此外，介入性冠状动脉造影(ica)具有价格昂贵和有创性的特点，限制了该方法成为常规临床检测手段的可能性。因此，基于冠状动脉血流储备分数(ffr)的心肌缺血性损伤诊断方法被广泛关注，并成为诊断冠心病的新“金标准”。

获取冠状动脉血流储备分数(ffr)的技术方案分为两类。第一类是利用血管内导管和压力导丝技术，实测患者冠脉狭窄处的血流储备分数(ffr)。第二类是利用计算流体力学(cfd)数值仿真，模拟计算患者全部冠状动脉的血流储备分数(ffr)。与方案一相比，方案二具有价格低、创伤小、时间短、灵活度高等显著优势。

基于计算流体力学(cfd)数值仿真的血流储备分数(ffr)计算方法，其计算精度主要受两方面因素制约。一是获取冠状动脉血管所包络成的血管腔体的三维空间几何形态数据，目前的技术手段包括冠状动脉cta、介入性冠状动脉造影(ica)、血管内超声(ivus)和血管内光学相干断层成相(oct)。二是对冠脉血管内血液动力学问题进行数学抽象和理论建模，特别是位于冠脉血管远端各分支出口的血液流量和微循环阻力分配问题，是制约数值仿真精度的关键。

为解决上述冠脉血管远端各分支出口的血液流量和微循环阻力分配问题，本专利提出一种基于孔隙介质理论的冠状动脉血流储备分数(ffr)计算方法。将冠状动脉血管腔体内的血液视为连续孔隙介质，血液运动规律遵守流体力学连续方程和动量方程。并引入任意朗格朗日欧拉方法(ale)和浸没物体法(iom)，分别模拟血管管壁柔性属性和血管腔体里存在静止或移动物体的因素。

技术实现要素：

本发明的目的是解决基于计算流体力学(cfd)数值仿真的冠脉血管血流储备分数(ffr)计算时，冠脉血管远端各分支出口的血液流量和微循环阻力分配不符合人体血液动力学分布的问题，提出一种更加精确的基于孔隙介质理论的冠状动脉血流储备分数(ffr)计算方法，并考虑血管管壁柔性属性和血管腔体里存在静止或移动物体的因素。

为实现上述发明目的，本发明采用的技术方案是：一种基于孔隙介质理论的冠状动脉血流储备分数计算方法，将冠状动脉血管腔体内的血液视为连续孔隙介质、考虑血管管壁弹性、考虑血管腔体里存在静止或移动物体，血液运动规律遵守连续方程和动量方程；所述冠状动脉血管腔体内的血液视为连续孔隙介质，介质参数包括孔隙率、粘性损失系数与惯性损失系数通过实验获得；所述血管管壁既是刚体属性，也是随血液流动而收缩、舒张的柔性属性；所述血管腔体内部包括代表血液的流体介质，代表血栓、纤维组织增生、钙质沉着、动脉粥样硬化、血管支架和介入性医疗器械的静止或移动固体介质；所述计算方法包括以下步骤：

步骤一，通过冠状动脉cta技术获取患者的冠状动脉血管所包络成的血管腔体的三维空间几何形态数据；通过测量舒张压与心肌质量获得患者冠脉血管主动脉入口压力pinlet及流量qinlet数据；通过血液测试结果获得孔隙介质的孔隙率n；

步骤二，创建冠脉血管血液流动分析计算域，并在无弹性变形的冠状动脉血管管壁、无弹性变形的管壁施加无滑移刚性壁面条件；在由于收缩、舒张引起管壁位移的冠状动脉血管管壁、有弹性变形的管壁施加收缩、舒张的无滑移柔性壁面条件；在冠脉血管入口、左冠状动脉入口、血管入口施加流量入口边界条件，流量与患者心肌质量成正比；在冠脉血管血液流动分析计算域出口、左前降第一出口、左前降第二出口、左前降第三出口、左前降第四出口、左前降第五出口、左前降第六出口、左前降第七出口、左回旋第一出口、左回旋第二出口、左回旋第三出口、左回旋第四出口、血管出口施加阻力出口边界条件，为考虑泊肃叶定律和莫里法则的远端阻力；以及初始条件，冠状动脉血管腔体内的孔隙介质化的血液流速和压力、依附于冠状动脉血管管壁的静止固体介质和冠状动脉血管腔体内随血液流动而移动的固体介质的位置分布情况；

步骤三，采用数值模型离散求解所述计算域的连续方程和动量方程，得到所述血管腔体内部的血液流动速度和压力；考虑血管管壁弹性，采用任意朗格朗日欧拉方法(ale)模拟由于管壁变形对血液的影响；考虑血液内部存在静止或移动的固体介质，采用浸没物体法(iom)模拟由于固体的存在对血液流动的影响；所述连续方程为式(1)和动量方程为式(2)：

其中为血管内血液的渗流速度分量，为血管内血液的流速分量，n为孔隙介质的孔隙率，ujm为冠脉血管由于管壁变形产生的计算域网格节点速度分量，ρ为患者血液的密度，p血管内血液压力，ν为患者血液的运动粘度，si为i方向的动量损失源项，fi为i方向的浸入物体动量损失源项，i为三个方向分量；

所述冠脉血管内血液视为连续孔隙介质，由于孔隙的存在，i方向的动量损失源项为式(3)：

其中为血管血液的粘性损失系数分量，为血管血液的惯性损失系数分量；

所述冠脉血管内由于代表血液的流体介质，代表血栓、纤维组织增生、钙质沉着、动脉粥样硬化、血管支架和介入性医疗器械的静止或移动固体介质的存在，i方向的浸入物体动量损失源项为式(4)：

其中ω为血管内部固体介质表示的计算域，r和a为数值离散参数，vo为固体介质的运动速度；

步骤四，采用通过上述步骤得到的血液压力p，计算冠脉血管的血流储备分数。

本发明的有益效果是：这种基于孔隙介质理论的冠状动脉血流储备分数计算方法，将冠状动脉血管腔体内的血液视为连续孔隙介质、考虑血管管壁弹性、考虑血管腔体里存在静止或移动物体，血液运动规律遵守连续方程和动量方程；冠状动脉血管腔体内的血液视为连续孔隙介质，介质参数包括孔隙率、粘性损失系数与惯性损失系数通过实验获得；考虑了血管管壁柔性属性和血管腔体里存在静止或移动物体的因素，相比于血管管壁为刚性边界属性且血液中只存在流体介质的计算方法，具有更良好的精确性和更广泛的适用性。

附图说明

图1是一种基于孔隙介质理论的冠状动脉血流储备分数(ffr)计算方法的工作流程图。

图2是冠脉血管血液流动分析计算域结构示意图。

图3是冠状动脉血管所包络成的血管腔体的三维空间几何形态示意图。

图4是考虑冠脉血管管壁弹性变形的计算算例。

图5是考虑冠脉血管内存在移动异物的计算算例。

图中：1、无弹性变形的冠状动脉血管管壁，2、冠状动脉血管腔体内的孔隙介质化的血液流速和压力，3、由于收缩、舒张引起管壁位移的冠状动脉血管管壁，4、依附于冠状动脉血管管壁的静止固体介质，5、冠状动脉血管腔体内随血液流动而移动的固体介质，6、所述冠脉血管入口，7、所述冠脉血管血液流动分析计算域出口，8、左冠状动脉入口，9、左前降第一出口,10、左前降第二出口,11、左前降第三出口,12、左前降第四出口，13、左前降第五出口,14、左前降第六出口，15、左前降第七出口，16、左回旋第一出口，17、左回旋第二出口，18、左回旋第三出口，19、左回旋第四出口,20、血管入口，21、血管出口，22、无弹性变形的管壁，23、有弹性变形的管壁。

具体实施方式

下面将结合本发明实施案例冠脉血管管壁弹性变形的计算算例(图4)和冠脉血管内存在移动异物的计算算例(图5)，对本发明实施案例中的技术方案进行清楚、完整地描述。所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施案例，本领域普通技术人员所获得的所有其他实施例，都属于本发明保护的范围。

图1是一种基于孔隙介质理论的冠状动脉血流储备分数(ffr)计算方法的工作流程图，首先(步骤一)，通过冠状动脉cta技术获取患者的冠状动脉血管所包络成的血管腔体的三维空间几何形态数据；通过测量舒张压与心肌质量获得患者冠脉血管主动脉入口压力pinlet及流量qinlet数据；通过血液测试结果获得孔隙介质的孔隙率n；

然后(步骤二)，创建冠脉血管血液流动分析计算域，并在无弹性变形的冠状动脉血管管壁、无弹性变形的管壁施加无滑移刚性壁面条件；在由于收缩、舒张引起管壁位移的冠状动脉血管管壁、有弹性变形的管壁施加收缩、舒张的无滑移柔性壁面条件；在冠脉血管入口、左冠状动脉入口、血管入口施加流量入口边界条件，流量与患者心肌质量成正比；在冠脉血管血液流动分析计算域出口、左前降第一出口、左前降第二出口、左前降第三出口、左前降第四出口、左前降第五出口、左前降第六出口、左前降第七出口、左回旋第一出口、左回旋第二出口、左回旋第三出口、左回旋第四出口、血管出口施加阻力出口边界条件，为考虑泊肃叶定律和莫里法则的远端阻力；以及初始条件，冠状动脉血管腔体内的孔隙介质化的血液流速和压力、依附于冠状动脉血管管壁的静止固体介质和冠状动脉血管腔体内随血液流动而移动的固体介质的位置分布情况；

接下来(步骤三)，采用数值模型离散求解所述计算域的连续方程和动量方程，得到所述血管腔体内部的血液流动速度和压力；考虑血管管壁弹性，采用任意朗格朗日欧拉方法(ale)模拟由于管壁变形对血液的影响；考虑血液内部存在静止或移动的固体介质，采用浸没物体法(iom)模拟由于固体的存在对血液流动的影响；

最后(步骤四)，采用通过上述步骤得到的血液压力p，计算冠脉血管的血流储备分数。

图2是冠脉血管血液流动分析计算域结构示意。图中3示意表示考虑冠脉血管管壁弹性情况下，产生收缩与扩张的血管管壁位置；图中4示意表示血管管壁粘附异物情况，图中5示意表示血管管腔存在移动异物情况。

图3是人体心脏左冠脉血管结构示意图。左冠脉血管包括两个分支，分别为左前降和左回旋，其中左前降共有7个出口，为9-15，左回旋共有4个出口，为16-19。

图4是冠脉血管管壁弹性变形的计算算例示意图。该算例模拟计算了左冠脉血管管段发生血管壁塌陷的算例。图中23为塌陷段。

冠脉血管管壁弹性变形的计算算例(图4)的具体实施方式为：

步骤一，通过冠状动脉cta技术获取患者的冠状动脉血管所包络成的血管腔体的三维空间几何形态数据(如图2)；以其中靠近冠状动脉血管入口处某一断面为基本断面，通过拉伸形成均匀管段，构成计算域；通过测量舒张压与心肌质量获得患者冠脉血管主动脉入口压力pinlet及流量qinlet数据。

步骤二，创建冠脉血管管壁弹性变形的计算算例计算域，并在无弹性变形的管壁22施加无滑移刚性壁面条件；在有弹性变形的管壁23施加收缩、舒张的无滑移柔性壁面条件，本算例管壁视为柔性薄膜，其归一化无量纲弹性模量e^*为0.1，归一化无量纲初始应力为5.964；在血管入口20施加流量入口边界条件流量qinlet；在血管出口21施加阻力出口边界条件，为考虑泊肃叶定律和莫里法则的远端阻力；共计算10秒。

步骤三，采用数值模型离散求解所述计算域的式(1)连续方程和式(2)动量方程，得到所述血管腔体内部的血液流动速度和压力；考虑血管管壁弹性，采用任意朗格朗日欧拉方法(ale)模拟由于管壁变形对血液的影响。

其中为所述血管内血液的渗流速度分量，为所述血管内血液的流速分量，n为所述孔隙介质的孔隙率，ujm为所述冠脉血管由于管壁变形产生的计算域网格节点速度分量，ρ为患者血液的密度，p血管内血液压力，ν为患者血液的运动粘度，i为三个方向分量。

步骤四，采用通过上述步骤得到的血液压力p，计算冠脉血管的血流储备分数(ffr)。

图5是考虑冠脉血管内存在移动异物的计算算例示意图。该算例模拟计算了血管横断面中心处放置一上下简谐运动的异物(该实例异物为单位长度的圆)时，管内血液的变化。图5(a)为相位角为0°时异物周围涡量场(左)和压力场(右)；图5(b)为相位角为36°时异物周围涡量场(左)和压力场(右)；图5(c)为相位角为96°时异物周围涡量场(左)和压力场(右)；图5(d)为相位角为108°时异物周围涡量场(左)和压力场(右)；图5(e)为相位角为192°时异物周围涡量场(左)和压力场(右)；图5(f)为相位角为216°时异物周围涡量场(左)和压力场(右)；图5(g)为相位角为288°时异物周围涡量场(左)和压力场(右)；图5(h)为相位角为324°时异物周围涡量场(左)和压力场(右)。

冠脉血管内存在移动异物的计算算例(图5)的具体实施方式为：

步骤一，本实例为简化算例，模拟在血管横断面中心处放置一上下简谐运动的圆，其无量纲化半径1。简谐运动幅值为无量纲尺寸2.38，频率为1.33赫兹。

步骤二，创建冠脉血管内存在移动异物的算例计算域，本算例计算域为无量纲化半径30的圆；共计算40周期。

步骤三，采用数值模型离散求解所述计算域的式(1)连续方程和式(2)动量方程，由于是二维计算，因此仅有两个方向。得到所述血管腔体内部的血液流动速度和压力；考虑血液内部存在静止或移动的固体介质，采用浸没物体法(iom)模拟由于固体的存在对血液流动的影响。

其中为所述血管内血液的渗流速度分量，为所述血管内血液的流速分量，n为所述孔隙介质的孔隙率，ujm为所述冠脉血管由于管壁变形产生的计算域网格节点速度分量，ρ为患者血液的密度，p血管内血液压力，ν为患者血液的运动粘度，i为两个个方向分量。

i方向的浸入物体动量损失源项为式(4)：

其中ω为血管内部固体介质表示的计算域，r和a为数值离散参数，vo为固体介质的运动速度；

步骤四，采用通过上述步骤得到的血液压力p，计算冠脉血管的血流储备分数(ffr)。