Higher-dimensional generalizations [of substitution systems]

The state of a d-dimensional substitution system can be represented by a nested list of depth d. The evolution of the system for t steps can be obtained from

SSEvolve[rule_, init_, t_, d_Integer] := Nest[FlattenArray[# /. rule, d]&, init, t]

FlattenArray[list_, d_] := Fold[Function[{a,n}, Map[MapThread[Join,#,n]&,a,-{d+2}]], list, Reverse[Range[d]-1]]

The analog in 3D of the 2D rule on page 187 is

{1->Array[If[LessEqual[##],0,1]&,{2,2,2}], 0->Array[0&,{2,2,2}]}

Note that in d dimensions, each black cell must be replaced by at least d+1 black cells at each step in order to obtain an object that is not restricted to a dimension d-1 hyperplane.