Actually this image gives us a hint that the implementation hasn't reached optimal time, because it isn't exactly a causal graph. We don't really need to calculate the squares at intermediary locations along straight lines between events, but we do anyways. We can probably refine our state update procedure, but what changes need to be made to the state structure? Does each square need to have its own local clock?

The answer is "yes", but due to incorporation of the Skip12 DataStructure, the code base has grown to considerable size. Here we will give a brief preview what to expect once the code is deployed through WFR. Assume we have a list of valid states in a variable outputData. The space time plot can be produced as follows:

PathGraphs[data_] := Module[
  {len = Length[data], actives, edges, ind},
  ind = 2;
  edges = Reap[While[ind <= len,
     actives = Select[data[[ind]]["PhaseSpace"],
       #["Time"] == data[[ind]]["Time"] &];
     Sow[MapThread[DirectedEdge, {
        Lookup[data[[ind - 1]]["PhaseSpace"],
         Keys[actives]], Values[actives]}]];
     ind++
     ]];
  Graph[Flatten[edges[[2]]]]
  ]

g1 = PathGraphs[outputData];

Show[Graph3D[g1,
  VertexShapeFunction -> "Cuboid",
  VertexSize -> 440,
  EdgeStyle -> Red,
  VertexStyle -> Blend[{Yellow, Orange}],
  VertexCoordinates -> (# -> Append[#["Position"], -#["Time"]
        ] & /@ VertexList[g1]),
  PerformanceGoal -> "Quality"
  ],
 ViewVertical -> {0, 0, 1},
 ViewPoint -> {10000, -2 20000, 2 10000},
 Lighting -> "ThreePoint",
 PlotRange -> {{0, 15}, {0, 15}, Automatic},
 Boxed -> True
 ]

This should be similar to one of the simple images above, except that we've removed some intermediary states along straight trajectories. This is made possible by the following line of code:

actives =  Select[outputData[[ind]][
   "PhaseSpace"], #["Time"] == outputData[[ind]]["Time"] &]

Which looks for masses whose local clock has the same time as the global clock. These masses then get new nodes added to their spacetime path graphs. The spacetime plot is not exactly the causal graph, but obviously very similar. Here's another spacetime graph, which is not so easy to decipher:

This is made just by changing outputData and sending it through graphics functions. Obtaining the causal graph is a simple one-liner:

Graph[Apply[GraphUnion,
  PathGraph[Position[#["Events"] & /@ outputData,
       #][[All, 1 ;; 2]],
     DirectedEdges -> True
     ] & /@ Range[16]],
 GraphLayout -> "LayeredDigraphEmbedding"]

The takeaway point is that this latest, optimized algorithm outputs data that is more fundamentally like the causal graph than previous algorithms. This is only possible because each mass keeps its own local time, which is then updated only at a relevant event.