砂山崩しモデルをEXCEL VBAでシミュレーションする方法

自然は放っておくと、自然に平衡状態（自己組織化臨界状態と同義語？）になるらしいです。その状態では、大小さまざまな出来事が、べき乗則に従って予測不可能なタイミングで起きます。例えば、地震の発生や、砂山の雪崩現象や、自然林の林冠ギャップの分布なども、この法則に従っています。

そこで問題なのは、自然環境の理想の姿は、この臨界状態なのか？はたまた、そうならないようにコントロールした方がよいのか？どっちなんだろうということです。これから、考えて行きたいテーマです。（世間では、もう結論でてるんでしょうか？世間知らずなもので。。。）

※余談ながら、現在の資本主義も、この臨界状態に達したんじゃなかろうかと思っています。そうなると、些細な出来事が原因で、大きな事件が予測不可能に発生しちゃう。

そのために、とりあえずやるべきことは、

自己組織化臨界状態を理解するためシミュレーションプログラムを実行してみる。
日本の森林は、臨界状態なのかどうかを確かめる。その方法を調べる。

今回は、まず自己組織化臨界現象の例として、砂山崩しモデルをEXCEL VBAでプログラミングして実行する方法を紹介します。

砂山崩しモデルEXCEL VBAプログラミング

いきなりですが、プログラムを張っておきます。解説は、また今度。
参考文献は、こちらなので、私の解説よりは、こっちを読んだ方が早いと思います。

Option Explicit
Const L = 64 'エリアの範囲
Const STEP = 60000 '砂を落とす回数
Dim d(L, L) As Integer '砂の高さ(厳密には高さの差分)
Dim df(L, L) As Integer '崩れたかどうかのフラグを記録
Dim nsize As Integer '崩れた箇所の総数（雪崩のサイズ）
Dim count(STEP) As Integer '砂を落としたステップごとの雪崩サイズを記録


Sub BTW()
Dim k, i As Long
Dim myrnd1, myrnd2 As Integer
Dim x, y As Integer

Randomize

For y = 1 To L 'セルの初期化
 For x = 1 To L
  Cells(x, y).Interior.Color = &H0
 Next
Next

k = 0
While k < STEP
 nsize = 0
 Erase df
 'エリア内の砂を落とす場所を乱数で決める
 myrnd1 = Int(Rnd() * L + 1)
 myrnd2 = Int(Rnd() * L + 1)
 '砂を落とす
 d(myrnd1, myrnd2) = d(myrnd1, myrnd2) + 1
 '崩れるかどうかを再帰的にチェックする
 Call check(myrnd1, myrnd2)
 
 If k >= L * L * 2 Then '定常状態になるまでまつ
  For y = 1 To L
   For x = 1 To L

   '雪崩が起きそうなところを色つけの場合
   '   If d(x, y) = 3 Then
   '    Cells(x, y).Interior.Color = RGB(255, 255, 255)
   '   Else
   '    Cells(x, y).Interior.Color = RGB(0, 0, 0)
   '   End If
   
   '砂の高さごとに色つけの場合
   '   If d(x, y) = 0 Then
   '    Cells(x, y).Interior.Color = RGB(0, 0, 0)
   '   ElseIf d(x, y) = 1 Then
   '    Cells(x, y).Interior.Color = RGB(85, 85, 85)
   '   ElseIf d(x, y) = 2 Then
   '    Cells(x, y).Interior.Color = RGB(170, 170, 170)
   '   ElseIf d(x, y) = 3 Then
   '    Cells(x, y).Interior.Color = RGB(255, 255, 255)
   '   End If
   
   '雪崩が起きたところを色付けの場合
    If df(x, y) = 1 Then
     Cells(x, y).Interior.Color = RGB(255, 255, 255)
    Else
     Cells(x, y).Interior.Color = RGB(0, 0, 0)
    End If
   
   Next
  Next
 End If

 count(k) = nsize
 k = k + 1
Wend

'Sheets2にステップごとの雪崩サイズを書き出し
 For i = 1 To STEP
  Worksheets("Sheet2").Cells(i, 1) = count(i - 1)
 Next
End Sub

Sub check(x, y)
 If d(x, y) >= 4 Then '砂の高さが4以上だと崩れる
  If df(x, y) = 0 Then '雪崩サイズをカウントする
   nsize = nsize + 1
   df(x, y) = 1
  End If
  'その場所が崩れ、上下左右に散らばる。
  '散らばった場所が、崩れないか再びチェック
  d(x, y) = d(x, y) - 4
  If x + 1 <= L Then d(x + 1, y) = d(x + 1, y) + 1
  If x - 1 >= 1 Then d(x - 1, y) = d(x - 1, y) + 1
  If y + 1 <= L Then d(x, y + 1) = d(x, y + 1) + 1
  If y - 1 >= 1 Then d(x, y - 1) = d(x, y - 1) + 1
  If x + 1 <= L Then Call check(x + 1, y)
  If x - 1 >= 1 Then Call check(x - 1, y)
  If y + 1 <= L Then Call check(x, y + 1)
  If y - 1 >= 1 Then Call check(x, y - 1)
 End If
End Sub

こちらは、臨界状態を理解するのに読みやすくて面白くて分かりやすいです。