Oregon Curriculum Network
Discovering Math with Python
By AlexanderAlUS  Own work, CC BYSA 3.0, Link
"Graphene" refers to an hexagonal grid of cells, the vertexes being carbon atoms. However any hexagonal mesh, such as for game boards, might be referred to as a "graphene pattern".
Quadrays are explained in other Notebooks. Four basis vectors emanate to the corners of a volume 1 tetrahedron of edges 2R or 1D, in the canonical version, where R and D refer respectively to the Radius and Diameter of imaginary spheres packed together, giving this home base tetrahedron.
The Quadrays {2, 1, 1, 0}, meaning all 12 permutations of those numbers, fan out from (0,0,0,0) to the corners of a cuboctahedron.
In [1]:
from itertools import permutations
g = permutations((2,1,1,0))
unique = {p for p in g} # set comprehension
print(unique)
I have elsewhere used this fact to algorithmically generate consecutive shells of 12, 42, 92, 162... spheres (balls) respectively; a growing cuboctahedron of $10 S^{2} + 2$ balls per shell S = 1,2,3... (1 when S=0).
However suppose we don't want to grow the grid omnidirectionally, but only in a plane. Each ball will be surrounded by six neighbors meaning at the center of a hexagon.
The cuboctahedron supplies four such hexagons i.e. its 24 edges comprise four hexagons orbiting the center. We may use any one of them.
The algorithm begins with a planar subset of the vectors {2, 1, 1, 0} used to compute the six vertexes surrounding (0,0,0,0). We may call these six vertexes "carbons".
Then hop to neighboring hexagon centers (where no carbon is located) using an additional set of vectors.
From these new centers, compute the six surrounding carbons again, some of which will have already been found, as neighbors share fences, with three faces (centers) sharing each fence post (carbon).
Using the Python set object, the algorithm filters to keep only unique carbons.
Keep track of hexagon centers, a dual mesh, in a separate set.
(0,0,0,0) will be the first center (ring0).
If qrays r, s are 60 degrees apart on the same hexagon, pointing to neighboring carbons, then r + s will be the "hop" vector over the fence (edge) to the neighboring "yard" (face), or center.
Once we have six vertex vectors from a center, computing the six hop vectors (for jumping over the fences) will be a matter of summing pairs of adjacent (60 degree separated) vectors. We only keep new centers i.e. those of the next ring (see below).
What about edges?
As we go around a hexagon in 60 degree increments, say in a clockwise direction, we will be finding edges in terms of adjacent ball pairs.
To avoid redundancy, (ball_a, ball_b)  any edge  will be sorted. Any two quadrays may be ordered as 4tuples e.g. (2, 1, 1, 0) is "greater than" (2, 1, 0, 1).
With unique representations of any edge, in the form of sorted tuples of qray namedtuples, we will be able to employ the same general technique employed with vertexes (carbons) and face centers: check the existing database for uniqueness and throw away (filter) anything already in the database. Sets will not allow duplicates.
The first step is to isolate six of the twelve from {2, 1, 1, 0} that define a hexagon.

One of the hexagons is TZOQXV. Do you see it in the above graphic? Another one is TYRQWS. If we regenerate all of the vectors AZ mentioned above, we'll have a vocabulary suitable for graphene grid development, and then some.
In [2]:
from qrays import Qvector, IVM
A, B, C, D = Qvector((1,0,0,0)), Qvector((0,1,0,0)), Qvector((0,0,1,0)), Qvector((0,0,0,1))
E,F,G,H = B+C+D, A+C+D, A+B+D, A+B+C
I,J,K,L,M,N = A+B, A+C, A+D, B+C, B+D, C+D
O,P,Q,R,S,T = I+J, I+K, I+L, I+M, N+J, N+K; U,V,W,X,Y,Z = N+L, N+M, J+L, L+M, M+K, K+J
In [3]:
# two "beacons" of six spokes
hexrays = [T, Z, O, Q, X, V] # to surrounding carbon atoms
hoprays = [T+Z, Z+O, O+Q, Q+X, X+V, V+T] # to neighboring (vacant) hex centers
In [4]:
(T.angle(Z), Z.angle(O), O.angle(Q), Q.angle(X), X.angle(V), V.angle(T))
Out[4]:
Lets verify that, going around the hexagon, each pair of consecutive hexrays is 60 degree apart. And ditto for hoprays, the vectors we'll use to jump over the fence to neighboring hexagon centers.
In [5]:
(hoprays[0].angle(hoprays[1]),
hoprays[1].angle(hoprays[2]),
hoprays[2].angle(hoprays[3]),
hoprays[3].angle(hoprays[4]),
hoprays[4].angle(hoprays[5]),
hoprays[5].angle(hoprays[0]))
Out[5]:
Looks like we're in business!
As with the growing cuboctahedron and the CCP packing, it makes sense to think in terms of consecutive rings.
The hexagonal coordination sequence is generated by:
In [6]:
def A008458(n):
# OEIS number
if n == 0:
return 1
return 6 * n
[A008458(x) for x in range(10)]
Out[6]:
I will use this as a check as the algorithm generates multiple rings.
In [7]:
centers = {IVM(0,0,0,0)} # center face
edges = set() # no duplicates permitted
carbons = set()
ring0 = [Qvector((0,0,0,0))]
def next_ring(ring):
"""
Use only the most recently added hexagonal ring
of face centers to compute the next ring, moving
outward: 1, 6, 12, 18, 24...
"""
new_faces = []
for face in ring:
verts = []
# CARBONS
for spoke in hexrays:
v = face + spoke
carbons.add(v.coords) # just the namedtuple is added to the set
verts.append(v)
# EDGES
for bond in zip(verts, verts[1:] + [verts[0]]):
# adding carbontocarbon bonds if not already in the set
edge = tuple(sorted([bond[0].coords, bond[1].coords]))
edges.add(edge)
# CENTERS
for jump in hoprays:
neighbor = face + jump
previous = len(centers)
centers.add(neighbor.coords)
if len(centers) > previous: # if True, face is new
new_faces.append(neighbor)
return new_faces
In [8]:
def rings(n):
prev = ring0
for ring in range(n):
print("Ring: {:3} Number: {:4}".format(ring, len(prev)))
nxt = next_ring(prev)
prev = nxt
rings(12)
Note these are the expected numbers for consecutive rings.
Now that we have our database, it's time to generate some graphical output. As with the FCC, I'll use POVRay's scene description language and then render in POVRay. We just want to look at the edges and carbon atom vertexes.
In [9]:
sph = """sphere { %s 0.1 texture { pigment { color rgb <1,0,0> } } }"""
cyl = """cylinder { %s %s 0.05 texture { pigment { color rgb <1.0, 0.65, 0.0> } } }"""
def make_graphene(fname="../c6xty/graphene.pov", append=True):
"""
Scan through carbons, edges, converting to XYZ and embedding
in POVRay Scene Description Language
"""
if append:
pov = open(fname, "a")
else:
pov = open(fname, "w")
# graphene will be included as a single object in the
# parent povray script, where lighting, camera position,
# and background are defined
print("#declare graphene = union{", file=pov)
for atom in carbons:
v = Qvector(atom).xyz()
s = sph % "<{0.x}, {0.y}, {0.z}>".format(v)
print(s, file=pov)
for bond in edges:
v0, v1 = bond
v0 = Qvector(v0).xyz()
v1 = Qvector(v1).xyz()
c = cyl % ("<{0.x}, {0.y}, {0.z}>".format(v0), "<{0.x}, {0.y}, {0.z}>".format(v1))
print(c, file=pov)
print("}\n", file=pov)
make_graphene(append=False)
Success!
The graphrene grid needs twelve of its hexagons to become pentagons in order to lose the 720 degrees necessary to turn it into a concave / convex "global matrix" or "hexapent" grid.
This 720 degrees is characteristic of polyhedrons more generally and is sometimes referred to as "Descartes' Deficit". We may visualize the subtration of 720 degrees as the removal of one tetrahedron, thereby likewise begetting a tetrahedron (something with an inside and outside, not flat).