Main Content

数値データの積分

この例では、一連の離散的な速度データを数値積分し、移動距離を概算する方法を示します。integral 関数群は関数ハンドルのみを入力として受け入れるので、離散データ セットには使用できません。関数式を積分に使用できないときは、trapz または cumtrapz を使用します。

速度データの表示

以下の速度データと、対応する時間データについて考えます。

vel = [0 .45 1.79 4.02 7.15 11.18 16.09 21.90 29.05 29.05 ...
29.05 29.05 29.05 22.42 17.9 17.9 17.9 17.9 14.34 11.01 ...
8.9 6.54 2.03 0.55 0];
time = 0:24;

このデータは、24 秒間、1 秒間隔で測定した自動車の速度 (m/秒) を表します。

速度のデータ点をプロットし、各点を直線で結びます。

figure
plot(time,vel,'-*')
grid on
title('Automobile Velocity')
xlabel('Time (s)')
ylabel('Velocity (m/s)')

加速中の勾配は正、速度が一定の期間の勾配は 1、減速中の勾配は負です。時間 t = 0 において、車両は速度 vel(1) = 0 m/秒で停止しています。車両は t = 8 秒で最高速度 vel(9) = 29.05 m/秒に達するまで加速し、その速度を 4 秒間維持します。次に vel(14) = 17.9 m/秒まで減速し、その速度を 3 秒間維持し、最終的には停車します。この速度曲線には不連続点が複数あるので、1 つの連続関数でこの曲線を表すことはできません。

合計移動距離の計算

trapz はデータ点を使用して台形を作成することにより離散積分を実行するので、不連続性をもつデータ セットの扱いに適しています。この方法ではデータ点間の動きが線形であると仮定しており、データ点間の動きが非線形である場合には正確性が低下することがあります。グラフのデータ点を頂点とする台形を作図して、これを図解することができます。

xverts = [time(1:end-1); time(1:end-1); time(2:end); time(2:end)];
yverts = [zeros(1,24); vel(1:end-1); vel(2:end); zeros(1,24)];
p = patch(xverts,yverts,'b','LineWidth',1.5);

trapz は一連の離散データの下側の領域を複数の台形に分割し、その台形の面積を計算します。次に、個々の台形の面積を加算して面積の合計を計算します。

trapz を使用して速度データを数値積分することにより、自動車の合計移動距離 (影付き部分の面積) を計算します。既定では、trapz(Y) 構文を使用する場合、点の間隔は 1 と仮定されます。ただし、trapz(X,Y) 構文を使用すると、別の等間隔または非等間隔 X を指定できます。この場合、time ベクトルの読み取り間隔は 1 なので、既定の間隔が使用可能です。

distance = trapz(vel)
distance = 345.2200

自動車が t = 24 秒までに移動した距離は約 345.22 m です。

累積移動距離のプロット

関数 cumtrapztrapz と密接に関連しています。trapz は最終積分値のみを返しますが、cumtrapz はベクトルの中間値も返します。

累積移動距離を計算して結果をプロットします。

cdistance = cumtrapz(vel);
T = table(time',cdistance','VariableNames',{'Time','CumulativeDistance'})
T=25×2 table
    Time    CumulativeDistance
    ____    __________________

      0                0      
      1            0.225      
      2            1.345      
      3             4.25      
      4            9.835      
      5               19      
      6           32.635      
      7            51.63      
      8           77.105      
      9           106.15      
     10            135.2      
     11           164.25      
     12           193.31      
     13           219.04      
     14            239.2      
     15            257.1      
      ⋮

plot(cdistance)
title('Cumulative Distance Traveled Per Second')
xlabel('Time (s)')
ylabel('Distance (m)')

参考

| |

関連するトピック