はじめに

はじめまして、りょーつと申します。大学院でメカの研究をしています。
先日トランスミッションの操作ができるラジコンを作りたくなったので、それについて書いていこうと思います。実機の動作まで確認できたのですが、今回はプログラムを書くために必要な数式をひとまずメモしていきます。そのうち図などを追加して分かりやすくしていきます！
実装とかについてはまた別の記事で紹介します。

自動車のモデル化

今回はMT車のモデルを作りたいと思います。ブリッピングとかエンジンブレーキも再現したいので、エンジンの回転数と車体の移動速度をそれぞれ計算しないといけません。ここではギアがつながってる場合とそうでない場合について場合分けしてモデルを作っていきます。ある程度ノリと勢いで作りましたが、それっぽくは動いているので大きく間違ってはいないと思います。
実機の動作についてはこちらのYouTubeチャンネルで公開中です。(https://www.youtube.com/@user-jm7lf1wn9o/featured)

ギアがつながってない場合のモデル

ギアがつながっていない場合は、エンジンの回転と車体はそれぞれ独立に動作するため、式を別々に立てることができます。

エンジンのモデル化

ここではエンジンのモデル化をやっていきます。とは言ってもエンジンの運動は吸気とか排気とかあって難しそうなので、モータに置き換えて考えていきます。つまり、エンジンが図1に示すような等価回路で書くことができるという仮定を置きます。

図1 モータの等価回路

まずは回路方程式を導出していきます。アクセルの入力によって生じる電源電圧を$u$ [V]、エンジンの内部抵抗を$R$ [$\Omega$]、逆起電力を$e$ [V]、回路を流れる電流を$i$ [A]とするとキルヒホッフの第二法則より、
$$u = Ri + e \tag{1}$$
がなりたちます。ここで逆起電力定数$k_e$、エンジンの回転数を$\omega$ [rad/s]とすれば
$$e = k_e \omega \tag{2}$$
となる。逆起電力定数$k_e$はモータにおける定格電圧$u_0$[V]と最高回転数$\omega_0$[rad/s]の比として与えられる定数です。ここではアクセルを踏み込んだ場合の入力を$u_0 = 1$、エンジンの最高回転数[rad/s]を$\omega_0$として計算します。
エンジンが発生させるトルク$\tau$ [Nm]はモータと同様、電流$i$に比例するものとして
$$\tau = k_t i \tag{3}$$
とします。ここで$k_t$はトルク定数とよびます(トルク定数の計算方法はあとで示します)。エンジンの力学的なモデルは、慣性モーメント$I$ [$kgm^2$]の回転軸がエンジンのトルク$\tau$によって回転されるだけのシンプルなものを考えます。この時運動方程式は
$$I\dfrac{d\omega}{dt} = \tau \tag{4}$$
というシンプルな形で記述されます。

図2　エンジンのモデル

(1)(2)式より$e$を消去し、$i$について解くと
$$i = \dfrac{u}{R} - \dfrac{k_e}{R}\omega \tag{5}$$
となり、(3)(4)(5)式から$i$と$\tau$を消去することで
$$\dfrac{d\omega}{dt} = -\dfrac{k_e k_t}{RI}\omega + \dfrac{k_t}{RI}u \tag{6}$$
が得られます。(6)式が求めたかったエンジンのモデルです。しかし未知の定数$k_t$と$R$と$I$があるため、それらを求める方法を示しておきます。モータのモデルにおいて、トルク定数$k_t$はモータの最大電流$i_{max}$と最大トルク$\tau_{max}$の比で表されます。また、内部抵抗$R$は最大電流$i_{max}$と定格電圧$u_0$の比で表されます。そのため、トルク定位数$k_t$と内部抵抗$R$の比は
$$\dfrac{k_t}{R} = \dfrac{\dfrac{\tau_{max}}{i_{max}}}{\dfrac{u_0}{i_{max}}} = \dfrac{\tau_{max}}{u_0} \tag{7}$$
で表されます。$\tau_{max}$をエンジンの最高トルク、アクセルを踏み込んだ場合の入力を$u_0 = 1$とすれば$k_t$と$R$をエンジンの諸定数で表すことができます。エンジンの慣性モーメント$I$については、計算方法を思いつかなかったので、フライホイールの慣性モーメントで表現することにしました。フライホイールを質量$m_f$ [kg]、直径$D_f$ [kg]の円板として近似すると、エンジンの慣性モーメント$I$は
$$I = \dfrac{1}{8}m_f {D_f}^2 \tag{8}$$
で計算できます。

図3　円柱で置き換えたフライホイール

車体のモデル化

ここでは、タイヤと車の運動方程式を考えていきます。まずはタイヤについて考えます。タイヤの慣性モーメントを$J_w$ [$kg m^2$]、地面を転がることで生じる転がり摩擦トルクを$T_f$ [Nm]、ブレーキをかけることで生じる動摩擦トルクを$T_b$ [Nm]、タイヤの回転速度を$\Omega$ [rad/s]とします。タイヤの慣性モーメント$J_w$はタイヤを質量が$m_w$ [kg]、直径$D_w$ [m]の円柱であると仮定すれば、タイヤが4本ついていることを考慮して
$$J_w = 4\times \dfrac{1}{8}m_w {D_w}^2 \tag{9}$$
とします。ほんとはゴムの部分とホイールの部分で密度が違いますが、そんなことはいったん無視してます。

図4　タイヤの慣性モーメント
今はギアがつながっていない状態を考えているため、タイヤ単体で考えた場合の運動方程式は
$$J_w \dfrac{d\Omega}{dt} = -T_b - T_f \tag{10}$$
となります。

図5　タイヤの運動方程式

次に車体の運動方程式を求めます。地面に角度$\theta$ [rad]の傾斜が付いているとすると、質量$M$ [kg]の車体が重力加速度$g$ [$kgm/s^2$]の空間内を速度$V$ [m/s]で走行する場合、
$$M \dfrac{dV}{dt} = -Mg\sin{\theta} \tag{11}$$
になります。タイヤと地面の間に滑りがないものとすれば
$$V = r\Omega \tag{12}$$
が成り立つので、(11)(12)式から$V$を消去することで
$$Mr \dfrac{d\Omega}{dt} = -Mg\sin{\theta} \tag{13}$$
が得られます。

図6　車体の運動方程式
さらに並進移動の等価慣性モーメント(参考: http://www.mekatoro.net/digianaecatalog/fujih-sougou/Book/fujih-sougou-P0230.pdf) (1)というものを作りたいので、さらに(13)式の両辺を$r$倍してやると
$$Mr^2 \dfrac{d\Omega}{dt} = -Mgr\sin{\theta} \tag{14}$$
となり、(10)式と(14)式の次元が一致しました。したがって(10)式と(14)式を足してやることで
$$(Mr^2 + J_w) \dfrac{d\Omega}{dt} = -Mgr\sin{\theta} -T_b - T_f \tag{15}$$
という式が得られます。なお、等価慣性モーメント$J$ [$kg m^2$]というのは(15)式の左辺の係数に当たる部分で、
$$J = Mr^2 + J_w \tag{16}$$
です。最後に(15)式へ(16)式を代入し、$J$で除してやると
$$\dfrac{d\Omega}{dt} = -\dfrac{Mgr\sin{\theta} + T_b + T_f}{J} \tag{17}$$
が得られます。(17)式が求めたかった車体(厳密にはタイヤ)のモデルです。

ギアがつながっている場合のモデル

ギアがつながっている場合は、エンジンの回転と車体の運動が連動するため、それらを考慮した式を作る必要があります。基本的にはギアがつながっていない場合のモデルと同様の手順で導出していきます。まずエンジンはこれまでと同様にモータの等価回路を用いてモデル化します。そのため(1)(2)(3)(5)の回路方程式はここでも使うことができます。また車体単体の運動方程式は(14)式になります。そのためここではエンジンとタイヤの運動方程式の導出をメインに書いていきます。
エンジンのトルク$\tau$ [Nm]とそれによりタイヤに発生する駆動トルク$T$ [Nm]の間には、トランスミッションのギア比を$k_G$としたとき、
$$T = k_G \tau \tag{18}$$
の関係があります。これより、タイヤは駆動トルク$T$と摩擦トルク$T_b$、$T_f$を受けて加減速するため、その運動方程式は
$$J_w \dfrac{d\Omega}{dt} = T -T_b - T_f \tag{19}$$
になります。実際にはエンジンの慣性モーメント$I$ [$kg m^2$]も考慮したほうがいいですが、等価慣性モーメントに比べて十分小さいので無視しました。(19)式に(14)式の両辺を足して、(16)(18)式を代入すると
$$J \dfrac{d\Omega}{dt} = k_G \tau + Mgr\sin{\theta} -T_b - T_f \tag{20}$$
が得られます。ここで(3)(5)の回路方程式から$i$を消去すると
$$\tau = \dfrac{k_t}{R}(u - k_e \omega) \tag{21}$$
になるので、(20)(21)式から$\tau$を消去し、両辺を等価慣性モーメント$J$で割ることで
$$\dfrac{d\Omega}{dt} = \dfrac{\dfrac{k_t}{R}(u - k_e \omega) k_G + Mgr\sin{\theta} -T_b - T_f}{J} \tag{22}$$
が求まります。この(22)式がギアのつながった状態を表すモデルとなります。なお、(22)式からエンジンの回転数を求めるためにはギア比$k_G$を考慮して
$$\omega = k_G \Omega \tag{23}$$
を用いればよいです。また、車速は(12)式を使って計算できます。

ギアチェンジのモデル

ギアがつながっていないときは、エンジン回転数$\omega$とタイヤ回転数$\Omega$は(6)式と(17)式で別々に計算されます。そのためギアをつなげるとき、(23)式をみたすように、別々に計算されている回転速度の辻褄を合わせるための計算式を作る必要があります。(ギアが離れるときは特に問題にはなりません)
ここでは角運動量保存則を用いて計算します。また、速度の辻褄合わせはすべてクラッチで発生しているものと考えます。角運動量保存則は、外力が発生しない場合、慣性モーメントを$J_i$、角速度を$\omega_i$としたときに
$$\sum J_i \omega_i = const \tag{24}$$
が満たされることを示した法則です。今回の例ではエンジンの角運動量とタイヤの角運動量の総和が一定になるようにすればいいのですが、トランスミッションのギア比$k_G$を考慮する必要があることに注意します。今回はクラッチに近いエンジン側で角運動量保存則を考えるので、タイヤ側の慣性モーメントはギアがある場合の等価慣性モーメント(参考: http://www.mekatoro.net/digianaecatalog/fujih-sougou/Book/fujih-sougou-P0230.pdf) (1)を使用する必要があります。タイヤのエンジン側からみた等価慣性モーメント$J'$は
$$J' = \dfrac{J}{{k_G}^2} \tag{25}$$
です。またタイヤの回転速度$\Omega$はエンジン側ではギア比$k_G$を考慮して
$$\Omega' = k_G \Omega \tag{26}$$
に変換されます。ギアがつながった後のエンジン回転速度を$\omega'$とすると、ギアがつながる前後の角運動量保存則により、
$$J' \Omega' + I \omega = (I + J')\omega' \tag{27}$$
がなりたちます。

図7　角運動量保存則

(27)式に(25)(26)式を代入し、$\omega'$について解くとギアチェンジ直後のエンジン回転速度$\omega'$が
$$\omega' = \dfrac{\dfrac{J}{{k_G}^2} k_G \Omega + I \omega}{I + \dfrac{J}{{k_G}^2}} = \dfrac{(I k_G \omega + J\Omega)k_G}{I {k_G}^2 + J} \tag{27}$$
によって得られます。実際には半クラッチを使って徐々に速度が変化していくのですが、摩擦のモデル化がややこしかったので、今回は半クラッチを使わずに走行すると仮定しました。

モデルのまとめ

長々と式を書いてしまいましたがここらでいったんまとめておこうと思います。
まずギアがニュートラルの場合、運動を記述する微分方程式は
$$\dfrac{d\omega}{dt} = -\dfrac{k_e k_t}{RI}\omega + \dfrac{k_t}{RI}u \tag{28}$$
$$\dfrac{d\Omega}{dt} = -\dfrac{Mgr\sin{\theta} + T_b + T_f}{J} \tag{29}$$
で与えられます。ギアがつながっている場合は
$$\dfrac{d\Omega}{dt} = \dfrac{\dfrac{k_t}{R}(u - k_e \omega) k_G + Mgr\sin{\theta} -T_b - T_f}{J} \tag{30}$$
で与えられます。実際の処理ではこれらを数値積分することで車っぽい動作ができるようになります。さらにギアがつながった直後の速度は
$$\omega' = \dfrac{(I k_G \omega + J\Omega)k_G}{I {k_G}^2 + J} \tag{31}$$
で計算できます。

計算例

車のパラメータを(28)(29)(30)式に代入する方法について具体的な例をもとに説明していきます。今回はマツダのRX-8の特性値をシミュレーションに反映していきます。そのためにネットのいろんなサイトをめぐって、RX-8の特性値を集めました。

エンジン最高回転数：9000 rpm
エンジン最高トルク：22 kgf m
フライホイール質量：約7.6 kg
フライホイール直径：約350 mm
タイヤ質量：約21 kg
タイヤ直径：約19 インチ
車体質量：約1350 kg

これらをもとに必要な定数を計算していきます。有効数字はいったんめんどうなので考慮しません。まずエンジン最高トルク$\tau_{max}$ [Nm]は単位変換により
$$\tau_{max} = 22 kgf m = 22\times 9.81 Nm = 215 Nm \tag{32}$$
になります。またエンジン最高回転数$\omega_0$ [rad/s]は単位変換により
$$\omega_0 = 9000rpm = 9000 \times \dfrac{2\pi}{60} rad/s = 942 rad/s \tag{33}$$
になります。アクセルの最大入力を$u_0 = 1$にしたので、逆起電力定数$k_e$は
$$k_e = \dfrac{1}{\omega_{max}} = \dfrac{1}{942} = 0.00106 \tag{34}$$
になります。また(7)式に代入することで
$$\dfrac{k_t}{R} = \dfrac{\tau_{max}}{u_0} = \dfrac{215}{1} = 215 \tag{35}$$
が得られます。あとはフライホイールとタイヤの慣性モーメントを求めるだけです。まずは単位変換をします。フライホイールの直径$D_f$とタイヤの直径$D_w$ [m]とタイヤの半径$r$ [m]は
$$D_f = 350 mm = 0.35 m \tag{36}$$$$D_w = 19 インチ = 19 \times 0.0254 m = 0.483 m \tag{37}$$$$r = \dfrac{D_w}{2} = 0.241 m \tag{38}$$
になります。あとは(8)(9)(16)式に必要な値を代入するだけでオッケーです。
$$I = \dfrac{1}{8}m_f {D_f}^2 = \dfrac{1}{8}\times 7.6 \times {0.35}^2 = 0.116 \tag{39}$$$$J_w = 4\times \dfrac{1}{8}m_w {D_w}^2 = 4\times \dfrac{1}{8}\times 21 \times {0.483}^2 = 2.45 \tag{40}$$$$J = Mr^2 + J_w = 1320\times {0.241}^2 + 2.45 = 79.1 \tag{41}$$
これで(28)~(31)式に出てくるパラメータの多くが明らかになりました。まだあと分かっていないのはブレーキの摩擦トルク$T_b$ [Nm]と転がり摩擦トルク$T_f$ [Nm]ですが、これらは摩擦なので値を見つけることができませんでした。そのため、これら2つはいい感じに操作しながら試行錯誤で設定していく必要があります。
ギア比$k_G$は適宜切り替わるのでその都度その都度代入する必要があります。ちなみにRX-8のギア比はこんな感じだそうです。

1速：15.5
2速：9.0
3速：6.2
4速：4.4
5速：3.6
R　：14.6

さいごに

この記事では理論に関することが中心となりましたが、今後は実装とか実機の作り方についても紹介していきたいと思います。

参考文献

(1) $GD^2$と慣性モーメント、http://www.mekatoro.net/digianaecatalog/fujih-sougou/Book/fujih-sougou-P0230.pdf