HTML Canvasで物理シミュレーション その2

前回(→ HTML Canvasで物理シミュレーション)はCanvasを使って斜方投射の物理シミュレーションを作成しました。今回は前回挙げた課題点を解決していきます。その上で、最後に空気抵抗がある場合の運動をシミュレートしたいと思います。

複数の物体への対応

前回は質点1個のシミュレーションでした。今回は複数同時にシミュレートできるように改良します。具体的には質点をクラス化して扱いやすくします。

// Ball クラス
class Ball {
  constructor() {}
}

Ball クラスには位置(x, y)と速度(vx, vy)をプロパティとして持たせます。コンストラクターの引数で初期パラメーターをオブジェクトとして渡せるようにしましょう。その際、初速度は速さと角度(v, deg)で指定できるようにします。全てデフォルト値は0とし、負の数も0とします。それと発射後の着地判定用のプロパティも追加しておきます。

// Ball クラス
class Ball {
  constructor(attr) {
    // 位置
    this.x = Math.max(attr.x || 0, 0);
    this.y = Math.max(attr.y || 0, 0);
    
    // 速度
    attr.v = Math.max(attr.v || 0, 0),
    attr.deg = Math.min(Math.max(attr.deg || 0, 0), 90); // degの範囲は [0, 90]

    this.vx = attr.v * Math.cos(attr.deg * Math.PI / 180); // 速度のx成分
    this.vy = attr.v * Math.sin(attr.deg * Math.PI / 180); // 速度のy成分

    // 着地したかどうか
    this.landed = false;
  }
}

位置の更新と描画のメソッドはこの Ball クラスにもたせます。

// Ball クラス
class Ball {
  〜
  // 更新
  update(dt) {
    if (this.landed) return; // 着地しているなら処理しない
    this.x += this.vx * dt;
    this.y += this.vy * dt;
    this.y = Math.max(this.y, 0);
    this.vy += g * dt;

    if (dt && this.y === 0) {
      this.landed = true; // 着地した
    }
  }

  // 描画
  draw() {
    ctx.beginPath();
    ctx.arc(this.x / scale, this.y / scale, 1, 0, Math.PI * 2);
    ctx.fill();
  }
}

それからこの Ball のインスタンスを格納しておく配列も準備しましょう。コンストラクター内でこの配列に自分自身を格納します。

const ballList = []; // Ball を管理する配列

// Ball クラス
class Ball {
  constructor(attr) {
    〜

    ballList.push(this); // 自身を配列に格納
  }
  〜
}

あとは requestAnimationFrame のコールバックで、ballList から一個ずつボールを取り出して、更新と描画を行うだけです。

const tick = (time) => {
  if (!lastTime) lastTime = time;
  const dt = (time - lastTime) / 1000; // 経過時間 [s]
  
  // 画面の消去
  ctx.clearRect(0, 0, canvas.width, canvas.height);
  
  // ボールの更新・描画
  ballList.forEach(ball => {
    ball.update(dt);
    ball.draw();
  });
  
  lastTime = time;
  requestAnimationFrame(tick);
};

実際に Ball をいくつか作ってアニメーションさせてみましょう。

new Ball({
  x: 8,
  y: 0,
  v: 15,
  deg: 45
});

new Ball({
  /* x, y を省略した場合はそれぞれ 0 になる */
  v: 12,
  deg: 54
});

new Ball({
  x: 0,
  y: 12,
  v: 12
  /** deg = 0 になるので水平投射 */
});

requestAnimationFrame(tick);

See
the Pen
Projectile Motion 2
by SHIN Inc. (@shin-inc)
on CodePen.0

これで好きなだけボールを飛ばせるようになりました(もちろん極端に数が多くなると処理落ちしたりメモリが不足します)。

以下は角度を0〜90度の範囲でちょっとずつ変えたものを200個一斉に飛ばす例です。

See the Pen
Projectile Motion 2 / 200 Balls
by SHIN Inc. (@shin-inc)
on CodePen.0


フレームレートの固定

これまで状態の更新は、前回の呼び出しからの経過時間を使用していました。ですが requestAnimationFrame
はブラウザの描画タイミングで呼び出されるため、フレームレート(FPS)を固定することができません。実行環境や処理の負荷によってその間隔が変わるため、同じ初期条件でも違うシミュレート結果が出る可能性があるのです。

そこで、開始時から現在までに何フレーム経過したかを計算し、前回呼び出しからのフレーム増分だけ更新処理を行うことで、FPSを固定したいと思います。

const fps = 30, // FPS
  dt = 1 / fps; // 1フレーム当たりの秒数

let startTime, // 開始時間
  frameCount = 0; // 経過フレーム数

const tick = (time) => {
  if (!startTime) startTime = time; // 初回呼び出し時に開始時間を設定
  const currentFrame = Math.floor((time - startTime) * fps / 1000); // 現在の経過フレーム数
  let df = currentFrame - frameCount // 前回からのフレーム増分

  // フレーム増分だけ更新処理を繰り返す
  while (df) {
    // ボールの更新
    ballList.forEach(ball => {
      ball.update(dt);
    });
    df--;
  }

  // 画面の消去
  ctx.clearRect(0, 0, canvas.width, canvas.height);

  // 再描画
  ballList.forEach(ball => {
    ball.draw();
  });

  frameCount = currentFrame; // フレーム数を更新
  requestAnimationFrame(tick);
};

これで経過時間がフレーム単位に満たない場合は更新処理がスキップされ、次回に持ち越されます。逆に間隔が伸びた場合は、複数回更新処理を行うことで辻褄を合わせるため、フレームレートが一定となります。

フレームレートを上げれば、それだけ刻み幅が小さくなり、誤差が少なくなります。でもその分、一回の処理の負荷が上がるため、処理落ちしやすくなります。処理速度と精度がトレードオフの関係になっているわけです。

ここで扱っている斜方投射のシミュレーションは、着地判定に誤差が生じるので、フレームレートによって到達距離が若干変わります。

空気抵抗がある場合の斜方投射

ここまではシミュレーション手法の改良でした。ここからはシミュレートする減少を変更します。今までの単純な斜方投射に空気抵抗を追加してみましょう。

流体の抵抗(抗力)は物体の速度が小さいときは速度に比例し、速度が大きいときは速度の二乗に比例することが知られていますが、空気中の運動では概ね速度の二乗に比例すると考えてよさそうです(日常的なスケールでは空気の粘性は小さく、流速の寄与の方が大きい)。

新たに Ball クラスを継承して、空気抵抗を受ける BallAirResistance クラスを作ります。空気抵抗の大きさ(係数)は物体の形状によって変わるものなので、任意の値を設定できるようにしています。また、空気抵抗が加わることでx軸方向の速度変化が加わります。
※抵抗が0のとき Ball と同じになるので、Ball を直接拡張しても構いません。

ここで、速度の二乗に比例する空気抵抗のx成分、y成分は次のように計算できます。
F は空気抵抗、k は空気抵抗係数、v は速さ、θ は速度が水平面と成す角度)

F = - k v 2 F x = F cos θ = - k v 2 cos θ = - k v 2 v x v = - k v v x F y = F sin θ = - k v 2 sin θ = - k v 2 v y v = - k v v y
// BallAirResistance クラス
class BallAirResistance extends Ball {
  constructor(attr) {
    super(attr);
    this.k = Math.max(attr.k || 0, 0); // 空気抵抗係数
  }
  update(dt) {
    if (this.landed) return; // 着地しているなら処理しない
    this.x += this.vx * dt;
    this.y += this.vy * dt;
    this.y = Math.max(this.y, 0);
    
    const v = Math.sqrt(this.vx ** 2 + this.vy ** 2);
    
    this.vx += -this.k * v * this.vx * dt;      // 空気抵抗を受けて減速
    this.vx = Math.max(this.vx, 0);             // vxの下限は0
    this.vy += (g - this.k * v * this.vy) * dt; // 重力と空気抵抗の両方を受ける
    
    if (dt && this.y === 0) {
      this.landed = true; // 着地した
    }
  }
}

補足:空気抵抗による加速度の係数は、本来質量に反比例する(k/m)ものですが、ここでは質量のパラメーターを使用していないので、質量込みの係数となっています。

以下に、同じ初速度で空気抵抗があるものとないものを同時に飛ばせるシミュレーターを用意しました。空気抵抗があると急激に速度が落ちる様子を確認してみてください。また、空気抵抗がある方は赤く、ない方は黄色に色付けしています(色付の仕組みはソースをご覧ください)。

See the Pen
Projectile Motion Simulator 2
by SHIN Inc. (@shin-inc)
on CodePen.0


次なる課題

ここまで重力加速度が一定の場での運動のシミュレーションを行ってきました。それなりにちゃんとしたものができたはずです。

が! これらは解析的に解ける(積分が計算できる)問題ばかりなので、シミュレーションするまでもなく正確な予測が可能です。もっと計算が難しい問題や、そもそも解析的に解けない問題こそシミュレーションのしがいがあるというものです。

例えば重力が位置によって変化するような系……というわけで、次回は天体の運動をシミュレートしてみたいと思います。