Page 310 - 《软件学报》2020年第10期

3286 Journal of Software 软件学报 Vol.31, No.10, October 2020

其次,以预测位置 p={p 0 ,p 1 ,…,p n }更新顶点位置 x={x 0 ,x 1 ,…,x n },即顶点位置 x=p.构造 E(x+Δx)的一阶泰勒
展开式 E(x+Δx)≈E(x)+∇ x E(x)·Δx=0,其中,E(·)为能量函数,包括弹性变形能 E S 、弯曲变形能 E B 以及各向异性变
形能 E A ,总能量为 3 种能量之和.计算各顶点位移修正量Δx=[Δx 1 ,Δx 2 ,…,Δx n ],使 E(x+Δx)=0.顶点位移修正量应
[5]
满足所有能量函数(设共有 N 个能量函数),需求解线性系统 :
∇ x E 1 ()⋅ Δ x =x − E 1 ()x
⋅⋅⋅ (5)
()⋅
∇ x E x Δ x = − E N ()x
N
[5]
本文基于位置动力学方法,以非线性高斯-赛德尔方法 求解式(5)约束函数,计算Δx,并以Δx 更新顶点位置.
最后,碰撞检测.检测发生碰撞的图元对,计算顶点位移修正量Δx,使 C vf (x+Δx)≥0,C ee (x+Δx)≥0,更新顶点位
[5]
置 .
t
设模型初始体积为 V 0 ,变形后体积为 V′,则体积比 r=V′/V 0 .定义阈值ε,当前时刻 t 顶点位置的集合为 x ,时间
步长为Δt,则每一帧薄壳变形算法流程描述如下.
t
输入:当前时刻顶点位置集合 x ;
t+1
输出:下一时刻顶点位置集合 x .
t+1
Step 1. 仅考虑外部压力,以 Sympletic 欧拉法计算顶点预测位置 p ;
t+1
t+1
Step 2. 更新顶点位置 x ←p ,基于公式(2)~公式(4)给出的 E S 、E B 和 E A 构建线性系统(5);
t+1
t+1
Step 3. 以非线性高斯-赛德尔方法迭代求解线性系统,计算顶点位移修正量Δx,更新顶点位置 x ←x +
Δx,迭代 5 次可得到较为精确的结果;
Step 4. 检测碰撞,若发生自碰撞,则计算新的顶点位置,使 C vf ≥0,C ee ≥0;
Step 5. 更新顶点位置,计算体积比 r,若 r 小于ε,则变形终止;
t
t+1
t+1
Step 6. 更新顶点速度 v ←(x -x )/Δt;
Step 7. 更新时间 t←t+Δt,跳转 Step 1,计算下一时刻顶点位置.
6 实验结果与分析
本文的实验环境为 Inter(R) Core(TM) i7 处理器,主频 3.4GHz,内存 32GB,PC,操作系统为 64 位 Windows.
编程语言环境为 Microsoft Visual Studio 2010 C++.

6.1 部分收缩变形实验结果
–2
图 8 为变形至 100 帧时,分别添加弹性变形能 E S 、E S +弯曲变形能 E B (k B =1)、E S +E B (k B =10 )、E S +E B +各向
异性变形能 E A 的情况下对应变形结果、平均曲率分布和能量分布.其中,红色区域能量较高,蓝色区域能量较低;
红色区域曲率值较大,蓝色区域曲率值较小.
首先,从图 8(a)可以看出,显然,在弹性变形能 E S 较高的区域,模型变形显著;弹性变形能 E S 较低,模型变形较
微弱.
其次,对比图 8(a)和图 8(c)可以发现,图 8(a)无弯曲变形能 E B ,红框区域被过度拉伸,平均曲率分布不均匀,
说明模型变形效果较粗糙;图 8(c)添加弯曲变形能 E B 后,红框区域平均曲率分布较均匀,模型表面较光滑,收缩
变形效果较好,无过度拉伸.对比图 8(b)和图 8(c),当弯曲变形能 E B 的系数 k B =1 时,图 8(b)中模型无明显收缩,抗
−2
弯曲现象较强,当弯曲变形能 E B 的系数 k B =10 时,图 8(c)中模型发生较好的收缩变形.
最后,对比图 8(c)和图 8(d),显然,在添加各向异性变形能 E A 之后,红框区域内局部类球结构变形较显著;对
比图 8(a)、图 8(c)和图 8(d)红框区域内平均曲率分布可见,在加入 E A 后,红框区域内产生了显著的凹陷变形;对
比图 8(a)、图 8(c)和图 8(d)红框区域内能量分布可见,在加入 E A 后,红框区域内能量分布发生了明显变化,该区
域内局部类球结构消失.显然,在添加 E A 后,红框内的局部类球结构收缩变形较添加 E A 前更为显著,局部类球结
构变形幅度更为明显.针对各向同性能量作用下薄壳局部类球结构收缩变形缓慢且不显著的不足,各向异性变
形能量 E A 能够使局部类球结构及邻近区域快速产生显著变形.