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

孙晓鹏 等:局部各向异性的薄壳收缩变形 3283

×
为 F = s − m 1 ∈ D D R 3 2 , 其应变能可定义为 [13]
E=A j Ψ(F),
其中,A j 为三角面片 f j 的面积,Ψ(F)为能量密度场,则面 f j 这 3 个顶点的应变能梯度分别记为
⎡ E ∂ E ⎤ ∂ ∂ E 2 ∂ E
()D
⎢ ⎥ = A PF − m T , = − ∑ (1)
j
⎢ ∂x 1 j ∂ 2 j ⎥ ⎣ ⎦ ∂x x 3 j i= 1 ∂ x ji
其中,P(F)为Ψ(F)的 1st Piola-Kirchhoff 应力张量.
基于变形梯度 F 计算面 f j 的二维变形梯度 =FpF ∈ R 22× , 其中,p 为投影矩阵 [23] ,则格林应变张量为 =ε
1 ( T − ) FF I Saint Venant-Kirchhoff 模型应变能量密度场 Ψ 和 1st Piola-Kirchhoff 应力张量 () 分别记为
PF
,
2 S
Ψ = μ 2 λ tr ( ), PF = ⎡ 2μ ε ε + tr ε λ + () ⎤ I ,
2
ε
() F
S F 2 ⎣ ⎦
其中,μ和λ为 Lame 系数,||·|| F 为 Frobenius 范数,tr(·)表示矩阵的迹,则薄壳模型 M 的弹性变形能 E S 可记为
= Ψ
E S ∫ Ω S dX (2)
根据公式(1)计算每个顶点的二维空间能量梯度,并基于 p 将能量梯度映射到三维空间,得到顶点三维空间
∂ E
)i ∈
能量梯度 S ( 1,2,3 .
∂x
ji
3.2 弯曲变形能
本节基于 Bender 等离散等距弯曲模型 [13,24] ,定义薄壳模型 M 变形过程中的弯曲变形能 E B ,以阻止模型收缩
变形过程中产生过度的面外弯曲,驱动模型恢复初始状态.

n1
n0
xs1 xs1
e1 e4
xs2 e0 xs3
xs0
xs3
xs2 e2 e3
xs0
θ
Fig.3 The geometry structure of a stencil
图 3 模板 s 的结构

对于边 e i ∈E 及其两个邻接三角面片,其 4 个顶点的集合记为 x s ={x s0 ,x s1 ,x s2 ,x s3 }、5 条边的集合记为 e s =
{x s0 x s1 ,x s1 x s2 ,x s2 x s0 ,x s0 x s3 ,x s3 x s1 },则模板 s 的定义如图 3 所示,其中,两个邻接三角面片的法向量分别记为 n 0 和 n 1 ,
其夹角为θ,则模板 s 的弯曲变形能定义为
1
s
T
E ( B s )x = ∑ Q x x (3)
si sj
, i j
2 , ij
3
T
这里, Q s = K K , A 0 和 A 1 分别为边 e i 的两个邻接三角面片的面积,其中,行向量 K 定义为(c 01 +c 04 ,c 02 +c 03 ,
A + A
0 1
–c 01 –c 02 ,–c 03 –c 04 ),其中,c ab =cot∠(e a ,e b ),则顶点 si 上的弯曲变形能梯度可定义为
∂ E
s
.
B =∑ Qx
∂x , ij sj
si j
图 4 展示模板变形过程中呈现的不同状态,设模板初始结构为黑色,无弯曲变形能 E B 实施变形后结构为紫
0
色,定义初始面片夹角为β .添加 E B 后,紫色模板受 E B 影响具有恢复初始结构的运动趋势,若 E B 过大,模板会发生
过度抗弯曲现象以抵抗收缩(如蓝色模板所示),使模型产生抖动.理想模板运动范围应位于紫色模板和黑色模
板之间(如绿色模板所示).