文摘

我们提出一个完全耦合的机电模型的心脏。模型集成了心脏电生理和心脏力学通过excitation-induced收缩和deformation-induced电流。基于有限元的数值方案实现在一台超级计算机。算例给出了使用薄心脏组织和一只狗心室与现实的几何。并行仿真方案的性能进行了研究。模型提供了一个有用的工具了解心血管动力学。

1。介绍

心血管疾病死亡的主要原因。复杂动力学的计算机仿真和可视化心脏有很大的潜力来提供定量指导心脏疾病的诊断和治疗。有开发准确的计算机模型深入研究促进心血管动力学机制的理解(1]。

灵感来自霍奇金和赫胥黎的开创性工作2),已研制出许多数学模型(3]。与此同时,各种各样的数学模型提出了机电模拟。纳什和Panfilov [4)提出了一个计算框架夫妇有三FitzHugh-Nagumo-type [5]excitation-tension模型控制方程的非线性应力平衡采用机电和mechanoelectric反馈。倪德勒et al。6]定量特征Ca的绑定^{2 +}过渡委员会,原肌球蛋白的动力学,结合位点的可用性,和crossbridge绑定的动力学扰动后肌节长度。Gurev et al。7]说明方法构造有限元机电模型的心脏和开发性器官的心室网基于磁共振扩散张量磁共振成像的心。的工作(7专注于心室网格的建设和没有考虑心脏电生理学上的机械收缩的影响。Goktepe和库尔8)提出了一个隐式和完全有限元素方法双向耦合的兴奋收缩的问题。电生理学被FitzHugh-Nagumo-type (FHN)模型(8]。柯南道尔et al。9)应用心脏力学的模拟的并行计算。他们使用非结构化网格评估模型的性能,并达到了最大加速因子为15.9当使用32个线程。Lafortune et al。10)建立了一个平行的机电模型的心脏。他们的模型能高效运行数百个处理器使用心室网的实际几何。Lafortune等人所描述的简单的电生理学有三FitzHugh-Nagumo-type FHN模型(5)或有三Fenton-Karma(颗)模型(11),采用“单向耦合”的位移不影响电生理学。此外,心脏力学行为的影响在其电子的行为,也称为mechanoelectric反馈,引起了研究者的关注(12- - - - - -15]。Mechanoelectric反馈可能是由于stretch-activated渠道或拉伸电气信号传播的影响。

这项工作的目标是开发一个双向耦合的机电模型的并行模拟复杂的心血管动力学。朝着这个目标,我们已经开发出一种完全耦合的机电模型的心脏,心脏电生理学一体,心脏产生的力学和双向耦合excitation-induced收缩和deformation-induced一代的电流。所描述的心脏电生理学Beeler-Reuter (BR)模型(16]。机电耦合问题是解决隐式使用有限元方法。并行计算算法使用消息传递接口(MPI)。模型是测试通过模拟薄心脏组织和狗心室与现实的几何。

本文组织如下。节2,我们介绍了生理模型。介绍了耦合问题的数值计算方法3。部分4给出了数值结果。我们将讨论论文的结果和结论部分5。

2。生理模型

心跳动的结果序列的电化学激发波从窦房结发起的。电脉冲诱导细胞内钙循环,进而导致心肌收缩。这一过程,称为兴奋收缩偶联(ECC),心脏的理解至关重要。另一方面,机械变化反应神经和激素的影响也对电气性能的影响。这种互补的概念被称为mechanoelectric反馈。参见图1电激活的关系、化学体内平衡和机械收缩。

2.1。心脏电生理学

几十个模型已经提出多年来模拟心脏电生理学(3]。这些模型大多来自霍奇金和赫胥黎的开创性工作2]。在这部作品中,Beeler-Reuter (BR)模型(16)采用数值插图。BR模型描述了单个细胞跨膜电压如下 在那里,代表跨膜电压,代表膜能力,总目前被描述为: 在这里,代表了电压门控钠电流,代表了长期有效的输出电流,代表了time-activated外电流,代表了一种缓慢的内向电流,和代表的外部刺激。注意,原始BR模型不包括stretch-activated通道的详细信息将在稍后讨论。刺激电流被选中的方波脉冲−80一个/F为1毫秒。我们参考读者16]BR模型的细节。

在心脏组织,(1)扩展到一个反应扩散形式包括空间扩散电流: 在哪里代表每个材料的空间坐标点的心;扩散张量,控制电的传导方向和速度波的激发心脏组织;膜电容,并设置为1F /厘米。

2.2。心脏力学

我们表示初始配置(舒张)Ω的心脏₀并通过Ω变形配置(收缩)。材料的位置向量的初始配置,在那里是一个矩形的单位基本向量笛卡尔坐标系统。表示质点的位置在时间通过。然后,物质点的空间坐标可以表示为。这个函数可以被看作是一个映射之间的初始配置和配置在时间吗。这两个措施,和由变形梯度有关,所示(4)。考虑以下: 有两种方法在描述一个连续体的变形:拉格朗日描述使用材料坐标作为独立变量和欧拉描述使用空间坐标作为独立的变量(17]。欧拉描述通常采用流体动力学。在这部作品中,利用拉格朗日描述方法。有两个配方的拉格朗日方法:拉格朗日公式和更新的拉格朗日公式。在总拉格朗日公式,方程离散对原始配置。相比之下,更新的拉格朗日公式是基于当前的配置和通常用于非线性、大变形。因此,我们在这项工作中使用更新的拉格朗日公式。

对于每一个时间步质点的位移,用定义的区别是,它的当前位置和以前的位置(5)。位移是由线性动量平衡,和方程是描述为(6)。在(6),基尔霍夫应力张量;占的体力或外部应用的压力。基尔霍夫应力由一个被动组件和一个活跃的组件。活跃的力量生成的电子激发后,将详细解释。被动组件是由基本力学方程(7)。在(7),和是瘸子常数控制各向同性的应激反应;代表肌纤维的被动的刚度。左Cauchy-Green张量表示和被定义为(8)。参数值被称为纳什和Panfilov (2004) (4)和Goktepe和库尔(2010)8]。考虑 让我们表示当地取向肌纤维的初始配置由一个单位向量在变形的配置,向量。在(7),代表了变形结构张量和被定义为(9)。在(9),在初始配置结构张量。结构张量在初始配置,在变形配置代表主导方向指定节点的社区(18,19]。考虑以下: 在(7),象征表示的系数决定肌纤维的刚度是否有效。它表明,当拉伸材料点,将1,否则为0。在数学上,它被定义为(10)。在(10),代表了心脏的拉伸材料点。在初始配置,在变形的配置,可能会导致一些材料点。此外,标量和在(7)被定义为(11)。考虑

2.3。机电耦合

总结部分2。1和2。2,耦合问题是由以下方程: 方程(12显示完整的心脏电生理和心脏力学的耦合。膜电位,同时解决了每个节点的空间坐标(12)。第一个方程(13)代表了无需通量边界条件对心脏的面域用。符号外表面正常吗。第二个方程(13)定义了强加的自然边界条件。第三个方程(13)显示必要的边界条件施加点是固定的,以确保机械的问题是定义良好的。的领域是对用本质边界条件。

2.3.1。Excitation-Induced收缩

节2。2基尔霍夫应力,我们介绍了活跃由excitation-induced生成收缩。从几何的角度看,积极的基尔霍夫应力的方向应由结构张量,它的大小是由跨膜电位控制的v。让积极压力的大小,我们有。

提出了相当多的模型来模拟压敏活性纤维张力(4,18]。在这项工作中,我们采用提出的简化方程(4]: 在(14),象征MPa / mV的最大活性纤维张力,和是静态电位−94.7 mV溴离子的心肌细胞模型。用的开关函数它决定了活性纤维张力速度将会改变对跨膜电位。参数的值,,,。当变化从−94.7 mV 20 mV,函数是如图2。

注意钙波动不是研究工作和活性纤维膜电位直接张力控制的简化,如图所示(14)。

2.3.2。扩散张量

在(3),细胞之间的互连是受扩散张量。它控制电的传导速度波激发的心脏组织。由于心脏组织的各向异性性质,实验中观察到的是传导显然是肌纤维方向的速度比其他方向。需要考虑额外的速度沿着纤维方向,扩散张量分为两个部分:。符号表示一个单位矩阵。的系数控制的速度各向同性转导四面八方和系数沿着纤维方向表示额外的速度。由于结构张量依赖于空间坐标,每个质点的扩散张量将会改变与心脏的重塑。

2.3.3。Deformation-Induced当前一代的

在(3),总离子跨膜电流由一个组件。拉伸激活通道的离子通道开放毛孔在回应机械变形的细胞膜19]。根据机制当前是什么引起的,有不同的配方激活渠道。在这项工作中,我们采用提出的公式(20.]: 在(15),最大电导;的静态电位stretch-activated渠道;的系数决定了肌纤维的刚度是否实际上,定义在(10);就是质点的心脏。

3所示。数值计算方法

控制方程(12)是解决采用算子分裂法(21]。首先,我们解决以下使用向前欧拉方法非线性常微分方程(22]: 那么,解决方案(16)是用来解决以下偏微分方程(17使用隐式欧拉方法)。这两个方程为每个时间步迭代得到解决。考虑以下: 弱形式的(17)构建遵循古典加勒金过程。疲软的形式获得通过的产品(17与测试函数)和域和集成它们。所需的时间独立测试函数C⁰并满足必要的边界条件。把测试函数相乘和的两个方程(18)和执行集成部分收益 应用自然边界条件(18)导致 在每个时间步,(19)线性化如下: 我们可以解决和从线性化方程。

传统的等参的伽辽金过程离散化连续弱式方程。域的心分解为子域,每个子域名是一种元素。然后字段变量和,两个相关的测试函数在每个子区域内插 在(21),每个元素和节点数量吗是C⁰interpolants,通常被称为有限元形状函数文献。隐式欧拉方法是利用离散化时时间导数项(18)。最后,我们取得了一个线性系统方程的形式: 线性系统方程的程度,在那里在网格的节点的数量。对于每个节点,是解决。他们是用于更新膜电位和每个节点的空间坐标。

模型是在c++中实现,并行使用消息传递接口(MPI) [23]。一个开源软件包称为美逖斯(24)是用于分区核心网,这样计算cpu之间的负载平衡。梅蒂斯人的算法是基于多级recursive-bisection,多级k路、多约束分区方案。

两个平行动力学被用来解决最后的线性系统。两个解决层次迭代并行解算器(臀部)25和Trilinos的解决方案26]。他们实现了广义极小剩余法(27]。当使用Trilinos,线性系统雅可比预调节器的先决条件。模拟运行在超级计算机,巨妖(28]。开源软件访问(29日)是用于可视化。

4所示。数值结果

4.1。一层薄薄的心脏组织

我们第一次在一层薄薄的广场进行了模拟心脏组织的尺寸0.4 * 0.4 * 0.001厘米³。在模拟中,左上的和右下方的角落是固定的,沿着垂直方向和纤维取向。我们使用向前欧拉方法解决ODE (16在0.005毫秒的时间步。在pde (17使用隐式欧拉方法)是解决时间步长为0.1 ms。

动作电位的一个节点(−−0.2厘米,0.2厘米,0.0005厘米),如图3获得了使用不同的网格尺寸。数值结果表明,一致的动作电位反应了使用不同的网格尺寸。因为心脏组织很薄,这是当作一个2 d组织。网格大小对和方向是每个0.4厘米长。

然后我们的机电模拟执行组织(0.4厘米* 0.4厘米* 0.001厘米)与网格大小= 0.004厘米。刺激的中心组织实施。图4显示了电波传播而不考虑收缩。电波传播的对称中心对整个组织。这个测试验证的反应扩散模型的一部分复制电波传播的基本现象在一个方形的组织。注意,系数在这个模拟设置为0,这样的组织是各向同性的。

图5显示了电子传播以及excitation-induced收缩在一个组织。在该测试中,垂直的纤维组织是一致的。因此,组织在垂直方向变形,清楚地显示在时间女士在图5。此外,与之前的测试,我们认为是额外的速度沿着纤维方向,如所示。的系数被设置为。它很容易观察到的从图5在垂直方向上的传播显然是更快。

4.2。狗心室与现实的几何

我们还模拟狗心室的收缩与现实的几何。检查两个网格。第一个880年由六面体网格元素。第二个元素有190080六面体网格是第一个使用软件中精炼出来的肘桑迪亚国家实验室开发(30.]。参见图6对原始网格和精制。心室网获得的(31日]。

本构和耦合模型,各向异性导电性和其他参数使用的一样。的刺激对上级部分心室如图7。一些节点上表面是克制的,问题是定义良好的。我们假设正常的纤维在任何时候是指向几何中心。在这种假设下,纤维形成的肌肉层的心。我们注意到,这种假设可能不接近现实的纤维布局的心脏。然而,假定结构允许我们测试的效率和鲁棒性的计算算法的心,从点对点的纤维取向的变化。

图8显示了一个狗心室机电模拟与现实的几何。电刺激对隔膜的上表面。收缩状态保持了100 ms,然后慢慢恢复到静息状态。

我们还模拟心室外表面附近的一个疤痕。我们采用了疤痕的简化描述。具体来说,我们认为疤痕没有传导能力和能保持被动机械收缩与其他细胞。疤痕的大小,如图9。约1厘米的疤痕有半径的表面积和心室壁厚度相似。

图9显示心室收缩的形状和膜电位在六个时刻的空间分布。比较图9与图8,可以观察到明显的区别。疤痕时,膜电位明显小,这可能导致更小的活性纤维张力,从而削弱了泵送能力。虽然没有严重的生理以来,本研究的结论可以给我们模拟器皿初步和缺乏实验验证,模拟演示的能力模型来研究真正的心。

4.3。性能分析

为我们的模型并行效率是至关重要的。使用网格的需求与成千上万的元素来实现高准确度和分辨率使计算效率一个巨大的挑战。并行效率测量分析的可伸缩性。使用156万六面体网格元素。在这个分析中,我们使用120、240和480年的核心。使用120 cpu时所花费的时间作为参考价值。图10显示了很强的可伸缩性480芯的巨妖的超级计算机。在这种情况下,240核的可伸缩性是大约90%的理想和480核的可伸缩性是大约70%的理想。可伸缩性的减少是由于增加的比例数量的通讯每个迭代和鬼的数量的比例的增加节点。扩展能力有限,PDE解算器的性能从Trilinos包26]。在未来的研究中,将受益于选择性能更好的并行求解算法。同样的性能得到改善,如果我们可以使用CUDA (32)这是一个并行计算平台和编程模型NVIDIA发明的。

5。讨论和结论

这是一个常见的方法在文献中以迭代的方式解决电机机械问题。在每个步骤中,电气问题已经解决了,然后结果从电气解决方案被提交到机械的问题,他的解决方案是使用在下一步解决电的问题。从电气和机械问题分别解决问题,这也是一个常见的做法采取细网格的电场和机械领域更粗的网格。

这项工作的贡献是开发一个完全耦合方案,准确解决机电的心脏。我们已经开发出一种心脏机电模型于一体的心脏电生理学,机械收缩,以及它们的相互作用。现实的生理模型用来描述的电气和机械功能的心。模型已使用隐式数值求解,有限元素方法。数值模拟已进行了使用并行模拟组织不同的几何图形。心脏力学所描述的更新拉格朗日方法,从当前配置视图问题,微分和积分的空间坐标。在网的角度描述,更新的拉格朗日描述的特征是使材料点保持一致与网格点。因此,拉格朗日描述简化了边界条件的实施以来,元素的边界网格保持与材料边界重合。开发模型和计算机代码验证每一步使用简单的测试例子来确保数值计算的准确性。多个模拟进行了使用各种网格和参数,以确保开发模型的数值鲁棒性。

由于计算涉及数以百万计的节点,当前框架为实时应用有一定的局限性。在未来的研究中,我们将进一步提高使用效率更高的性能模型PDE解决或在CUDA实现模型(32]。本文采用了一种简化的纤维配置中,在正常的纤维被认为是指向心脏的几何中心。更现实的心形和纤维配置应该利用在未来为生理上忠实的工作参数的研究。与材料,此外,由于拉格朗日网格变形网格可能会变得扭曲,如果心太大的变形。

未来的工作也可以考虑更详细的模型等积极强调混合模型(18]。在混合模型,依赖于活跃的力量,而不是跨膜电位,这是更多的生物物理合理。此外,细胞内的混合模型还考虑绑定^{2 +}肌钙蛋白C、配置更改的原肌球蛋白和肌动蛋白和肌凝蛋白的相互作用,这是一个更准确的描述在兴奋收缩交互。数值模拟进行了验证模型实现本文的目的,但更多的数值试验将使用平台执行调查机电功能之间的相互作用在心脏和心律失常的影响。

确认

这项工作是在NSF在格兰特没有支持的一部分。cmmi - 0845753。这项研究受到了分配通过TeraGrid先进的支持项目,在批准号BCS090012。一些模拟进行牛顿在UTK高性能集群。h·夏支持通过研究生研究助理职务国家数学研究所和生物合成,由美国国家科学基金会赞助的一个研究所,美国国土安全部和美国农业部通过NSF奖。ef - 0832858,额外的田纳西大学的支持下,诺克斯维尔。