Robust Quasistatic Finite Elements and Flesh Simulation

文章来源:
Robust Quasistatic Finite Elements and Flesh Simulation

3.Quasistatic Formulation

{
公式xt→=v⃗ 以及vt→=M−1f⃗ (t,x⃗ ,t⃗ )中下标t表示对t求导
}

4.Strain Energy

That is, the global stiffness matrix −∂f⃗ /∂x⃗  is always symmetric, as a result of the hyperelastic energy having continuous second derivatives with respect to the spatial configuration.

{
[引用于https://en.wikipedia.org/wiki/Symmetry_of_second_derivatives#Schwarz.27_theorem]
Schwarz’ theorem
In mathematical analysis, Schwarz’ theorem (or Clairaut’s theorem) named after Alexis Clairaut and Hermann Schwarz, states that if
f:Rn→R
has continuous second partial derivatives at any given point in Rn , say, (a1,…,an), then ∀i,j∈{1,2,…,n},

∂2f∂xi∂xj(a1,…,an)=∂2f∂xj∂xi(a1,…,an),

The partial derivations of this function are commutative at that point.
}
{
所以对于Ψ的黑塞H矩阵中的元素来说Hij=∂2Ψ∂xi∂xj=∂2Ψ∂xj∂xi=Hji, 所以这里说the global stiffness matrix是对称的.
}

Furthermore, a steady state corresponds to a local minimum of the hyperelastic energy indicating that the energy Hessian,∂2Ψ/∂x⃗ 2, (or equivalently the global stiffness matrix)is positive definite in the vicinity of an isolated steady state.
{
[引用于Hessian Matrix]
也就是在临界点时,黑塞矩阵正定则对应函数值的最小点
}
{
Critical points
If the gradient (the vector of the partial derivatives) of a function f is zero at some point x, then f has a critical point (or stationary point) at x.
}

7.Diagonalization

公式(4)的推导

δP=∂UP(UTFV)VT∂(UTFV)∣∣∣F:δ(UTFV)=U{∂P(F)∂F∣∣∣UTFV:UTδFV}VT

{
这里要证明两点:
(1)δ(UTFV)=UTδFV
(2)∂UP(UTFV)VT∂(UTFV)∣∣∣F:δ(UTFV)=U{∂P(F)∂F∣∣∣UTFV:UTδFV}VT
证明:
(1)

δ(UTFV)=δ(Uije^j⊗e^iFmne^m⊗e^nVpqe^p⊗e^q)=Uije^j⊗e^iδFmne^m⊗e^nVpqe^p⊗e^q=UTδFmnV

(2)

∂UP(UTFV)VT∂(UTFV)∣∣∣F:δ(UTFV)=∂UP(F)VT∂(F)∣∣∣UTFV:UTδ(F)V=∂(Uije^i⊗e^jPmne^m⊗e^nVpqe^q⊗e^p)∂(Fxye^x⊗e^y):ϵabe^a⊗e^b=Uij∂Pmn∂FxyVpq(e^i⊗e^j)⋅(e^m⊗e^n)⋅(e^q⊗e^p)⊗(e^x⊗e^y):ϵabe^a⊗e^b=Uij(e^i⊗e^j){∂Pmn∂Fxy(e^m⊗e^n)⋅(e^x⊗e^y:ϵabe^a⊗e^b)}Vpq(e^q⊗e^p)=Uij(e^i⊗e^j){∂Pmn∂Fxye^m⊗e^n⊗e^x⊗e^y:ϵabe^a⊗e^b}Vpq(e^q⊗e^p)=U{∂P(F)∂F∣∣∣UTFV:UTδFV}VT

注意这里运用到的一些公式:
引用来源:tensors.pdf p31,p83

}

8. Enforcing Positive Definiteness

{
我们先看下P与F的关系
引用于The classical FEM method and discretization methodology: Course Notes (Eftychios Sifakis) p23

将F写成F=UF^VT，则P=U(ΨI⋅2F^+ΨII⋅4F^3+ΨIII⋅2III⋅F^−1)VT
写成索引形式:
Pij=f(σ1,σ2,σ3)kUikVjk
将F写成索引形式:
Fij=σkUikVjk
那么
∂Pij∂Fij∣∣∣F^=∂Pij∂Fij∣∣∣σ1,σ2,σ3
注意这里自变量为F^,即把变量限制为σ1,σ2,σ3,当然此时Uik,Vjk此时就变成了δik,δjk
另外约定: ∂c∂a,∂a∂c为0,∂c1∂c2是无意义的,其中a为变量，c,c1,c2为变量
所以通过这种方式可以确定9*9矩阵里面的前三行和前三列，交叉部分的偏导可以把变量限制为σ1,σ2,σ3,其余部分为0,那么求A的时候我们会用这种形式来求.B12,B13,B23我们用另外一种方式来求.

先来求B
另外注意一点σ1,σ2,σ3与F的非对角线元素关于变量σ1,σ2,σ3相互独立的
即
∂f(σ1,σ2,σ3)∂Fij∣∣∣F^=∂f(σ1,σ2,σ3)∂σkUikVjk∣∣∣F^=∂f(σ1,σ2,σ3)∂σkδij∣∣∣(σ1,σ2,σ3)=∂f(σ1,σ2,σ3)∂0∣∣∣(σ1,σ2,σ3)=0
注意由于是非对角线元素所以，i≠j,即δij=0
下面开始求B,注意下面中的i≠j,p≠q
∂Fij∂Fpq=δipδjq

∂(FFTF)ij∂Fpq=∂FikFlkFlj∂Fpq=∂Fik∂FpqFlkFlj+Fik∂Flk∂FpqFlj+FikFlk∂Flj∂Fpq=δipδkqFlkFlj+δlpδkqFikFlj+δlpδjqFikFlk=δipδkqFlkFlj+FiqFpj+δlpδjqFikFlk

这里强调下δ的化简方式通常是往常量去化,比如说:δliδkjFikFlj这里i,j为常量,l,k为变量 ,所以说应该化FijFij

下面引用来源:tensors.pdf p50
A−1=adj(A)|A|=cof(A)T|A|
PDF上其实有点错误,修正因为
C1i=ϵijkA2jA3k,C2i=ϵijkA3jA1k,C3i=ϵijkA1jA2k
这里C代表A的余子阵cofacter

所以

F−Tij=cof(F)ij|F|=ϵjmnF(i⊕1)mF(i⊕2)n|F|

注意这里\oplus记成模3循环加法,我也不知道怎么定义=.=,自创的,即2⊕2=1

∂F−Tij∂Fpq=∂ϵjmnF(i⊕1)mF(i⊕2)n∂|F|Fpq

现在小结一下:

P1=∂Fij∂Fpq=δipδjq

P2=∂(FFTF)ij∂Fpq=δipδkqFlkFlj+FiqFpj+δlpδjqFikFlk

P3=∂F−Tij∂Fpq=∂ϵjmnF(i⊕1)mF(i⊕2)n∂|F|Fpq

这样就把B中的α,β求出来了

下面来求A
Pii=UikVik(ΨI⋅2σk+ΨII⋅4σ3k+ΨIII⋅2III⋅1σk)
Fjj=UjkVjkσk
如果i=j

∂Pii∂Fii∣∣∣F^=∂(ΨI⋅2σi+ΨII⋅4σ3i+ΨIII⋅2III⋅1σi)∂σi=∂([2σi,4σ3i,2IIIσi]⎡⎣⎢ΨIΨIIΨIII⎤⎦⎥)∂σi=∂([2σi,4σ3i,2IIIσi])∂σi⋅⎡⎣⎢ΨIΨIIΨIII⎤⎦⎥+[2σi,4σ3i,2IIIσi]⋅∂⎡⎣⎢ΨIΨIIΨIII⎤⎦⎥∂σi=αii+βii+4IIIΨIIIσ2i+[2σi,4σ3i,2IIIσi]⋅⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢ΨI,I⋅∂I∂σiΨI,II⋅∂II∂σiΨI,III⋅∂III∂σiΨII,I⋅∂I∂σiΨII,II⋅∂II∂σiΨI,III⋅∂III∂σiΨIII,I⋅∂I∂σiΨIII,II⋅∂II∂σiΨIII,III⋅∂III∂σi⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥=αii+βii+4IIIΨIIIσ2i+[2σi,4σ3i,2IIIσi]⋅⎡⎣⎢ΨI,IΨI,IIΨI,IIIΨII,IΨII,IIΨI,IIIΨIII,IΨIII,IIΨIII,III⎤⎦⎥⋅⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢∂I∂σi∂II∂σi∂III∂σi⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥

引用于The classical FEM method and discretization methodology: Course Notes (Eftychios Sifakis) p23

∂Pii∂Fii∣∣∣F^=αii+βii+4IIIΨIIIσ2i+[2σi,4σ3i,2IIIσi]⋅⎡⎣⎢ΨI,IΨI,IIΨI,IIIΨII,IΨII,IIΨI,IIIΨIII,IΨIII,IIΨIII,III⎤⎦⎥⋅⎡⎣⎢⎢⎢2σi4σ3i2IIIσi⎤⎦⎥⎥⎥=αii+βii+γii

当i≠j

∂Pii∂Fjj∣∣∣F^=∂([2σi,4σ3i,2IIIσi]⎡⎣⎢ΨIΨIIΨIII⎤⎦⎥)∂σj=∂([2σi,4σ3i,2IIIσi])∂σj⋅⎡⎣⎢ΨIΨIIΨIII⎤⎦⎥+[2σi,4σ3i,2IIIσi]⋅∂⎡⎣⎢ΨIΨIIΨIII⎤⎦⎥∂σj=4IIIΨIIIσiσj+[2σi,4σ3i,2IIIσi]⋅⎡⎣⎢ΨI,IΨI,IIΨI,IIIΨII,IΨII,IIΨI,IIIΨIII,IΨIII,IIΨIII,III⎤⎦⎥⋅⎡⎣⎢⎢⎢⎢2σj4σ3j2IIIσj⎤⎦⎥⎥⎥⎥=γij

}