Hi everyone,

A couple of days ago I came across this fascinating interview of Stephen Wolfram where he addressed the different concepts arising from his latest and most promising project: the Wolfram Physics Project. If you are curious about how his Theory of Everything relates to physics, quantum computation or the origin of life I highly recommend having a look at the short clips available on Lex Fridman’s channel, they are much easier to digest than the whole +3h interview.

Among the different topics covered I was particularly interested in the central notion of hypergraphs and their visual aesthetics in 3D space.

From what I understand these graphs are similar to Cellular Automaton in that they follow a discrete model of computation but are also reminiscent of the Lindenmayer L-Systems because of the rewriting system they are built on. Those who are familiar with these types of algorithms (and I know there are many inside the Processing community) will easily catch up with this new model as it mainly consists in applying a rule over and over again until it produces something that looks really complicated.

What I find particulariy fascinating is not only the complexity that emerges from simple rules but also the remarkable diversity of hypergraphs that it can produce. Without any notion of geometry the rewriting system can yet lead to various looking shapes in 2D and 3D space: trees, open/closed networks, meshes, reticular systems, Sierpiński structures,…etc.

Like for Cellular Automata, the Wolfram community is keeping a registry of the notable “universe” models (as they call them) found so far with their corresponding rules. Here is a naive interpretation in Python mode of the first and most popular example used by Stephen Wolfram in his presentations.

## Hypergraph Example (Python Mode)

```
from collections import defaultdict
from itertools import chain
add_library('toxiclibs')
add_library('peasycam')
W, H = 1400, 800
####### Multiple relations example (see part 2.7)
####### -> https://www.wolframphysics.org/technical-introduction/basic-form-of-models/rules-depending-on-more-than-one-relation/
# RULE (if edges corresponding to pattern 'a' -> replace by new edges corresponding to pattern 'b')
a = [[0, 1], [0, 2]]
b = [[0, 2], [0, 3], [1, 3], [2, 3]]
# Starting axiom
s = [[0, 0], [0, 0]]
# Number of steps/generations
G = 13
def setup():
size(W, H, P3D)
stroke(105, 140, 200)
smooth(8)
global physics, uniques, a, b, s
cam = PeasyCam(this, 1000)
physics = VerletPhysics3D()
physics.setDrag(.2)
#Transform / Translate axiom according to rule
a_,b_,s_ = [list(chain(*e)) for e in (a,b,s)]
ws = set(k for k in b_ if k > max(a_))
while ws:
a_ += [ws.pop()]
s_ += [max(s_)+1]
d = {k:i for i, k in enumerate(a_)}
t = [map(lambda x: s_[d[x]], e) for e in b]
print "Translated axiom: ", t
#Finding matching relations in the array of edges
for g in xrange(G-1):
temp = []
ids = set()
w = max(chain(*t))+1
#This part (until line 66 + line 70) should be modified depending on the input rule
for i, (xa, ya) in enumerate(t):
for j, (xb, yb) in enumerate(t):
if i != j:
if xa == xb and yb > ya: # reminder: [[xa, xb], [ya, yb]] = [[0, 1], [0, 2]]
if {i, j}.isdisjoint(ids):
temp.append((xa, ya, xb, yb))
ids.update([i, j])
for id in sorted(ids, reverse=True): del t[id]
for xa, ya, xb, yb in temp[::-1]:
t.extend([[xa, yb], [xa, w], [ya, w], [yb, w]])
w += 1
#Get id of each node
uniques = list(set(chain(*t)))
print "Hypergraph after %.i steps: %.i nodes // %.i edges" % (G, len(uniques), len(t))
#Store neighboring nodes for each node in a dictionary (handle both binary and ternary edges)
d = defaultdict(set)
for edge in t:
n = len(edge)
if n == 2:
n1, n2 = edge
d[n1].add(n2)
d[n2].add(n1)
elif n == 3:
n1, n2, n3 = edge
d[n1].add(n2)
d[n3].add(n2)
d[n2].update([n1, n3])
else:
pass
#Find smallest and largest number of neighbors (not mandatory, using this for computing edge length later)
nn = sorted(len(d[k]) for k in d)
lo, hi = nn[0], nn[-1]+1
#Create a particle for each node + give it a repulsion force
for n in uniques:
p = VerletParticle3D(Vec3D.randomVector().scale(W).add(Vec3D(W>>1, H>>1, 0)))
physics.addParticle(p)
physics.addBehavior(AttractionBehavior3D(p, 200, -.5))
#Create spring between each pair of nodes (handle both binary and ternary edges)
for edge in t:
n = len(edge)
if n == 2:
n1, n2 = edge
p1, p2 = [physics.particles.get(uniques.index(n)) for n in edge]
l = ( len(d[n1]) + len(d[n2]) ) * .5 # make length of edge depend on the average number of neighbors of its end vertices
f = map(l, lo, hi, 1, .3) #the higher the degree of the end vertices, the lower the factor
s = VerletSpring3D(p1, p2, l*l*f, .3)
physics.addSpring(s)
elif n == 3:
n1, n2, n3 = edge
p1, p2, p3 = [physics.particles.get(uniques.index(n)) for n in edge]
l1 = ( len(d[n1]) + len(d[n2]) ) * .5
l2 = ( len(d[n2]) + len(d[n3]) ) * .5
f1 = map(l1, lo, hi, 1, .2)
f2 = map(l2, lo, hi, 1, .2)
s1 = VerletSpring3D(p1, p2, l1*l1*f1, .2)
s2 = VerletSpring3D(p2, p3, l2*l2*f2, .2)
physics.addSpring(s1)
physics.addSpring(s2)
else:
pass
def draw():
background('#000000')
translate(-W>>1, -H>>1)
#Update physics
physics.update()
#Draw springs
for s in physics.springs:
line(s.a.x(),s.a.y(), s.a.z(), s.b.x(), s.b.y(), s.b.z())
#Draw nodes
pushStyle()
strokeWeight(3.5)
stroke('#BEF0FF')
for p in physics.particles:
point(p.x(), p.y(), p.z())
popStyle()
```

It works with other rules (pictures above), but you will have to modify lines 59 to 66 according to the pattern you are trying to find.

My hope is that someone will come up with a better implementation in Java mode and share it here for everyone to play with. Who knows, some day one of us might come across an amazing pattern that will be worth including in Wolfram registry !