shigureni.dev

入り組み具合を測る — フラクタル次元

「どの縮尺でも似ている」形は、長さや面積では捉えきれない。海岸線は測る物差しを細かくするほど長くなって、ふつうの長さに落ち着かない。

長さの代わりに次元を測る。格子の目を細かくしながら、形が通過するマスの数を数える。マスを細かくしたとき、通過数が何乗で増えるか。その指数が次元になる。直線なら 1、面を塗り潰せば 2。その間の半端な値が出る形がフラクタル。

箱数え次元。マスの一辺を s、形が通過したマス数を N(s) とすると、自己相似な形では N ∝ s^(-D) が成り立つ。両辺の対数を取れば log N が log(1/s) の一次関数に乗り、その傾きが次元 D。シェルピンスキー三角形なら D = log3 / log2 ≈ 1.585 で、線でも面でもない値になる。海岸線の長さが物差し依存で発散するという問いが、この一本の傾きに畳まれる。

シェルピンスキーの三角形を、三隅への中点写像をくり返すカオスゲームで点群にして溜めておく。この点群に格子をかぶせ、点が一つでも入ったマスを塗り、塗ったマスの数を数える。下は格子の目を 4・8・16・32 と段々細かくして、通過マスを数え直していく。

const buildPoints = (verts, count) => {
  const points = []
  let x = verts[0][0]
  let y = verts[0][1]
  for (let i = 0; i < count; i++) {
    const v = verts[(Math.random() * 3) | 0]
    x = (x + v[0]) / 2
    y = (y + v[1]) / 2
    if (i > 20) points.push(x, y)
  }
  return points
}

const countCells = (points, ox, oy, side, n, hit) => {
  hit.fill(0)
  const cellSize = side / n
  for (let i = 0; i < points.length; i += 2) {
    const cx = ((points[i] - ox) / cellSize) | 0
    const cy = ((points[i + 1] - oy) / cellSize) | 0
    if (cx < 0 || cx >= n || cy < 0 || cy >= n) continue
    hit[cy * n + cx] = 1
  }
  let total = 0
  for (let i = 0; i < hit.length; i++) total += hit[i]
  return total
}

export default (c) => {
  const w = c.canvas.width
  const h = c.canvas.height
  const side = Math.min(w, h) - 32
  const ox = (w - side) / 2
  const oy = (h - side) / 2

const verts = [
    [ox + side / 2, oy],
    [ox, oy + side],
    [ox + side, oy + side],
  ]
  const points = buildPoints(verts, 60000)
  const grids = [4, 8, 16, 32]
  const maxN = grids[grids.length - 1]
  const hit = new Uint8Array(maxN * maxN)
  let frame = 0

return () => {
    frame += 1
    const n = grids[((frame / 60) | 0) % grids.length]
    countCells(points, ox, oy, side, n, hit)
    const cellSize = side / n

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

c.fillStyle = colors.tintAlpha(80, 0.5)
    for (let cy = 0; cy < n; cy++) {
      for (let cx = 0; cx < n; cx++) {
        if (hit[cy * n + cx] === 0) continue
        c.fillRect(ox + cx * cellSize, oy + cy * cellSize, cellSize, cellSize)
      }
    }

c.fillStyle = colors.tint(20)
    for (let i = 0; i < points.length; i += 2) {
      c.fillRect(points[i], points[i + 1], 1, 1)
    }

c.strokeStyle = colors.tintAlpha(130, 0.6)
    c.lineWidth = 1
    c.beginPath()
    for (let i = 0; i <= n; i++) {
      const p = i * cellSize
      c.moveTo(ox + p, oy)
      c.lineTo(ox + p, oy + side)
      c.moveTo(ox, oy + p)
      c.lineTo(ox + side, oy + p)
    }
    c.stroke()
  }
}

buildPoints は中点写像で点群を作り、最初の 20 個を捨てて配列に溜める。countCells は点が落ちたマスに 1 を立て、立った数を合計する。目を細かくすると塗られるマスは増える。ただし面を塗り潰すほどには増えない。隙間が同じ割合で空き続けるので、増え方が s² より緩い。

増え方を読むには、各格子サイズの (log(1/s), log N) を平面に打つ。自己相似なら点は一直線に乗り、その傾きが次元。下は格子を 2・4・8・16・32・64 と倍々に細かくして通過マスを数え、log-log 平面に点を打って、最小二乗で直線を当てる。傾きの数値を縦に伸びる棒で示し、1 と 2 の目盛りに添える。

const fitSlope = (xs, ys) => {
  const len = xs.length
  let sx = 0
  let sy = 0
  let sxx = 0
  let sxy = 0
  for (let i = 0; i < len; i++) {
    sx += xs[i]
    sy += ys[i]
    sxx += xs[i] * xs[i]
    sxy += xs[i] * ys[i]
  }
  return (len * sxy - sx * sy) / (len * sxx - sx * sx)
}

export default (c) => {
  const w = c.canvas.width
  const h = c.canvas.height
  const side = Math.min(w, h)

const verts = [
    [side / 2, 0],
    [0, side],
    [side, side],
  ]
  const points = buildPoints(verts, 80000)
  const grids = [2, 4, 8, 16, 32, 64]
  const maxN = grids[grids.length - 1]
  const hit = new Uint8Array(maxN * maxN)

const logInvS = []
  const logN = []
  for (let g = 0; g < grids.length; g++) {
    const n = grids[g]
    const total = countCells(points, 0, 0, side, n, hit)
    logInvS.push(Math.log(n))
    logN.push(Math.log(total))
  }
  const slope = fitSlope(logInvS, logN)

const pad = 56
  const plotW = w - pad * 2
  const plotH = h - pad * 2
  const xMin = logInvS[0]
  const xMax = logInvS[logInvS.length - 1]
  const yMax = logN[logN.length - 1]
  const PX = (lx) => pad + ((lx - xMin) / (xMax - xMin)) * plotW
  const PY = (ly) => h - pad - (ly / yMax) * plotH
  let reveal = 0

return () => {
    if (reveal < grids.length) reveal += 0.04

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

c.strokeStyle = colors.tintAlpha(150, 0.7)
    c.lineWidth = 1
    c.beginPath()
    c.moveTo(pad, h - pad)
    c.lineTo(w - pad, h - pad)
    c.moveTo(pad, h - pad)
    c.lineTo(pad, pad)
    c.stroke()

c.strokeStyle = colors.tint(40)
    c.lineWidth = 2
    c.beginPath()
    c.moveTo(PX(xMin), PY(slope * xMin + (logN[0] - slope * logInvS[0])))
    c.lineTo(PX(xMax), PY(slope * xMax + (logN[0] - slope * logInvS[0])))
    c.stroke()

const shown = Math.min(grids.length, Math.floor(reveal))
    c.fillStyle = colors.tint(20)
    for (let i = 0; i < shown; i++) {
      const cx = PX(logInvS[i])
      const cy = PY(logN[i])
      c.beginPath()
      c.arc(cx, cy, 4, 0, Math.PI * 2)
      c.fill()
    }

const barX = w - pad + 12
    c.strokeStyle = colors.tintAlpha(150, 0.5)
    c.lineWidth = 1
    c.beginPath()
    c.moveTo(barX, PY((1 * yMax) / 2))
    c.lineTo(barX + 10, PY((1 * yMax) / 2))
    c.moveTo(barX, PY(yMax))
    c.lineTo(barX + 10, PY(yMax))
    c.stroke()

c.strokeStyle = colors.tint(20)
    c.lineWidth = 6
    c.beginPath()
    c.moveTo(barX + 5, h - pad)
    c.lineTo(barX + 5, h - pad - (slope / 2) * plotH)
    c.stroke()
  }
}

fitSlope が六点に最小二乗で直線を当て、その傾きを返す。六点はほぼ一直線に乗り、傾きは 1.585 あたりに落ち着く。三角形の隙間が同じ割合で空き続けるので、面（2）には届かず、線（1）よりは込み入る。verts を一直線上に並べれば傾きは 1 へ、点群を面いっぱいに撒けば 2 へ寄る。

buildPoints の中点写像を別の自己相似な形の点列に取り替えても、同じ countCells と fitSlope でその形の次元が出る。形は反復で立て、その細かさを変えながら測ると、込み入り具合が一本の傾きに畳まれる。