发明名称 软组织形变仿真方法
摘要 本发明涉及一种软组织形变仿真方法,包括以下步骤:建立软组织的生物力学模型、并对其中的各个质点进行初始化;力反馈设备对软组织施加作用力,进行碰撞检测;用改进的欧拉算法计算运动状态信息;将模型每个时间步长的状态输出到显示屏上,动态显示软组织形变过程;计算反馈力并输出触觉反馈。通过采用上述步骤,能够有效解决虚拟手术仿真中的实时性、精确性以及反馈力光滑性问题,提高了软组织形变仿真的精度和实时性,从而满足虚拟手术仿真的需要。
申请公布号 CN103400023B 申请公布日期 2016.11.02
申请号 CN201310267557.7 申请日期 2013.06.28
申请人 华北水利水电大学 发明人 刘雪梅;王瑞艺;孙新娟;姚建斌;宋连公;冯飞*;郭松;高阳;李海瑞;朱韶楠
分类号 G06F19/00(2011.01)I 主分类号 G06F19/00(2011.01)I
代理机构 郑州大通专利商标代理有限公司 41111 代理人 陈大通
主权项 一种软组织形变仿真方法,其特征在于包含以下步骤:步骤1):采集软组织的数据信息,并采用基于四面体的质点‑弹簧体模型建立软组织的生物力学模型,该生物力学模型由n个质点组成,并且对于其中任意一个质点i,满足以下方程:<maths num="0001"><math><![CDATA[<mrow><msub><mi>m</mi><mi>i</mi></msub><mfrac><mrow><msup><mo>&part;</mo><mn>2</mn></msup><msub><mi>x</mi><mi>i</mi></msub></mrow><mrow><mo>&part;</mo><msup><mi>t</mi><mn>2</mn></msup></mrow></mfrac><mo>+</mo><msubsup><mi>F</mi><mrow><mi>i</mi><mi>n</mi></mrow><mi>i</mi></msubsup><mo>=</mo><msubsup><mi>F</mi><mrow><mi>e</mi><mi>x</mi><mi>t</mi></mrow><mi>i</mi></msubsup></mrow>]]></math><img file="FDA0001010635540000011.GIF" wi="362" he="103" /></maths><maths num="0002"><math><![CDATA[<mrow><msubsup><mi>F</mi><mrow><mi>i</mi><mi>n</mi></mrow><mi>i</mi></msubsup><mo>=</mo><msubsup><mi>&Sigma;</mi><mrow><mi>j</mi><mo>=</mo><mi>j</mi><mn>1</mn></mrow><mrow><mi>j</mi><mi>q</mi></mrow></msubsup><mi>&mu;</mi><mrow><mo>(</mo><mo>|</mo><msub><mi>x</mi><mi>j</mi></msub><mo>-</mo><msub><mi>x</mi><mi>i</mi></msub><mo>|</mo><mo>-</mo><mo>|</mo><msub><mi>x</mi><mi>j</mi></msub><mo>-</mo><msub><mi>x</mi><mi>i</mi></msub><msup><mo>|</mo><mn>0</mn></msup><mo>)</mo></mrow><mfrac><mrow><mo>(</mo><msub><mi>x</mi><mi>j</mi></msub><mo>-</mo><msub><mi>x</mi><mi>i</mi></msub><mo>)</mo></mrow><mrow><mo>|</mo><msub><mi>x</mi><mi>j</mi></msub><mo>-</mo><msub><mi>x</mi><mi>i</mi></msub><mo>|</mo></mrow></mfrac><mo>+</mo><msubsup><mi>&Sigma;</mi><mrow><mi>j</mi><mo>=</mo><mi>j</mi><mn>1</mn></mrow><mrow><mi>j</mi><mi>q</mi></mrow></msubsup><mi>&eta;</mi><mrow><mo>(</mo><msub><mi>v</mi><mi>j</mi></msub><mo>-</mo><msub><mi>v</mi><mi>i</mi></msub><mo>)</mo></mrow></mrow>]]></math><img file="FDA0001010635540000012.GIF" wi="1380" he="111" /></maths><maths num="0003"><math><![CDATA[<mrow><msub><mi>a</mi><mi>i</mi></msub><mo>=</mo><mfrac><mrow><msubsup><mi>F</mi><mrow><mi>e</mi><mi>x</mi><mi>t</mi></mrow><mi>i</mi></msubsup><mo>-</mo><msubsup><mi>F</mi><mrow><mi>i</mi><mi>n</mi></mrow><mi>i</mi></msubsup></mrow><msub><mi>m</mi><mi>i</mi></msub></mfrac></mrow>]]></math><img file="FDA0001010635540000013.GIF" wi="219" he="111" /></maths>其中,m<sub>i</sub>为质点i的质量,x<sub>i</sub>表示质点i的位置矢量,<img file="FDA0001010635540000014.GIF" wi="67" he="68" />表示与质点i相连的所有质点j施加于质点i的内力,与质点i相连的所有质点j包括j1至jq,x<sub>j</sub>表示质点j的位置矢量,<img file="FDA0001010635540000015.GIF" wi="82" he="63" />表示质点i所受外力,|x<sub>j</sub>‑x<sub>i</sub>|表示两个质点间位置矢量之差,|x<sub>j</sub>‑x<sub>i</sub>|<sup>0</sup>是发生形变前弹簧的长度,v<sub>i</sub>为质点i的速度,v<sub>j</sub>为质点j的速度,a<sub>i</sub>为质点i的加速度,μ为弹簧的弹性系数,η为阻尼器的阻尼系数;步骤2):依据步骤1)所建的生物力学模型,对其中的各个质点进行初始化,包括初始化各个质点的位置、质量、速度、加速度和受力信息,构建模型的初始状态,并且计算该模型中的每个质点与相连质点间的弹簧初始长度;步骤3):外接的力反馈设备对软组织施加作用力,进行碰撞检测,确定软组织上发生碰撞的质点和受力发生形变的区域,以及软组织被按压或者拉伸的长度;步骤4):用改进的欧拉算法计算发生形变区域内的各个质点的运动状态信息,其中运动状态信息包括发生形变区域内的各个质点的位置、速度和受力随时间变化的信息;计算过程具体为:步骤4.1):经力反馈设备对软组织施加外作用力、并经碰撞检测,返回发生碰撞的质点的序列号r及所受外力<img file="FDA0001010635540000016.GIF" wi="106" he="61" />并且当时刻k=0时,<img file="FDA0001010635540000017.GIF" wi="246" he="62" />步骤4.2):计算质点r的加速度:<maths num="0004"><math><![CDATA[<mrow><msub><mi>a</mi><mi>r</mi></msub><mrow><mo>(</mo><mi>k</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><msubsup><mi>F</mi><mrow><mi>e</mi><mi>x</mi><mi>t</mi></mrow><mi>r</mi></msubsup><mo>-</mo><msubsup><mi>F</mi><mrow><mi>i</mi><mi>n</mi></mrow><mi>r</mi></msubsup><mrow><mo>(</mo><mi>k</mi><mo>)</mo></mrow></mrow><msub><mi>m</mi><mi>r</mi></msub></mfrac></mrow>]]></math><img file="FDA0001010635540000018.GIF" wi="390" he="100" /></maths>步骤4.3):用改进的欧拉算法求解质点r的速度和位移,其中用显示欧拉法对速度v<sub>r</sub>进行迭代求解,用隐式欧拉法对位置矢量x<sub>r</sub>进行求解:<maths num="0005"><math><![CDATA[<mrow><msubsup><mi>v</mi><mi>r</mi><mrow><mi>k</mi><mo>+</mo><mn>1</mn></mrow></msubsup><mo>=</mo><msubsup><mi>v</mi><mi>r</mi><mi>k</mi></msubsup><mo>+</mo><mi>&Delta;</mi><mi>t</mi><mo>&CenterDot;</mo><msub><mi>a</mi><mi>r</mi></msub><mrow><mo>(</mo><mi>k</mi><mo>)</mo></mrow></mrow>]]></math><img file="FDA0001010635540000019.GIF" wi="523" he="68" /></maths><maths num="0006"><math><![CDATA[<mrow><msubsup><mi>x</mi><mi>r</mi><mrow><mi>k</mi><mo>+</mo><mn>1</mn></mrow></msubsup><mo>=</mo><msubsup><mi>x</mi><mi>r</mi><mi>k</mi></msubsup><mo>+</mo><mfrac><mrow><mi>&Delta;</mi><mi>t</mi></mrow><mn>2</mn></mfrac><mo>&CenterDot;</mo><mrow><mo>(</mo><msubsup><mi>v</mi><mi>r</mi><mi>k</mi></msubsup><mo>+</mo><msubsup><mi>v</mi><mi>r</mi><mrow><mi>k</mi><mo>+</mo><mn>1</mn></mrow></msubsup><mo>)</mo></mrow></mrow>]]></math><img file="FDA0001010635540000021.GIF" wi="614" he="95" /></maths>其中,<img file="FDA0001010635540000022.GIF" wi="56" he="63" />为质点r在时刻k的速度向量,<img file="FDA0001010635540000023.GIF" wi="102" he="61" />为质点r在时刻k+1的速度向量,<img file="FDA0001010635540000024.GIF" wi="54" he="62" />为质点r在时刻k的位置向量,<img file="FDA0001010635540000025.GIF" wi="102" he="61" />为质点r在时刻k+1的位置向量;步骤4.4):计算质点r开始运动后与其相连质点间弹簧的长度及弹簧的形变量,并计算质点r所受内力<img file="FDA0001010635540000026.GIF" wi="83" he="63" />步骤4.5):循环执行步骤4.1)至4.4),计算受力区域内其他质点的位置、速度及受力信息;步骤5):依据步骤3)和步骤4)计算反馈力,并将该反馈力输出至力反馈设备;步骤6):循环执行步骤4)至步骤5),计算每个时间步长中各质点的运动状态,并在显示器上动态显示软组织的形变过程。
地址 450011 河南省郑州市北环路36号