shigureni.dev

二重振り子 — 初期値への鋭敏さ

腕の先にもう一本腕を吊るす。上の腕の角度 a1、下の腕の角度 a2 の二つだけで姿勢が決まる。この二角度を運動方程式で進めると、振り方が周期に収まらず、二度と同じ軌道を描かない。

二重振り子。支点に吊った腕の先にもう一本の腕を吊った系で、状態は二つの角度 a1, a2 とその角速度 v1, v2 の四つ。各角の角加速度は両方の角度・角速度が絡む連立式で決まり、解析的に解けない。一本の振り子は小振幅なら単振動に収まるが、腕が一本増えるだけで非線形項が効き、カオス系の標準例になる。

下の腕の先 tip の位置は a1, a2 から三角関数で出る。角加速度 acc1, acc2 を連立式で計算し、角速度を足して角度を進める。一回の刻みを大きく取ると式が発散するので、DT を小さく刻んで毎フレーム数回まわす。

const G = 1
const DT = 0.06
const SUBSTEPS = 8

const accel = (a1, a2, v1, v2, l1, l2) => {
  const s12 = Math.sin(a1 - a2)
  const c12 = Math.cos(a1 - a2)
  const denom = 3 - Math.cos(2 * a1 - 2 * a2)
  const acc1 =
    (-G * 3 * Math.sin(a1) -
      G * Math.sin(a1 - 2 * a2) -
      2 * s12 * (v2 * v2 * l2 + v1 * v1 * l1 * c12)) /
    (l1 * denom)
  const acc2 =
    (2 * s12 * (v1 * v1 * l1 * 2 + G * 2 * Math.cos(a1) + v2 * v2 * l2 * c12)) / (l2 * denom)
  return [acc1, acc2]
}

export default (c) => {
  const w = c.canvas.width
  const h = c.canvas.height
  const px = w / 2
  const py = h * 0.42
  const l1 = h * 0.2
  const l2 = h * 0.2
  let a1 = Math.PI * 0.8
  let a2 = Math.PI * 1.1
  let v1 = 0
  let v2 = 0

c.fillStyle = colors.paper
  c.fillRect(0, 0, w, h)

return () => {
    c.fillStyle = colors.tintAlpha(238, 0.12)
    c.fillRect(0, 0, w, h)

for (let i = 0; i < SUBSTEPS; i++) {
      const acc = accel(a1, a2, v1, v2, l1, l2)
      const dt = DT / SUBSTEPS
      v1 += acc[0] * dt
      v2 += acc[1] * dt
      a1 += v1 * dt
      a2 += v2 * dt
    }

const midX = px + l1 * Math.sin(a1)
    const midY = py + l1 * Math.cos(a1)
    const tipX = midX + l2 * Math.sin(a2)
    const tipY = midY + l2 * Math.cos(a2)

c.strokeStyle = colors.tint(28)
    c.lineWidth = 3
    c.beginPath()
    c.moveTo(px, py)
    c.lineTo(midX, midY)
    c.lineTo(tipX, tipY)
    c.stroke()

c.fillStyle = colors.tint(28)
    c.beginPath()
    c.arc(midX, midY, 5, 0, Math.PI * 2)
    c.arc(tipX, tipY, 6, 0, Math.PI * 2)
    c.fill()
  }
}

腕が振れて、下の先が同じ道を二度と通らずに揺れ回る。半透明の paper を毎フレーム上から塗っているので、腕の残像が薄く尾を引く。denom の 3 - cos(2a1-2a2) がゼロに近づくと角加速度が跳ねるので、DT を小さく刻んで SUBSTEPS 回まわし、一回ぶんの角度変化を小さく抑えている。

同じ式の振り子をもう一本、a2 だけ 0.001 ずらして重ねる。最初は二本がぴたりと重なって動く。

const G = 1
const DT = 0.06
const SUBSTEPS = 8

const advance = (p, l1, l2) => {
  for (let i = 0; i < SUBSTEPS; i++) {
    const acc = accel(p.a1, p.a2, p.v1, p.v2, l1, l2)
    const dt = DT / SUBSTEPS
    p.v1 += acc[0] * dt
    p.v2 += acc[1] * dt
    p.a1 += p.v1 * dt
    p.a2 += p.v2 * dt
  }
}

export default (c) => {
  const w = c.canvas.width
  const h = c.canvas.height
  const px = w / 2
  const py = h * 0.42
  const l1 = h * 0.2
  const l2 = h * 0.2

const arms = [
    { a1: Math.PI * 0.8, a2: Math.PI * 1.1, v1: 0, v2: 0, g: 200 },
    { a1: Math.PI * 0.8, a2: Math.PI * 1.1 + 0.001, v1: 0, v2: 0, g: 40 },
  ]

c.fillStyle = colors.paper
  c.fillRect(0, 0, w, h)

return () => {
    c.fillStyle = colors.paper
    c.fillRect(0, 0, w, h)

for (const p of arms) {
      advance(p, l1, l2)
      const midX = px + l1 * Math.sin(p.a1)
      const midY = py + l1 * Math.cos(p.a1)
      const tipX = midX + l2 * Math.sin(p.a2)
      const tipY = midY + l2 * Math.cos(p.a2)

c.strokeStyle = colors.tint(p.g)
      c.lineWidth = 3
      c.beginPath()
      c.moveTo(px, py)
      c.lineTo(midX, midY)
      c.lineTo(tipX, tipY)
      c.stroke()

c.fillStyle = colors.tint(p.g)
      c.beginPath()
      c.arc(tipX, tipY, 6, 0, Math.PI * 2)
      c.fill()
    }
  }
}

二本が重なって一本に見える時間がしばらく続き、ある瞬間から急に離れて、片方が上に跳ねたとき片方は下に垂れる。初期の差は 0.001 しかないのに、進めるほど差が拡大して別の軌道になる。

初期値への鋭敏な依存（カオスの定義の一つ）。近い二つの初期値の差が時間とともに指数的に広がる性質。差が e 倍になる速さがリアプノフ指数で、正なら系はカオス。差が小さいうちは指数の肩がまだ小さく二本は揃って見えるが、肩がある大きさを越えると一気に分かれる。初期値を完全に測れない限り長期の予測ができない、というのがこの系の含意。

初期値をわずかずつ変えた振り子を何本も重ね、腕を消して下の先の点だけを打つ。最初は一点に固まり、進むと一本の弧に伸び、その先で画面に散る。

const G = 1
const DT = 0.06
const SUBSTEPS = 8
const COUNT = 28

export default (c) => {
  const w = c.canvas.width
  const h = c.canvas.height
  const px = w / 2
  const py = h * 0.42
  const l1 = h * 0.2
  const l2 = h * 0.2

const arms = Array.from({ length: COUNT }, (_v, i) => ({
    a1: Math.PI * 0.8,
    a2: Math.PI * 1.1 + i * 0.0004,
    v1: 0,
    v2: 0,
    g: 20 + (i / COUNT) * 200,
  }))

c.fillStyle = colors.paper
  c.fillRect(0, 0, w, h)

return () => {
    c.fillStyle = colors.paper
    c.fillRect(0, 0, w, h)

c.fillStyle = colors.tint(p.g)
      c.beginPath()
      c.arc(tipX, tipY, 4, 0, Math.PI * 2)
      c.fill()
    }
  }
}

28 個の点が一塊から動きだし、最初の一秒は団子のまま揺れ、そこから先頭と後尾がほどけて弧に伸び、最後は支点のまわり全域にばらける。初期値の差 0.0004 刻みの並びが、進むほど引き伸ばされて軌道全体を埋める。