MATLAB代码（沙漏模型）

%sand pile
clear all
clf
nx=52; %must be divisible by 4
ny=100;
Pbridge = .05;
z=zeros（nxny）;
o=ones（nxny）;
sand = z ;
sandNew = z;
gnd = z ;
gnd（1:nxny-3）=1 ; % the ground line
gnd（nx/4:nx/2+4ny-15）=1; %the hole line
gnd（nx/2+6:nxny-15）=1; %the hole line
gnd（nx/4 ny-15:ny） = 1; %side line
gnd（3*nx/4 1:ny） = 1 ;
imh = image（cat（3z‘sand‘gnd‘））;
set（imh ‘erasemode‘ ‘none‘）
axis equal
axis tight
for i=1:1000
p=mod（i2）; %margolis neighborhood
sand（nx/2ny/2） = 1; %add a grain at the top
%upper left cell update
xind = [1+p:2:nx-2+p];
yind = [1+p:2:ny-2+p];
vary = rand（nxny）<.95 ;
vary1 = 1-vary;
sandNew（xindyind） = ...
gnd（xindyind）.*sand（xindyind） + ...
（1-gnd（xindyind））.*sand（xindyind）.*sand（xindyind+1） .* ...
（sand（xind+1yind+1）+（1-sand（xin