砂山モデルと自己組織化臨界 — BTWモデルをブラウザで動かす
このトピックを深掘りする
クリエイティブコーディング
ブラウザでの数学的シミュレーション・自然現象の可視化・生成アートの実装技法
全12本中 8 本目。
次に読む記事
二重振り子のカオスを数値で見る — RK4積分と初期値鋭敏性の実装
二重振り子シミュレーションの実装を通じて、ラグランジュ力学による運動方程式の導出、RK4積分の安定性、NaN対策、初期値鋭敏性の可視化までを解説する。
Fourier級数で音を作る — 倍音構成と波形合成の実装
矩形波・鋸歯状波・三角波をFourier級数で合成し、エピサイクルで可視化する。倍音の加算が音色をどう変えるか、Web Audio APIによる音合成の実装を通じて体感する。
アリの巣シミュレーション — フェロモンで創発する群知能
フェロモンベースの経路探索をエージェントベースモデルで実装し、アリの群知能をブラウザで再現する。蒸発・拡散パラメータの調整で創発パターンがどう変わるかも検証した。
シリーズ全体の流れを見ながら、次に読む記事へ進めます。 初めての方はホームへ →
自己組織化臨界とは何か
砂山に砂粒を1つずつ落としていく。最初は何も起きない。砂が積もり、傾斜が急になると、あるとき雪崩が起きる。小さな雪崩もあれば、山全体を巻き込む大きな雪崩もある。
この「小さな入力が予測不能な規模の連鎖反応を引き起こす」現象を自己組織化臨界(Self-Organized Criticality, SOC)と呼ぶ。1987年にPer Bak、Chao Tang、Kurt Wiesenfeldが提唱した概念で、地震、森林火災、株価変動など、自然界の至るところに現れるべき乗則の背後にあるメカニズムとして注目された。
重要なのは「自己組織化」の部分。外部からパラメータを調整しなくても、系が自然に臨界状態へ向かう。砂山モデル(BTWモデル)はこの現象を最もシンプルに表現する数理モデルだ。
BTWモデルのルール
ルールは驚くほど単純:
- 2Dグリッドの各セルに「砂粒数」を持たせる
- 中央に砂粒を1つ落とす
- 砂粒数が閾値(4)以上のセルは崩壊(topple): 自身から4を引き、上下左右の4近傍に1ずつ配る
- 崩壊が連鎖しなくなるまで繰り返す
- 端から落ちた砂粒は消失する(開放境界条件)
崩壊前: 崩壊後:
0 0 0 0 1 0
0 4 0 → 1 0 1
0 0 0 0 1 0たった4行のルールから、べき乗分布に従う雪崩が生まれる。
コアロジック: topple関数
実装の核心は崩壊処理。Uint8Array で1Dフラット配列としてグリッドを管理する:
const TOPPLE_THRESHOLD = 4
function topple(grid: Uint8Array, cols: number, rows: number): number {
let totalTopples = 0
let unstable = true
while (unstable) {
unstable = false
for (let y = 0; y < rows; y++) {
for (let x = 0; x < cols; x++) {
const idx = y * cols + x
if (grid[idx] >= TOPPLE_THRESHOLD) {
unstable = true
const distribute = Math.floor(grid[idx] / 4)
grid[idx] -= distribute * 4
totalTopples += distribute
// Von Neumann近傍(上下左右)
if (y > 0) grid[(y - 1) * cols + x] += distribute
if (y < rows - 1) grid[(y + 1) * cols + x] += distribute
if (x > 0) grid[y * cols + (x - 1)] += distribute
if (x < cols - 1) grid[y * cols + (x + 1)] += distribute
}
}
}
}
return totalTopples
}なぜ while (unstable) でループするのか。1回のスキャンでは崩壊が収まらない場合があるからだ。あるセルが崩壊して近傍に砂を配ると、その近傍も閾値を超える。この連鎖がアバランシェ(雪崩)。ループは連鎖が完全に収まるまで回り続ける。
Math.floor(grid[idx] / 4) としているのは、閾値を大幅に超えた場合への対応。砂粒を64個まとめて落とすケースでは、1セルに大量の砂が載る。grid[idx] -= 4 だと何度もループが必要になるが、一度に floor(n/4) 回分の崩壊を処理することで収束を速める。
開放境界条件の意味
境界の扱いが物理的に重要な意味を持つ:
// 端のセルが崩壊すると、外側への分配分は消失
if (y > 0) grid[(y - 1) * cols + x] += distribute
if (y < rows - 1) grid[(y + 1) * cols + x] += distribute
if (x > 0) grid[y * cols + (x - 1)] += distribute
if (x < cols - 1) grid[y * cols + (x + 1)] += distributeif 文でガードしているだけ。条件を満たさない方向への砂粒は単に捨てられる。これが開放境界条件(open boundary)。砂山の端から砂が落ちていくイメージだ。
もし境界をトーラス(周期境界条件)にすると、砂粒が系から出ていかない。エネルギーが散逸しないので臨界状態に到達せず、SOCは観察できない。開放境界は散逸系であることの表現であり、SOCの必要条件。
ImageDataによる高速描画
グリッドを毎フレーム描画するには、セルごとに fillRect を呼ぶのではなく Canvas API の ImageData を直接操作する:
const CELL_SIZE = 6
const GRAIN_COLORS: [number, number, number][] = [
[17, 24, 39], // 0粒 - 暗い背景
[59, 130, 246], // 1粒 - 青
[34, 197, 94], // 2粒 - 緑
[245, 158, 11], // 3粒 - 琥珀
]
function render(
ctx: CanvasRenderingContext2D,
grid: Uint8Array,
cols: number,
rows: number,
ratio: number,
) {
const w = cols * CELL_SIZE * ratio
const h = rows * CELL_SIZE * ratio
const imgData = ctx.createImageData(w, h)
const data = imgData.data
const scaledCell = CELL_SIZE * ratio
for (let y = 0; y < rows; y++) {
for (let x = 0; x < cols; x++) {
const val = grid[y * cols + x]
const color = GRAIN_COLORS[Math.min(val, 3)]
const pxStartX = Math.floor(x * scaledCell)
const pxStartY = Math.floor(y * scaledCell)
const pxEndX = Math.floor((x + 1) * scaledCell)
const pxEndY = Math.floor((y + 1) * scaledCell)
for (let py = pxStartY; py < pxEndY; py++) {
for (let px = pxStartX; px < pxEndX; px++) {
const idx = (py * w + px) * 4
data[idx] = color[0]
data[idx + 1] = color[1]
data[idx + 2] = color[2]
data[idx + 3] = 255
}
}
}
}
ctx.putImageData(imgData, 0, 0)
}色は砂粒数に対応する4色のみ。Math.min(val, 3) で3以上はすべて同じ色にクランプしている。安定状態では各セルが0〜3のいずれかなので、4色で十分。崩壊中の一瞬だけ4以上になるが、描画タイミングではほぼ見えない。
ratio は devicePixelRatio 。Retinaディスプレイでは2倍のピクセルを描く必要がある。Canvas の内部解像度と CSS サイズを分離するのは Canvas 描画の基本テクニック。
べき乗則の観察
砂山モデルの醍醐味は、雪崩のサイズ分布がべき乗則に従うこと:
P(s) ~ s^(-tau)小さな雪崩は頻繁に、大きな雪崩はまれに起きる。しかしその分布は指数分布ではなく、べき乗分布。つまり「ありえないほど大きな雪崩」が指数分布の予想よりもずっと高い確率で起きる。
実装では topple の戻り値(totalTopples)がアバランシェサイズ。これをヒストグラムに記録して対数プロットすれば、直線的な分布が見えてくる。SandPileコンポーネントでは簡略化して、アバランシェ回数と最大サイズだけを表示している:
const topples = topple()
if (topples > 0) {
avalanchesRef.current++
if (topples > maxAvalancheRef.current) {
maxAvalancheRef.current = topples
}
}Auto Dropモードで砂粒を落とし続けると、系が臨界状態に達した後は最大サイズが突発的にジャンプする。これがSOCの本質: 系は自発的に「ギリギリ不安定」な状態を維持する。
自動投下のタイミング制御
Auto Dropモードでは requestAnimationFrame のフレームカウントでタイミングを制御する:
let frameCount = 0
const loop = () => {
frameCount++
if (autoDropRef.current && frameCount % Math.max(1, 60 - speed * 5) === 0) {
const cx = Math.floor(cols / 2)
const cy = Math.floor(rows / 2)
dropGrains(cx, cy, dropAmount)
}
render()
rafRef.current = requestAnimationFrame(loop)
}speed が大きいほどモジュロの値が小さくなり、投下頻度が上がる。Math.max(1, ...) は0除算防止。setInterval ではなく rAF ベースにしているのは、描画と同期させるため。タブが非アクティブ時に rAF が停止する副作用も都合が良い。
なぜ砂山モデルが面白いのか
BTWモデルが示す根本的な洞察は、複雑さはルールの複雑さから生まれるのではないということ。4つのルール、1つの閾値、1つの近傍定義。これだけで、自然界のべき乗則を再現できる。
セルオートマトン全般に言えることだが、局所的なルールがグローバルな振る舞いを生む。ライフゲームは「生命っぽさ」を、SandFallは「物質っぽさ」を、そして砂山モデルは「臨界っぽさ」を生む。どれもグリッドと近傍ルールという同じ骨格の上に成り立っている。
使った道具たち
自己組織化臨界と複雑系の全体像を掴むなら、道具棚の書籍が参考になる。単純なルールから複雑な振る舞いが創発するメカニズムは、ソフトウェア設計における「疎結合なモジュールが複雑なシステムを構成する」という思想と地続きだ。