## 2 Dynamic Model of Quadrotor Aircrafts  
The coordinate system and structure of the quadrotor are shown in Figure 1. We use six variables to describe the motion of the quadrotor, namely the quadrotor position \(x\), \(y\), \(z\), the attitude angle \(\phi\), \(\theta\) and \(\psi\)[23]. The centroid of quadrotor moves upward continuously from below the origin of the fuselage mount when the load is gradually reduced. Referring to the body coordinate system B and the ground fixed coordinate system E, the dynamic model of the variable mass quadrotor aircraft is established [24]. As shown in Figure 1, the rotation matrix from coordinate system B to E is \(R_{B\to E}\)[25; 26]. The positional motion of the quadrotor satisfies [PERSON]'s second law,  
\[\left\{ \begin{array}{l}\dot{X}^{E}=V^{E}\\ (M-\Delta m)\dot{V}^{E}=L^{E}+G_{M}^{E}+G_{\Delta m}^{E}+f_{1}^{E}\end{array}  \right., \tag{1}\]  
where \(X^{E}=[x\quad y\quad z]^{T}\); \(V^{E}\) is the linear velocity; \(M\) is the mass of the quadrotor; \(L^{E}\) is the total lift; \(G_{\Delta m}^{E}= \begin{bmatrix}0&0\quad\Delta mg\end{bmatrix} ^{T}\); \(G_{M}^{E}= \begin{bmatrix}0&0\quad-Mg\end{bmatrix} ^{T}\); \(f_{1}^{E}= \begin{bmatrix}-k\dot{x}&-k\dot{y}&-k\dot{z}\end{bmatrix} ^{T}\) is the matrix of air resistance. The lift force of the quadrotor is expressed as  
\[L^{B}= \begin{bmatrix}L_{x}\quad L_{y}\quad L_{z}\end{bmatrix} ^{T}= \begin{bmatrix} 0\quad 0\end{bmatrix} ^{T}, \tag{2}\]  
where \(L=\sum\limits_{i=1}^{4}L_{i}=L_{1}+L_{2}+L_{3}+L_{4}\); from \(R_{B\to E}\), we can obtain \(L^{E}\):  
\[L^{E}=R_{B\to E}L^{B}=L \begin{bmatrix}\sin\psi\sin\phi+\cos\phi\cos\psi\sin \theta\\ -\sin\phi\cos\psi+\cos\phi\sin\psi\sin\theta\\ \cos\phi\cos\theta\end{bmatrix} . \tag{3}\]  
In summary, we can obtain the positional dynamics equation of the quadrotor. \[\left\{ \begin{array}{l}\ddot{x}=\frac{-\dot{k}_{f}\dot{x}}{(M-\Delta m)}+ \frac{L(\sin\psi\sin\phi+\cos\psi\cos\phi\sin\theta)}{(M-\Delta m)}\\ \ddot{y}=\frac{-\dot{k}_{f}\dot{y}}{(M-\Delta m)}+\frac{L(-\cos\psi\sin\phi+ \sin\psi\cos\phi\sin\theta)}{(M-\Delta m)}\\ \ddot{z}=\frac{-\dot{k}_{f}\dot{z}}{(M-\Delta m)}+\frac{L\cos\phi\cos\theta}{ (M-\Delta m)}-g\end{array} \right.. \tag{4}\]  
As the load mass decreases, the centroid of quadrotor gradually moves upward, thereby affecting the attitude angle of the quadrotor. Because the quadrotor is loaded by the mounting platform, the two can be considered as a whole. The initial coordinates of the centroid of quadrotor are \((x_{1},y_{1},z_{1})\), the coordinates of the centroid of the part of the load  
Figure 1: Structure diagram of a quadrotor. that is reduced during operation are \((x_{2},y_{2},z_{2})\). When the load mass decreases, the change in the coordinates of the centroid of quadrotor is shown in the following formula:  
\[ \begin{cases}x_{a}=(x_{1}M-x_{2}\Delta m)/(M-\Delta m)\\ y_{a}=(y_{1}M-y_{2}\Delta m)/(M-\Delta m)\\ z_{a}=(z_{1}M-z_{2}\Delta m)/(M-\Delta m)\\ \end{cases} , \tag{5}\]  
where the coordinates of the centroid of the quadrotor are \((x_{a},y_{a},z_{a})\). According to Figure 1, initial coordinate of the centroid of the quadrotor is \((0,0,0)\), the coordinate of the centroid of the reduced load mass during the working process is \((0,0,z_{c})\). Substituting it into Formula (6), the coordinate of the centroid of the quadrotor can be obtained as \((0,0,z_{a})\). \[z_{a}=z_{c}m/(M-\Delta m). \tag{6}\]  
According to the Newton-Euler equation, the attitude dynamics equation can be expressed as [27]. The attitude angle motion of the quadrotor satisfies the Newton-Euler equation, and its attitude dynamic model can be expressed as [28; 29]  
\[\left\{ \begin{array}{l}I_{3\times 3}\dot{\tau}^{B}+\tau^{B}\times I_{3\times 3 }\tau^{B}=\varsigma^{B}\\ \varsigma^{B}=Q^{B}+f_{2}^{B}+P^{B}\end{array}