1. 二维沙堆模型
算法规则：
- \(L\times L\) 的正方格子上，随机地分布着沙子，每一格上的沙子数不超过 4 。
- 有一个小孩，在 \([L/2, L/2]\) 这个格子上，一粒粒地增加沙子。
- 任何一个格子的沙子数大于等于 4 的时候，就会发生一次沙崩，这个格子上的沙子数 - 4，四个相邻的格子上的沙子数 +1，即沙子流向四个邻居。
- 如果边缘上的沙子发生沙崩，那么一部分沙子流出整个沙盘，例如，如它有 3 个邻居在格子上，则这次沙崩这些邻居的沙子数分别 +1，还有一粒沙子流失。
2. 意义
- 这个沙堆是一个耗散系统，沙子源源不断地从外界流入，又源源不断地流出，系统内的物质并不守恒。
- 如果沙堆的初态确定，那么后续的所有动作就已经完全确定。所以是确定性的系统。
- 随着沙子的累积，小的沙崩是随时会发生的，大的沙崩需要较多的沙子累积，所以较少发生。
- 大的沙崩是很多小的沙崩的叠加。如何从解析上描述自相似性？是否存在标度不变，如何解析描述？
- 如何进一步得到 power law？
3. 代码
以下是一段 python 代码，可以实现小规模的沙堆模型。
import numpy as np
import matplotlib.pyplot as plt
size = 10
nsand = np.zeros([size, size], dtype='uint8')
for i in range(size):
for j in range(size):
nsand[i,j] = np.random.randint(0,4)
def avalanch( nsand, size ):
for i in range( size ):
for j in range( size ):
if nsand[i,j] >=4 :
nsand[i,j] -= 4;
if i-1 >= 0: nsand[i-1,j] += 1
if i+1 < size: nsand[i+1,j] += 1
if j-1 >= 0: nsand[i,j-1] += 1
if j+1 < size: nsand[i,j+1] += 1
return 1
return 0
nstep = 10000
magnitudes = []
for i in range(nstep):
print("i=",i)
nsand[(int)(size/2),(int)(size/2)] += 1; magnitude=0
while avalanch( nsand, size ) > 0 :
magnitude +=1
#print("--------------------------------- \n", nsand)
if magnitude > 0:
#print("magnitude = ", magnitude)
magnitudes.append( magnitude )
#print("magnitudes = ", magnitudes)
x = np.bincount(magnitudes)
data = []
for i in range( x.size ): # 去除次数为０的点
if x[i] > 0: data.append( [ i, x[i] ] )
xdata = []; ydata= [] #
for t in data:
xdata.append( np.log10(t[0]) ); ydata.append( np.log10(t[1]) )
plt.plot(xdata, ydata, "o")
plt.xlabel(r"$log_{10}$ magnitude")
plt.ylabel(r"$log_{10}$ (number of avalanches)")
plt.savefig("PowerLaw.png")
plt.show()
4. 结果
把沙崩的震级定义为一次大沙崩中包括多少次小沙崩，然后把震级、发生的次数取了 log 以后画在图上，得到下面的图片。
- 曲线非常紧凑，这说明存在一些规律。
- 曲线并不完全是直的，这说明什么？这个弯曲的意义是什么？
