欧拉方法

显式欧拉方法

$\begin{cases}\vec{v}\left(t_{1}\right)=\vec{v}\left(t_{0}\right)+M^{-1} \vec{F}\left(t_{0}\right)\Delta t \\ \vec{x}\left(t_{1}\right)=\vec{x}\left(t_{0}\right)+\vec{v}\left(t_{0}\right) \Delta t\end{cases}$

用当前的状态计算速度和加速度$(\text{Forward} )$

能量增加，开放的系统

$\int_{t^{[0]}}^{t^{[1]}}\boldsymbol v(t)\text dt=\boldsymbol v(t^{[0]})\Delta t$

隐式欧拉方法

$\vec{v}\left(t_{1}\right)=\lambda(\vec{v}\left(t_{0}\right)+M^{-1} \vec{F}\left(t_{1}\right)\Delta t) \\ \vec{x}\left(t_{1}\right)=\vec{x}\left(t_{0}\right)+\vec{v}\left(t_{1}\right) \Delta t$

用未来的状态去计算速度和加速度再来更新$(\text{Backward} )$

$\int_{t^{[0]}}^{t^{[1]}}\boldsymbol v(t)\text dt=\boldsymbol v(t^{[1]})\Delta t$

中点方法

$\int_{t^{[0]}}^{t^{[1]}}\boldsymbol v(t)\text dt=\boldsymbol v(t^{[0.5]})\Delta t$

二阶准确

半隐式欧拉方法

$\vec{v}\left(t_{1}\right)=\vec{v}\left(t_{0}\right)+M^{-1} \vec{F}\left(t_{0}\right)\Delta t \\ \vec{x}\left(t_{1}\right)=\vec{x}\left(t_{0}\right)+\vec{v}\left(t_{1}\right) \Delta t$

本质上为显式方法

旋转

欧拉角

$\mathbf{R}_{x y z}(\alpha, \beta, \gamma)=\mathbf{R}_{x}(\alpha) \mathbf{R}_{y}(\beta) \mathbf{R}_{z}(\gamma)$

任意旋转都可以分解为绕$x,y,z$三个轴的旋转

顺序非常重要

$\mathbf{R}_{x}(\alpha) \mathbf{R}_{y}(\beta) \mathbf{R}_{z}(\gamma)\ne\mathbf{R}_{x}(\alpha)\mathbf{R}_{z}(\gamma) \mathbf{R}_{y}(\beta)$

飞机滚转（Roll）；飞机俯仰（Pitch）；飞机偏航（Yaw）

罗德里格斯旋转公式

$\mathbf{R}(\mathbf{n}, \alpha)=\cos (\alpha) \mathbf{I}+(1-\cos (\alpha)) \mathbf{n} \mathbf{n}^{T}+\sin (\alpha) \underbrace{\left(\begin{array}{ccc} 0 & -n_{z} & n_{y} \\ n_{z} & 0 & -n_{x} \\ -n_{y} & n_{x} & 0 \end{array}\right)}_{\mathbf{N}}$

四元数

$q=(q_1, q_ 2, q _3, q _4)=\left(\cos \left(\frac{\theta}{2}\right), \sin \left(\frac{\theta}{2}\right) n_ x, \sin \left(\frac{\theta}{2}\right) n_ y, \sin \left(\frac{\theta}{2}\right)n_ z\right)$

对应于绕$\mathbf{n}$的旋转

$p'=qpq^{-1}$

质点弹簧系统

$f_{a\to b}$为弹簧作用在$a$上的力

能量耗散

引入阻力

$f=-k_d\dot{\boldsymbol x}$

但表示合速度，表示不了弹簧内部的损耗

$f_a=-k_s\frac{\boldsymbol a-\boldsymbol b}{||\boldsymbol a-\boldsymbol b||}(\dot{\boldsymbol a}-\dot{\boldsymbol b})\cdot\frac{\boldsymbol a-\boldsymbol b}{||\boldsymbol a-\boldsymbol b||}$

流体

不可压缩流体

$h_{i}\left(t_{0}+\Delta t\right)=2 h_{i}\left(t_{0}\right)-h_{i}\left(t_{0}-\Delta t\right)+\frac{\Delta t^{2}}{\Delta x^{2} \rho}\left(\left(\frac{h_{i-1}+h_{i}}{2}\right)\left(P_{i-1}-P_{i}\right)+\left(\frac{h_{i+1}+h_{i}}{2}\right)\left(P_{i+1}-P_{i}\right)\right)\\ h_{i}\left(t_{0}+\Delta t\right)=h_{i}\left(t_{0}\right)+\beta\left(h_{i}\left(t_{0}\right)-h_{i}\left(t_{0}-\Delta t\right)\right)+\alpha\left(h_{i+1}\left(t_{0}\right)+h_{i-1}\left(t_{0}\right)-2 h_{i}\left(t_{0}\right)\right)$

有限元

形函数

$\begin{pmatrix}X\\Y\\Z\end{pmatrix}=\sum_{i=0}^{i=3}\begin{pmatrix}X_i\\Y_i\\Z_i\end{pmatrix}\phi_i\left(\mathbf{X}\right) \\ \begin{pmatrix}X-X_0\\Y-Y_0\\Z-Z_0\end{pmatrix}=\begin{pmatrix}\Delta X_1&\Delta X_2&\Delta X_3\\\Delta Y_1&\Delta Y_2&\Delta X_3\\\Delta Z_1&\Delta Z_2&\Delta X_3\end{pmatrix}\begin{pmatrix}\phi_1\left(\mathbf{X}\right)\\\phi_2\left(\mathbf{X}\right)\\\phi_3\left(\mathbf{X}\right)\end{pmatrix}$

所以

$\begin{pmatrix}\phi_1\left(\mathbf{X}\right)\\\phi_2\left(\mathbf{X}\right)\\\phi_3\left(\mathbf{X}\right)\end{pmatrix}=T^{-1}(\mathbf{X}-\mathbf{X}_0)$

齐次坐标

变换

$\mathbf{x}\left(\mathbf{X}\right)=\sum_{i=0}^{i=3}\mathbf{x}_{i}\left(t\right){\phi_{i}\left(\mathbf{X}\right)}\\ \mathbf{x}\left(\mathbf{X}\right)=\begin{pmatrix}\phi_0\mathbf{I}&\phi_1\mathbf{I}&\phi_2\mathbf{I}&\phi_3\mathbf{I}\end{pmatrix}\begin{pmatrix}\mathbf{x}_0\\\mathbf{x}_1\\\mathbf{x}_2\\\mathbf{x}_3\end{pmatrix}$

求动能

$\mathbf{v}\left(\mathbf{X}\right)=\dot{\mathbf{x}}\left(\mathbf{X}\right)=\mathbf{N}\left(\mathbf{X}\right)\dot{\mathbf{q}}\\ \frac12\int_{\Omega}\rho\mathbf{v}^T\mathbf{v}d\Omega =\frac12\int_{\Omega}\rho\left(\dot{\mathbf{q}}^T\mathrm{N}\left(\mathbf{X}\right)^T\mathrm{N}\left(\mathbf{X}\right)\dot{\mathbf{q}}\right)d\Omega$

重组

$\frac12\dot{\mathbf{q}}^T\left(\int_{\Omega}\rho\mathrm{N}\left(\mathbf{X}\right)^T\mathrm{N}\left(\mathbf{X}\right)d\Omega\right)\dot{\mathbf{q}}=\frac12\dot{\mathbf{q}}^T\mathbf{M}_0\dot{\mathbf{q}}$

$\mathbf{M}_0$为单元的质量矩阵

$6\rho\cdot\boldsymbol{vol}\cdot\int_0^1\int_0^{1-\phi_1}\int_0^{1-\phi_1-\phi_2}\left(\phi_r\phi_s\right)d\phi_3d\phi_2d\phi_1$