摘要
目的 构建包含静脉窦以及不包含静脉窦的两种不同结构的静脉模型,利用流固耦合仿真研究静脉窦结构对于下肢静脉系统血液流动的影响,分析静脉窦结构在静脉瓣开闭过程中的作用。方法 基于超声图像和动物解剖图像,使用建模软件建立两种不同结构的静脉系统三维几何模型,通过浸没边界法完成静脉瓣膜的流固耦合仿真,以获取瓣膜形态变化和流线分布、壁面剪切力分布、Q准则强度分布等血流动力学参数。结果 带窦模型的瓣膜几何开口面积大于无窦模型。血液流线的分布证明了窦区涡流与血管内螺旋流的存在。同端的静脉瓣对上呈现不对称的壁面剪切力分布,带窦模型静脉瓣背面的高壁面剪切力区域和远、近心端瓣膜间的Q准则强度均小于无窦模型。结论 静脉窦结构有助于静脉瓣膜的打开,对螺旋流动有一定的削弱作用,但螺旋流动模式有利于静脉窦区涡流的产生。静脉窦对静脉系统中的血流动力学有着重要影响,在静脉瓣的研究中应将静脉窦结构作为不可或缺的部分。
关键词:
静脉
正交静脉瓣膜
静脉窦
流固耦合仿真
血流动力学
静脉是人体心血管循环系统中不可或缺的构成部分,它不仅储存着人体80%的血液,还有效保证了缺氧血液回流到心脏这一过程的顺利进行。慢性静脉功能不全(chronic venous insufficiency,CVI)是一种常见的下肢疾病,我国下肢静脉疾病的患病率为8.89%,年发病率为0.5%~3.0%,患者主要临床表现为下肢静脉曲张、腿部胀痛,严重者会出现水肿、溃疡等症状,极大影响了患者的身体健康与经济状况 。而静脉瓣膜功能不全是慢性静脉功能不全的主要原因,当静脉瓣膜受损或静脉系统出现其他病变,导致静脉瓣膜功能出现异常,血液便无法回流到心脏从而在静脉内滞留,使得静脉内压力持续升高,进而引起静脉血管壁的病理性扩张、炎症和内皮的促血栓表型变化,最终导致慢性静脉功能不全 。因此,开展关于静脉瓣膜问题的研究对静脉疾病的诊断和治疗具有重要意义。
近年来,国内外很多学者对静脉瓣的动力学行为进行了研究。Buxton等 最早对静脉进行仿真研究,可视化展示了静脉瓣膜防止反流的功能特点。Zervides等 建立了静脉瓣膜的三维模型,但并没有将静脉窦这一重要结构考虑在模型中。Tien等 通过粒子图像测速(particle image velocimetry,PIV)实验对血液流经静脉窦区时的血流动力学进行了研究,其结果显示,静脉窦可以通过在瓣膜小叶后产生涡流来避免血流停滞。研究发现,大隐静脉和股静脉中两对相邻瓣之间呈正交分布,导致血液螺旋流动,从而提高静脉回流的效率 。Chen等 通过对不同配置的相邻瓣膜进行模拟仿真发现,正交瓣膜对之间存在较强的螺旋流,有效降低了静脉血栓形成的风险,但没有考虑静脉窦对静脉血液流动的影响。汪波等 通过建立三维模型再现了静脉瓣运动与血液流动的生理过程,但研究对象为单个瓣膜,没有考虑相邻瓣膜对血液流动的协同作用。
针对上述研究中存在的不足,本文考虑加入静脉窦结构,建立两种结构不同且包含两对相邻正交静脉瓣的静脉模型。对这两对瓣膜在血液流动中的过程进行三维流固耦合仿真,研究静脉窦结构对于静脉血液流动和静脉瓣功能的影响。
首先,使用SolidWorks 2022(Dassault System公司,法国)建立两种不同结构的静脉血管三维几何模型,然后在Hypermesh 2020(Altair公司,美国)中进行流体域和固体域网格划分,将网格文件导入LS-PrePost 4.6(ANSYS公司,美国)进行前处理,赋予其材料属性和边界条件,再使用LS-Run 2022 R1(ANSYS公司,美国)对流固耦合仿真模型进行求解,最后借助LS-PrePost 4.6、Fluent 2022 R1(ANSYS公司,美国)和CFD-Post 2022 R1(ANSYS公司,美国)等软件对求解结果进行后处理与分析。
根据华中科技大学同济医学院附属协和医院提供的人体下肢静脉系统的超声图像[见 图1 (a)],确定下肢静脉的基本形状,再结合在动物实验中得到的静脉解剖图像 [见 图1 (b)],发现静脉瓣膜是一种静脉内壁向管腔突出而形成的半月形囊袋状薄膜结构,其瓣叶的叶尖朝向近心端,瓣膜根部的静脉血管向外膨出呈现窦状 。下肢静脉瓣膜绝大部分以成对的形式出现在静脉系统中,并且相邻两对静脉瓣膜之间的朝向通常呈现近乎垂直的大角度偏转[见 图1 (c)] 。
基于静脉瓣以上的结构特点和从解剖结果中得到的各部分参数 ,本文使用SolidWorks 2022构建带窦静脉血管三维几何模型,该模型由静脉血管壁、静脉窦与两对静脉瓣膜组成,带窦模型如 图2 (a)所示。为了进行对比分析,研究静脉窦在血液流经两对相邻正交静脉瓣膜时对血流产生的影响,同时为确保实验条件的一致性,在同一参数设置下,构建了一个不包含静脉窦的对照模型[见 图2 (b)]。静脉瓣膜模型的三视图如 图2 (c)所示。
将两种不同结构的静脉血管三维几何模型导入Hypermesh中,对模型进行网格划分。静脉瓣瓣叶被划分为单元尺寸为0.5 mm、单元厚度为0.05 mm四边形壳单元,静脉血管壁被划分为单元尺寸为1 mm四边形壳单元,血液流域模型被划分为六面体实体单元,单元尺寸设置为0.5 mm,这一尺寸与瓣叶尺寸一致,从而提高流固耦合计算的稳定性。以上述网格为基网格,使用更小的最大单元尺寸生成加密网格,通过比较静脉出口的流速验证网格敏感性。最终两种网格的结果差异小于5%,故采用基网格进行仿真分析。两种模型的计算网格参数如 表1 所示。
将静脉瓣膜设置为弹性模量为0.15 MPa的线弹性材料 。由于本文重点关注瓣膜与血液之间的相互作用,故将静脉血管壁设置为刚性壁,且流体在管壁面无滑移。血液为不可压缩的牛顿流体,密度为1 060 kg/m 3 ,黏度为3.5 mPa·s 。
为了模拟人体在静息仰卧姿态下,由心脏搏动作用驱动的静脉血液回流运动,本文以人体在仰卧时的生理条件为参考,在静脉系统入口处即远心端设置如 图3 所示的压力曲线,模型周期时长为1 s。前0.5 s采用正弦函数 A sin (2π t ) 的前半段,其中 A =2 mmHg(1 mmHg=0.133 kPa,下同),旨在以脉冲函数形式的压力模拟心脏搏动对静脉血液的动力作用;后0.5 s对瓣膜施加方向为由近心端向远心端的1 mmHg负压,模拟人体的呼吸作用以促进瓣膜关闭,检验瓣膜闭合功能是否正常。同时在模型的出入口设置4 mmHg初始静水压强,模拟人体仰卧姿态下静脉血管内的静水压。
图3 压力边界条件(1 mmHg=0.133 kPa)
单个计算周期为1 s,为消除初始条件对仿真带来的误差,共计算3个周期,第2、3周期的计算结果仅相差0.23%。因此,为节约计算成本,取第2个周期数据进行后续分析。
图4 (a)展示了通过流固耦合仿真获得的两种模型中静脉瓣膜处于4个阶段的不同形态。在开启阶段,此时静脉血液由于心脏搏动作用从静脉远心端向近心端流动,静脉瓣膜受到血液冲击而逐步打开。当瓣膜几何开口面积(geometric orifice area, GOA)达到较大水平时,瓣叶形态维持稳定,此时称为平衡阶段,由带窦静脉模型可以明显地看到瓣叶向窦区内凸出,形成正常生理情况下静脉瓣膜的“漏斗型”造型。在关闭阶段,瓣叶腹部由于跨瓣压差向远心端凹进。在闭合阶段,两种模型的瓣膜对都几乎完全闭合,说明了本文构建的静脉瓣膜模型开闭功能正常。
瓣膜在开闭过程中的开口面积是衡量瓣膜是否功能正常的一项关键指标,过小的开口面积会阻碍血液回流,引起静脉高压。本文计算了两种模型随时间变化的GOA [见 图4 (b)],对比发现带窦静脉模型的GOA最大可达到0.14 cm 2 ,而无窦静脉模型的GOA最大只有0.08 cm 2 ,并且前者的GOA在绝大部分时间都要大于后者。
当静脉瓣膜处于平衡阶段稳定打开时,流经瓣膜的血液会分为三部分,其具体的血液流线分布情况如 图5 所示。主体血液由于文丘里效应形成射流,加速流向下一对瓣膜,冲击到相邻瓣膜的部分血液回落,流回远心端瓣膜的瓣叶背部区域,此时血液便会在瓣叶背部与血管壁之间形成涡流,它可以在防止血液滞留的同时阻止瓣叶持续打开与静脉内壁相贴,在图中表现为流线在瓣膜背部区域弯曲旋转形成漩涡状。结果表明,带窦静脉模型中的血液最大流速为0.689 m/s,小于无窦静脉模型的血液最大流速(0.793 m/s)。此外,观察两种模型的流线形状发现,在流线分布正视图中,流线在通过近心端瓣膜后,都呈现出明显向左倾斜弯曲的趋势;在流线分布左视图中,远近心端瓣膜之间的流线均是先向左侧血管壁偏移,再向右偏移通过远心端瓣膜。
通过将平衡阶段的瓣膜几何形态导出到仿真软件Fluent中进行单纯流体仿真,得到两种静脉模型处于平衡阶段的壁面剪切力(wall shear stress, WSS)分布图。对比 图6 (a)中两种模型远、近心端4对瓣膜的瓣叶背部WSS分布发现,不论是远心端还是近心端,带窦静脉模型瓣叶的背部,高剪切力即WSS>0.3 Pa的区域都要明显小于无窦静脉模型。并且在两种模型的剪切力分布图中,同一端两片瓣叶背部的WSS分布都呈现出不对称的一侧大一侧小的现象。 图6 (b)展示了两种模型血管壁上的WSS分布,带窦静脉模型高剪切力区域明显小于无窦静脉模型。
Q准则(Q criterion)是一种常用的涡识别方法,Q的正值部分表示流场中涡量占主导的区域,负值部分表示应变速率或黏性应力占主导的区域。因此,当Q为正值时代表该处易产生涡流,为负值则相反,且数值大小在一定程度上反映了该位置涡流的强弱。 图7 展示了两种模型在两种正交视角下的Q准则强度分布。血液在流经两种模型远心端瓣膜后,立刻在瓣膜开口处呈现出较强的朝向侧边血管壁的涡流[见 图7 (a)]。血液在远心端与近心端瓣膜之间的区域两侧呈现出极强的涡流分布,无窦模型的强涡流区域要明显大于带窦模型[见 图7 (b)]。但通过纵向对比带窦模型近心端与远心端瓣膜的背部Q准则强度分布,发现近心端处涡流强度较远心端有明显增大 [见 图7 (c)]。
本文基于华中科技大学同济医学院附属协和医院提供的人体下肢静脉系统的超声图像,建立两种不同结构的三维静脉系统模型,研究静脉窦在静脉系统中对于静脉血流的影响,借助浸没边界法完成了静脉血液回流过程的流固耦合仿真,获得了静脉瓣膜的形态学、流线分布、WSS分布、Q准则强度分布等血流动力学结果。
本文发现,在静脉系统中,静脉窦结构对含有两对正交静脉瓣的静脉血管内的血液流动产生了显著影响。形态学结果显示,无静脉窦模型中的静脉瓣膜GOA在绝大部分时刻都小于带静脉窦模型,说明静脉窦结构对于瓣膜的打开功能具有显著的促进作用,这一方面是由于静脉窦的外凸结构为瓣膜的打开提供了更加充足的背部空间;另一方面无窦模型由于窦区管径小于带窦模型,故无窦模型瓣膜间血液流速更大,瓣间压强更小,使得瓣叶正面与背面压差增大,从而导致无窦模型的开口面积进一步缩小。
平衡阶段在整个瓣膜运动周期中占据了接近50%的时间,故对该阶段的血流动力学研究十分重要。从流线分布来看,两种模型中流线在不同视角下均表现出大幅度倾斜、偏移等极不对称的分布规律,与常见泊肃叶流的对称流场不同,这些都是显著的螺旋流特征。该特征与Lee等 获得的超声影像相符,说明了本研究所建立两种仿真模型的可信度。
静脉瓣膜背面区域上的WSS分布结果显示,不论是近心端瓣膜还是远心端,两种模型的两侧静脉瓣膜剪切力分布都呈现出极不对称的现象,这是由于血液在流经静脉瓣后产生了螺旋流动,大量血液向一侧血管壁偏移,造成两侧血液流速急缓程度不一,流速较大的血流在同侧瓣膜背部产生较大的WSS,另一侧则较小。与此同时,横向比较两种模型同一端两对瓣膜上的WSS分布,带静脉窦模型的高WSS分布区域要显著小于无静脉窦模型,这一现象与Tien等 得出的结论相符。静脉窦结构对血液流动起到一定的调节作用,它使得血液在窦区的流速减缓,并易形成涡流,减小WSS,这从另一个角度验证了本研究的可靠性。
前人对于相邻静脉瓣呈正交分布这一更符合生理情况的研究只停留在不加入静脉窦结构的静脉系统仿真 ,但静脉窦结构是静脉血管中显著的生理特征,它对相邻正交静脉瓣导致的螺旋流动模式会产生何种影响具有重要的研究价值。由Q准则强度分布可见,正视图中远心端瓣膜开口处存在两股朝向两侧血管壁的强涡量区,这是由于原先沿血管直线均匀流动的血液在此处发生扭转,形成螺旋流动方式,导致开口处出现局部涡量增大的情况。而左视图中,带窦模型主体的涡量强度分布显著小于无窦模型,提示其螺旋流动的强度小于无窦模型,静脉窦结构会削弱正交静脉瓣导致的螺旋流动。这一结果也可由静脉血管壁上的WSS分布图佐证,较弱的螺旋流动使得带窦模型血管壁上的WSS也较小。血液流经远心端瓣膜前是泊肃叶流动状态,而流经近心端瓣膜前是螺旋流动状态,结合带窦模型远、近心端瓣膜背部Q准则强度分布的差异推测,螺旋流动模式与静脉窦结构的协同作用使得窦区更易产生涡流,防止血液停滞,减小血栓形成风险。由此可见,静脉窦对下肢静脉系统中的血流动力学有着非常重要的作用,且静脉窦的加入更加贴合生理情况。因此,本研究对未来人工静脉瓣膜的设计工作有一定指导意义。
本研究的局限性如下:① 对静脉血液流动的研究还停留在仰卧体位下,没有考虑直立情况下重力对血液回流的阻碍作用;② 研究中给予静脉血液的压力边界条件来源于心脏搏动,但在生理情况下,人在直立行走时,小腿肌肉泵的收缩作用才是驱动静脉血液回流到心脏的主要动力来源;③ 静脉瓣膜被设置为线弹性材料,但根据相关文献,静脉瓣膜实际是各向异性的超弹性材料 ;④ 带窦模型与无窦模型管内流速不同,可能会影响WSS分布、Q准则分布等血流动力学参数。在后续工作中,需要获取直立体态下的静脉血管内压力信息,将静脉瓣膜的线弹性材料替换为更合理的超弹性材料,将小腿肌肉泵的收缩作用转化为动态压力边界条件对模型进行优化。
本文利用流固耦合仿真方法,研究静脉窦结构对包含两对正交静脉瓣的静脉血管中血液流动的影响,得出以下主要结论:静脉窦结构有助于静脉瓣膜的打开,促进静脉血液的回流;血液在相邻两对正交静脉瓣之间呈现螺旋流动的特点;静脉窦结构使得瓣膜背部区域较易产生涡流,降低WSS,减小对内皮细胞的损伤;螺旋流动模式有助于血液在窦区形成涡流,但静脉窦结构会对正交静脉瓣导致的螺旋流动有一定的削弱作用,为临床上设计和优化静脉瓣膜置换手术提供了参考依据。
作者贡献声明: 李杰负责建模仿真、数据处理和论文撰写;黎健明负责建模和仿真方法指导;王志红、党一平、王盛章负责论文写作指导及修改。