二輪車の運動特性(5) 旋回挙動の固有モードとの関係

二輪車の直進から旋回への移行は、操舵トルクのステップ入力によって実現可能であり、その挙動がウイーブモードとキャプサイズモードの安定性によって支配されていることを前回の記事で示しました。

全体の挙動は固有モードごとの応答の重ね合わせによって表されることから、モード展開することによって各々の固有モードの応答に分解することで寄与率を求めることができます。

今回は各固有モードがどのように影響を及ぼしているのか、ステップ応答計算をモード展開することによって明らかにしてみたいと思います。

モード展開によるステップ応答の計算

これまで運動方程式の解を求めるのに、外力がない場合について微分方程式の解法を紹介してきましたが、今回はステップ応答を求めたいため、外力を考慮した微分方程式の解法が必要になります。

一般に状態方程式は

$$\dot{\bm{x}}(t) = \bold{A}\bm{x}(t) + \bold{B}\bm{u}(t)$$

で表され、このシステムの自由応答はこれまでにご紹介したように、$\bm{u}(t)=\bm{0}$としてシステム行列$\bold{A}$の固有値・固有ベクトルから求めることができます。このような形式の微分方程式は同次方程式と呼ばれています。

一方、今回求めたいステップ応答のように$\bm{u}(t) \ne \bm{0}$の場合は非同次方程式と呼ばれ、解法は複雑となって、新たに定数変化法などの手法を適用する必要があります。その詳細は教科書に譲ってここではその結果だけを使わせてもらうことにします^[1][2]。

状態方程式のステップ応答

状態方程式の一般解は

$$\bm{x}(t) = e^{\bold{A}t}\bm{x}(0) + \int_0^t e^{\bold{A}(t-\tau)}\bold{B}\bm{u}(\tau)d\tau$$

で表されますが、入力をステップとして時間によらず一定とし、$\bm{u}(t) = \bold{U_c} = [u_{\phi c}\, u_{\delta c}]^T$、また初期条件として$\bm{x}(0) = [0~~0]^T$として式展開していきます。

$s=t-\tau$とおくと、$ds=-d\tau,\; \tau=0$のとき$s=t,\; \tau=t$のとき$s=0$を代入すると、

$$\begin{align*} \bm{x}(t) &= \int_t^0 -e^{\bold{A}s}\bold{BU_c}ds \\ &= \int_0^t e^{\bold{A}s}\bold{BU_c}ds \end{align*}$$

行列$\bold{A}$の固有値$\lambda_1, \lambda_2, \ldots, \lambda_n$, 固有ベクトル$\bm{v}_1, \bm{v}_2, \ldots, \bm{v}_n$を既知とし、

$$\begin{align*} \bold{\Lambda} &= \begin{pmatrix} \lambda_1 & 0 & \ldots & 0 \\ 0 & \lambda_2 & \ldots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \dots & \lambda_n \end{pmatrix} \\ \bold{V} &= [\bm{v}_1\; \bm{v}_2\; \ldots\; \bm{v}_n] \end{align*}$$

とすると、

$$\begin{align*} \bold{A} &= \bold{V \Lambda V}^{-1} \\ e^{\bold{A}t} &= \bold{V}e^{\bold{\Lambda}t}\bold{V}^{-1} = \bold{V} \begin{pmatrix} e^{\lambda_1 t} & 0 & \ldots & 0 \\ 0 & e^{\lambda_2 t} & \ldots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \ldots & e^{\lambda_n t} \end{pmatrix} \bold{V}^{-1} \end{align*}$$

となり、状態方程式の解は次のように表せます。

$$\begin{align*} \bm{x}(t) &=\bold{V}\left( \int_0^t e^{\bold{\Lambda}s} ds \right) \bold{V}^{-1}\bold{BU_c} \\ &= \bold{V} \begin{pmatrix} \int_0^t e^{\lambda_1 s}ds & 0 & \ldots & 0 \\ 0 & \int_0^t e^{\lambda_2 s}ds & \ldots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \ldots & \int_0^t e^{\lambda_n s}ds \end{pmatrix} \bold{V}^{-1}\bold{BU_c} \end{align*}$$

ここで、$\lambda_i \ne 0$とすると、

$$\int_0^t e^{\lambda_i s}ds = \left[ \frac{e^{\lambda_i s}}{\lambda_i} \right]_0^t = \frac{e^{\lambda_i t} - 1}{\lambda_i}$$

$$\therefore \bm{x}(t) = \bold{V} \begin{pmatrix} \frac{e^{\lambda_1 t} - 1}{\lambda_1} & 0 & \ldots & 0 \\ 0 & \frac{e^{\lambda_2 t} - 1}{\lambda_2} & \ldots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \ldots & \frac{e^{\lambda_n t} - 1}{\lambda_n} \end{pmatrix} \bold{V}^{-1}\bold{BU_c}$$

となり、モードごとに独立な時間関数の重ね合わせとしてステップ応答についてのモード展開が求まります。

※ $\lambda_i=0$の場合は、$\int_0^t e^{\lambda_i s}ds = t$

ちなみに、$\lambda_i$の実部が負であれば、$t\rightarrow \infin$のとき$e^{\lambda_i} \rightarrow 0$となり、$x$は定常値に収束して、

$$\begin{align*} \bm{x}(\infin) &= \bold{V} \begin{pmatrix} -\frac{1}{\lambda_1} & 0 & \ldots & 0 \\ 0 & -\frac{1}{\lambda_2} & \ldots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \ldots & -\frac{1}{\lambda_n} \end{pmatrix} \bold{V}^{-1}\bm{BU_c} \\ &= -\bold{V\Lambda}^{-1}\bold{V}^{-1}\bold{BU_c} \\ &= -\bold{A}^{-1}\bold{BU_c} \end{align*}$$

となりますが、これは状態方程式の$\dot{\bm{x}}=0$とした定常時の釣合式$\bold{A}\bm{x}+\bold{B}\bm{u}=0$と矛盾しないことがわかります。

ベンチマーク自転車のステップ応答モード展開

車速$v=5$m/sで走行中のベンチマーク自転車に、時刻$t=0$で操舵トルクとして$T_\delta=-0.2$Nmのステップ入力された時のロール角と舵角の応答を、固有モードごとにプロットすると次のようになります。

青線がキャスタリングモード、緑線がウィーブモード、赤線がキャプサイズモードの時間応答を表し、水色は全ての合計です。すべてのモードの時間応答を合計した水色の応答は前回の記事で示した数値積分の結果と一致しています。

入力開始直後に表れる逆操舵（舵角が旋回とは逆を向く）挙動はウイーブモードによって引き起こされていますが、全体の挙動はキャプサイズモードによって支配されていることがわかります。

また、定常状態に関しては、ロール角にはウイーブモードとキャスタリングモードがほとんど寄与しないのに対し、舵角についてはウイーブモードは操舵トルクと同じ向きの挙動となるのにたいし、キャスタリングモードが逆向きとなり、互いに打ち消しあって最終的にほぼキャプサイズモードの応答となっています。

さらに車速を変化させて、ウイーブモードとキャプサイズモードの固有値を変化させたときの応答を見てみます。

まず$v=4.5$m/sではウイーブモードの減衰が小さくなり、初期の応答が振動的になりますが、全体の挙動はキャプサイズモードの応答に支配されていて応答は早く、定常値は小さくなることがわかります。

逆に$v=5.5$m/sではウイーブモードの減衰が大きくなり、初期応答の振動は早く収束する一方、キャプサイズモードの応答は緩慢になり、定常値が大きくなっていることがわかります。

計算に使用したScilabスクリプトを末尾に載せておきます。これまでと同様、固有値解析のスクリプト（benchmark_bicycleEIG.sce）を実行後に、車速（5行目のv）や入力（6行目のu）を変更して計算を実行してみてください。

まとめ

ステップ応答をモード展開することにより、旋回開始時の逆操舵挙動はウイーブモードによって引き起こされ、旋回挙動全体はキャプサイズモードによって支配されていることを示しました。特にキャプサイズモードが直進から旋回に移行するときの応答を支配しており、操縦性を評価する上で重要な固有モードであるといえるでしょう。

これまで安定性の観点から、キャプサイズモードは重要でないという扱いをされており、あまり注目されてきませんでした。また、操縦性についてもライダーの感性と密接に関わるため定量的な評価が困難とされてきましたが、キャプサイズモードと操縦性の関連については研究する価値があるのではないかと思います。

以上のように、二輪車の安定性と操縦性に関する基礎的な知見について、解析手法とともに紹介してきました。これをもって第一部終了としたいと思います。

今後はこの解析のベースとして用いた運動方程式の導出方法とその限界や、より高速走行が可能なオートバイに特有な挙動への拡張などの紹介を計画しています。

benchmark_bicycleME.sce

//
// Step Response Analysis with modal expansion
//
clear y;
v=5.0;
u=[0; -0.2];
sl = BicycleSys(v,bicycle.g,M,K0,K2,C1);
[V D]=spec(sl.A);
[D k]=gsort(diag(D));
V=V(:,k);
alpha=inv(V)*sl.B*u;
beta=sl.C*V;
//
t=0:0.01:30;
for i=1:4
    E=zeros(4,4);
    E(i,i)=1/D(i);
    y(:,:,i)=beta*E*alpha*(exp(D(i)*t)-1)*180/%pi;
end
ya=sum(y,3);
//
// Data Plot
a=gda();
a.font_size=2;
a.title.font_size=4;
a.x_label.font_size=3;
a.y_label.font_size=3;
f=scf();
f.axes_size=[1200 400];
//
subplot(1,2,1);
plot(t,y(1,:,1),t,2*real(y(1,:,2)),t,y(1,:,4),t,ya(1,:));
replot([0 -5 30 30]);
title('Step Response of Roll Angle at speed ' + string(v) +' [m/s]');
xlabel('Time [sec]');
ylabel('Angle [deg]');
xgrid();
legend(['Castering'; 'Weave'; 'Capsize'; 'Total']);
//
subplot(1,2,2);
plot(t,y(2,:,1),t,2*real(y(2,:,2)),t,y(2,:,4),t,ya(2,:));
replot([0 -2 30 14]);
title("Step Response of Steering Angle at speed " + string(v) + " [m/s]");
xlabel("Time [sec]");
ylabel("Angle [deg]");
xgrid();
legend(['Castering'; 'Weave'; 'Capsize'; 'Total']);

脚注

1. 荒井 迅, 「常微分方程式の解法」, 2020, p.109, ISBN978-4-320-11188-2

2. 野波健蔵, 西村秀和, 「MATLABによる制御理論の基礎」, 1998, p.85, ISBN4-501-31930-5