CN107330972B - 模拟生物力学特性的实时软组织形变方法和系统 - Google Patents

模拟生物力学特性的实时软组织形变方法和系统 Download PDF

Info

Publication number
CN107330972B
CN107330972B CN201710508385.6A CN201710508385A CN107330972B CN 107330972 B CN107330972 B CN 107330972B CN 201710508385 A CN201710508385 A CN 201710508385A CN 107330972 B CN107330972 B CN 107330972B
Authority
CN
China
Prior art keywords
soft tissue
constraint
collision detection
force
spring
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201710508385.6A
Other languages
English (en)
Other versions
CN107330972A (zh
Inventor
刘谦
徐浪
卢昱华
李文成
郭琛琛
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Huazhong University of Science and Technology
Ezhou Institute of Industrial Technology Huazhong University of Science and Technology
Original Assignee
Huazhong University of Science and Technology
Ezhou Institute of Industrial Technology Huazhong University of Science and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Huazhong University of Science and Technology, Ezhou Institute of Industrial Technology Huazhong University of Science and Technology filed Critical Huazhong University of Science and Technology
Priority to CN201710508385.6A priority Critical patent/CN107330972B/zh
Publication of CN107330972A publication Critical patent/CN107330972A/zh
Application granted granted Critical
Publication of CN107330972B publication Critical patent/CN107330972B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00Three dimensional [3D] modelling, e.g. data description of 3D objects

Abstract

一种模拟生物力学特性的实时软组织形变方法,其包括如下步骤:S1、基于软组织器官三维可视化数据生成软组织物理模型数据;S2、基于软组织物理模型数据生成碰撞检测模型数据;S3、通过游戏引擎加载物理模型数据和碰撞检测模型数据进行碰撞检测;S4、通过约束的最优化求解方法来计算软组织变形后的状态。

Description

模拟生物力学特性的实时软组织形变方法和系统
技术领域
本发明涉及计算机图形学技术领域,特别涉及一种模拟生物力学特性的实时软组织形变方法和系统。
背景技术
手术模拟器能够无损耗的帮助训练者熟悉真实手术流程和操作。它通过软件上的可视化渲染和硬件机械手上的交互反馈实现了科学和可重复的训练。相较于基于生物器官的手术训练,它的学习成本低,而且训练有效性也得到保证。软组织变形的模拟涉及到手术模拟器的三个关键部分:手术环境的渲染,器械器官的交互以及力反馈。因此它是手术模拟器开发中最重要的部分。通常一个优秀的形变算法使得模拟器更接近临床表现。同时准确的形变,撕裂,流血及力反馈让训练更有效。为了达到流畅和真实的操作。软组织变形模拟算法实现时需要达到真实性和实时性的兼容。然而软组织作为粘弹性体,表现出更复杂的非线性形变行为,比较于完全弹性和完全黏性体。包括蠕变,应力松弛或者不可压缩等几何特征。因此开发一个稳定而快速的实时模拟软组织变形的非线性行为的方法很重要。实时形变算法的发展伴随着硬件计算能力的提升。在早期工作中,由于较低的计算表现,几何技术被用于模拟形变。在这些情况中,物理的准确性被牺牲给计算效率,系统也不了解对象的材料特性。随着计算能力的进步,基于物理的方法发展了起来。在这个领域也有着两个不同的方向。一个方向是求解基于连续物理模型的偏微分方程。比较典型的有基于网格划分的有限元方法,该方法使用同样的策略比如有限体积法和有限差分法。同时另一类基于邻近支持域离散的无网格方法最近也得到应用,如SPH,RPIM等。基于连续力学的物理模型的优势是准确性足够到能够反映真实对象的物理属性。他们被广泛应用于结构分析和电磁辐射模拟计算。但是在手术模拟的应用中的目标是达到能够反映软组织的重要物理特征,比如粘弹性不可压缩性等。因此,在这些实时的领域,选择另一类所谓看起来像物理的方法是合适的,比如弹簧振子模型,链式模型和基于位置的动态模型。通过抽象出典型的特征到这类模型中,他们能够达到很高的计算效率。
发明内容
有鉴于此,本发明提出一种模拟生物力学特性的实时软组织形变方法和系统。
一种模拟生物力学特性的实时软组织形变方法,其包括如下步骤:
S1、基于软组织器官三维可视化数据生成软组织物理模型数据;
S2、基于软组织物理模型数据生成碰撞检测模型数据;
S3、通过游戏引擎加载物理模型数据和碰撞检测模型数据进行碰撞检测;
S4、通过约束的最优化求解方法来计算软组织变形后的状态。
在本发明所述的模拟生物力学特性的实时软组织形变方法中,
所述步骤S1包括:
分析软组织的拓扑结构,引入相应的约束和质点分布实现表征器官的变形变现与拓扑结构相关的效果。
在本发明所述的模拟生物力学特性的实时软组织形变方法中,
所述步骤S2包括:
通过三角形和球形的基本单元来构成软组织物理模型;
通过改进的空间hash算法对软组织物理模型中基本单元进行AABB包围盒的计算;通过空间hash算法将AABB包围盒索引存储于一个hash数组中;遍历所有的hash数组对数组内的元素进行一一碰撞检测;AABB包围盒用于投影基本单元对象到坐标轴的三个不同方向,同时维持一个二维的数组来存储每一个基元与基元之间的碰撞检测状态来避免重复检测碰撞;
在碰撞检测结束后,进行碰撞响应的计算实现达到移动基本单元到不碰撞的状态。
在本发明所述的模拟生物力学特性的实时软组织形变方法中,
所述步骤S3包括:
通过Unity3D游戏引擎进行虚拟手术场景的渲染构建;
渲染构建包括场景的碰撞检测和物理计算,场景的渲染和机械手设备数据采集处理功能被置于Unity3D运行时环境中。
在本发明所述的模拟生物力学特性的实时软组织形变方法中,
所述步骤S4包括:
将需要模拟的对象以质点的集合,约束的集合以及弹簧的集合来进行离散,在每一个固定时间步长中,需要通过计算得到的数据为质点集合的位置;
对于生物力学相关的非线性和粘弹性通过弹簧振子模型来描述,对于拓扑结构相关的特性通过基于位置的动态模拟来描述。
本发明还提供一种模拟生物力学特性的实时软组织形变系统,其包括如下单元:
软组织物理模型生成单元,用于基于软组织器官三维可视化数据生成软组织物理模型数据;
碰撞检测模型生成单元,用于基于软组织物理模型数据生成碰撞检测模型数据;
碰撞检测单元,用于通过游戏引擎加载物理模型数据和碰撞检测模型数据进行碰撞检测;
软组织变形计算单元,用于通过约束的最优化求解方法来计算软组织变形后的状态。
在本发明所述的模拟生物力学特性的实时软组织形变系统中,
所述软组织物理模型生成单元包括:
分析软组织的拓扑结构,引入相应的约束和质点分布实现表征器官的变形变现与拓扑结构相关的效果。
在本发明所述的模拟生物力学特性的实时软组织形变系统中,
所述碰撞检测模型生成单元包括:
通过三角形和球形的基本单元来构成软组织物理模型;
通过改进的空间hash算法对软组织物理模型中基本单元进行AABB包围盒的计算;通过空间hash算法将AABB包围盒索引存储于一个hash数组中;遍历所有的hash数组对数组内的元素进行一一碰撞检测;AABB包围盒用于投影基本单元对象到坐标轴的三个不同方向,同时维持一个二维的数组来存储每一个基元与基元之间的碰撞检测状态来避免重复检测碰撞;
在碰撞检测结束后,进行碰撞响应的计算实现达到移动基本单元到不碰撞的状态。
在本发明所述的模拟生物力学特性的实时软组织形变系统中,
所述碰撞检测单元包括:
通过Unity3D游戏引擎进行虚拟手术场景的渲染构建;
渲染构建包括场景的碰撞检测和物理计算,场景的渲染和机械手设备数据采集处理功能被置于Unity3D运行时环境中。
在本发明所述的模拟生物力学特性的实时软组织形变系统中,
所述软组织变形计算单元包括:
将需要模拟的对象以质点的集合,约束的集合以及弹簧的集合来进行离散,在每一个固定时间步长中,需要通过计算得到的数据为质点集合的位置;
对于生物力学相关的非线性和粘弹性通过弹簧振子模型来描述,对于拓扑结构相关的特性通过基于位置的动态模拟来描述。
实施本发明提供的模拟生物力学特性的实时软组织形变方法和系统与现有技术相比具有以下有益效果:首先针对于实体模型和空腔模型采用了不同的特征约束来表现拓扑结构上的不同,而对于非线性,粘弹性等生物力学特征通过粘弹簧模型来表现。在最后质点位置求解时,通过基于高斯赛德尔迭代的位置动态模拟实现无条件稳定的求解。此外本发明的算法简单,速度快,可有效地用于虚拟手术系统的软组织形变仿真中。
附图说明
图1是本发明实施例的模拟生物力学特性的实时软组织形变方法流程图;
图2是实体肝脏物理模型以及空腔胆囊物理模型构建过程示意图;
图3是构建腹腔镜胆囊切除场景图;
图4是a部分为原始状态,b部分为器械提拉肝脏的效果,c部分为器械同时夹持肝脏和胆囊的效果;
图5是点与点碰撞约束示意图;
图6是点与三角形碰撞约束示意图;
图7是变形测试图;
图8是体约束效果图;
图9是全局约束效果图;
图10是方法稳定性测试图。
具体实施方式
如图1所示,本发明实施例提供一种模拟生物力学特性的实时软组织形变方法,包括如下步骤:
S1、基于软组织器官三维可视化数据生成软组织物理模型数据;
S2、基于软组织物理模型数据生成碰撞检测模型数据;
S3、通过游戏引擎加载物理模型数据和碰撞检测模型数据进行碰撞检测;
S4、通过约束的最优化求解方法来计算软组织变形后的状态。
在本发明一实施例中,所述步骤S1的具体实现方式为:首先分析软组织的拓扑结构,因为器官的变形变现与组织结构特征具有重要联系。因此表征这些与拓扑结构相关的效果,通过引入相应的约束和质点分布能够实现。举例实现的两类典型器官,肝脏和胆囊。在组织结构特点上,肝脏内部由紧密排列的组织细胞构成,而胆囊是由表皮细胞成囊状排列构成。对于肝脏的物理模型,使用四面体的拓扑结构来描述。通过四面体的边点面体信息来生成约束。而对于胆囊,使用三角面片表面约束来描述。通过三角面片边点面来生成约束。仔细来说,在模型构建上,首先是基于中国数字人的彩色切片数据进行三维重建得到器官的表面几何模型。而对于肝脏使用Tetgen进行有效的四面体网格划分,分别得到点边面体信息。除了基本生物组织模拟结构部分。对每一个四面体施加体积不变约束,如下,来表示实体模型的体积不变效果。而对于胆囊模型,直接使用表面几何特征,点表面信息。除了基本生物组织模拟部分。对于封闭表面施加了一个全局的体积约束,如下,来表示囊状模型的体积不变效果。
在本发明一实施例中,所述步骤S2的具体实现方式为:为了进行准确而快速的碰撞检测,使用了改进的空间hash算法。物理模型的基本单元为三角形和球形来构成。然后对基元进行AABB包围盒的计算。通过空间hash算法将AABB包围盒索引存储于一个hash数组中。遍历所有的hash数组对数组内的元素进行一一碰撞检测。AABB包围盒只是被用来投影基元对象到坐标轴的三个不同方向。同时维持了一个二维的数组来存储每一个基元与基元之间的碰撞检测状态来避免重复检测碰撞。在碰撞检测结束后,需要进行碰撞响应的计算来达到移动基元到不碰撞的状态。碰撞响应的过程,借助于传统的PBD,通过三种类型的距离约束来进行碰撞的移动。
在本发明一实施例中,所述步骤S3的具体实现方式为:通过Unity3D游戏引擎进行虚拟手术场景的渲染构建。包括场景的碰撞检测和物理计算,场景的渲染和机械手设备数据采集处理等功能被至于Unity3D运行时环境中。
在本发明一实施例中,所述步骤S4的具体实现方式为:将需要模拟的对象以质点的集合,约束的集合以及弹簧的集合来进行离散,在每一个固定时间步长中,需要通过计算得到的数据为质点集合的位置。质点的位置由作用于其上的力载荷所控制,力载荷ftatal分为约束力fconstraint,弹簧弹力fspringand外力fout,因此,基于牛顿第二定律有.
Figure GDA0003024831760000081
在上面的式子中,x=[x1 x2 … xn]T为质点位置列向量,M是质点质量对角矩阵,fconstraint=[fc1 fc2 … fcm]T为约束力列向量,fspring和fout是弹力和外力列向量。为了导出fconstraint,约束力势能被定义为
Figure GDA0003024831760000082
其中,C(x)=[C1C2…Cm]T是约束函数列向量,α是与约束强度相关的顺应程度,对应于重力势能和弹性势能等保守力势能力方向是平行于能量变化最大的方向.因此约束力被定义为势能的梯度方向
Figure GDA0003024831760000091
为了求解质点的位置,将等式(1)以时间步长Δt进行离散,n+1代表需要迭代到的时间步长位置。
Figure GDA0003024831760000092
通过引入拉格朗日乘子来分解力为标量和向量部分。通过定义
Figure GDA0003024831760000093
来简化后面的描述
Figure GDA0003024831760000094
通过将等式(5)(3)代入等式(4)中然后重新分布得到
Figure GDA0003024831760000095
等式右边的部分可以化简为
Figure GDA0003024831760000096
这是在约束投影前需要进行计算的预测位置。为了求解时间步长为n+1的质点位置。可以得到接下来的约束最优化求解方程组。
Figure GDA0003024831760000097
Figure GDA0003024831760000098
分析上面的方程,的目标是获取λ和x在时间步长为n+1的时候.上述方程组的变量个数是等于方程的个数,因此该方程组具有唯一解.为了求解该方程组,通过牛顿拉夫逊方法来线性化然后通过迭代的数值解法来求解.将迭代次数简化为i然后省略掉n+1.得到
Figure GDA0003024831760000101
Figure GDA0003024831760000102
的目标是在每一帧的情况下获取Δx和Δλ当迭代完成时通过λi+1=λi+Δλ和xi+1=xi+Δx来更新点的位置.应用了两个近似的方法来简化求解的过程.首先,考虑积分迭代次数为0的时候i=0而g(xi,λi)=0,当迭代进行的过程中g(xi,λi)保持一个很小的值,因此省略掉g(xi,λi).第二
Figure GDA0003024831760000103
省略掉第一部分并不会对结果的准确性有什么影响.因此让
Figure GDA0003024831760000104
经过这两个步骤之后,求解线性化后的方程(9)(10)为
Figure GDA0003024831760000105
Figure GDA0003024831760000106
计算等式(11)获取
Figure GDA0003024831760000107
将其代入到等式(12)获得
Figure GDA0003024831760000108
当Δλ被解出时,Δx也被求解出。
现在回过头来导出弹力fspring,非线性弹力方程如下,对于质点i有
Figure GDA0003024831760000111
衰减力fdemper
Figure GDA0003024831760000112
弹簧弹力为
Figure GDA0003024831760000113
Figure GDA0003024831760000114
Figure GDA0003024831760000115
kij
Figure GDA0003024831760000116
弹簧的伸长量直接相关.因此在伸长量较小时为伸长量的三阶相关,当伸长量较长时为线性相关。
从上面的过程可以看出,的模型是通过两部分来分别反应软组织力学特性。对于生物力学相关的非线性和粘弹性就通过弹簧振子模型来描述,对于拓扑结构相关的特性就通过基于位置的动态模拟来描述。
本部分结合附图对发明的具体实施步骤流程进行解释。
本发明提供了一种模拟生物力学特性的实时软组织形变方法,主要步骤介绍如下:
1、特征器官模型处理
首先依据需要模拟的器官组织结构,将器官分为实体模型和空腔囊状模型。对于这两类模型,首先根据形状特征构建三维表面面片数据。该数据的定义格式为fbx文件格式。基本的数据单元为,顶点集合,三角形集合和贴图材质等渲染数据。如图5,对于实体模型采集三维表面面片数据中的顶点信息和三角形信息利用TetGen四面体生成算法生成空间四面体信息。对于表面模型采集顶点信息和三角形信息判断表面模型的封闭情况,消除重复顶点和边,实现完整封闭的表面网格模型。
在预处理步骤结束后,对于实体模型遍历每一个空间四面体,生成对应的四面体体积不变约束。体积不变约束对应于顶点坐标为P1,P2,P3,P4的四面体对应的约束函数如下:
Figure GDA0003024831760000121
对于空腔模型直接对所有的表面三角形,生成对应的全局体积不变约束。对应的全局体积不变约束函数如下:
Figure GDA0003024831760000122
在特征模型约束生成完毕后,对模型的通用约束进行生成。包括特定的距离约束,过拉伸约束。为了表现出生物力学特性中的非线性和粘弹性,需要生成对应的弹簧振子模型,综合模型的控制参数如下:
表1.物理控制参数描述
Figure GDA0003024831760000123
Figure GDA0003024831760000131
依据以上的控制参数能够调试模型的模拟效果。
2、碰撞检测数据生成
在物理模型数据生成完毕后,为了实现器械和器官变形之间的交互功能,需要对模型进行碰撞检测数据的生成,如图2,图2中其中从a到c为实体肝脏物理模型构建过程,从d到f为空腔胆囊物理模型构建过程。如图3所示,首先对于器械使用空间球结构来进行模拟,对于器官组织使用空间三角形结构来进行模拟。通过遍历所有的碰撞检测基元单元(分别包括,空间三角形和空间球)直接生成AABB碰撞盒集合。AABB碰撞盒集合作为欧拉空间坐标系下的轴向六面体,通过判断空间基元在x轴y轴和z轴三个方向的最大最小值来构建。构建完成所有的AABB碰撞盒集合后。通过单元长度的大小,将碰撞检测空间离散为等大小的立方体集合。该集合为空间分布的三维数组,分别具有三个方向的参数。使用一个三维到一维的空间hash映射函数能够实现将立方体单元映射为哈希数组的一个元素。因此在该步骤中,能够获取一个一维分布的空间哈希数组。哈希数组里面存储的结构为空间立方体的索引。而在同一个索引中的碰撞基元才可能发生碰撞,因此只检测同一索引中的基元的碰撞情况。图4构建腹腔镜胆囊切除场景图。a部分为原始状态。b部分为器械提拉肝脏的效果。c部分为器械同时夹持肝脏和胆囊的效果
而对于基元的碰撞情况分为两类进行求解。一类碰撞情况是空间球和空间球之间的碰撞响应,使用点与点之间的距离约束求解,如图4,图4中a部分为原始状态。b部分为器械提拉肝脏的效果。c部分为器械同时夹持肝脏和胆囊的效果。二类碰撞情况是空间球和空间三角形之间的碰撞响应,使用点与三角形之间的距离约束求解,如图5。
3、场景渲染和物理计算
在物理模型数据生成和碰撞检测数据生成结束后,为了渲染出用户操作可见的虚拟现实场景。直接将场景数据导入到Unity3D游戏引擎中,该游戏引擎能够实现次时代的三维场景渲染。而的模拟代码直接运行在Unity3D支持的Mono运行环境时,代码以C#语言结构写成.整个方法的流程被总结如下。首先非线性和粘弹性等力学特性在弹簧受力计算步骤被表现。然后通过扩展后的弹力作为外力来估计质点位置。在这一步骤结束之后,约束投影步骤开始来修正质点的位置来满足约束函数的方程。尽管使用扩展后的约束动态来避免迭代次数和形变效果相关。但是,得益于高斯塞德尔方法,最终求得的质点位置会保持稳定。
Figure GDA0003024831760000141
Figure GDA0003024831760000151
Figure GDA0003024831760000161
在每帧迭代求解完成之后,顶点的信息数据被输出到Unity3D的Mesh面片结构中,用于修改场景的模型几何,达到场景中的器官变形效果。在应用本发明实现的腹腔镜胆囊切除手术模拟场景中,如图6,能够看到实际的效果,图6为以垂直于组织上表面的外力拉伸组织变形的效果。
图7至图10为生物力学特性和组织结构特性测试效果图。
图7表示以垂直于组织上表面的外力拉伸组织变形的效果。
图8表示在重力作用下下垂过程中,悬梁臂模型保持局部体积不变效果。
图9表示在平板压缩过程中,气球模型保持全局体积不变效果。
图10表示从a到c显示了本发明的方法从模型压缩后恢复原状的过程,从d到f显示了传统的弹簧振子模型从模型压缩后无法恢复原状。
本发明的技术关键点和欲保护点在于:
(1)技术关键点一:将生物力学特性和拓扑结构特性分别使用非线性粘弹簧和基于位置的动态方法来模拟,实现了快速稳定的求解。
(2)技术关键点二:对于不同组织结构的器官,采用具有拓扑结构特征的不同约束来表现,实现了组织特性的针对性表现。
可以理解的是,对于本领域的普通技术人员来说,可以根据本发明的技术构思做出其它各种相应的改变与变形,而所有这些改变与变形都应属于本发明权利要求的保护范围。

Claims (4)

1.一种模拟生物力学特性的实时软组织形变方法,其特征在于,其包括如下步骤:
S1、基于软组织器官三维可视化数据生成软组织物理模型数据,其中,所述步骤S1包括:分析软组织的拓扑结构,引入相应的约束和质点分布实现表征器官的变形变现与拓扑结构相关的效果,且所述器官分别为肝脏和胆囊,肝脏的物理模型使用四面体的拓扑结构来描述,胆囊的物理模型使用三角面片表面约束来描述;
S2、基于软组织物理模型数据生成碰撞检测模型数据,其中,所述步骤S2包括:通过三角形和球形的基本单元来构成软组织物理模型;通过改进的空间hash算法对软组织物理模型中基本单元进行AABB包围盒的计算;通过空间hash算法将AABB包围盒索引存储于一个hash数组中;遍历所有的hash数组对数组内的元素进行一一碰撞检测;AABB包围盒用于投影基本单元对象到坐标轴的三个不同方向,同时维持一个二维的数组来存储每一个基元与基元之间的碰撞检测状态来避免重复检测碰撞;在碰撞检测结束后,进行碰撞响应的计算实现达到移动基本单元到不碰撞的状态;
S3、通过游戏引擎加载物理模型数据和碰撞检测模型数据进行碰撞检测;
S4、通过约束的最优化求解方法来计算软组织变形后的状态,其中,所述步骤S4的具体实现方式为:将需要模拟的对象以质点的集合,约束的集合以及弹簧的集合来进行离散,在每一个固定时间步长中,需要通过计算得到的数据为质点集合的位置;质点的位置由作用于其上的力载荷所控制,力载荷ftotal分为约束力fconstraint,弹簧弹力fspring和外力fout,基于牛顿第二定律有:
Figure FDA0003024831750000028
其中,x=[x1 x2…xn]T为质点位置列向量,M是质点质量对角矩阵,fconstraint=[fc1fc2…fcm]T为约束力列向量,fspring和fout是弹簧的弹力和外力列向量:
为了导出fconstraint,约束力势能定义为:
Figure FDA0003024831750000021
其中,C(x)=[C1 C2…Cm]T为约束函数列向量,α是与约束强度相关的顺应程度,则约束力定义为势能的梯度方向:
Figure FDA0003024831750000022
对等式(1)以时间步长Δt进行离散,n+1代表迭代到的时间步长位置,即
Figure FDA0003024831750000023
定义
Figure FDA0003024831750000024
简化描述:
Figure FDA0003024831750000025
将等式(3)(5)代入等式(4)中重新分布得到:
Figure FDA0003024831750000026
求解约束最优化求解方程组:
Figure FDA0003024831750000027
Figure FDA0003024831750000031
其中,公式(7)为公式(6)的整理形式,公式(8)为公式(5)的离散形式;
通过牛顿拉夫逊方法来线性化然后通过迭代的数值解法来求解,将迭代次数简化为i然后省略掉n+1得到:
Figure FDA0003024831750000032
Figure FDA0003024831750000033
在每一帧的情况下获取Δx和Δλ当迭代完成时通过λi+1=λi+Δλ和xi+1=xi+Δx更新点的位置;
当积分迭代次数为0的时候i=0而g(xi,λi)=0;
Figure FDA0003024831750000034
求解线性化后的方程(9)(10)为:
Figure FDA0003024831750000035
Figure FDA0003024831750000036
计算等式(11)获取
Figure FDA0003024831750000037
将其代入到等式(12)获得:
Figure FDA0003024831750000038
当Δλ被解出时,Δx也被求解出;
导出弹力fspring,非线性弹力方程如下,对于质点i有:
Figure FDA0003024831750000039
衰减力fdamper为:
Figure FDA00030248317500000310
弹簧弹力为:
Figure FDA0003024831750000041
Figure FDA0003024831750000042
Figure FDA0003024831750000043
kij为弹性系数与
Figure FDA0003024831750000044
弹簧的伸长量直接相关,当伸长量较长时为线性相关。
2.如权利要求1所述的模拟生物力学特性的实时软组织形变方法,其特征在于,所述步骤S3包括:
通过Unity3D游戏引擎进行虚拟手术场景的渲染构建;
渲染构建包括场景的碰撞检测和物理计算,场景的渲染和机械手设备数据采集处理功能被置于Unity3D运行时环境中。
3.一种模拟生物力学特性的实时软组织形变系统,其特征在于,其包括如下单元:
软组织物理模型生成单元,用于基于软组织器官三维可视化数据生成软组织物理模型数据,其中,所述软组织物理模型生成单元包括:分析软组织的拓扑结构,引入相应的约束和质点分布实现表征器官的变形变现与拓扑结构相关的效果,且所述器官分别为肝脏和胆囊,肝脏的物理模型使用四面体的拓扑结构来描述,胆囊的物理模型使用三角面片表面约束来描述;
碰撞检测模型生成单元,用于基于软组织物理模型数据生成碰撞检测模型数据,其中,所述碰撞检测模型生成单元包括:通过三角形和球形的基本单元来构成软组织物理模型;通过改进的空间hash算法对软组织物理模型中基本单元进行AABB包围盒的计算;通过空间hash算法将AABB包围盒索引存储于一个hash数组中;遍历所有的hash数组对数组内的元素进行一一碰撞检测;AABB包围盒用于投影基本单元对象到坐标轴的三个不同方向,同时维持一个二维的数组来存储每一个基元与基元之间的碰撞检测状态来避免重复检测碰撞;在碰撞检测结束后,进行碰撞响应的计算实现达到移动基本单元到不碰撞的状态;
碰撞检测单元,用于通过游戏引擎加载物理模型数据和碰撞检测模型数据进行碰撞检测:
软组织变形计算单元,用于通过约束的最优化求解方法来计算软组织变形后的状态,其中,所述步骤S4的具体实现方式为:将需要模拟的对象以质点的集合,约束的集合以及弹簧的集合来进行离散,在每一个固定时间步长中,需要通过计算得到的数据为质点集合的位置;质点的位置由作用于其上的力载荷所控制,力载荷ftotal分为约束力fconstraint,弹簧弹力fspring和外力fout,基于牛顿第二定律有:
Figure FDA0003024831750000052
其中,x=[x1 x2…xn]T为质点位置列向量,M是质点质量对角矩阵,fconstraint=[fc1fc2…fcm]T为约束力列向量,fspring和fout是弹簧的弹力和外力列向量;
为了导出fconstraint,约束力势能定义为:
Figure FDA0003024831750000051
其中,C(x=[C1 C2…Cm]T为约束函数列向量,α是与约束强度相关的顺应程度,则约束力定义为势能的梯度方向:
Figure FDA0003024831750000061
对等式(1)以时间步长Δt进行离散,n+1代表迭代到的时间步长位置,即
Figure FDA0003024831750000062
定义
Figure FDA0003024831750000063
简化描述:
Figure FDA0003024831750000064
将等式(3)(5)代入等式(4)中重新分布得到:
Figure FDA0003024831750000065
求解约束最优化求解方程组:
Figure FDA0003024831750000066
Figure FDA0003024831750000067
其中,公式(7)为公式(6)的整理形式,公式(8)为公式(5)的离散形式;
通过牛顿拉夫逊方法来线性化然后通过迭代的数值解法来求解,将迭代次数简化为i然后省略掉n+1得到:
Figure FDA0003024831750000068
Figure FDA0003024831750000069
在每一帧的情况下获取Δx和Δλ当迭代完成时通过λi+1=λi+Δλ和xi+1=xi+Δx更新点的位置;
当积分迭代次数为0的时候i=0而g(xi,λi)=0:
Figure FDA0003024831750000071
求解线性化后的方程(9)(10)为:
Figure FDA0003024831750000072
Figure FDA0003024831750000073
计算等式(11)获取
Figure FDA0003024831750000074
将其代入到等式(12)获得:
Figure FDA0003024831750000075
当Δλ被解出时,Δx也被求解出;
导出弹力fspring,非线性弹力方程如下,对于质点i有:
Figure FDA0003024831750000076
衰减力fdamper为:
Figure FDA0003024831750000077
弹簧弹力为:
Figure FDA0003024831750000078
Figure FDA0003024831750000079
Figure FDA00030248317500000710
kij为弹性系数与
Figure FDA00030248317500000711
弹簧的伸长量直接相关,当伸长量较长时为线性相关。
4.如权利要求3所述的模拟生物力学特性的实时软组织形变系统,其特征在于,所述碰撞检测单元包括:
通过Unity3D游戏引擎进行虚拟手术场景的渲染构建;
渲染构建包括场景的碰撞检测和物理计算,场景的渲染和机械手设备数据采集处理功能被置于Unity3D运行时环境中。
CN201710508385.6A 2017-06-28 2017-06-28 模拟生物力学特性的实时软组织形变方法和系统 Active CN107330972B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710508385.6A CN107330972B (zh) 2017-06-28 2017-06-28 模拟生物力学特性的实时软组织形变方法和系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710508385.6A CN107330972B (zh) 2017-06-28 2017-06-28 模拟生物力学特性的实时软组织形变方法和系统

Publications (2)

Publication Number Publication Date
CN107330972A CN107330972A (zh) 2017-11-07
CN107330972B true CN107330972B (zh) 2021-06-25

Family

ID=60198187

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710508385.6A Active CN107330972B (zh) 2017-06-28 2017-06-28 模拟生物力学特性的实时软组织形变方法和系统

Country Status (1)

Country Link
CN (1) CN107330972B (zh)

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108256246B (zh) * 2018-01-25 2021-12-07 深圳意动航空科技有限公司 一种基于Unity3D的多孔材料的设计方法和装置
CN108577806B (zh) * 2018-03-12 2021-01-19 华中科技大学鄂州工业技术研究院 基于低帧率激光散斑衬比成像的粘弹性检测系统及方法
CN108511074B (zh) * 2018-03-26 2021-11-09 福建师范大学福清分校 一种基于空间核映射和子空间聚集的软组织形变方法
CN109993830B (zh) * 2019-04-10 2022-05-03 福建师范大学福清分校 一种软组织表面刺穿形变仿真方法及终端
CN110074866B (zh) * 2019-05-07 2022-02-15 广东工业大学 一种软组织手术机器人的控制方法、装置及设备
CN111467796B (zh) * 2020-04-02 2023-05-02 北京代码乾坤科技有限公司 分布式物理结算方法和装置
CN111569423B (zh) * 2020-05-14 2023-06-13 北京代码乾坤科技有限公司 碰撞形态的修正方法和装置
CN112802203B (zh) * 2021-01-18 2023-02-28 华中科技大学鄂州工业技术研究院 一种基于特征的空间哈希连续碰撞检测方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101587372A (zh) * 2009-07-03 2009-11-25 东南大学 一种用于虚拟现实人机交互的增强力触觉建模方法
JP2012202707A (ja) * 2011-03-23 2012-10-22 Toyota Central R&D Labs Inc ダミー計測装置
CN104679958A (zh) * 2015-03-12 2015-06-03 北京师范大学 基于弹簧模型的球b样条编针织物形变仿真的方法
CN104778894A (zh) * 2015-04-28 2015-07-15 关宏刚 一种虚拟仿真正骨手法培训系统及其建立方法
CN106528993A (zh) * 2016-11-03 2017-03-22 南京信息工程大学 基于由碟形弹簧片构成的组合弹簧虚拟模型建模方法

Family Cites Families (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
AU2007233560A1 (en) * 2006-03-30 2007-10-11 Monash University Method of modelling the interaction between deformable objects
CN101404039B (zh) * 2008-03-28 2010-06-16 华南师范大学 一种虚拟手术的方法及其装置
CN101441831A (zh) * 2008-11-17 2009-05-27 江苏科技大学 基于力反馈的虚拟手术仿真系统
CN103793552B (zh) * 2013-12-18 2017-07-11 哈尔滨工程大学 一种软组织形变的局部质点弹簧模型的实时动态生成方法
US10296707B2 (en) * 2014-04-10 2019-05-21 Siemens Healthcare Gmbh System and method for patient-specific image-based guidance of cardiac arrhythmia therapies
CN104200730B (zh) * 2014-09-09 2017-05-10 华中科技大学 一种用于虚拟腹腔镜手术的装置、方法及系统
US10140745B2 (en) * 2015-01-09 2018-11-27 Vital Mechanics Research Inc. Methods and systems for computer-based animation of musculoskeletal systems
CN105513130B (zh) * 2016-02-01 2018-02-16 福建师范大学福清分校 一种基于网格与无网格混合的软组织形变方法
CN106781941A (zh) * 2016-11-24 2017-05-31 北京理工大学 一种用于模拟微创穿刺手术的方法及其系统

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101587372A (zh) * 2009-07-03 2009-11-25 东南大学 一种用于虚拟现实人机交互的增强力触觉建模方法
JP2012202707A (ja) * 2011-03-23 2012-10-22 Toyota Central R&D Labs Inc ダミー計測装置
CN104679958A (zh) * 2015-03-12 2015-06-03 北京师范大学 基于弹簧模型的球b样条编针织物形变仿真的方法
CN104778894A (zh) * 2015-04-28 2015-07-15 关宏刚 一种虚拟仿真正骨手法培训系统及其建立方法
CN106528993A (zh) * 2016-11-03 2017-03-22 南京信息工程大学 基于由碟形弹簧片构成的组合弹簧虚拟模型建模方法

Also Published As

Publication number Publication date
CN107330972A (zh) 2017-11-07

Similar Documents

Publication Publication Date Title
CN107330972B (zh) 模拟生物力学特性的实时软组织形变方法和系统
Teran et al. Finite volume methods for the simulation of skeletal muscle
Bender et al. Position-based Methods for the Simulation of Solid Objects in Computer Graphics.
Müller et al. Real-time simulation of deformation and fracture of stiff materials
Faure et al. Sparse meshless models of complex deformable solids
San-Vicente et al. Cubical mass-spring model design based on a tensile deformation test and nonlinear material model
Lim et al. Real time simulation of nonlinear tissue response in virtual surgery using the point collocation-based method of finite spheres
CN106875462B (zh) 一种基于元球模型和混合驱动方法的实时数字器官切割方法
Li et al. A surface mass-spring model with new flexion springs and collision detection algorithms based on volume structure for real-time soft-tissue deformation interaction
Lee et al. Volumetric Object Modeling Using Internal Shape Preserving Constraint in Unity 3D.
Nesme et al. Physically realistic interactive simulation for biological soft tissues
CA3182690A1 (en) Method and system for real-time simulation of elastic body
Liu et al. Deformation of soft tissue and force feedback using the smoothed particle hydrodynamics
Zou et al. A new deformation simulation algorithm for elastic-plastic objects based on splat primitives
Li et al. Soft articulated characters in projective dynamics
Aras et al. An analytic meshless enrichment function for handling discontinuities in interactive surgical simulation
Peng et al. Bi-potential and co-rotational formulations applied for real time simulation involving friction and large deformation
Faure et al. An Implicit Tensor-Mass solver on the GPU for soft bodies simulation
Tagawa et al. Adaptive and embedded deformation model: An approach to haptic interaction with complex inhomogeneous elastic objects
Va et al. Real-time volume preserving constraints for volumetric model on GPU
Halic et al. Soft tissue deformation and optimized data structures for mass spring methods
Mendoza et al. Issues in deformable virtual objects simulation with force feedback
Huang et al. Virtual surgery deformable modelling employing GPU based computation
Zhang et al. Real time simulation of tissue cutting based on gpu and cuda for surgical training
Shi et al. GPU in haptic rendering of deformable objects

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant