元胞自动机-沙堆模型代码

function[]=testCA(n)z=zeros(n,n);cells=z;cells(n/2,.25*n:.75*n)=1;cells(.25*n:.75*n,n/2)=1;imh=image(cat(3,cells,z,z));set(imh,'erasemode','none')axisequalaxistight%Ö÷º¯Êý²ÎÊýnx=52;%mustbedivisibleby4ny=100;Pbridge=.05;z=zeros(nx,ny);o=ones(nx,ny);sand=z;sandNew=z;gnd=z;gnd(1:nx,ny-3)=1;%thegroundlinegnd(nx/4:nx/2+4,ny-15)=1;%theholelinegnd(nx/2+6:nx,ny-15)=1;%theholelinegnd(nx/4,ny-15:ny)=1;%sidelinegnd(3*nx/4,1:ny)=1;%Ö÷º¯Êýfori=1:1000p=mod(i,2);%margolisneighborhoodsand(nx/2,ny/2)=1;%addagrainatthetop%upperleftcellupdatexind=[1+p:2:nx-2+p];yind=[1+p:2:ny-2+p];%randomizetheflow--10%ofthetimevary=rand(nx,ny).9;vary1=1-vary;sandNew(xind,yind)=...gnd(xind,yind).*sand(xind,yind)+...(1-gnd(xind,yind)).*sand(xind,yind).*sand(xind,yind+1).*...(sand(xind+1,yind+1)+(1-sand(xind+1,yind+1)).*sand(xind+1,yind));sandNew(xind+1,yind)=...gnd(xind+1,yind).*sand(xind+1,yind)+...(1-gnd(xind+1,yind)).*sand(xind+1,yind).*sand(xind+1,yind+1).*...(sand(xind,yind+1)+(1-sand(xind,yind+1)).*sand(xind,yind));sandNew(xind,yind+1)=...sand(xind,yind+1)+...(1-sand(xind,yind+1)).*...(sand(xind,yind).*(1-gnd(xind,yind))+...(1-sand(xind,yind)).*sand(xind+1,yind).*(1-gnd(xind+1,yind)).*sand(xind+1,yind+1));sandNew(xind+1,yind+1)=...sand(xind+1,yind+1)+...(1-sand(xind+1,yind+1)).*...(sand(xind+1,yind).*(1-gnd(xind+1,yind))+...(1-sand(xind+1,yind)).*sand(xind,yind).*(1-gnd(xind,yind)).*sand(xind,yind+1));%scramblethesitestomakeitlookbettertemp1=sandNew(xind,yind+1).*vary(xind,yind+1)+...sandNew(xind+1,yind+1).*vary1(xind,yind+1);temp2=sandNew(xind+1,yind+1).*vary(xind,yind+1)+...sandNew(xind,yind+1).*vary1(xind,yind+1);sandNew(xind,yind+1)=temp1;sandNew(xind+1,yind+1)=temp2;sand=sandNew;set(imh,'cdata',cat(3,z',sand',gnd'))drawnowend%buildtheGUI%definetheplotbuttonplotbutton=uicontrol('style','pushbutton','string','Run','fontsize',12,'position',[100,400,50,20],'callback','run=1;');%definethestopbuttonerasebutton=uicontrol('style','pushbutton','string','Stop','fontsize',12,'position',[200,400,50,20],'callback','freeze=1;');%definetheQuitbuttonquitbutton=uicontrol('style','pushbutton','string','Quit','fontsize',12,'position',[300,400,50,20],'callback','stop=1;close;');number=uicontrol('style','text','string','1','fontsize',12,'position',[20,400,50,20]);stop=0;%waitforaquitbuttonpushrun=0;%waitforadrawfreeze=0;%waitforafreezewhile(stop==0)if(run==1)%nearestneighborsumsum(x,y)=cells(x,y-1)+cells(x,y+1)+cells(x-1,y)+cells(x+1,y)+cells(x-1,y-1)+cells(x-1,y+1)+cells(3:n,y-1)+cells(x+1,y+1);%TheCArulecells=(sum==3)|(sum==2&cells);%drawthenewimageset(imh,'cdata',cat(3,cells,z,z))%updatethestepnumberdiaplaystepnumber=1+str2num(get(number,'string'));set(number,'string',num2str(stepnumber))endif(freeze==1)run=0;freeze=0;enddrawnow%needthisintheloopforcontrolstoworkend