Your data structure needs to be a set of grids like so, doubling in their resolution each layer.
so the resolution of layer n is 1 << n by 1 << n.
Each square just holds a boolean value. 1 or 0.
The new value at each square is a function f of its neighbours values (including its own present value), for example here are the neighbours of the red square:
You might say that the orange squares are its siblings, the yellows are its children and parent. I also sometimes include the grey squares as neighbours of the red square.
Each square is then updated independently according to function f. Like other cellular automata, you need to have a separate copy of the grid hierarchy so you can 'double buffer' and the change in one cell doesn't affect the others on that frame.
The layers must be updated in a specific order, which appears like so:
So on frame 1 you update layer 3 based on its siblings and its parent(s). Since there is no layer below you can treat the layer below as being zeros, just as we do for the neighbours outside of the edges of the grid. On the next frame, update layer 2, then the next frame update layer 3 and so on.
Despite appearances the order is quite easy to calculate. Take the frame number in binary and find the position of the first 1 in it-
In our example we have four layers, so our depth is 4. If p is the position of the first 1, then:
layerToUpdate = depth - 1 - p.
To prevent overflow, add:
if layerToUpdate < 0
frame = 1;
layerToUpdate = depth-1;
The function f
The function f is the set of rules determining whether the centre square 'lives or dies', or you can call it the mapping from the neighbour values to 0,1. The number of mappings is so large that it helps to try out different simpler rule sets. For example a mapping based on the number of 'on' neighbours without worrying about where they are placed.