Approximation of Euler angles with small rotation hypothesis Sorry for boring you during summer vacation, my friends. I am haunted by the approximated expression of Euler angle rotation matrix found in this textbook. In the appendix, the author declares that the conventional Euler (Tait–Bryan) angle matrix could be approximated into the following second order form as
$$
\begin{bmatrix}
1 & & \\ & 1 & \\ & & 1
\end{bmatrix}
+
\begin{bmatrix}
0 & -\theta_z & \theta_y \\ \theta_z & 0 & -\theta_x\\ -\theta_y& \theta_x& 0
\end{bmatrix}
-\frac{1}{2}
\begin{bmatrix}
\theta_y^2+\theta_z^2 & & \\ & \theta_x^2+\theta_z^2 & \\ & & \theta_x^2+\theta_y^2
\end{bmatrix}
+
\begin{bmatrix}
0 & \theta_x\theta_y & \theta_x\theta_z \\ \theta_x\theta_y & 0 & \theta_y\theta_z\\ \theta_x\theta_z& \theta_y\theta_z& 0
\end{bmatrix}
$$
Normally, the first two matrices are used for approximating the Euler angles following small rotation hypothesis. And in my mind, the second matrix represents physically the rotational motion of small rotations. However, no further explication and analysis are presented in the textbook for all the other matrices. Thus, I wonder if someone could share some insight with this form of approximation (I cannot find a similar approximation in anywhere else but in the textbook). Thanks a lot in advance.
 A: This follows from Rodrigues' rotation formula.
Define $\theta$ as the magnitude of the rotation, $\theta=\sqrt{\theta_x^2+\theta_y^2+\theta_z^2}$, and the matrix $\mathbf K$ as
$$\mathbf K = \frac1{\theta}\,\begin{bmatrix}
0&-\theta_z &\phantom{-}\theta_y \\
\phantom{-}\theta_z & 0 & -\theta_x \\
-\theta_y & \phantom{-}\theta_x & 0
\end{bmatrix}$$
Note that the square of this matrix is
$$\mathbf K^2 = \frac1{\theta^2}\,\begin{bmatrix}
-\theta_y^2-\theta_z^2&\theta_x\theta_y &\theta_x\theta_z \\
\theta_x\theta_y & -\theta_x^2-\theta_z^2 & \theta_y\theta_z \\
\theta_x\theta_z & \theta_y\theta_z & -\theta_x^2-\theta_y^2
\end{bmatrix}$$
By Rodriques' rotation formula, the rotation matrix is
$$\mathbf R = \mathbf I + \sin\theta\,\mathbf K + (1-\cos\theta)\,\mathbf K^2$$
To second order, this becomes
$$\begin{aligned}
\mathbf R &\approx \mathbf I + \theta\,\mathbf K + \frac12\theta^2\,\mathbf K^2 \\
&= \begin{bmatrix}1&0&0\\0&1&0\\0&0&1\end{bmatrix}
+ \begin{bmatrix}
0&-\theta_z &\phantom{-}\theta_y \\
\phantom{-}\theta_z & 0 & -\theta_x \\
-\theta_y & \phantom{-}\theta_x & 0
\end{bmatrix}
+ \frac12\begin{bmatrix}
-\theta_y^2-\theta_z^2&\theta_x\theta_y &\theta_x\theta_z \\
\theta_x\theta_y & -\theta_x^2-\theta_z^2 & \theta_y\theta_z \\
\theta_x\theta_z & \theta_y\theta_z & -\theta_x^2-\theta_y^2
\end{bmatrix}
\end{aligned}$$
