浅谈牛顿-欧拉方程及其推导

前文

最近阅读课题论文时看到采用牛顿-欧拉方程来建立动力学方程，尝试系统学习一下这一知识，发现网上资源推导说明不清晰不利于新手入门理解，此文用于记录推导过程以加深理解。
若文中出现错误或误导内容，欢迎指正！

---

---

为何物（背景）

描述刚体的力学变化和运动过程有很多种方法，例如基于能量的拉格朗日公式、牛顿力学法以及现在讨论的牛顿-欧拉方程。不同的方法以不同的基本定理和角度来描述刚体运动过程（能量守恒，动量守恒等）。牛顿-欧拉方程则是从刚体的转动和平动两个运动类型来描述力学与运动之间的关系。
刚体的平动描述的是描述对象的运动轨迹，刚体的转动描述的是描述对象自身的转动过程。在实际生活中，刚体的运动往往是平动和转动的复合过程。并且，刚体运动的描述是基于参考系下的，此处要引入惯性坐标系和非惯性坐标系的概念。

惯性坐标系即绝对坐标系，常说的地球坐标系就是典型的惯性坐标系。
非惯性坐标系就是基于运动描述对象建立的坐标系。
举个栗子（故意不小心打错的）
一辆火车上，你坐在座位上静止不动，在火车的坐标系上，你的运动速度就是为0。但是从地面坐标系来看，你跟着火车以一个非常快的速度运动着。因此基于不同的参考系，你的运动描述是不同的。
这里要埋个伏笔，由于非惯性系的与惯性系之间关于运动描述的差别，在进行力学分析时，引入了惯性力的概念。后续在进行推导时会提到这一概念
再举一个栗子（肥肠自信.jpg）
根据牛顿第一定理，如果物体上没有受力，它将一直保持当前的运动状态，这种运动状态也称作惯性。
假设游乐园的旋转木马匀速的旋转着。在旋转木马的坐标系中，它相对地面坐标系进行旋转运动。当你坐到木马上之后，在旋转木马的坐标系中，你就是静止不动的，但是在地面坐标系中，你就是在转圈。
你的转圈运动不是惯性运动，你的运动速度一直被改变，存在一个看不见的力持续的改变你的运动状态。这个就是因为非惯性坐标系的缘故产生的一个改变运动状态的力，称作惯性力。

从何来（公式推导）

线运动推导

根据一个图来讲解，惯性坐标系及非惯性坐标系的描述。

此图摘自https://blog.csdn.net/huangjunsheng123/article/details/110249073，侵删！
对于旋转过程的描述，通常采用的是旋转矩阵的形式来描述，关于坐标系旋转变换的知识，可以参考我之前写过的这篇博客 坐标系旋转变换。
图纸$U$为惯性坐标系，$B$为非惯性坐标系，$A$为空间中一点，${\vec{r}}_A$表示该点在惯性坐标系下的位置矢量，${\vec{r}}_{A|B}$表示非惯性坐标系下该点的位置矢量。${\vec{r}}_B$为非惯性坐标系在惯性坐标系中的位置矢量。
从简单的矢量运算开始
$$\overrightarrow{r_A}=\overrightarrow{r_B}+\overrightarrow{r_{A|B}}$$
这里就涉及到一个问题，${\vec{r}}_{A|B}$是基于$B$坐标系的向量，它的向量值不是基于惯性坐标系的，因此做矢量相加时，需要将这个向量变换到惯性坐标系中，才可以做相加。这里把${\vec{r}}_{A|B}$平移到惯性坐标系原点，设旋转矩阵为$C_B^U$表示从非惯性坐标系变换到惯性坐标系。
则对上面的公式进行修正
$$\overrightarrow{r_A}=\overrightarrow{r_B}+C_B^U\cdot \overrightarrow{r_{A|B}}$$
这里的修正很重要，我看到很多帖子上都是按照未修正的公式进行推导的，看的我一脸懵。我觉得这些细节是需要明确指出，否则无法快速理解后续的推导。
——
这里要做一个强调，后面的推导是统一在惯性坐标系下进行的。
——
对两边求导，得到
$$\frac{d\overrightarrow{r_A}}{dt}=\frac{d\overrightarrow{r_B}}{dt}+\frac{d{C}_B^U\cdot \overrightarrow{r_{A|B}}}{dt}$$

这里要补充一个内容点——旋转矩阵同样是关于时间的函数，因此对时间求导时需要对其一并计算。

• 对于旋转矩阵的求导公式，可以参考这个博客 https://blog.csdn.net/weixin_38632538/article/details/106085426
• 对向量求导存在以下公理
  $$\frac{d\vec{r}}{dt}=\vec{\omega }\times \vec{r}$$
  这个公式很好理解，即坐标系只进行旋转运动时，向量的大小不会改变，只会改变其角度，改变的趋势与坐标系的旋转角速度有关。
• 假设向量$\vec{r}$在非惯性坐标系内，随着非惯性坐标系做旋转运动，其中 $\overrightarrow{{\omega }_U}$表示该向量在惯性坐标系下的角速度,$\overrightarrow{{\omega }_B}$表示该向量在非惯性坐标系下的角速度。设旋转矩阵是惯性坐标系下任意三个线性无关列向量组合成的矩阵，则有
  $$\frac{d{C}_B^U}{dt}=\frac{d\left[\begin{array}{ccc}\overrightarrow{b_1^U}& \overrightarrow{b_2^U}& \overrightarrow{b_3^U}\end{array}\right]}{dt}=\left[\begin{array}{ccc}\overrightarrow{{\omega }_U}\times \overrightarrow{b_1^U}& \overrightarrow{{\omega }_U}\times \overrightarrow{b_2^U}& \overrightarrow{{\omega }_U}\times \overrightarrow{b_3^U}\end{array}\right]$$
  又根据坐标系旋转变换规律
  $$\overrightarrow{{\omega }_U}=C_B^U\cdot \overrightarrow{{\omega }_B}$$
  并且满足矩阵计算满足以下规律
  $$\begin{gathered} \overrightarrow{b_1^U}=C_B^U\cdot \left[\begin{array}{c}1\\ 0\\ 0\end{array}\right]=C_B^U\cdot \overrightarrow{e_1} \\ \overrightarrow{b_2^U}=C_B^U\cdot \left[\begin{array}{c}0\\ 1\\ 0\end{array}\right]=C_B^U\cdot \overrightarrow{e_2} \\ \overrightarrow{b_3^U}=C_B^U\cdot \left[\begin{array}{c}0\\ 0\\ 1\end{array}\right]=C_B^U\cdot \overrightarrow{e_3} \end{gathered}$$
  联立上述三个式子可以得到
  $$\begin{gathered} \frac{d{C}_B^U}{dt}=\left[\begin{array}{ccc}(C_B^U\cdot \overrightarrow{{\omega }_B})\times (C_B^U\cdot \overrightarrow{e_1})& (C_B^U\cdot \overrightarrow{{\omega }_B})\times (C_B^U\cdot \overrightarrow{e_2})& (C_B^U\cdot \overrightarrow{{\omega }_B})\times (C_B^U\cdot \overrightarrow{e_3})\end{array}\right] \\ \Rightarrow \frac{d{C}_B^U}{dt}=C_B^U\left[\begin{array}{ccc}\overrightarrow{{\omega }_B}\times \overrightarrow{e_1}& \overrightarrow{{\omega }_B}\times \overrightarrow{e_2}& \overrightarrow{{\omega }_B}\times \overrightarrow{e_3}\end{array}\right] \end{gathered}$$
• 这里就是涉及到一个如何将叉乘变化为点乘的知识点——反对称矩阵
  设任意三维列向量 $\vec{A}\vec{B}$，则满足以下关系
  $$\vec{A}\times \vec{B}=\left[\begin{array}{ccc}{a}_2{b}_3-a_3{b}_2& a_3{b}_1-a_1{b}_3& a_1{b}_2-a_2{b}_1\end{array}\right]=S(\vec{A})\cdot \vec{B}$$
  其中$S(\vec{A})$为反对称矩阵，其表达式为
  $$S(\vec{A})=\left[\begin{array}{ccc}0& -a_3& a_2\\ a3& 0& -a_1\\ -a_2& a_1& 0\end{array}\right]$$
  将上式重新代入求导公式可得
  $$\begin{gathered} \frac{d{C}_B^U}{dt}=C_B^U\left[\begin{array}{ccc}S(\overrightarrow{{\omega }_B})\cdot \overrightarrow{e_1}& S(\overrightarrow{{\omega }_B})\cdot \overrightarrow{e_2}& S(\overrightarrow{{\omega }_B})\cdot \overrightarrow{e_3}\end{array}\right] \\ \Rightarrow \frac{d{C}_B^U}{dt}=C_B^{U}S(\overrightarrow{{\omega }_B})\cdot I=C_B^{U}S(\overrightarrow{{\omega }_B}) \end{gathered}$$

按照上面的公式结论来继续求导。
$$\begin{gathered} \frac{d\overrightarrow{r_A}}{dt}=\frac{d\overrightarrow{r_B}}{dt}+\frac{d{C}_B^U\cdot \overrightarrow{r_{A|B}}}{dt}=\overrightarrow{v_B}+C_B^U\cdot \overrightarrow{v_{A|B}}+\frac{d{C}_B^U}{dt}\overrightarrow{r_{A|B}} \\ \overrightarrow{v_A}=\overrightarrow{v_B}+C_B^U\cdot \overrightarrow{v_{A|B}}+C_B^{U}S(\overrightarrow{{\omega }_B})\overrightarrow{r_{A|B}} \\ \overrightarrow{v_A}=\overrightarrow{v_B}+C_B^U\cdot (\overrightarrow{v_{A|B}}+S(\overrightarrow{{\omega }_B})\overrightarrow{r_{A|B}}) \\ \Rightarrow \overrightarrow{v_A}=\overrightarrow{v_B}+C_B^U\cdot (\overrightarrow{v_{A|B}}+\overrightarrow{{\omega }_B}\times \overrightarrow{r_{A|B}}) \end{gathered}$$

注意这里括号内的向量$(\overrightarrow{v_{A|B}}+\overrightarrow{{\omega }_B}\times \overrightarrow{r_{A|B}})$，是基于非惯性坐标系的，经过旋转矩阵变换后才变换到惯性坐标系的。
这里是一个伏笔，后面会考！

接下来就是通过速度求导得到加速度方程，也就可以得到线运动的动力学方程。
如果前面看懂了的话，后面这一步就是高数里简单的复合函数求导题目，有手就行。
$$\begin{gathered} \frac{d\overrightarrow{v_A}}{dt}=\frac{d\overrightarrow{v_B}}{dt}+\frac{d{C}_B^U}{dt}(\overrightarrow{v_{A|B}}+\overrightarrow{{\omega }_B}\times \overrightarrow{r_{A|B}})+C_B^U\frac{d(\overrightarrow{v_{A|B}}+\overrightarrow{{\omega }_B}\times \overrightarrow{r_{A|B}})}{dt} \\ \Rightarrow \overrightarrow{a_A}=\overrightarrow{a_B}+C_B^U(\overrightarrow{{\omega }_B}\times (\overrightarrow{v_{A|B}}+\overrightarrow{{\omega }_B}\times \overrightarrow{r_{A|B}}))+C_B^U(\overrightarrow{a_{A|B}}+\overrightarrow{{\omega }_B^{\prime }}\times \overrightarrow{r_{A|B}}+\overrightarrow{{\omega }_B}\times \overrightarrow{v_{A|B}}) \\ \Rightarrow \overrightarrow{a_A}=\overrightarrow{a_B}+C_B^U(\overrightarrow{{\omega }_B}\times \overrightarrow{v_{A|B}})+C_B^U(\overrightarrow{{\omega }_B}\times \overrightarrow{{\omega }_B}\times \overrightarrow{r_{A|B}})+C_B^U\overrightarrow{a_{A|B}}+C_B^U(\overrightarrow{{\omega }_B^{\prime }}\times \overrightarrow{r_{A|B}})+C_B^U(\overrightarrow{{\omega }_B}\times \overrightarrow{v_{A|B}}) \\ \Rightarrow \overrightarrow{a_A}=\overrightarrow{a_B}+C_B^U\cdot 2(\overrightarrow{{\omega }_B}\times \overrightarrow{v_{A|B}})+C_B^U(\overrightarrow{{\omega }_B}\times \overrightarrow{{\omega }_B}\times \overrightarrow{r_{A|B}})+C_B^U\overrightarrow{a_{A|B}}+C_B^U(\overrightarrow{{\omega }_B^{\prime }}\times \overrightarrow{r_{A|B}}) \end{gathered}$$

到这里就已经推导完成，下面简单讨论上面各项的物理意义
$\overrightarrow{a_B}$ 是非惯性坐标系的加速度；
$2(\overrightarrow{{\omega }_B}\times \overrightarrow{v_{A|B}})$ 是由于非惯性系旋转产生的惯性力，称作科氏力 （Coriolis Force）
$\overrightarrow{{\omega }_B}\times \overrightarrow{{\omega }_B}\times \overrightarrow{r_{A|B}}$ 是法向加速度，可以把向量看成可以拉伸的弹力绳，在非惯性系旋转过程中，这跟绳子会被甩出去变长，使其变长的力。
$\overrightarrow{a_{A|B}}$ 是矢量$A$在非惯性坐标系下的加速度，再经过坐标变换到惯性坐标系下。
$\overrightarrow{{\omega }_B^{\prime }}\times \overrightarrow{r_{A|B}}$ 是由于非惯性的变速旋转运动，导致矢量的切向方向运动变化的力。

角运动推导

角运动的过程与线运动本质上相同。直接给公式
$$\overrightarrow{{\omega }_A}=\overrightarrow{{\omega }_B}+C_B^U\cdot \overrightarrow{{\omega }_{A|B}}$$
两边求导得到
$$\overrightarrow{{\alpha }_A}=\frac{d\overrightarrow{{\omega }_A}}{dt}=\overrightarrow{{\alpha }_B}+C_B^U(\overrightarrow{{\omega }_B}\times \overrightarrow{{\omega }_{A|B}})+C_B^U\cdot \overrightarrow{{\alpha }_{A|B}}$$

总结

总体来说，整个推导过程相对简单，但是在我学习过程中查阅的资料里并没有对旋转变换过程进行描述，导致我自己推导过程遇到许多阻碍，希望这篇文章可以帮助到入门学习的同学。
若文中出现错误或者描述有误的地方，欢迎指正！