（计算机科学与技术专业论文）基于物理的弹性变形模拟技术研究.pdf

浙江大学博士学位论文摘要 摘要 模型弹性变形的模拟在计算机图形学中有广泛的应用，与之相关的研究也已 成为图形学研究的一个重要方向。迄今为止，该领域已经发展了非物理方法和基 于物理的模拟两大类技术。近年来，随着计算机硬件的飞速发展和虚拟现实应用 需求的不断提高，基于物理的模拟日益受到重视。一方面，此类方法因遵循物体 的客观运动规律，而能够生成更加逼真合理的形变和动态效果；另一方面，对于 复杂的变形过程，容易对变形对象进行简易合理的控制；同时，由于方法上和工 程力学仿真具有高度一致性，而更加容易融合于工程领域的可视化辅助分析应 用。 本文重点研究了基于物理的弹性变形模拟中的三方面问题：基于骨架驱动的 弹性变形的有限元模拟、无网格模拟和在弹性变形中存在拓扑变化情况的处理。 在应用有限元方法对骨架驱动变形的模拟中，提出一种基于骨架旋转扩散的 技术求取模型的旋转场，并用以补偿模型大变形模拟中因采用线性计算而引起的 变形效果失真。该技术通过扩散方程的构建、旋转特征向量的求取以及模型中单 元旋转的计算等步骤，实现了模型旋转场的非迭代求解方法。与已有技术相比， 该方法在静力学变形中不需要划分子步，在动力学模拟中减少了误差累积，提高 了模拟的稳定性和速度。 无网格方法是近年来受到广泛关注的微分方程求解技术。本文在对骨架驱动 弹性变形模拟中进行了无网格模拟方面的探索，采用了伽辽金型移动最小二乘近 似，在模型的空间坐标域内进行了模拟。同时，提出了一种无网格框架内更为直 接的基于骨架旋转扩散的旋转场求取技术，并将有限元框架内发展起来s t i f f n e s s w a r p i n g 技术成功地运用到无网格框架内，解决了无网格框架内对大变形线性模 拟所导致的失真问题。 模型在变形过程中的几何拓扑变化，如遭受切割和断裂等是计算机图形学中 模拟的难点。本文研究扩展有限元方法，通过改进标准有限元的单元逼近函数， 使之能够表达单元内的几何不连续，从而实现切割和断裂的模拟。同时，针对计 算机图形学领域的应用，对其中多裂纹和相交裂纹处理等有关问题做了深入的研 究，提出了改进方法或解决策略。本文的方法在实现过程中，不需要对模型有限 浙江大学博士学位论文 摘要 元网格进行重新剖分，避免了由此带来的计算量大、稳定性差等问题；并且能够 较为方便地与在标准有限元模拟发展起来的技术，如非线性失真的补偿技术等有 效地融合。 关键词：基于物理的方法；弹性变形模拟；非线性；有限元；无网格；几 何拓扑改变 浙江大学博士学位论文a b s t r a c t a b s t r a c t s i m u l a t i o no fe l a s t i cd e f o r m a t i o nh a se x t e n s i v ea p p l i c a t i o na n dr e s e a r c ho nl th a sb e e n a ni m p o r t a n ta r e ai nc o m p u t e rg r a p h i c s t h e r ea r et w oc a t e g o r i e so ft e c h n i q u e sf o r t h i sp r o b l e m ：n o n p h y s i c a lm e t h o d sa n dp h y s i c a l l yb a s e ds i m u l a t i o n w i t ht h e o v e r g r o w t ho ft h ec o m p u t e rh a r d w a r ea n dt h ei m p r o v e m e n to ft h ea p p l i c a t i o n ， i n c r e a s i n g l yi m p o r t a n c eh a sb e e na t t a c h e dt op h y s i c a l l yb a s e ds i m u l a t i o n i nf a c t ， p h y s i c a l l yb a s e ds i m u l a t i o nh a s a tl e a s tt h r e eo u t s t a n d i n ga d v a n t a g e s f i r s t l y , p l a u s i b l e r e s u l tc a nb e o b t a i n e db e c a u s ep h y s i c a l l yb a s e dm e t h o di si n c o m p l i a n c e w i t h o b j e c t i v em o v e m e n tl a w s s e c o n d l y , t h ec o m p l e xd e f o r m a t i o nc a nb es i m p l ya n d r a t i o n a l l yc o n t r o l l e d f i n a l l y , t h es i m u l a t i o nr e s u l t sa r em o r ee a s i l ya p p l i e db yt h e v i s u a l a i da n a l y s i so ft h ee n g i n e e r i n ga r e a t h r e ep r o b l e m sa r es t u d i e di nt h i st h e s i s ：f i n i t ee l e m e n ts i m u l a t i o no fe l a s t i c d e f o r m a t i o nd r i v e nb ys k e l e t o n ，m e s h l e s ss i m u l a t i o no fe l a s t i cd e f o r m a t i o nd r i v e nb y s k e l e t o na n dt h et r e a t m e n to fc h a n g i n go f g e o m e t r yt o p o l o g yi ne l a s t i cd e f o r m a t i o n i nt h ef i n i t ee l e m e n ts i m u l a t i o no fe l a s t i cd e f o r m a t i o nd r i v e nb ys k e l e t o n ，an o v e l m e t h o di s p r o p o s e dw h i c he s t i m a t e s f i n i t ee l e m e n t sr o t a t i o nb yd i f f u s i n gf r o m s a m p l i n gp o i n t si na t t a c h e d s k e l e t o ns t r u c t u r ea n du s i n ga ne m p i r i c a lf o r m u l ao f r o t a t i o nq u a t e m i o n i nc a l c u l a t i n gd e f o r m a t i o nb a s e do nt h es t a t i ce q u a t i o n ，o u r m e t h o dc a ng r e a t l ya c c e l e r a t ec o m p u t i n go fd e f o r m a t i o nw i t h o u tf u s s ys u b s t e p c a l c u l a t i o n ，w h i l ei tc a nm a k et h ew h o l es o l u t i o nm o r es t a b l ei nk i n e t i c sc a l c u l a t i o n m o r e o v e r , o u rm e t h o dc a nr e d u c et h eu s e ri n t e r a c t i o na n do b t a i nr a t i o n a lr e s u l t s am e s h l e s ss i m u l a t i o ns y s t e mi s p r e s e n t e df o re l a s t i cd e f o r m a t i o nd r i v e nb y s k e l e t o ni nc h a p t e r5 i nt h i ss y s t e m ，g a l e r k i nm e s h l e s sm e t h o db a s e dm o v i n gl e a s t s q u a r ea p p r o x i m a t i o ni sa p p l i e di ns p a c ec o o r d i n a t e d o m a i n i na d d i t i o n ，an e w m e t h o df o rc a l c u l a t i n gn o d er o t a t i o ni sp r o p o s e dw h i l ea p p l y i n gas i m i l a rt e c h n i q u e w i t hs t i f f n e s sw a r p i n gt ot a c k l et h en o n l i n e a rl a r g ed e f o r m a t i o n i no u rm e t h o d ，a l l n o d er o t a t i o n sa r ee v a l u a t e df r o ms a m p l i n gp o i n t si na t t a c h e ds k e l e t o nb yc o n s t r u c t i n g a n ds o l v i n gt h ed i f f u s i o np a r t i a ld i f f e r e n t i a le q u a t i o n t h ee x p e r i m e n t si n d i c a t e dt h a t t h em e t h o dc a ne n h a n c et h es t a b i l i t yo ft h ed y n a m i c sa n da v o i df u s s ys u b - s t e p c a l c u l a t i o ni ns t a t i cd e f o r m a t i o ne d i t i o n a sf o rt h ep r o b l e mo ft h ec h a n g i n go fg e o m e t r yt o p o l o g y , t h a ti sc u t t i n ga n d 浙江大学博主兰位笙奎 竺! 皇竺! _ _ _ - _ _ - - _ _ _ - - i - _ - _ - _ _ - - _ _ _ - - - _ _ i 一一一 f r a c t u r i n go fm o d e l ，af i n i t ee l e m e n ts i m u l a t i o nf r a m e w o r k w i t h o u tr e m e s h i n gi s p r e s e n t e d t h e k e r n e li sa d d i n gad i s c o n t i n u o u s f u n c t i o nf o r t h es t a n d a r d a p p r o x i m a t i o nt oa c c o u n tf o rt h ec r a c k af e a s i b l et e c h n i q u ei sa d o p t e df o rd e a l i n g w i t hm u l t i p l ec r a c k sa n di n t e r s e c t i n gc r a c k s s e v e r a li n v o l v e dp r o b l e m si n c l u d i n g e x t e n d e df r e e d o m so ff i n i t ee l e m e n tn o d e sa n dm a s sm a t r i xc a l c u l a t i o n e t ca r e d i s c u s s e d t h et e c h n i q u ei se a s yt os i m u l a t ec h a n g i n gg e o m e t r i ct o p o l o g y m o r e o v e r , p r e v i o u sm e t h o d sd e v e l o p e di ns t a n d a r df i n i t ee l e m e n tf r a m e w o r k ，s u c ha ss t i f f n e s s w a r p i n gm e t h o d ，c a n b ee x t e n d e da n du t i l i z e d k e y w o r d s ： p h y s i c a l l yb a s e ds i m u l a t i o n ；e l a s t i cd e f o r m a t i o n ；f i n i t ee l e m e n t ， m e s h l e s s ；g e o m e t r yt o p o l o g y 浙江大学博上学位论文 图表目录 图表目录 图2 1 一个质点弹簧系统一9 图2 2 基于质点弹簧系统的人的面部表情模拟1 0 图2 3b u n n y 模型被离散化为四面体单元构成的有限元网格1 1 图2 4 用有限元方法模拟的球撞击模型变形1 4 图2 5 用粒子系统实现的表面模型变形和拓扑重构1 5 图2 6 无网格方法实现的变形模拟序列1 7 图2 7 边界元方法对物体变形的模拟1 9 图2 8 基于模态分析技术的变形模拟2 3 图2 9 线性模拟大变形的失真2 5 图2 1 0 骨架驱动的弹性变形2 6 图2 1 1 肘关节处的变形失真2 7 图2 1 2 碗的破碎模拟2 8 图3 1 应用不同方法得到的关节区变形比较一3 1 图3 2 平板变形模拟的失真31 图3 3 光线求交法交互生成骨架3 4 图3 4 骨架的局部标架3 4 图3 54 节点四面体单元3 6 图3 6 基于扩散方程的旋转场一4 1 图3 7 旋转特征向量4 2 图3 8 关节区域变形结果对比4 4 图3 9 杆的扭转变形结果4 5 图3 。1 0c a c t u s 模型的动力学变形4 5 图3 1 1 马模型的动力学变形4 6 图4 1 模型的点离散5 0 图4 2 二维移动最小二乘近似5 2 图4 3 无网格方法计算的大变形失真与补偿5 6 图4 4 骨架旋转扩散5 8 图4 5 骨架关节区的无网格模拟一6 1 图4 6d i n o 模型的动力学变形模拟6 1 图4 7a m o d i l i o 模型的动力变形6 2 图5 1 虚节点方法的基本原理一6 5 图5 2 断裂单元节点的附加自由度6 8 图5 3 断裂单元上的h 函数6 9 图5 4 运动的独立性7 0 图5 5 多条裂纹情况7 4 i i i 浙江人学博上学位论文 图表目录 图5 6 单元内两条裂纹情况7 5 图5 7 模型的切割7 6 图5 8 具有切割面的b u n n y 模型的运动序列7 7 图5 9 平板受拉断裂8 0 图5 1 0 板下落过程中遭受撞击断裂8 0 图5 1 1 头部遭受撞击 f l y m a n 模型81 表3 1 本文方法与m i i l l e r 方法所得旋转场数值比较4 4 表5 1 实例物理参数及运行统计8 1 i v 浙江大学研究生学位论文独创性声明 本人声明所呈交的学位论文是本人在导师指导下进行的研究工作及取得的 研究成果。除了文中特别加以标注和致谢的地方外，论文中不包含其他人已经发 表或撰写过的研究成果，也不包含为获得逝垄盘堂或其他教育机构的学位或 证书而使用过的材料。与我一同工作的同志对本研究所做的任何贡献均己在论文 中作了明确的说明并表示谢意。 学位论文作者签名：旁、矗之签字日期：j 帕矛年月7 日 学位论文版权使用授权书 本学位论文作者完全了解逝婆盘堂有权保留并向国家有关部门或机 构送交本论文的复印件和磁盘，允许论文被查阅和借阅。本人授权逝鎏盘堂 可以将学位论文的全部或部分内容编入有关数据库进行检索和传播，可以采用影 印、缩印或扫描等复制手段保存、汇编学位论文。 ( 保密的学位论文在解密后适用本授权书) 学位论文作者签名：每墨乙 签字日期：加孑年月7 日 导师签名： 签字日期：1 d 目 浙江大学博上学位论文第1 章绪论 第1 章绪论 在计算机图形学中，基于物理的弹性变形模拟因其广泛应用而成为研究热点 之一，本文的研究内容围绕该领域的关键问题展开。作为开始，本章对研究的应 用背景、问题的描述和全文的内容及组织结构做简单介绍。 1 1 弹性变形模拟在计算机图形学中的应用 计算机图形学中，弹性变形模拟综合运用了计算机科学、数学、物理学和其 他相关学科知识，在计算机上快速地生成具有真实感的模型弹性变形效果。它在 数字娱乐和虚拟现实领域具有广泛的应用，概括地说，主要在于以下三个方面： ( 1 ) 影视产品。9 0 年代末期开始出现了三维动画片和利用计算机动画生成 特技画面的大制作电影，如怪物史莱克、冰河世纪、侏罗纪公园和终 结者i i 等。在这些影视产品中，有相当部分的角色动作、特技场景如汽车的剧 烈碰撞等都是用计算机模拟生成的。此外，在这些影视作品中的部分自然景物运 动，如动画片中树木随风运动等的展现都应用了弹性变形模拟技术。应用包括弹 性变形模拟在内的当代虚拟特效技术可以实现真实情况下人体或其他角色的高 难度或夸张动作；能够视觉真实地再现灾难性场面和降低制作成本。 ( 2 ) 数字娱乐。在大量网络游戏中，模型的弹性变形场景都广泛存在。如 各类打斗类游戏中，角色的肌肉变形以及武器与周围物体的碰撞等都需要游戏开 发者对相互接触的物体的变形进行模拟。因此，变形和运动的真实程度是影响此 类游戏产品质量的重要因素之一，其中所采用的弹性物体变形模拟技术往往是至 关重要的。 ( 3 ) 计算机辅助工程及辅助医学诊断。弹性变形模拟在工程中的应用历史 悠久，十分广泛，变形程度已经成为各个工程行业设计的一个重要衡量标准。尽 管本文所涉及的弹性变形模拟与严格的工程模拟有所区别，但是图形学中的弹性 变形模拟仍可在一定范围内为工程设计，如虚拟装配等提供参考。此外，在医学 辅助诊断和治疗中，弹性变形的模拟也发挥着重要作用，如近年来研究广泛的虚 拟手术，一方面可以用来作为医师虚拟训练系统，另一方面，在辅助诊断和治疗 的本身也发挥着越来越重要的作用。 浙江大学博士学位论文第l 章绪论 除上述三个主要方面外，弹性变形模拟还在计算机图形学中有很多其他方面 的应用，在此不一一列举。因此，弹性变形模拟在计算机图形学研究中越来越受 到关注，已经成为一个重要的研究领域。 1 2 弹性变形模拟的问题描述 计算机图形学中弹性变形模拟问题可数学描述为：一般地，对任意变形模型 区域表示成qcr 3 ，对任意一点p eq ，变形前后p ：f 2 x r - - r 3 ：( x ，0 ) t - - - x ( x ，t ) ， 其中x 为点p 的初始位置坐标( 材料坐标) ，x ( x ，t ) = x + u ( x ，t ) 为点p 在t 时刻的 位置坐标，u ( x ，t ) 是相对于初始时刻的位移，设在r t 。面上存在位移约束，r ，面上 存在力约束，并且位移约束和力约束随着时间变化可以发生变化。变形模拟的基 本目标即是求解q 内所有点在任意时刻的位置坐标。 同时，图形学中变形模拟一般还需满足如下要求： 1 较强的真实感。计算机图形学的真实感包括形象的真实感和运动的真实 感，简单的说，形象的真实感就是要求模拟的对象( 如场景、角色动画的角色等) 具有视觉的真实性；运动的真实感指模拟对象的变形和运动具有运动的合理性。 前者对建模技术，真实感绘制技术方面较高的要求，而后者则需要符合物理学上 的物体运动规律，因为这些规律描述了物体在真实世界中如何运动和在运动过程 中相互的影响。运动的真实感是本文重点关注和研究的内容。 2 足够快的模拟速度。图形学中的很多应用要求变形模拟速度足够快，如电 脑游戏中需要实时的模拟速度等，这是变形模拟的显著特征之一。 3 良好的稳定性。因为一般需要在场景中交互的应用而要求模拟系统要有良 好的稳定性。 根据模拟过程中的不同方面，弹性变形模拟可以分为几类。如根据模拟计算 是否为线性求解可以分为线性模拟和非线性模拟；根据模型本身是否具有骨架或 者变形运动是否由骨架驱动来完成可以分为骨架驱动的弹性变形和非骨架驱动 的弹性变形；根据在变形过程中是否存在结构拓扑的改变分为弹性变形模拟和切 割、断裂模拟；根据是否根据物理学定律来进行模拟分为非物理的弹性变形模拟 和基于物理的模拟等。 1 3 弹性变形模拟的方法 迄今为止，计算机图形学对弹性变形模拟已经发展了两大类方法：非物理的 2 浙江大学博七学位论文 第1 章绪论 方法和基于物理的方法，两类方法在都有广泛的应用和深入的研究。 1 3 1 非物理方法 非物理方法主要通过几何操作方法等非物理学技术对模型的变形加以模拟 或编辑。本文的重点在基于物理的方法，但不可忽视的是，两大类方法不是截然 分开的，某些非物理的方法在基于物理的模拟框架也得到了应用，如自由曲面变 形技术等。以下简要介绍几大类非物理方法。 ( 1 ) 曲线曲面技术 早期对变形模拟的研究主要在计算机辅助几何设计( c a g d ) 领域，应用曲 线曲面技术如b e z i e r 曲线、b 样条、有理b 样条、非均匀有理b 样条( n u r b s ) 等技术对模型( 曲线、曲面) 的变形进行控制。这些技术主要使模型逼近相应的 曲线曲面，通过控制点的运动对模型变形加以控制和模拟，具体方法在文献【l 】中 有详细讨论。曲线曲面技术基于模型的参数化表达，对于变形的控制和运动交互 较为方便。但是对变形的细节控制需要很大工作量，甚至即使很小的变形有时也 需要变动很多的控制点，需要有足够经验的用户才能够完成。 ( 2 ) 自由曲面变形技术 自由曲面变形技术( f r e e f o r md e f o r m a t i o n ，f f d ) 最早由a b a r r 等在19 8 4 年 【2 】提出。此技术把变形模型映射到包含模型的另一更容易模拟的模型或空间中去， 通过模拟该映射模型或空间的变形来实现对原模型变形的模拟。与曲线曲面通过 控制顶点来控制模型变形相比，f f d 技术提供了更方便、有力的体控制或空间控 制方法。因此，近年来得到了很大的发展和广泛的应用 3 - 6 。 ( 3 ) 几何微分域技术 基于几何微分域的变形模拟技术 7 - 9 】，通过使用模型网格的微分坐标构造网格 几何能量等范函，求解该范函问题，进而得到模型网格变形后的构形，这类方法 能够很好地保持模型变形中的细节，因而近年来得到了蓬勃的发展【1 0 - 1 2 】。 此外，还有许多非物理方法，如蒙皮技术【13 1 ，多分辨率编辑技术【1 4 1 等，都针 对骨架驱动变形、复杂模型的变形等方面做了研究。 1 3 2 基于物理的方法 基于物理的方法主要是根据物理学、力学中的基本规律和计算方法，对物体 浙江大学博士学位论文第1 章绪论 弹性变形及其运动进行模拟。具体而言，根据场景特点建立弹性变形问题的力学 模型，选择包括静力平衡方程、牛顿第二定律、能量守恒原理或弹性力学变分原 理等刻画变形和运动基本规律的定理、原理，为模拟系统建立求解方程并确定边 界条件，然后运用方程的数值求解方法如有限元方法、无网格方法以及时间积分 技术求解，进而确定变形模型中任意点在模拟时间区间内任一时刻的位置坐标。 此类方法能较好保持运动真实感，对非物理方法较难解决的问题，如变形过程中 模型几何细节的保持等不需要专门考虑。 近年来，基于物理的方法受到研究者的广泛关注。一方面，相对于非物理方 法，基于物理的方法因其遵循物体的客观运动规律，而能够生成更加逼真合理的 形变和动态效果；另一方面，对于复杂的变形和运动过程，此类方法往往有直观 明晰的参数设定，更容易对变形对象进行简易合理的控制；同时，由于方法上和 工程仿真具有一致性，从而更加容易融合于面向工程领域的可视化辅助分析应 用，如虚拟医学辅助诊断等。 必须指出的是，尽管图形学中模拟弹性变形问题与物理学和力学中应用的基 本运动规律是一致的，但是由于与后者面向的应有领域不同，所以存在一些显著 区别： 第一，变形研究范围的区别。计算机图形学领域的变形问题，大部分是大变 形问题，而物理学和力学领域的问题很多是小变形问题。 第二，模拟的精度区别。图形学中很多要求视觉的真实感，而物理学和力学 许多是面向工程应用领域，有着更为严格的精度要求。 第三，模拟的速度需要区别。计算机图形中很多是用以满足交互模拟的需要， 对模拟速度有较严格的要求，所以需要结合很多加速技术和方法。事实上，图形 学中追求的是变形模拟速度与真实感的平衡。 1 4 本文的内容及结构 本文在前人工作的基础上，对基于物理弹性变形模拟的两方面问题做深入探 讨，改进或提出相应的解决方法。具体地说，针对应用广泛的骨架驱动的弹性变 形的非线性模拟和模型的切割、断裂的模拟展开研究。 在本文的第二章，对基于物理的弹性变形模拟的主要技术进行综述，侧重于 讨论与本文研究内容相关方法，同时兼顾该方向研究的其他重要方面。 4 浙江大学博士学位论文第l 章绪论 第三章，将在模型有限元离散框架内，对在角色动画等应用广泛的骨架驱动 的弹性变形模拟做相关研究。为能实现快速稳定的变形，此章针对已有对大变形 模拟进行旋转刚度补偿所使用迭代方法的不足，提出了一种通过构造和求解扩散 微分方程，从骨架采样点旋转得到单元的旋转特征向量，进而给出了旋转四元数 经验预测公式，直接构造计算域内单元旋转场的新方法。该方法是非迭代的，实 现了在静力学模拟中，避免了繁琐的变形子步划分，提高了计算速度；在动力学 模拟中，允许使用更大的求解时间步长，增强了变形系统的求解稳定性。同时， 对于骨架关节点区域变形模拟，较大程度降低了用户干预，取得了合理的效果。 第四章，同样针对骨架驱动弹性变形，采用无网格伽辽金方法对弹性变形进 行模拟。尽管伽辽金型无网格方法在模拟速度方面与有限元方法相当，但是具有 独特的优点：( 1 ) 与配点型无网格方法相比更加稳定；( 2 ) 更容易扩展；( 3 ) 计 算量较大的逼近函数可以通过预计算完成；( 4 ) 更适合于点采样模型的模拟。因 此基于此类型无网格计算框架的弹性模拟有较高的研究价值。在本章，还实现了 无网格框架内的大变形线性模拟失真进行旋转补偿，同时也提出了更为直接的的 模型旋转场的扩散求取方法。 在第五章，重点研究模型的切割和断裂模拟。与此前研究者应用标准有限元 方法、无网格方法等进行断裂模拟不同，本文首次在图形学中应用了扩展有限元 法。采用改进有限单元近似函数使之能够表达单元内不连续，并对相关重要问题 做了深入的讨论与研究。此章方法能够避免标准有限元法模拟此类问题所需要的 网格重新剖分及其带来的物理参量的重新映射等，保证了模拟系统具有良好的稳 定性。同时，许多标准有限元方法发展起来的技术，如非线性变形补偿等，能与 之成功结合，应用于切割和断裂的模拟。 最后，将对本文的研究工作做总结，并提出未来的研究工作设想。 浙江大学博士学位论文第2 章基于物理的弹性变形模拟技术综述 第2 章基于物理的弹性变形模拟技术综述 本章介绍基于物理的弹性变形模拟主要技术，将从模拟所采用的控制方程、 模拟的主要数值方法、时间积分和加速技术等方面阐述该方向的主要研究进展。 同时，重点讨论本文研究所涉及的难点和相关工作，如骨架驱动弹性变形、切割 和断裂等。 2 1 控制方程 大多数文献中使用的变形控制方程可归结为如下三类： 1 牛顿第二定律。 m i + c 童+ k x = 匕 ( 2 1 ) 其中x ，童和置分别表示质点的空间坐标，速度和加速度，m i 为表征模型质点系 的惯性力，c 支是运动中所受的阻尼力，k x 是质点系之间相互作用的内力，匕是 质点系所受外力。 2 能量原理。如最小势能原理等，将变形问题转化为如下的能量范函问题 兀= t + u 一缈 ( 2 2 ) 其中丁为模型所具有的动能，u 为模拟系统内能，通常指应变能、热能等，形为 外力所做功。由此方程能够推导出最小势能原理、虚功原理、弹性力学的各种变 分原理。 3 欧拉拉格朗日方程。在连续介质力学中，求解( 2 2 ) 的范函问题最终可以 通过求解该范函的欧拉拉格朗日方程 要，要 + 罢蛾： ( 2 3 ) 磊l 面j + 瓦+ c p 鲥 【2 3 ) 一般地，对求解区域应用数值离散后，( 2 3 ) 式可以写成 m i + c i + 婴：f e z , ( 2 4 ) 在方程( 2 3 ) 、( 2 4 ) 中，m 、c 分别为质量矩阵和阻尼矩阵，在具体数值方法中求 取方法可能有所不同。 可以把方程( 2 3 ) 和( 2 4 ) 看做是公式( 2 2 ) 的导出结果，由于图形学中大量文献 直接应用( 2 3 ) 和( 2 4 ) ，所以这里也单独列出。事实上，方程( 2 1 ) ( 2 4 ) 在本质上是 6 浙江大学博士学位论文第2 章基于物理的弹性变形模拟技术综述 统一的，只是基于不同的力学原理推导出来的结果，例如方程( 2 4 ) 的孚为内能 o x 对空间坐标的偏导，该项的物理意义可以理解为模型内力，和方程( 2 1 ) 的k x 的 物是一致的。在以下的讨论中，本文将不再区分物理规律的来源，而直接应用上 述方程。 在实际应用中，往往把位移作为基本求解未知量，在控制方程中，将质点空 间坐标x 用位移i l l 替换，将不会改变公式意义。即 m i + c 宝+ k ( x x ) = 匕或m i i + c u + k u = f e x t ( 2 5 ) 和 m i i + c f i + 罢- u = 匕 ( 2 6 ) 在对模型弹性变形的编辑等应用中，由于变形速度没有改变，或者改变不是 很大，可以认为加速度为零，从而采用静力学求解。此时方程( 2 3 ) 和( 2 4 ) 可以写 成静力学形式 k ( x - x ) = 匕或k u = ( 2 7 ) 和 i 6 u ：匕 ( 2 8 ) 次 翻 、7 在求解模型内力时，通常需要求解应力、应变以及应变能。在力学中，应力 用来度量受力程度，应变用来度量连续介质变形程度。根据不同的度量标准和不 同的应用需要，有很多种应力与应变。以下将对计算机图形学中常用的应力和应 变做简要介绍。 1 柯西应力 柯西应力6 的意义是：若作用于面积鲋上的力为p ，削的法向为n ，则 n 6 a a = p 。在文献中”5 ，1 6 1 出现第- - p i o l a k i r c h h o f f e 力和第- - p i o l a k i r c h h o f f 应 力，均为对特殊材料的应用，此处不做介绍。 2 格林应变和柯西应变 应变在力学中是描述物体变形程度的物理量。它与变形梯度密切相关，变形 梯度是表征变形特征的一个重要参量，定义为f ：o x ，其中x 为质点的空间坐标， d a 浙江大学博士学位论文第2 章基于物理的弹性变形模拟技术综述 x 为材料坐标，数学上变形梯度f 是模型运动p ：q r - - - r 3 ：( x ，0 ) i - - x ( x ，t ) 的 j a c o b i a n 矩阵。 图形学中常用的应变有格林应变和柯西应变两种。格林应变的定义为 g = 去( f r f i ) ( 2 9 ) 格林应变是一个二阶张量，其分量可以表示成白= ( 瑶一4 ) 2 ，指标满足爱因 斯坦求和约定，磊为k r o n e c k e r 算子，即当江时磊= 1 否则磊= 0 。格林应变的 常用位移表示形式 g = 去( ( v u ) 7 + v u + v u ( v u ) r ) ( 2 1 0 ) 容易证明，在刚体运动中格林应变为零。如果略去二阶导数项，在形式上格 林应变即变为柯西应变 c = 去( ( v u ) 7 + v u ) ( 2 1 1 ) 当面向的应用为小变形问题，应变的二阶导数项对问题的影响较小，可以采 用柯西应变。在图形学中应用中，基本上都是大的变形，此时若采用小变形的柯 西应变，因在刚体运动中，模型的应变可能不为零，导致误差过大而模拟效果视 觉失真。此问题本文将在后面章节展开讨论。 在三维问题中，可以用向量表示应力和应变为6 = 吒，吒) r 和= 岛，岛，乞，岛，乞) 2 。同时用矩阵表示应变与位移u 的关系 = b u 或= b ( x x )( 2 1 2 ) 其中，b 称为应变矩阵，显然采用不同种类的应变时b 是不同的。 3 应变能 本文研究所采用的材料均为线弹性材料，即应力应变满足用矩阵表示的方程 式 o = d f 2 1 3 ) 其中d 材料属性矩阵。 模型的应变能( 也称弹性势能) 可以由下式计算 u = 土e r 6 d f l = 吉r d 耐q(214)9 垃 ，) 比 、 浙江大学博十学位论文第2 章基于物理的弹性变形模拟技术综述 2 2 弹性变形模拟的主要数值方法 基于物理的弹性变形模拟的关键是根据问题的实际，建立和求解控制方程。 动力学方程( 2 1 ) 和( 2 3 ) 一2 6 ) 均为对时间的常微分方程，静力学方程( 2 7 ) 和( 2 8 ) 为 普通的代数方程。对于计算机图形学应用而言，控制方程中常微分方程一般是没 有解析解的，需要进行数值求解。根据不同模型特点和具体模拟要求，建立数值 计算过程的主要方法有质点弹簧系统、有限元方法、无网格方法、有限差分方法 以及边界元方法等。 2 2 1 质点一弹簧系统 稚m 翻口l 零葫g 图2 1 一个质点弹簧系统( 摘自文献 1 7 】) 2 2 1 1 质点弹簧系统模型 质点弹簧系统把整个模拟对象模型离散成有限个质点，质点间按一定规则采 用弹簧连接( 如图2 1 ) 。模型的运动状态由所有质点的空间坐标x ，、速度文，和加 速度量，确定，此处i = 1 ，z 。每个质点所受力f f 由其连接的所有弹簧计算，对任一 质点f ，由牛顿第二运动定律有 鳘，= + ( 2 1 5 ) j 其中他是质点的质量，为作用在质点上除弹簧力外的外力，如重力等；易是 与质点z 之间存在弹簧连接的质点，对i 施加的弹簧作用力，对于线弹性弹簧 f ：【，毡( | x 小勺) 南咆x ( 2 1 6 ) 其中x ，- - - - x ，一x ，赶是弹簧的刚度，f l f 是弹簧的初始长度，屯为弹簧速度阻尼系 数。对于整个模型，刀个质点最终组成一个系统，即形式如方程( 2 1 ) 的方程组， 浙江大学博上学位论文第2 章基于物理的弹性变形模拟技术综述 其中系数矩阵m 、c 和k ，均为3 n x 3 n 的矩阵，分别称为质量矩阵、阻尼矩阵和 刚度矩阵，通过求解该方程组得到不同时刻所有质点的位置，最终确定模型该时 刻的构形。 事实上，质点弹簧模型在图形学中得到了飞速发展【17 1 ，已经超越了传统意 义上的图2 1 的结构模式，由质点和质点间相互作用力为系统主要力学特征 1 8 - 2 0 的大多数粒子系统模型可以包含在质点一弹簧模型中。另外，在具体的应用中弹簧 的类型也可以扩展为铰弹簧、面弹簧、体弹簧等。 图2 2 基于质点- 弹簧系统的人的面部表情模拟( 摘自文献 2 4 ) 2 2 1 - 2 质点一弹簧系统在图形学中的应用 质点弹簧系统最早被应用静力学方程( 2 7 ) 或( 2 8 ) 模拟人脸表情【2 1 ，2 2 1 ，随后应 用动力学方程模拟皮肤、肌肉等工作相继出现 2 3 - 2 5 】，如图2 2 即为用质点弹簧系 统模拟人的面部表情实例。m i l l e r 等质点弹簧系统模拟蛇、蚕和鱼等模型的运动 2 6 ，2 7 1 ，c h r i s t e n s e n 等结合f f d 技术实现复杂模型的弹性变形模拟2 8 1 等。另外， 质点弹簧系统也用来模拟在变形过程中出现断裂的情况【2 9 ，3 0 1 。近年来，质点弹 簧系统被广泛用来模拟布料的运动 1 8 ，1 9 , 3 1 , 3 2 1 ，这里不再详述。 质点弹簧模型的最大优点是简单、直观、求解计算量相对于有限元等方法较 小，所以模拟速度较快。但是也有一些不足，主要有：( 1 ) 对于模拟一般的弹性 模型，由于不是建立在连续介质相关理论上的，所以模拟系统的稳定性较差，模 拟误差较大；( 2 ) 对质点弹簧系统的拓扑依赖较重；( 3 ) 对于约束的处理不是很 自然。( 4 ) 不容易对模拟模型施加材料属性。 浙江大学博士学位论文第2 章基于物理的弹性变形模拟技术综述 2 2 2 有限元方法 2 2 2 1 有限元方法的基本思想 有限元法方法( f e m ) 是求解微分方程的最常用数值方法之一。该方法分别 在应用数学别3 3 1 和工程技术界刚被独立提出，1 9 6 3 年，m e l o s h 3 5 1 使数学界和工 程界对同一问题的不同描述得以统一。有限元方法的基本思想如下【3 6 j ： ( 1 ) 从方法论角度，有限元方法是一种分析综合法。它先将模型整体结构分 解为单元，再将单元合成为整体结构。 分解 ( 2 ) 从数学角度，有限元方法是连续问题离散化的近似解法，将无限自由度 的问题转化为有限自由度问题；将微分方程问题转化为代数方程问题。 l 微分方程问题离散化离散问题 l 无限自由度问题有限自由度问题 有限元方法的求解变形问题基本过程为：首先将模型的区域离散化，划分为 有限个基本单元的集合体，单元与单元的连接点称为节点；然后选择位移近似函 数，将单元内任意点的位移表示为以单元节点位移为自变量的函数；对每个单元 建立方程；然后将所有单元集合起来，形成整体方程，得到一组以节点位移为未 知量的线性方程组；最后求解线性方程组，得到各节点的位移值。 ( a ) 原始模型( b ) 有限元网格 图2 3b u n n y 模型被离散化为四面体单元构成的有限元网格 浙江大学博上学位论文第2 章基于物理的弹性变形模拟技术综述 具体地，有限元方法的步骤如下： ( 1 ) 将模拟对象模型离散化。即选择某种类型的基本单元，将区域划分为 由该单元组成的集合体。对于二维区域，常见的单元类型有三角形单元和矩形单 元；对于三维区域，常见的单元类型有四面体单元和六面体单元，如图2 3 即为 b t m w 模型的四面体单元离散实例。 ( 2 ) 选择位移近似函数。表示单元内任意点以该点坐标为自变量的位移函 数。单元内任意点的位移u 可用单元节点位移i l l 8 表示： i l l = n u 9 ( 2 1 7 ) 其中，n 称为形函数矩阵。 ( 3 ) 求得单元的各种类型的矩阵和向量。根据能量原理或变分原理，得到各 单元的质量矩阵、刚度矩阵、阻尼矩阵和节点载荷向量 m 8 = i p n 7 n d q 娅l b r d b d q c 8 2 j l , y n 7 n d q = l n 礓