技术领域
[0001] 本
发明属于
人机交互领域，尤其涉及一种模拟血流与血管壁作用的方法和装置。
背景技术
[0002] 在许多外科手术的训练和规划的模拟血管手术系统中，研究血流与血管壁相互作用，从而真实地模拟血流及其与周围组织的相互作用，可以为模拟血管手术系统提供实时、有效的交互信息，在遇到特殊情况时，可以更为有效的帮助医生了解如何采取正确的行动。因而，血流与血管壁的相互作用的研究工作，对评价系统好坏具有十分重要的作用。
[0003] 现有的模拟血流和血管壁相互作用的方法，包括基于粒子的方法，其具有独立性较强的复杂网状结构及较高的性能，因此常被应用在
对流体的实时模拟中，Muller等人把光滑粒子
流体动
力学方法(SPH)应用于模拟手术系统中模拟血液的流动，然而其并没有提供血流与
变形组织(血管壁)之间相互作用的解决方案，对于血流与血管壁之间存在相互作用的模拟血管手术系统中，不能真实的反应血管壁和血流作用的交互信息。
发明内容
[0004] 本发明的目的在于提供一种模拟血流与血管壁作用的方法，以解决
现有技术因不能血流与血管壁之间相互作用，在模拟血管手术系统中不能真实的反应血流与血管壁作用的交互信息问题。
[0005] 本发明是这样实现的，一种模拟血流与血管壁作用的方法，所述方法包括：
[0006] 通过光滑粒子流体动力学方法计算血管内的流体粒子的压力
密度、粘性力密度和外力密度；
[0007] 使用质点
弹簧模型计算位于血管壁与血流的交界处的边界粒子受到的弹回边界力；
[0008] 采用Runge-Kuntta模型整合所述计算的血管内的流体粒子的压力密度、粘性力密度、外力密度和位于血管壁与血流的交界处的边界粒子受到的弹回边界力，更新所述血流与血管壁的动态信息。
[0009] 本发明
实施例的另一目的在于提供一种模拟血流与血管壁作用的装置，所述装置包括：
[0010] 第一计算单元，用于通过光滑粒子流体动力学方法计算血管内的流体粒子的压力密度、粘性力密度和外力密度；
[0011] 第二计算单元，用于使用质点弹簧模型计算位于血管壁与血流的交界处的边界粒子受到的弹回边界力；
[0012] 更新单元，用于采用Runge-Kuntta模型整合所述计算的血管内的流体粒子的压力密度、粘性力密度、外力密度和位于血管壁与血流的交界处的边界粒子受到的弹回边界力，更新所述血流与血管壁的动态信息。
[0013] 在本发明中，利用光滑粒子流体动力学方法计算血管内的流体粒子的压力密度、粘性力密度和外力密度，使用质点弹簧模型计算位于血管壁与血流的交界处的边界粒子受到的弹回边界力，从而使得在模拟血流流动的同时，能够模拟血管壁的变形，处于边界上的粒子既可以看成是光滑粒子流体动力学方法中的边界颗粒，又可以看作是质点弹簧模型中的边界颗粒，通过弹回边界力的计算，能够为虚拟血管手术系统提供更真实、更实时、更有效的互动信息。
附图说明
[0014] 图1是本发明第一实施例提供的模拟血流与血管壁作用的方法的实现
流程图；
[0015] 图2是本发明第一实施例提供的采用光滑粒子动力学方法求解血管内的流体粒子的动态特征过程图；
[0016] 图3是本发明第一实施例提供的边界粒子所受的弹回边界力示意图；
[0017] 图4为本发明第一实施例提供的基于光滑粒子动力学方法获取的粒子的实时性能表；
[0019] 图5b为动脉瘤模型经过20个时间步长后大量粒子注入动脉瘤空腔的示意图；
[0020] 图5c为经过50个时间步长后，粒子从血管壁涌出的示意图；
[0021] 图5d为基于移动立体方法的动脉瘤破裂示意图；
[0022] 图6为本发明第二实施例提供的模拟血流与血管壁作用的装置的结构
框图。
具体实施方式
[0023] 为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
[0024] 当人们把大量时间和精力花费在人体组织与虚拟设备之间的软组织变形模型及实时交互的
力反馈模型开发工作时，却忽略了对流体与变形组织之间相互作用的模拟工作，如血流与血管壁之间的相互作用等。
[0025] 事实上，真实地模拟血流及其与周围组织的相互作用，可以帮助医生了解如何在受到干预时采取正确的行动，在许多外科训练和规划系统中，研究血流与血管壁相互作用对评价系统好坏具有十分重要的作用。
[0026] 本发明所提出的模拟血流及其与周围组织，即血管壁的相互作用的方法，包括：
[0027] 通过光滑粒子流体动力学方法计算血管内的流体粒子的压力密度、粘性力密度和外力密度；
[0028] 使用质点弹簧模型计算位于血管壁与血流的交界处的边界粒子受到的弹回边界力；
[0029] 采用Runge-Kuntta模型整合所述计算的血管内的流体粒子的压力密度、粘性力密度、外力密度和位于血管壁与血流的交界处的边界粒子受到的弹回边界力，更新所述血流与血管壁的动态信息。
[0030] 本发明利用光滑粒子流体动力学方法计算血管内的流体粒子的压力密度、粘性力密度和外力密度，使用质点弹簧模型计算位于血管壁与血流的交界处的边界粒子受到的弹回边界力，从而使得在模拟血流流动的同时，能够模拟血管壁的变形，处于边界上的粒子既可以看成是光滑粒子流体动力学方法中的边界颗粒，又可以看作是质点弹簧模型中的边界颗粒，通过弹回边界力的计算，能够为虚拟血管手术系统提供更真实、更实时、更有效的互动信息。
[0031] 实施例一：
[0032] 图1示出了本发明第一实施例提供的模拟血流与血管壁作用的方法的实现流程，详述如下：
[0033] 在步骤S101中，通过光滑粒子流体动力学方法计算血管内的流体粒子的压力密度、粘性力密度和外力密度。
[0034] 具体的，通常情况下，血流运动的基本方程是由流体的
质量守恒、动量守恒和
能量守恒推导出。我们假设血流的
温度几乎是恒定的且没有热量输入到系统中。因此，能量平衡方程与机械能平衡的方程是一样的，能量守恒导致没有新的独立方程。此外，粒子系统通过维持粒子数量与单个粒子质量不变而达到质量守恒。因此，我们可以将重点放在如何解决Navier-Stokes方程，其动量守恒公式为：
[0035]
[0036] 其中，ρ是流体粒子密度； 是速度矢量场；适量运算符 p是压力；是自身的力，如引力；μ是粘性系数。根据方程(1)可以计算出每个时间步长内，粘性力为 压力为 外力为 的力密度。
[0037] 与其他的模拟流体的方法相比，光滑粒子流体动力学SPH方法在运用及快速计算方面相对简单。光滑粒子流体动力学SPH的基本原理是：
[0038]
[0039] 其中，标量 在 处的标量位是邻近的粒子贡献的加权总和。mj和ρj分别是粒子j的质量和密度，且 Aj是在xj处的场量。函数 称为满足：及 条件的平滑
内核，其中h为光滑粒子流体动力学方法的支持半径，δ为狄拉克函数。
[0040] 采用光滑粒子动力学方法求解血管内的流体粒子的动态特征过程如图2所示，可以包括如下步骤：
[0041] 在步骤S201中，在单位时间步长中，在设定的支持半径范围内查找与目标粒子i相邻近的相邻粒子j。
[0042] 在通常情况下，寻找邻近粒子的方法有三种，即成对接近的方法、链接列表的方法和基于树的方法。采用链接列表法在速度方面远远快于成对接近法，而链接列表法与基于树的方法相比，具有更高的并行度。
[0043] 在步骤S202中，根据光滑粒子流体动力学原理，获取血管内的目标粒子i的密度ρi以及粒子i的压力密度fip。
[0044] 可以采用两种计算密度的方法，即求和密度法与继续密度法。求和密度法可直接应用于与光滑粒子流体动力学SPH近似的密度本身。继续密度法是利用SPH近似的概念，引入变化的连续性方程估计密度。虽然求和密度法可能会引起一些边界效应导致边界结果不准确，但它可以使总质量守恒。因此，在光滑粒子流体动力学SPH实际应用中很方便。在本方案中，在计算血管内的目标粒子的密度时，可以通过在边界周围设置一些虚拟流体粒子来降低边界效应。
[0045] 血管内的流体粒子有三个初始量：质量、
位置和速度。单位时间步长的评估中，光滑粒子流体动力学SPH是用来评估密度ρj和由粒子i受到的压力、粘性力和外力共同作用产生的力密度fip、fiex、fiv。我们通过基本方程，可以计算出了粒子i的密度和压力密度。
[0046]
[0047]
[0048] 其中，pj、pi分别表示粒子j和粒子i的压力，其它参数与上同。
[0049] 在步骤S203中，设定引力为唯一的外力，获取外力密度fiex。
[0050] 把引力作为唯一的外力时，可以由公式
[0051] fiex＝ρig (5)
[0052] 计算出外力密度。
[0053] 在步骤S204中，结合层
粘度与子粒子尺度计算粘度力密度fiv。
[0054] 具体的，血液粒子流经小口径的孔或形状不规则的容器(如动脉瘤)时，不能简单的将血液看作是
层流，应结合层粘度与子粒子尺度来计算粘度力密度：
[0055]
[0056] 可以理解，步骤S101中所采用的基于光滑粒子动力学方法获取血管内流体粒子的压力密度、粘性力密度和外力密度只是其中一种可选的实施方式，本领域一般技术人员可以明白，还可以通过其它方法，如采用欧拉方法计算流体粒子的状态、或者采用有限元法模拟流体的运动等。
[0057] 在步骤S102中，使用质点弹簧模型计算位于血管壁与血流的交界处的边界粒子受到的弹回边界力。
[0058] 具体的，本步骤中采用弹回边界条件模拟血流和血管壁之间的相互作用。在不降低速度的条件下，将位于血管壁和血流的交界处的边界粒子，可以血流的流体粒子，也可以血管壁上的粒子上的力，分解为与边界垂直的方向上的力N(x)和与边界相切方向上的力T(y)，如图3所示，那么边界粒子所受的弹回边界力可以表示为：
[0059]
[0060] 其中，N(x)、T(y)分别表示垂直方向的力和切线方向的力；代表力的方向的单位矢量；x和y分别为边界粒子与流体粒子在垂直方向和切线方向上的距离。
[0061] 其中所述与边界垂直的方向上的力N(x)具体为：
[0062] 当边界粒子与流体粒子间的距离大于预定的
阈值距离时，所述与边界垂直的方向上的力N(x)为零；
[0063] 当边界粒子与流体粒子间的距离小于或等于预定的阈值距离、且所述与边界垂直的方向上的力N(x)小于或等于边界粒子弹簧连接中拉伸力的最大值时，所述与边界垂直的方向上的力N(x)为 所述 h为光滑粒子动力学方法的支持半径，x为边界粒子与流体粒子间的距离，p为流体粒子a与边界交互的压力，ρ是粒子密度，Vab为流体粒子a相对于相对边界粒子b的速度；nb为边界粒子b的法向方向；
[0064] 当所述N(x)大于边界粒子弹簧连接中拉伸力的最大值时，所述与边界垂直的方向上的力N(x)为零。
[0065] 在设计函数T(y)时，需保证当流体粒子平行移动到边界曲线时，从这些沿切线方向的边界粒子产生的合力保持不变，使流体粒子可以沿着边界平稳的流动。
[0066] 在步骤S103中，采用Runge-Kuntta模型整合所述计算的血管内的流体粒子的压力密度、粘性力密度、外力密度和位于血管壁与血流的交界处的边界粒子受到的弹回边界力，更新所述血流与血管壁的动态信息。
[0067] 具体的，根据所计算的血管内的流体粒子的压力密度、粘性力密度以及与边界粒子交互产生的外力，结合粒子的当前位置和粒子密度，根据
牛顿第二定律计算各粒子的
加速度、速度，根据所述计算的速度与加速度值重新确定粒子的位置。
[0068] 为了评估本发明实施例的可行性，下面进行了相应的实验进行验证。
[0069] 本实验针对不同数量的粒子进行时间性能评估。首先，构造一个动脉瘤模型验证改进后的弹回边界条件，通过实验对本方案中关键参数的性能的影响进行调查。为了对本方案中关键参数具有可比性，所有的实验都是在同一
硬件配置的电脑上进行。
[0070] 基于光滑粒子动力学方法与基于网格的解决方案相比，其优势在于具有较高的性能。图4为本发明实施例所述实验的实时性能表，从图4中所示实时性能表可以看出，当粒子数量增加时，
帧速率下降。当流体粒子达到大约5000个、边界粒子达到大约600个时，
帧速率可以保持在大约20fps，这样的结果可以满足大部分虚拟交互手术系统的要求。当粒子数量增加至10000左右时，帧速率将小于10fps。在实验过程中，我们设置时间步长为0.1s、支持半径为0.01、搜索半径为0.02且阈值距离为0.02。
[0071] 图5a-5d示出了一个包含了6961个粒子的动脉瘤模型。图5a为初始状态下的动脉瘤。图5b为经过20个时间步长后越来越多的粒子流入动脉瘤空腔内的状态。图5c经过50个时间步长后，边界粒子和流体粒子之间的排斥力超过边界弹簧的最大拉伸力，粒子从血管壁涌出即动脉瘤破裂，图5d为基于移动立方体方法的动脉瘤破裂示意图。在这组实验中，我们将时间步长设置为0.1s，支持半径设置为0.01。
[0072] 需要注意的是，一些关键参数可能会在极大程度上影响本方法的效果。通常是根据支持半径的大小、设置搜索半径、距离阈值。
[0073] 实施例二：
[0074] 图6示出了本发明第二实施例提供的模拟血流与血管壁作用的装置的结构框图，详述如下：
[0075] 本发明实施例所述模拟血流与血管壁作用的装置，包括：
[0076] 第一计算单元601，用于通过光滑粒子流体动力学方法计算血管内的流体粒子的压力密度、粘性力密度和外力密度；
[0077] 第二计算单元602，用于使用质点弹簧模型计算位于血管壁与血流的交界处的边界粒子受到的弹回边界力；
[0078] 更新单元603，用于采用Runge-Kuntta模型整合所述计算的血管内的流体粒子的压力密度、粘性力密度、外力密度和位于血管壁与血流的交界处的边界粒子受到的弹回边界力，更新所述血流与血管壁的动态信息。
[0079] 进一步的，所述第一计算单元601包括：
[0080] 查找子单元，用于在单位时间步长中，在设定的支持半径范围内查找与目标粒子i相邻近的相邻粒子j；
[0081] 第一获取子单元，用于根据光滑粒子流体动力学原理，获取血管内的目标粒子i的密度ρi以及粒子i的压力密度fip；
[0082] 第二获取子单元，用于设定引力为唯一的外力，获取外力密度fiex；
[0083] 计算子单元，用于结合层粘度与子粒子尺度计算粘度力密度fiv。
[0084] 具体的，所述查找子单元具体用于在单位时间步长中，在设定的支持半径范围内，采用成对接近方法、链接列表方法或者基于树的方法，查找与目标粒子i相邻近的相邻粒子j。
[0085] 所述第一获取子单元具体用于根据光滑粒子流体动力学原理，采用求和密度法或者继续密度法，获取血管内的目标粒子i的密度ρi。
[0086] 进一步的，第二计算单元602具体用于：
[0087] 当边界粒子与流体粒子间的距离大于预定的阈值距离时，所述与边界垂直的方向上的力N(x)为零；
[0088] 当边界粒子与流体粒子间的距离小于或等于预定的阈值距离、且所述与边界垂直的方向上的力N(x)小于或等于边界粒子弹簧连接中拉伸力的最大值时，所述与边界垂直的方向上的力N(x)为 所述 h为光滑粒子动力学方法的支持半径，x为边界粒子与流体粒子间的距离，p为流体粒子a与边界交互的压力，ρ是粒子密度，Vab为流体粒子a相对于相对边界粒子b的速度；nb为边界粒子b的法向方向；
[0089] 当所述N(x)大于边界粒子弹簧连接中拉伸力的最大值时，所述与边界垂直的方向上的力N(x)为零。
[0090] 所述获取单元用于根据所计算的血管内的流体粒子的压力密度、粘性力密度以及与边界粒子交互产生的外力，结合粒子的当前位置和粒子密度，根据牛顿第二定律计算各粒子的加速度、速度，根据所述计算的速度与加速度值重新确定粒子的位置。
[0091] 本发明实施例所述模拟血流与血管壁作用的装置与实施例一中所述模拟血流与血管壁作用的方法相对应，在此不作重复赘述。
[0092] 以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何
修改、等同替换和改进等，均应包含在本发明的保护范围之内。
