航空機の運動方程式

何よりもまず結論から示すと、

並進運動　

(\begin{matrix} F_{x} \\ F_{y} \\ F_{z} \end{matrix}) = m (\begin{matrix} \dot{u} + q w - r v \\ \dot{v} + r u - p w \\ \dot{w} + p v - q u \end{matrix})

回転運動　

(\begin{matrix} \bar{L} \\ \bar{M} \\ \bar{N} \end{matrix}) = (\begin{matrix} I_{x} \dot{p} - (I_{y} - I_{z}) q r - I_{x y} (r p - \dot{q}) - I_{y z} (q^{2} - r^{2}) - I_{z x} (\dot{r} + p q) \\ I_{y} \dot{q} - (I_{z} - I_{x}) r p - I_{x y} (q r - \dot{p}) - I_{y z} (\dot{r} + p q) - I_{z x} (r^{2} - p^{2}) \\ I_{z} \dot{r} - (I_{x} - I_{y}) p q - I_{x y} (p^{2} - q^{2}) - I_{y z} (\dot{q} + r p) - I_{z x} (\dot{p} - q r) \end{matrix})

以下は、この式の導出から変形等を示す。

使用する記号等

機体固定座標系
　　x軸：機首方向　y軸：右翼方向　z軸：機体の腹方向　の右手系
　　重心：C

地面固定座標

x_{e}

軸 ,

y_{e}

軸 ,

z_{e}

軸　　　earth
　　原点：O
　　オイラー角　

(ϕ, θ, ψ)

ϕ

：

x_{e}

軸周りの角度 roll角
　　

θ

：

y_{e}

軸周りの角度 pitch角
　　

ψ

：

z_{e}

軸周りの角度 yaw角

機体の質量

m

機体の一部の微小質量

∆ m

重心の位置ベクトル

r_{c}

微小質量の位置ベクトル

r_{m}

微小質量の重心からの位置ベクトル

r

速度　

V = (u, v, w)

角速度（機体の重心C周りの回転ベクトル）　

ω = (p, q, r)

全外力　

F = (X, Y, Z)

モーメント　

(L, M, N)

角運動量　

H = (H_{x}, H_{y}, H_{z})

x軸,y軸,z軸に関する慣性能率　

(I_{x x}, I_{y y}, I_{z z})

x-y面,y-z面,z-x面に関する慣性乗積　

(I_{x y}, I_{y z}, I_{z x})

迎角

α

：

\tan α = \frac{w}{u}

横滑り角

β

：

\sin β = \frac{v}{V}

飛行経路角

γ

：

γ = θ - α

上反角

Γ

：

Γ = \frac{- z_{e}^{'}}{\sqrt{{x_{e}^{'}}^{2} + {y_{e}^{'}}^{2}}}

舵角

(δ_{a}, δ_{e}, δ_{r})

: エルロン角、エレベータ角、ラダー角

u^{'}, v^{'}, w^{'} = f (u, v, w, p, q, r, θ, ϕ, ψ, m, T, ρ, V, S, C_{x}, C_{y}, C_{z})

p^{'}, q^{'}, r^{'} = f (p, q, r, I_{x}, I_{y}, I_{z}, I_{x z}, ρ, V, S, c, b, C_{l}, C_{m}, C_{n})

回転座標系の時間微分の仕方

通常、固定座標系においては、
　

x_{e}

を微分すると、

\frac{d x_{e}}{d t}

となるが、

角速度　

ω

　で回転する回転座標系を基準に考える場合、
　

x_{e}

を微分すると、

\frac{d x_{e}}{d t} = \frac{d x}{d t} + ω \times x

となる。

航空機の運動方程式の導出（並進運動方程式）

重心の位置ベクトル

r_{c}

微小質量の位置ベクトル

r_{m}

重心から微小質量の位置ベクトル

r

　　の関係は、

　

r_{m} = r_{c} + r

微分すると
　

\frac{d r_{m}}{d t} = \frac{d r_{c}}{d t} + (\frac{d r}{d t} + ω \times r)

r

は時間によって変化することはない定数なので、

\frac{d r}{d t} = 0

\frac{d r_{m}}{d t} = \frac{d r_{c}}{d t} + ω \times r

\frac{d r_{c}}{d t}

は機体重心の速度ベクトルであり、機体座標系の成分で表すと　

\frac{d r_{c}}{d t} = v_{c} = (u, v, w)

なので、

　

\frac{d r_{m}}{d t} = v_{c} + ω \times r

　（速度の式）

もう一回微分すると

\begin{aligned} \frac{d^{2} r_{m}}{d t^{2}} & = (\frac{d r_{c}}{d t} + ω \times r) + (\frac{d}{d t} (ω \times r) + ω \times (ω \times r)) \\ = (\frac{d r_{c}}{d t} + ω \times r) + (\frac{d ω}{d t} \times r + ω \times \frac{d r}{d t}) + ω \times (ω \times r) \\ = \frac{d r_{c}}{d t} + ω \times v_{c} + \frac{d ω}{d t} \times r + ω \times (ω \times r) \end{aligned}

　（加速度の式）

微小質量に係る外力を

∆ F

とすると、
ニュートンの第2法則　

F = m α

　より、

\begin{aligned} ∆ F & = ∆ m \frac{d^{2} r_{m}}{d t^{2}} \\ = ∆ m (\frac{d r_{c}}{d t} + ω \times v_{c} + \frac{d ω}{d t} \times r + ω \times (ω \times r)) \end{aligned}

積分して、機体全体の質量にして考えると

\begin{aligned} F & = \int d m (\frac{d r_{c}}{d t} + ω \times v_{c} + \frac{d ω}{d t} \times r + ω \times (ω \times r)) \\ = \int d m (\frac{d r_{c}}{d t} + ω \times v_{c}) + \frac{d ω}{d t} \times \int r d m + ω \times (ω \times \int r d m) \end{aligned}

rは時間によって変化することはない定数なので、

\int r d m = 0

よって、

\begin{aligned} F & = \int d m (\frac{d r_{c}}{d t} + ω \times v_{c}) \\ = m (\frac{d r_{c}}{d t} + ω \times v_{c}) \end{aligned}

成分表示をすると

\begin{aligned} (\begin{array}{c} F_{x} \\ F_{y} \\ F_{z} \end{array}) & = m {\frac{d}{d t} (\begin{array}{c} u \\ v \\ w \end{array}) + (\begin{array}{c} p \\ q \\ r \end{array}) \times (\begin{array}{c} u \\ v \\ w \end{array})} \\ = m {(\begin{array}{c} \dot{u} \\ \dot{v} \\ \dot{w} \end{array}) + (\begin{array}{c} q w - r v \\ r u - p w \\ p v - q u \end{array})} \\ = m (\begin{array}{c} \dot{u} + q w - r v \\ \dot{v} + r u - p w \\ \dot{w} + p v - q u \end{array}) \end{aligned}

　

　（並進運動の方程式）

航空機の運動方程式の導出（回転運動方程式）

回転運動のモーメントを考える。
微小質量

∆ m

にかかるモーメントは、

\begin{aligned} r \times ∆ F & = r \times ∆ m (\frac{d v_{c}}{d t} + ω \times v_{c} + \frac{d ω}{d t} \times r + ω \times (ω \times r)) \\ = ∆ m (r \times (\frac{d v_{c}}{d t} + ω \times v_{c}) + r \times (\frac{d ω}{d t} \times r) + r \times (ω \times (ω \times r))) \end{aligned}

行列の外積の計算の仕方　

[A \times (B \times C) = (A ・ C) ・ B - (A ・ B) ・ C]

　なので

\begin{aligned} r \times ∆ F & = ∆ m (r \times (\frac{d v_{c}}{d t} + ω \times v_{c}) + (r ・ r) \frac{d ω}{d t} - (r ・ \frac{d ω}{d t}) r + r \times ((ω ・ r) ω - (ω ・ ω) r)) \\ = ∆ m (r \times (\frac{d v_{c}}{d t} + ω \times v_{c}) + (r ・ r) \frac{d ω}{d t} - (r ・ \frac{d ω}{d t}) r + (ω ・ r) ・ (ω \times r) - (ω ・ ω) ・ (r \times r) \end{aligned}

r \times r = 0

　なので、

r \times ∆ F = ∆ m (r \times (\frac{d v_{c}}{d t} + ω \times v_{c}) + (r ・ r) \frac{d ω}{d t} - (r ・ \frac{d ω}{d t}) r + (ω ・ r) (ω \times r))

積分して、機体全体の質量にして考えると

\begin{aligned} r \times F & = \int d m (r \times (\frac{d v_{c}}{d t} + ω \times v_{c}) + (r ・ r) \frac{d ω}{d t} - (r ・ \frac{d ω}{d t}) r + (ω ・ r) (ω \times r)) \\ = \int r d m \times (\frac{d v_{c}}{d t} + ω \times v_{c}) + \int d m ((r ・ r) \frac{d ω}{d t} - (r ・ \frac{d ω}{d t}) r + (ω ・ r) (ω \times r)) \\ = 0 + \int d m ((r ・ r) \frac{d ω}{d t} - (r ・ \frac{d ω}{d t}) r + (ω ・ r) (ω \times r)) \\ = \int d m ((r ・ r) \frac{d ω}{d t} - (r ・ \frac{d ω}{d t}) r + (ω ・ r) (ω \times r)) \end{aligned}

成分表示をすると

(\begin{matrix} \bar{L} \\ \bar{M} \\ \bar{N} \end{matrix}) = \int d m {(\begin{matrix} x \\ y \\ z \end{matrix}) ･ (\begin{matrix} x \\ y \\ z \end{matrix})} \frac{d}{d t} (\begin{matrix} p \\ q \\ r \end{matrix}) - {(\begin{matrix} x \\ y \\ z \end{matrix}) ･ \frac{d}{d t} (\begin{matrix} p \\ q \\ r \end{matrix})} (\begin{matrix} x \\ y \\ z \end{matrix}) + {(\begin{matrix} p \\ q \\ r \end{matrix}) ･ (\begin{matrix} x \\ y \\ z \end{matrix})} {(\begin{matrix} p \\ q \\ r \end{matrix}) \times (\begin{matrix} x \\ y \\ z \end{matrix})}

= \int d m {(x^{2} + y^{2} + z^{2}) (\begin{matrix} \dot{p} \\ \dot{q} \\ \dot{r} \end{matrix}) - (x \dot{p} + y \dot{q} + z \dot{r}) (\begin{matrix} x \\ y \\ z \end{matrix}) + (p x + q y + r z) (\begin{matrix} y r - z q \\ z p - x r \\ x q - y p \end{matrix})}

= \int d m {(\begin{matrix} (x^{2} + y^{2} + z^{2}) \dot{p} \\ (x^{2} + y^{2} + z^{2}) \dot{q} \\ (x^{2} + y^{2} + z^{2}) \dot{r} \end{matrix}) - (\begin{matrix} x^{2} \dot{p} + x y \dot{q} + z x \dot{r} \\ x y \dot{p} + y^{2} \dot{q} + y z \dot{r} \\ z x \dot{p} + y z \dot{q} + z^{2} \dot{r} \end{matrix}) + (\begin{matrix} x y r p - z x p q + y^{2} q r - y z q^{2} + y z r^{2} - z^{2} q r \\ z x p^{2} - x^{2} r p + y z p q - x y q r + z^{2} r p - z x r^{2} \\ x^{2} p q - x y p^{2} + x y q^{2} - y^{2} p q + z x q r - y z r p \end{matrix})}

= \int d m (\begin{matrix} (y^{2} + z^{2}) \dot{p} + (y^{2} - z^{2}) q r - x y (r p - \dot{q}) - y z (q^{2} - r^{2}) - z x (\dot{r} + p q) \\ (z^{2} + x^{2}) \dot{q} + (z^{2} - x^{2}) r p - x y (q r + \dot{p}) - y z (\dot{r} - p q) - z x (r^{2} - p^{2}) \\ (x^{2} + y^{2}) \dot{r} + (x^{2} - y^{2}) p q - x y (p^{2} - q^{2}) - y z (\dot{q} + r p) - z x (\dot{p} - q r) \end{matrix})

= (\begin{matrix} \int (y^{2} + z^{2}) d m \dot{p} + \int (y^{2} - z^{2}) d m q r - \int x y d m (r p - \dot{q}) - \int y z d m (q^{2} - r^{2}) - \int z x d m (\dot{r} + p q) \\ \int (z^{2} + x^{2}) d m \dot{q} + \int (z^{2} - x^{2}) d m r p - \int x y d m (q r + \dot{p}) - \int y z d m (\dot{r} - p q) - \int z x d m (r^{2} - p^{2}) \\ \int (x^{2} + y^{2}) d m \dot{r} + \int (x^{2} - y^{2}) d m p q - \int x y d m (p^{2} - q^{2}) - \int y z d m (\dot{q} + r p) - \int z x d m (\dot{p} - q r) \end{matrix})

ここで、　

{\begin{cases} \int (y^{2} + z^{2}) d m = I_{x} ， \int (z^{2} + x^{2}) d m = I_{y} ， \int (x^{2} + y^{2}) d m = I_{z} \\ \int x y d m = I_{x y} = I_{y x} ， \int y z d m = I_{y z} = I_{z y} ， \int z x d m = I_{z x} = I_{x z} \end{cases}

とすると、

(\begin{matrix} \bar{L} \\ \bar{M} \\ \bar{N} \end{matrix}) = (\begin{matrix} I_{x} \dot{p} + (I_{y} - I_{z}) q r - I_{x y} (r p - \dot{q}) - I_{y z} (q^{2} - r^{2}) - I_{z x} (\dot{r} + p q) \\ I_{y} \dot{q} + (I_{z} - I_{x}) r p - I_{x y} (q r + \dot{p}) - I_{y z} (\dot{r} - p q) - I_{z x} (r^{2} - p^{2}) \\ I_{z} \dot{r} + (I_{x} - I_{y}) p q - I_{x y} (p^{2} - q^{2}) - I_{y z} (\dot{q} + r p) - I_{z x} (\dot{p} - q r) \end{matrix})

　　　（回転運動の方程式）

左右対称な機体の場合の回転運動方程式

もしも機体が　左右対称形　であった場合　

I_{x y} = I_{y z} = 0

　となる。
そのため、回転運動の方程式は

(\begin{matrix} \bar{L} \\ \bar{M} \\ \bar{N} \end{matrix}) = (\begin{matrix} I_{x} \dot{p} + (I_{y} - I_{z}) q r - I_{z x} (\dot{r} + p q) \\ I_{y} \dot{q} + (I_{z} - I_{x}) r p - I_{z x} (r^{2} - p^{2}) \\ I_{z} \dot{r} + (I_{x} - I_{y}) p q_{z x} (\dot{p} - q r) \end{matrix})

この時、
　

\bar{M}

　の第２項の　

(I_{z} - I_{x}) r p

　を「ロール・ヨー・カップリング」と言い、
　

\bar{N}

　の第２項の　

(I_{x} - I_{y}) p q

　を「ロール・ピッチ・カップリング」と言う。

「ロール・ヨー・カップリング」
　ロール運動中にヨー運動が加わるとピッチングモーメントが発生する。

「ロール・ピッチ・カップリング」
　ロール運動中にピッチ運動が加わるとヨーイングモーメントが発生する

外力と（並進運動方程式の左辺）について

航空機に働く外力は　（重力）と（推力）と（揚力）と（抗力）です。
（揚力）と（抗力）は、合わせて（空気力）と言うことができます。

つまり、航空機に働く外力は　（重力）と（推力）と（空気力）です。

１　重力

W = (- m g \sin θ, m g \cos θ \sin ϕ, m g \cos θ \cos ϕ)

２　推力
　　

T = (T \cos i_{T}, 0, - T \sin i_{T})

３　空気力
　　

(\frac{1}{2} C_{x} ρ S V^{2}, \frac{1}{2} C_{y} ρ S V^{2}, \frac{1}{2} C_{z} ρ S V^{2})

よって、外力Fは

(\begin{matrix} F_{x} \\ F_{y} \\ F_{z} \end{matrix}) = (\begin{matrix} - m g \sin θ + T \cos i_{T} + \frac{1}{2} C_{x} ρ S V^{2} \\ m g \cos θ \sin ϕ + \frac{1}{2} C_{y} ρ S V^{2} \\ m g \cos θ \cos ϕ - T \sin i_{T} + \frac{1}{2} C_{z} ρ S V^{2} \end{matrix})

なので、航空機の並進運動の運動方程式
　

(\begin{matrix} F_{x} \\ F_{y} \\ F_{z} \end{matrix}) = m (\begin{matrix} \dot{u} + q w - r v \\ \dot{v} + r u - p w \\ \dot{w} + p v - q u \end{matrix})

は、

(\begin{matrix} - m g \sin θ + T \cos i_{T} + \frac{1}{2} C_{x} ρ S V^{2} \\ m g \cos θ \sin ϕ + \frac{1}{2} C_{y} ρ S V^{2} \\ m g \cos θ \cos ϕ - T \sin i_{T} + \frac{1}{2} C_{z} ρ S V^{2} \end{matrix}) = m (\begin{matrix} \dot{u} + q w - r v \\ \dot{v} + r u - p w \\ \dot{w} + p v - q u \end{matrix})

と書くことができ、これを整理すると

(\begin{matrix} \dot{u} \\ \dot{v} \\ \dot{w} \end{matrix}) = m (\begin{matrix} - q w + r v - g \sin θ + \frac{T}{m} \cos i_{T} + \frac{1}{2 m} C_{x} ρ S V^{2} \\ - r u + p w + g \cos θ \sin ϕ + \frac{1}{2 m} C_{y} ρ S V^{2} \\ - p v + q u + g \cos θ \cos ϕ - \frac{T}{m} \sin i_{T} + \frac{1}{2 m} C_{z} ρ S V^{2} \end{matrix})

　　　これも（並進運動方程式）

モーメント（回転運動方程式の左辺）について

航空機に働くモーメントは、
（空気力によるモーメント）と（エンジンジャイロモーメント）に分けられます。

１　空気力によるモーメント
　　平均空力翼弦を　

\bar{c}

　、翼幅を　

b

　、空力モーメント係数を

(C_{l}, C_{m}, C_{n})

　とすると、

（空気力によるモーメント）は　

(\frac{1}{2} C_{l} ρ S V^{2} b, \frac{1}{2} C_{m} ρ S V^{2} \bar{c}, \frac{1}{2} C_{n} ρ S V^{2} b)

２　エンジンジャイロモーメント
　　エンジンによる慣性モーメントを

I_{R}

、エンジンによる回転角速度を

ω_{R}

とすると、

（エンジンジャイロモーメント）は　

(0, - I_{R} ω_{R} ･ r, I_{R} ω_{R} ･ q)

航空機に働くモーメントは

(\begin{matrix} \bar{L} \\ \bar{M} \\ \bar{N} \end{matrix}) = (\begin{matrix} \frac{1}{2} C_{l} ρ S V^{2} b \\ \frac{1}{2} C_{m} ρ S V^{2} \bar{c} - I_{R} ω_{R} ･ r \\ \frac{1}{2} C_{n} ρ S V^{2} b + I_{R} ω_{R} ･ q \end{matrix})

なので、

航空機の回転運動の運動方程式

(\begin{matrix} \bar{L} \\ \bar{M} \\ \bar{N} \end{matrix}) = (\begin{matrix} I_{x} \dot{p} - (I_{y} - I_{z}) q r - I_{x y} (r p - \dot{q}) - I_{y z} (q^{2} - r^{2}) - I_{z x} (\dot{r} + p q) \\ I_{y} \dot{q} - (I_{z} - I_{x}) r p - I_{x y} (q r - \dot{p}) - I_{y z} (\dot{r} + p q) - I_{z x} (r^{2} - p^{2}) \\ I_{z} \dot{r} - (I_{x} - I_{y}) p q - I_{x y} (p^{2} - q^{2}) - I_{y z} (\dot{q} + r p) - I_{z x} (\dot{p} - q r) \end{matrix})

は、

(\begin{matrix} \frac{1}{2} C_{l} ρ S V^{2} b \\ \frac{1}{2} C_{m} ρ S V^{2} \bar{c} - I_{R} ω_{R} ･ r \\ \frac{1}{2} C_{n} ρ S V^{2} b + I_{R} ω_{R} ･ q \end{matrix}) = (\begin{matrix} I_{x} \dot{p} - (I_{y} - I_{z}) q r - I_{x y} (r p - \dot{q}) - I_{y z} (q^{2} - r^{2}) - I_{z x} (\dot{r} + p q) \\ I_{y} \dot{q} - (I_{z} - I_{x}) r p - I_{x y} (q r - \dot{p}) - I_{y z} (\dot{r} + p q) - I_{z x} (r^{2} - p^{2}) \\ I_{z} \dot{r} - (I_{x} - I_{y}) p q - I_{x y} (p^{2} - q^{2}) - I_{y z} (\dot{q} + r p) - I_{z x} (\dot{p} - q r) \end{matrix})

機体が左右対称であれば、

I_{x y} = I_{y z} = 0

　なので、

(\begin{matrix} \frac{1}{2} C_{l} ρ S V^{2} b \\ \frac{1}{2} C_{m} ρ S V^{2} \bar{c} - I_{R} ω_{R} ･ r \\ \frac{1}{2} C_{n} ρ S V^{2} b + I_{R} ω_{R} ･ q \end{matrix}) = (\begin{matrix} I_{x} \dot{p} + (I_{y} - I_{z}) q r - I_{z x} (\dot{r} + p q) \\ I_{y} \dot{q} + (I_{z} - I_{x}) r p - I_{z x} (r^{2} - p^{2}) \\ I_{z} \dot{r} + (I_{x} - I_{y}) p q_{z x} (\dot{p} - q r) \end{matrix})

と書くことができ、これを整理すると

(\begin{matrix} I_{x} \dot{p} \\ I_{y} \dot{q} \\ I_{z} \dot{r} \end{matrix}) = (\begin{matrix} (I_{y} - I_{z}) q r + I_{z x} (\dot{r} + p q) + \frac{1}{2} C_{l} ρ S V^{2} b \\ (I_{z} - I_{x}) r p + I_{z x} (r^{2} - p^{2}) + \frac{1}{2} C_{m} ρ S V^{2} \bar{c} - I_{R} ω_{R} ･ r \\ (I_{x} - I_{t}) p q + I_{z x} (\dot{p} - q r) + \frac{1}{2} C_{n} ρ S V^{2} b + I_{R} ω_{R} ･ q \end{matrix})

ここで、

(\begin{matrix} (I_{y} - I_{z}) q r + I_{z x} p q + \frac{1}{2} C_{l} ρ S V^{2} b \\ (I_{z} - I_{x}) r p + I_{z x} (r^{2} - p^{2}) + \frac{1}{2} C_{m} ρ S V^{2} \bar{c} - I_{R} ω_{R} ･ r \\ (I_{x} - I_{t}) p q - I_{z x} q r + \frac{1}{2} C_{n} ρ S V^{2} b + I_{R} ω_{R} ･ q \end{matrix}) = (\begin{matrix} L \\ M \\ N \end{matrix})

とおくと、

　

(\begin{matrix} I_{x} \dot{p} \\ I_{y} \dot{q} \\ I_{z} \dot{r} \end{matrix}) = (\begin{matrix} L + I_{z x} \dot{r} \\ M \\ N + I_{z x} \dot{p} \end{matrix})

さらに整理すると、

　

(\begin{matrix} \dot{p} \\ \dot{q} \\ \dot{r} \end{matrix}) = (\begin{matrix} \frac{\frac{L}{I_{x}} + \frac{I_{z x}}{I_{x}} ･ \frac{N}{I_{z}}}{1 - \frac{{I_{z x}}^{2}}{I_{x} I_{z}}} \\ \frac{M}{I_{y}} \\ \frac{\frac{N}{I_{z}} + \frac{I_{z x}}{I_{z}} ･ \frac{L}{I_{x}}}{1 - \frac{{I_{z x}}^{2}}{I_{x} I_{z}}} \end{matrix})

　これも（回転運動方程式）