Gfyuki

I'm trying to plot a directed graph of $2^{16}$ nodes and $2^{16}$ edges (but not simply a cycle). Ultimately, I need to be able to share the graph with others in a way that it can be "explored" (zooming, panning, labels).

Mathematica did a fine job of drawing this graph in a way that minimized edge overlaps:

What you're seeing is a mass of nodes (blue) mashed up, totally hiding all the edges. This isn't a viable solution because it (1) requires an installation of Mathematica, (2) takes several minutes to generate, and (3) can't be saved—exporting the plot as SVG crashed all SVG viewers I tried.

SigmaJS with random initial positions, then ForceAtlas2

It seems for large graphs, rendering with HTML5 Canvas is the way to go, and SigmaJS is a popular choice.

The first problem with SigmaJS though was that it doesn't automatically place nodes the way Mathematica did. So to apply any force-directed drawing algorithm, first I had to supply all nodes with initial positions.

Randomly dispersing the 65,536 nodes in a square caused such a hopeless tangle of edges that, even after several minutes of running ForceAtlas2, I could only see this:

SigmaJS with ring-adjacent placement, then ForceAtlas2

Well, no big deal. Instead of random placement, I decided to do a naive depth-first search and place nodes in a ring as I found them. This way the majority of nodes would start right next to a neighbor. Here's what the evolution of that looked like with ForceAtlas2 at start, a few minutes in, an hour later, and a few hours later:

But this made it really evident that the results of force-directed graphing algorithms depend heavily on their initial states. I can see each of those radial "islands" being stuck in local optima, never reaching the configuration Mathematica did.

About this particular graph, and investigating Mathematica's algorithm

The graph I'm studying is a pseudo-random number generator of the form

$$x_{next} = 5x_{current}+273 bmod 65536$$

for the most part. (A quirk in the actual assembly code implementation actually causes shift-by-1s for ~700 of the 65536 edges.) In other words, what I'm graphing is the succession of "random" numbers generated by that formula, e.g.

$$0 rightarrow 273 rightarrow 1365 rightarrow 7098 rightarrow 35763 rightarrow 48016 rightarrow ldots$$

Eventually this succession yields a number we've already seen, closing the loop and forming one of the 3 cycles of this graph. I know this isn't really about data science (much less "Big Data") but the technique I'm looking for would help similar visualization problems on large, sparsely connected graphs.

To see what Mathematica's doing, first I plotted just a single succession for the first 1,000 integers, i.e.

begin{align}0 &rightarrow 273\1 &rightarrow 279\2 &rightarrow 285\ &ldots \999 &rightarrow 5268end{align}

(There are some 3- and 4-length chains due to the quirk mentioned before.) And here's the same for 10,000 integers:

Clearly Mathematica organizes subgraphs in some order to do with the size of each subgraph.

"Life" begins to form around 40,000 nodes, and as edges connect subgraphs at varying midpoints to produce more and more interesting shapes, and we converge toward the graph we began with:

Okay, so...

a. Is there a library or known algorithm that can do this in Javascript (or any other reasonably shareable output format)? This would be ideal.

b. Or, does anyone have any insights about Mathematica's algorithm? It's clearly not just a force-directed drawing algorithm. I wonder if I can develop a SigmaJS plugin that does something similar.

Here's the Mathematica implementation if anyone is interested in playing around with it.

  Hex[exp_] := FromDigits[exp, 16];

  LByte[exp_] := BitAnd[exp, Hex@"00ff"];

  HByte[exp_] := BitAnd[exp, Hex@"ff00"]~BitShiftRight~8;

  PRNG[v_] := Module[{L5, H5, v1, v2, carry},

    L5 = LByte@v*5;

    H5 = HByte@v*5;

    v1 = LByte@H5 + HByte@L5 + 1;

    carry = HByte@v1~BitGet~0;

    v2 = BitShiftLeft[LByte@v1, 8] + LByte@L5;

    Mod[v2 + Hex@"0011" + carry, Hex@"ffff" + 1]

    ];

  mappings = # -> PRNG@# & /@ Range[0, Hex@"ffff"];

  (* WARNING! GraphPlot takes a long time to generate. *)

  (* GraphPlot[mappings] *)