CN101425188A - 带刚体核的广义弹簧振子形变仿真方法 - Google Patents

Abstract

本发明公开了一种带刚体核的广义弹簧振子形变仿真方法，目的是解决对表面网格模型进行形变仿真时存在不能保持初始形状和不能表现整体性形变行为的问题。技术方案是首先构建由质点和无质量的弹簧组成的带刚体核的广义弹簧振子模型，然后在实时形变仿真的视觉显示循环中对带刚体核的广义弹簧振子模型进行形变仿真，包括：采用对刚体核进行更新、获得模型中各质点所受的内力、结合模型所受到的外力得到模型中每个质点受到的合力，并根据运动方程获得模型中质点的运动状态、利用模型中每个质点在当前时刻的质点位置在显示器中绘制模型的当前形状四步。采用本发明解决了传统仿真方法不能够保持初始形状和不能对整体性形变效果进行仿真的问题。

Description

带刚体核的广义弹簧振子形变仿真方法

技术领域

本发明涉及一种计算机仿真领域对物体形变进行仿真的方法，尤其是一种带刚体核的广义弹簧振子形变仿真方法。

背景技术

计算机图形学中基于物理的形变仿真方法从物体的物理属性入手，对物体的形变进行仿真模拟，它主要应用于虚拟现实、电影、动画等领域。许多模型被用于基于物理的形变仿真，例如有限元方法(finite elementmethod)、弹簧振子模型(mass-spring model)以及基于模态分析的精简形变模型等。与其他模型相比，弹簧振子模型实现简单、求解高效，因此它是最早被应用于计算机图形学的实时形变模型之一。1987年，K.Waters在ACM SIGGRAPH会议上发表的“A Muscle Model for AnimationThree-Dimensional Facial Expression”(《一种用于三维面部表情动画的肌肉模型》)论文中采用弹簧振子模型对脸部进行建模。1994年，D.Breen等人在ACM SIGGRAPH会议上发表的“Predicting the Drape of Woven ClothUsing Interacting Particles”(《使用交互粒子预测织物窗帘行为的方法》)论文中首次将弹簧振子模型应用于布料仿真模拟中，并指出弹簧振子模型比其他连续模型能更好地模拟布料的行为。弹簧振子模型的另一个优点是易于实现模型拓扑结构的修改，可以方便地模拟物体的断裂、切割等行为，而不像有限元方法等需要使用一些复杂的剖分方法。2003年，刘卉等人在《软件学报》上发表的《基于改进的弹簧-质点模型的三维服装模拟》论文中也采用了弹簧振子模型对布料进行模拟。2002年，蔡及时等人在《软件学报》上发表的《应用于网上虚拟手术的自适应变形模型》论文中使用弹簧振子模型对人体软组织进行形变模拟。2006年，谌海新等在《***仿真学报》上发表的《质点-弹簧模型驱动的体数据模型快速变形模拟》论文中，使用弹簧振子模型对三维体数据模型进行快速变形模拟。2007年，王彦臻等人在《计算机辅助设计与图形学学报》上发表的《改进的采用表面网格的弹簧振子模型》中对传统的弹簧振子模型进行了改进，修正了线弹性弹簧引起的“超弹性”现象。

基本弹簧振子模型由质点和无质量的结构型(structural)弹簧组成。对于一个弹性物体的弹簧振子模型Ω，其质点集合为V，每个质点的初始位置为

质量为m_i(1≤i≤‖V‖)。结构型弹簧集合为S，每一条结构型弹簧连接相邻的两个质点，用于表示模型的结构信息以及在质点间进行力的传播，连接质点i和j之间的结构型弹簧初始长度为l_ij(1≤i，j≤‖V‖)。使用基本弹簧振子模型对弹性物体进行形变仿真的方法在每一个视觉显示循环中，根据弹性物体受到外力的情况，使用弹簧振子模型和运动方程获得模型中各个质点的位置和运动情况，用于对弹性物体的形变状态进行实时绘制。每一个视觉显示循环中，基本弹簧振子模型对弹性物体进行形变仿真的方法包括以下步骤。

第一步，根据模型中结构型弹簧的形变状态获得模型中质点所受到的内力。

基本弹簧振子模型使用线弹性的结构型弹簧。如果每一个质点的当前位置为x_i(1≤i≤‖V‖)，当前速度为v_i(1≤i≤‖V‖)，根据胡克定律，连接质点i和j之间的结构型弹簧施加给质点i的力为

$f_{\mathrm{ij}}=k_s(\left|x_{\mathrm{ij}}\right|-l_{\mathrm{ij}})\frac{x_{\mathrm{ij}}}{\left|x_{\mathrm{ij}}\right|}+k_d(v_j-v_i)$

其中，f_ij为质点i和质点j之间的结构型弹簧施加给质点i的力，k_s为弹簧的劲度系数，x_ij是两个质点位置向量之差(x_j-x_i)，l_ij是弹簧的原始长度，k_d为弹簧的阻尼系数。

质点i此时所受的模型内力为

$f_i^{\mathrm{int}}=\underset{j\∈N\left(i\right)}{\mathrm{\Σ}}{f}_{\mathrm{ij}}$

其中，N(i)为与质点i之间有结构型弹簧连接的质点编号集合。

第二步，结合模型所受到的外力(重力、操纵力等)得到模型中每个质点受到的合力，并根据运动方程获得模型中质点的运动状态(位置、速度)。

模型中的质点i所受的合力为其所受模型内力与所受外力之和：

$f_i=f_i^{\mathrm{int}}+f_i^{\mathrm{ext}}$

其中，

为质点i所受的外力。

根据拉格朗日运动方程 $\stackrel{\·\·}{x}=M^{-1}f(x,v),$ 结合模型上一时刻(t时刻)的已知运动状态，采用向前-向后欧拉法得到当前时刻(t+Δt时刻)模型中每个质点的运动状态：

$\left\{\begin{array}{c}{v}_i(t+\mathrm{\Δt})=v_i\left(t\right)+\mathrm{\Δt}{f}_i\left(t\right)/m_i\\ x_i(t+\mathrm{\Δt})=x_i\left(t\right)+{\mathrm{\Δtv}}_i(t+\mathrm{\Δt})\end{array}\right.$ (1≤i≤‖V‖)

其中，□(t)为t时刻时的物理量，Δt为时间步长。

第三步，利用模型中每个质点在当前时刻的位置向量x_i(t+Δt)在显示器中绘制模型的当前形状。

在仿真方法实施中，这三个步骤循环执行，实现时间的推移。在每一个视觉显示循环中得到模型当前的形变状态并进行绘制，对弹性物体模型的动态形变进行仿真。

在形变仿真模拟中，直接基于弹性物体的体网格(四面体网格、六面体网格等)建立的体弹簧振子模型能够体现物体的体行为，对弹性物体的形变进行逼真的仿真。然而在实际应用中，交互式三维建模工具一般都采用表面网格模型，因此物体的表面网格模型比体网格模型更容易获得。直接采用基本弹簧振子形变仿真方法对表面网格模型进行形变仿真模拟时，质点对应于网格顶点，结构型弹簧对应于网格边。由于表面网格模型的网格边只存在于模型表面，因此模型缺乏对弹性物体内部属性的描述，在形变仿真模拟时模型将缺乏体行为，模型在外力下的形变行为类似于一个空壳：模型的一边受到的外力难以影响模型的另一边，而且在外力作用下模型不能够保持体积，在外力作用消失后模型也不能够恢复初始的形状。为了解决上述问题，2006年，鲍春波等在《***仿真学报》上发表的《一种用于软组织变形仿真的动态质点弹簧模型》论文中向表面弹簧模型中添加了一种虚拟体弹簧，用于体现软体的体特征。但这种方法只适用于人体组织器官等位置固定的特定弹性物体，不能够对弹性物体的刚体运动、伸缩、切变等整体性形变行为进行仿真模拟。

发明内容

本发明要解决的技术问题在于针对背景技术在对表面网格模型进行形变仿真时存在不能保持初始形状和不能表现整体性形变行为的不足，提出一种带刚体核的广义弹簧振子形变仿真方法。

本发明的技术方案是：第一步，构建带刚体核的广义弹簧振子模型。

为了使得表面网格模型在受到外力作用时能够保持初始形状，本发明定义了表面网格模型对应的刚体核，并根据刚体核向基本弹簧振子模型中增加了“初始位置”广义弹簧。具体的，带刚体核的广义弹簧振子模型由质点和无质量的弹簧组成，质点集合对应于表面网格模型Ω的网格顶点集合V，弹簧包括两种类型：结构型弹簧和“初始位置”广义弹簧，它们都采用线弹性弹簧。结构型弹簧的集合对应于表面网格模型Ω的网格边集合E，每条结构型弹簧连接相邻的两个质点，初始长度为网格边的长度。

表面网格模型对应的刚体核的定义如下：

对于表面网格模型Ω，其顶点集合为V，初始状态下各个顶点的位置为

称满足下列条件的结点集合H为表面网格模型Ω的刚体核：

i.集合H中的结点个数与模型Ω中的顶点个数相同；

ii.模型Ω中的每一个顶点都对应集合H中的一个结点，第i个结点的位置为 $h_i=x_i^0,$ 称h_i为顶点i的初始位置；

iii.集合H中所有结点组成的整体结点位置向量h不会由于受力而改变。

所述的“初始位置”广义弹簧是指连接质点与其初始位置之间的无质量弹簧，“初始位置”广义弹簧的初始长度为0。“初始位置”广义弹簧的引入使得模型中的质点受到回复初始位置的弹性力。因此，带刚体核的广义弹簧振子模型在受到外力的作用时，具有保持初始形状的能力。

第二步，对带刚体核的广义弹簧振子模型进行实时形变仿真，即在实时形变仿真的视觉显示循环中，随着时间的推移，以Δt为时间步长循环执行以下操作，在每一个时间步中得到模型当前的形变状态并进行绘制，对弹性物体模型的动态形变进行仿真。Δt根据当前所使用的计算机软硬件条件和被模拟的物体属性而定，同时需满足 $\mathrm{\&Delta;t}<2\sqrt{\min \left(m_i\right)/k_s}(1\&le;i\&le;\left|\right|V\left|\right|),$ 其中m_i为质点i的质量，k_s为弹簧的劲度系数。具体的，在一个视觉循环时间步中的操作包括以下四个步骤：

步骤1，更新刚体核。

为了体现弹性物体的整体性形变行为，需要对刚体核进行更新。根据模型中质点的当前位置，采用形状匹配方法对表面网格模型对应的刚体核进行更新，以体现模型的整体性形变效果，具体步骤是：

(1)用最小二乘法，获取表面网格模型从初始形状到当前形状的最优线性变换矩阵A。

$A=\left({\mathrm{\&Sigma;}}_i{m}_i{p}_i{q}_i^T\right){\left({\mathrm{\&Sigma;}}_i{m}_i{q}_i{q}_i^T\right)}^{-1}$

其中，q_i＝h_i-h_cm，p_i＝x_i-x_cm。h_cm为模型初始刚体核的质心坐标，x_cm为模型当前形状的质心坐标，x_i为模型当前形状中质点i的位置坐标。

(2)根据A获取表面网格模型初始形状与当前形状间的最优旋转变换矩阵R。

首先对A进行奇异值分解(SVD，singular value decomposition)，得到A的对角化形式：

$A=U\hat{A}{V}^T$

其中，

为对角矩阵，U和V均为正交矩阵。则R为：

R＝UV^T

(3)更新表面网格模型对应的刚体核。

为了保证模型具有恢复初始形状的能力，首先结合最优线性变换矩阵A与最优旋转变换矩阵R，生成刚体核更新变换矩阵T：

T＝αA+βR

其中，α，β∈[0，1)且α+β＝1，用于控制模型的整体刚性，α越小，模型的整体刚性越强。

然后使用刚体核更新变换矩阵T对模型的刚体核进行更新：

g_i＝T(h_i-h_cm)+x_cm(1≤i≤‖V‖)

其中，g_i为质点i在更新后的刚体核中对应的顶点位置坐标，h_i为质点i的初始位置。

步骤2，根据模型中弹簧的形变状态获得带刚体核的广义弹簧振子模型中各质点所受内力。

模型中质点受到的内力包括结构型弹簧施加给质点的力和“初始位置”广义弹簧施加给质点的力。

与质点i相连的结构型弹簧施加给质点i的力为

其中N(i)为与质点i存在结构型弹簧连接的相邻质点编号集合。

“初始位置”广义弹簧施加给质点i的力：

f_ii＝k_sx_ii

其中，x_ii是更新后的刚体核中质点i所对应的结点位置与它的当前位置之差(g_i-x_i)。

带刚体核的广义弹簧振子模型中质点i所受的内力

为：

$f_i^{\mathrm{int}}\&prime;=f_{\mathrm{ii}}+\underset{j\&Element;N\left(i\right)}{\mathrm{\&Sigma;}}{f}_{\mathrm{ij}}$

步骤3，结合模型各质点所受到的外力(重力、操纵力等)得到模型中每个质点受到的合力，并根据运动方程获得模型中质点的运动状态(位

和速度)。

模型中的质点i所受的合力为其所受模型内力与所受外力之和：

$f_i\&prime;=f_i^{\mathrm{int}}\&prime;+f_i^{\mathrm{ext}}$

其中，

为质点i所受的外力。

根据拉格朗日运动方程 $\stackrel{\&CenterDot;\&CenterDot;}{x}=M^{-1}f(x,v),$ 结合模型上一时刻(t时刻)的已知运动状态，采用向前-向后欧拉法得到当前时刻(t+Δt时刻)模型中每个质点的运动状态，即位

x_i(t+Δt)和速度v_i(t+Δt)：

$\left\{\begin{array}{c}{v}_i(t+\mathrm{\&Delta;t})=v_i\left(t\right)+\mathrm{\&Delta;t}{f}_i\&prime;\left(t\right)/m_i\\ x_i(t+\mathrm{\&Delta;t})=x_i\left(t\right)+{\mathrm{\&Delta;tv}}_i(t+\mathrm{\&Delta;t})\end{array}\right.$ (1≤i≤‖V‖)

其中，

为t时刻时的物理量。

步骤4，利用模型中每个质点在当前时刻的质点位置向量x_i(t+Δt)在显示器中绘制模型的当前形状。

采用本发明可获得如下有益效果：

第一，本发明直接采用物体模型的表面网格数据作为输入，与背景技术中采用体网格的方法相比，简化了模型的获取和修改；

第二，本发明向基本弹簧振子模型中引入刚体核和“初始位置”广义弹簧，解决了基本弹簧振子形变仿真方法不能够在外力作用时保持表面网格模型初始形状的问题；

第三，本发明通过使用形状匹配方法对模型的刚体核进行更新，解决了背景技术中添加虚拟体弹簧的方法不能对物体的整体性形变效果进行仿真模拟的问题，使得表面网格模型表现出逼真的体行为。

附图说明

图1为本发明总体流程图；

图2为带刚体核的广义弹簧振子模型局部示意图；

具体实施方式

本实施例目的在于对表面网格模型描述的弹性物体进行实时动态的形变仿真模拟，采用3DS Max软件生成的VRML格式表面网格模型作为输入。弹性物体的表面网格模型表示为Ω＝(V，E)，其中V为网格顶点集合，E为网格边集合。如图1中的流程图所示，首先构建与表面网格模型Ω对应的带刚体核的广义弹簧振子模型，然后根据模型的受力信息对带刚体核的广义弹簧振子模型进行实时动态的形变仿真模拟，具体步骤为：

第一步，构建带刚体核的广义弹簧振子模型，带刚体核的广义弹簧振子模型局部示意如图2所示，其中实心圆点表示模型中的质点，空心圆点表示模型中质点对应的初始位置结点，实线表示结构型弹簧，虚线表示“初始位置”广义弹簧。

第二步，对带刚体核的广义弹簧振子模型进行形变仿真

在实时形变仿真的视觉显示循环中，随着时间的推移，以Δt为时间步长循环执行以下操作，在每一个时间步中得到模型当前的形变状态并进行绘制，对弹性物体模型的动态形变进行仿真。具体的，在一个视觉循环时间步中的操作包括以下四个步骤：

步骤1，更新刚体核，即根据模型中质点的当前位置，采用形状匹配方法对表面网格模型对应的刚体核进行更新。

步骤2，根据模型中弹簧的形变状态获得模型中各质点所受内力。

步骤3，结合质点所受外力，采用向前-向后欧拉法求解拉格朗日运动方程，得到模型中每个质点的运动状态。

步骤4，对模型的形变结果进行绘制，即利用模型中每个质点在当前时刻的整体质点位置向量在显示器中绘制模型的当前形状。

Claims (1)

1.一种带刚体核的广义弹簧振子形变仿真方法，其特征在于包括以下步骤：

第一步，构建带刚体核的广义弹簧振子模型，带刚体核的广义弹簧振子模型由质点和无质量的弹簧组成，质点集合对应于表面网格模型Ω的网格顶点集合V，弹簧包括结构型弹簧和“初始位置”广义弹簧，它们都采用线弹性弹簧；“初始位置”广义弹簧是指连接质点与其初始位置之间的无质量弹簧，“初始位置”广义弹簧的初始长度为0；结构型弹簧的集合对应于表面网格模型Ω的网格边集合E，每条结构型弹簧连接相邻的两个质点，初始长度为网格边的长度；表面网格模型对应的刚体核的定义如下：

对于表面网格模型Ω，其顶点集合为V，初始状态下各个顶点的位置为 $x_i^0(1\&le;i\&le;\left|\right|V\left|\right|),$ 称满足下列条件的结点集合H为表面网格模型Ω的刚体核：

i.集合H中的结点个数与模型Ω中的顶点个数相同；

ii.模型Ω中的每一个顶点都对应集合H中的一个结点，第i个结点的位置为 $h_i=x_i^0,$ 称h_i为顶点i的初始位置；

iii.集合H中所有结点组成的整体结点位置向量h不会由于受力而改变；

第二步，对带刚体核的广义弹簧振子模型进行实时形变仿真，即在实时形变仿真的视觉显示循环中，随着时间的推移，以Δt为时间步长循环执行以下操作，在每一个时间步中得到模型当前的形变状态并进行绘制，对弹性物体模型的动态形变进行仿真；Δt根据当前所使用的计算机软硬件条件和被模拟的物体属性而定，满足 $\mathrm{\&Delta;t}<2\sqrt{\min \left(m_i\right)/k_s}(1\&le;i\&le;\left|\right|V\left|\right|),$ 其中m_i为质点i的质量，k_s为弹簧的劲度系数。

步骤1，采用形状匹配方法对表面网格模型对应的刚体核进行更新，方法是：

(1)采用最小二乘法，获取表面网格模型从初始形状到当前形状的最优线性变换矩阵A，

$A=\left({\mathrm{\&Sigma;}}_i{m}_i{p}_i{q}_i^T\right){\left({\mathrm{\&Sigma;}}_i{m}_i{q}_i{q}_i^T\right)}^{-1}$

q_i＝h_i-h_cm，p_i＝x_i-x_cm，h_cm为模型初始刚体核的质心坐标，x_cm为模型当前形状的质心坐标，x_i为模型当前形状中质点i的位置坐标；

(2)根据A获取表面网格模型初始形状与当前形状间的最优旋转变换矩阵R，方法是：首先对A进行奇异值分解，得到A的对角化形式： $A=U\hat{A}{V}^T,$

为对角矩阵，U和V均为正交矩阵，R＝UV^T；

(3)更新表面网格模型对应的刚体核，方法是：首先结合最优线性变换矩阵A与最优旋转变换矩阵R，生成刚体核更新变换矩阵T：T＝αA+βR，其中，α，β∈[0，1)且α+β＝1，用于控制模型的整体刚性，α越小，模型的整体刚性越强；然后使用T对模型的刚体核进行更新：

g_i＝T(h_i-h_cm)+x_cm(1≤i≤‖V‖)，g_i为质点i在更新后的刚体核中对应的顶点位置坐标，h_i为质点i的初始位置；

步骤2，根据模型中弹簧的形变状态获得带刚体核的广义弹簧振子模型中各质点i所受的内力

$f_i^{\mathrm{int}}\&prime;=f_{\mathrm{ii}}+\underset{j\&Element;N\left(i\right)}{\mathrm{\&Sigma;}}{f}_{\mathrm{ij}},$

为与质点i相连的结构型弹簧施加给质点i的力，N(i)为与质点i存在结构型弹簧连接的相邻质点编号集合，f_ii为“初始位置”广义弹簧施加给质点i的力，f_ii＝k_sx_ii，x_ii是更新后的刚体核中质点i所对应的结点位置与它的当前位置之差(g_i-x_i)；

步骤3，结合模型各质点所受到的外力得到模型中每个质点受到的合力，并根据运动方程获得模型中质点的运动状态：

模型中的质点i所受的合力为其所受模型内力与所受外力之和：

$f_i\&prime;=f_i^{\mathrm{int}}\&prime;+f_i^{\mathrm{ext}},$ 其中，

为质点i所受的外力；

根据拉格朗日运动方程 $\stackrel{\&CenterDot;\&CenterDot;}{x}=M^{-1}f(x,v),$ 结合模型上一时刻即t时刻的已知运动状态，采用向前-向后欧拉法得到当前时刻即t+Δt时刻模型中每个质点的运动状态，即位置x_i(t+Δt)和速度v_i(t+Δt)：

$\left\{\begin{array}{c}{v}_i(t+\mathrm{\&Delta;t})=v_i\left(t\right)+\mathrm{\&Delta;t}{f}_i\&prime;\left(t\right)/m_i\\ x_i(t+\mathrm{\&Delta;t})=x_i\left(t\right)+\mathrm{\&Delta;t}{v}_i(t+\mathrm{\&Delta;t})\end{array}\right.(1\&le;i\&le;\left|\right|V\left|\right|)$

其中，

为t时刻时的物理量；

步骤4，利用模型中每个质点在当前时刻的质点位置向量x_i(t+Δt)在显示器中绘制模型的当前形状。

Images