A328506 - OEIS

%I #27 Oct 24 2019 10:52:24

%S 1,5,16,36,66,101,160,218,285,374,464,565,680,815,969,1124,1282,1467,

%T 1659,1863,2091,2346,2559,2824,3100,3411,3690,4043,4380,4697,5060,

%U 5468,5833,6266,6670,7132,7595,8006,8502,9004,9518,10039,10609,11155,11740,12304,12971,13603,14202,14861,15532,16217

%N Iteration of Abelian sandpile model where the n-th matrix expansions occurs. Begins with infinite sand in 1 X 1 matrix.

%C The Abelian sandpile model is a cellular automaton considering the behavior of sand grains on a square grid. Any square containing 4 or more grains will topple, sending a grain to each of its 4 neighbors and subtracting 4 grains from itself.

%C Here, expansion refers to the addition of a boundary layer to the outside of the existing matrix when the model reaches beyond the previous matrix boundary.

%e                                                       _ _ _ _ _

%e           _ _ _      _ _ _      _ _ _      _ _ _     |0|0|1|0|0|

%e    _     |0|1|0|    |0|2|0|    |0|3|0|    |0|4|0|    |0|2|1|2|0|

%e   |∞| -> |1|∞|1| -> |2|∞|2| -> |3|∞|3| -> |4|∞|4| -> |1|1|∞|1|1| -> ...

%e    ‾     |0|1|0|    |0|2|0|    |0|3|0|    |0|4|0|    |0|2|1|2|0|

%e           ‾ ‾ ‾      ‾ ‾ ‾      ‾ ‾ ‾      ‾ ‾ ‾     |0|0|1|0|0|

%e                                                       ‾ ‾ ‾ ‾ ‾

%e             ^                                             ^

%e      1st expansion on                              2nd expansion on

%e    1st iteration (a(1) = 1)                      5th iteration (a(2) = 5)

%o (MATLAB)

%o L = 3;

%o plane = zeros(3,3);

%o plane(2,2) = 99999999999999999999999999999999999999999999999;

%o listn = [];

%o for n = 1:50000

%o     plane2 = plane;

%o     for r = 1:L

%o         for c = 1:L

%o             if plane(r,c) > 3

%o                 plane2(r,c) = plane2(r,c) - 4;

%o                 plane2(r-1,c) = plane2(r-1,c)+1;

%o                 plane2(r+1,c) = plane2(r+1,c)+1;

%o                 plane2(r,c-1) = plane2(r,c-1)+1;

%o                 plane2(r,c+1) = plane2(r,c+1)+1;

%o             end

%o         end

%o     end

%o     if sum(plane2(:,1))+sum(plane2(1,:)) > 0

%o         plane2 = padarray(plane2,[1,1]);

%o         L = L+2;

%o         listn = [listn n];

%o     end

%o     plane = plane2;

%o end

%o fprintf('%s\n', sprintf('%d,', listn))

%o (PARI)

%o Step(M)={my(n=#M, R=matrix(n,n)); for(i=2, n-1, for(j=2, n-1, if(M[i,j]>=4, R[i,j]-=4; R[i,j+1]++; R[i,j-1]++; R[i-1,j]++; R[i+1,j]++))); M+R}

%o Expand(M)={my(n=#M, R=matrix(n+2, n+2)); for(i=1, n, for(j=1, n, R[i+1, j+1]=M[i,j])); R}

%o seq(n)={my(L=List(), M=matrix(3,3), k=0); while(#L<n, k++; my(o=#L+2); M[o,o]=4; M=Step(M); if(M[1,], M=Expand(M); listput(L,k))); Vec(L)} \\ _Andrew Howroyd_, Oct 23 2019

%Y Cf. A259013, A249872.

%K nonn

%O 1,2

%A _Parker Grootenhuis_, Oct 22 2019