3D Zonotiles

by
Russell Towle
version 1.9, January 4, 2004
rtowle@inreach.com

Preface

The extension of the 2D zonotile function to three dimensions has proved very difficult. For several years I grappled with the problem, only to turn away in frustration. Not only do the methods employed here cost much time, when many planes are involved (say, more than twenty), but it fails outright to create many zonotiles.

Nevertheless, some lovely local space-fillings, of zonohedra tiled by smaller zonohedra, can be obtained. It is a curious fact, though far be it from me to prove this assertion, that, given a zonohedron tiled by some suite of smaller zonohedra, the entire space outside the zonohedron can be tiled by this same suite, as well.

These tilings of zonohedra by smaller zonohedra—these 3D zonotiles—are based upon arrangements of planes. Wherever three or more planes intersect at one point, a zonohedron arises. Any convex polyhedron ipso facto yields an arrangement of planes—the planes of its faces. These faces, produced into infinite planes, intersect, perhaps, in points outside the original polyhedron; and these points are the vertices of various stellations of that polyhedron. The exact relationship between 3D zonotiles and stellations of convex polyhedra seems an especially interesting subject for future study.

—Russell Towle, December, 2003.

Remarks about version 1.9: There is so much which could be improved in this notebook. For instance, variable names are often opaque; no matter whether a single number, or a collection of one hundred zonohedra, I will call it a, or qaq, or some such thing. Then again, I will sometimes go back and change a name in one part of my code, but not another, which depends upon the original name; and then all hell breaks loose.

New in this version: stellation diagrams, and a brief look at multiple stellations, and their related zonotiles; solid zones of zonohedra, within a zonotile; intriguing new examples based upon Platonic polyhedra (see the “Platonic Zonotiles II” randomization subsubsections); and the addition of some few references.

I always appreciate hearing from people who enjoy my notebooks, and especially like to learn about errors, or anything which needs improvement.

—Russell Towle, January, 2004.

Introduction

In an expedition to discover the El Dorado of geometric beauty, I am only a scout. Thus there is much lacking in this introduction to 3D zonotiles. For instance, what is a zonohedron? One simple definition: any convex polyhedron, bounded entirely by centrally-symmetrical polygons, is a zonohedron. Famous and familiar examples include the Cube; any right regular hexagonal prism, or right regular 2n-gonal prism; and Kepler’s Rhombic Dodecahedron. And, generally, an n-dimensional zonotope is a convex polytope bounded by (n-1)-dimensional zonotopes. A two-dimensional zonotope is called a zonogon.

Very well.

Based upon some unattributed PostScript code I found on the internet, which created the line segments of a 2D zonotile (a zonogon tiled by smaller zonogons), I wrote Mathematica code to create the zonotiles as lists of polygons. The procedure is known as the generalized dual method. In extending the algorithm to the third dimension, further changes were made. To handle the cases where a 3D zonotile (a zonohedron tiled by smaller zonohedra) contains zonohedra bounded in part by “higher” zonogons (hexagons, octagons, etc.), I modified George Hart’s “zonohedrify” function slightly, and named it zedron. For the cases where a component zonohedron is bounded only by rhombs or parallelograms, my own Z function is used.

2D zonotiles are based upon arrangements of lines, 3D zonotiles, on arrangements of planes. It is easy for us to construct and examine an arrangement of lines; but very difficult for us to visualize an arrangement of planes. For, one plane hides another.

The procedure is somewhat as follows: a set of planes is defined, in the form {x, y, z, d}, where {x, y, z} is the unit normal vector, and d is the distance from the origin. To every point where the planes intersect corresponds a zonohedron (determined by the normals to the planes meeting at that point), and to every region of space defined by the planes corresponds a vertex of the zonotile. A set of sign vectors is found, each vector equal in length to the number of planes. These sign vectors make the encoding for creating the zonohedra which tile the larger zonohedron (the zonotile).

Some sets of plane normals {x, y, z} (as for instance with the vectors which determine a polar zonohedron, which are easily found in cyclic order) give results which are largely analogous to those of the 2D zonotile algorithm. Neither the 2D nor the 3D zonotile algorithm can admit lines or planes with opposite normals, say, in 3D, {x, y, z} and {-x, -y, -z}. To obtain admissible parallel planes one must use different distances d; for instance, to get parallel planes one unit to either side of the origin: {x, y, z, 1} and {x, y, z, -1}.  But, when vectors are chosen from the octahedral or icosahedral symmetry groups, the options become more limited. For instance, in such a case one cannot assign distances of 0 and 1 to each plane (vector); the distances must “balance” across the origin, so, for instance, the one distance 0 is permissible (but this just gives the one undivided zonohedron determined by the set of vectors (planes)), or the distances {-1, 1}, or the distances, {-1, 0, 1}, or the distances, {-2, -1, 0, 1, 2}, etc.

So, the algorithm is in a precarious state, and not at all robust. I feel somewhat as though I were a witch, casting—I don’t know what: crystals, toads, sacred relics, and flowers, into an ancient cauldron. The cauldron steams and shudders for a good long while; perhaps the toads eat the crystals; but whatever it is which happens in there, I hardly know and cannot say—and at long last, out pops a zonotile.

As implemented below, the main zonotile creation portion of the procedure returns two lists: allpoly, in which the individual zonohedra are sorted into sets, according to the number of faces which bound them; and polyflat, a list of all the zonohedra in the form {center of zonohedron, zonohedron}. This latter allows for easy “explosion” of a zonotile, by translating each zonohedral tile away from the origin by some positive multiple of its center. It also facilitates selecting regions of zonohedral tiles, say, all tiles within some distance of the z=0 plane, or all tiles within some distance from the center of the zonotile-at-large.

Thus, some methods to select regions of zonohedra within a zonotile are sketched in, as well as to explode the zonohedra apart slightly or drastically, and to sort zonohedra into subsets according to volume, and then assign a different color to each subset. It is difficult to visualize a close-packing of polyhedra. It helps to explode the polyhedra apart, or to use animation to build up larger and larger regions of close-packed polyhedra.

Meta-Geometry

The subject of zonotopes tiled by smaller zonotopes is rich with ramifications. They relate to quasicrystals. Then again, zonotopes are affine projections of hypercubes, and under certain conditions, orthogonal projections of hypercubes. When we dissect an n-polar zonohedron into rhombic hexahedra, each hexahedron is a projection of a discrete 3-cube bounding the “parent” n-cube. But only if we dissect the polar zonohedron into hexahedra in every possible way do we obtain projections of all the 3-cubes bounding the n-cube. Then the hexahedra interpenetrate in a horrible mish-mash. However, any one dissection may at least sometimes be regarded as a hidden-detail-removed projection of the n-cube into a 3-space. For instance, the 4-cube projects into an n=4 polar zonohedron, and if this zonohedron is dissected (in either of two ways) into four rhombic hexahedra, that is exactly a hidden-detail-removed projection from the 4-space into the 3-space—as “seen” from some one direction in the 4-space. That is, the 4-cube is bounded by eight 3-cubes, and were we in the 4-space, and were the 4-cube opaque, so we could not see right through it, we could see at most four of the eight 3-cubes at one time. Some of the “detail” of the 4-cube would be hidden.

So we may at least vaguely relate zonotiles to hidden-detail-removed projections of hypercubes. To construct unequivocal examples of hidden-detail-removed projections of 4-polytopes into polyhedra tiled by smaller polyhedra, see my Regular Polytopes notebook (there is a hyperlink in the References section of this notebook).

Tilings of zonotopes by smaller zonotopes are in bijection with single-element liftings of oriented matroids. I must stop exactly there, having already exceeded my own understanding of this relation. Works cited in the References section amplify upon zonotopes and oriented matroids.

Any zonotope is “determined by” a set of vectors. These are sometimes called the generators of the zonotope. For instance, a regular hexagon is determined by three vectors, a square, by two. A few years ago I devised an unusual construction I called “zonohedral completion” (or zonogonal completion, in two dimensions). Given a “seed” convex polyhedron, say, the Platonic icosahedron, bounded by twenty equilateral triangles, and centered, let’s say, upon the origin, for each face, take its vertices as the generators of a zonohedron, and construct that zonohedron, rooted at the origin, and passing through and beyond the face. In the case of the icosahedron, the twenty sets of three vectors give twenty equal rhombic hexahedra radiating from a common point (the origin). Now we can discard the “seed” icosahedron; its job is done. Between these radiating hexahedra are various hollows. Find each hollow, find the edges meeting in the hollow, take these as vectors, and construct the determined zonohedron, rooted in each hollow. Now a new set of hollows arises, and the process is continued, through however many stages, until one closes up a convex zonohedron. In the case of the icosahedron, one closes up Kepler’s Rhombic Triacontahedron, in three stages.

In two dimensions, given a regular n-gon, one closes up a regular 2n-gon, after however many stages. The procedure generalizes to higher space, and in all cases, given the “seed” convex polytope, one closes up, after however many stages, the very zonotope which is determined by the vectors from the center to the vertices of the seed. This zonotope is tiled by smaller zonotopes, hence, it is a zonotile.

Finally, in the zonotile construction employed here, and in my 2D_Zonotiles notebook (there is a hyperlink in the References section of this notebook), “arrangements of hyperplanes” are used to set up sign vectors which encode the exact construction of the zonotile, zonohedron by zonohedron, or zonogon by zonogon. Here I will only remark that every convex polygon gives an arrangement of lines, and that every convex polyhedron gives an arrangement of planes. Thus, in the general, n-dimensional case, every convex n-polytope corresponds to an arrangement of hyperplanes which determines a particular n-zonotile.

But these same arrangements of lines, planes, and hyperplanes, derived from convex polygons, polyhedra, and polytopes, also determined the stellations of the polygons, polyhedra, and polytopes.

The stellation of polyhedra is a devilishly difficult subject. For instance, there are millions of stellations of Kepler’s Rhombic Triacontahedron. All are obtained by merely producing the planes of the triacontahedron’s faces. Let those planes be produced to infinity (rather than “stopping them off” at various discrete stellations). Every vertex where three or more planes intersect corresponds to a zonohedron in the zonotile. Every open region enclosed at least partly by planes corresponds to a vertex of the zonotile. However, the reciprocity between plane arrangements and zonohedra tiled by smaller zonohedra is more complete, yet.

Imagine a set of planes, with unit normals {x, y, z}, at non-zero distances from the origin, which intersect pairwise in lines, and in some several places, three (say) planes intersect in single points. Draw a sphere around one such point, of a size small enough that no other points are within.  Fixing our attention within the sphere, we find a single 0-dimensional object (the point), six 1-dimensional line segments, twelve 2-dimensional plane regions, and eight 3-dimensional solid regions of non-zero volume. Reciprocally, these correspond to a rhombic hexahedron, its six faces, its twelve edges, and its eight vertices. Hence the generalized dual method.

I cannot doubt but that there is some extremely close and significant connection between stellations and zonotiles. It remains to properly characterize this connection, and discover whatever numerical and geometrical patterns actually come into play. However, it may at least be remarked that, if one studies the stellation diagram of the Platonic Dodecahedron, points of intersection are at three distances from the center; while in the stellation diagram of the Icosahedron, points of intersection are at eight distances from the center. Meanwhile, there are three sets, and eights sets, of zonohedra, concentrically arranged, in the corresponding zonotiles. And the correspondence, in detail, is even more exact.

Initialization

Off[General::spell]
Off[General::spell1]

<<Graphics`Polyhedra`

Hart

(*faster than Cross*)
cross[ {ax_, ay_, az_}, {bx_, by_, bz_} ] :=
   {ay bz - az by, az bx - ax bz, ax by - ay bx}

collinear[v1_, v2_]:= Chop[cross[v1, v2]]=={0,0,0}

setStar[vecs_]:=
Block[{n,b,d,a,t,v},
v=DeleteCases[vecs, {0,0,0}, Infinity];
n=Range[Length[v]];
d={};
While[n!={},
t=v[[n]];
a=Map[collinear[First@t,#]&, t];
b=n[[ Flatten[Position[a,True]] ]];
n=Complement[n,b];
d={d,b[[1]]}];(*end While*)
gStar = v[[ Flatten[d] ]] (*set global variable*)
]
    
pToXYZ[p_] := Plus @@ MapThread[Times, {gStar/2, p}]

(*where p is a sign vector; whatever length gon, first two vectors suffice*)
faceNormal[p_] := 
  With[{indices = Position[p, 0]}, 
   cross[gStar[[indices[[1,1]]]], gStar[[indices[[2,1]]]]]]
   
endPoints[p_] := {(If[#1 == 0, -1, #1] & ) /@ p, (If[#1 == 0, 1, #1] & ) /@ p}
   
hemi[n_Integer]:=
	Block[{i,j},Flatten[Table[{i,j},{i,n-1},{j,i+1,n}],1]]
   
makeFace[{i_, j_}] := 
  With[{normal = cross[gStar[[i]], gStar[[j]]]}, 
   (Sign[Chop[normal . #1]] &) /@ gStar]

allFaces := With[{halfTheFaces = makeFace /@ hemi[Length[gStar]]}, 
   Union[halfTheFaces, -halfTheFaces]]
   
makeEdge[i_, p_, faceNormal_] := 
  With[{normal = cross[gStar[[i]], faceNormal]}, 
   MapThread[If[#1 != 0, #1, Sign[Chop[#2 . normal]]] & , {p, gStar}]]
   
makeEdges[p_] := 
  With[{normal = faceNormal[p]}, 
   With[{halfEdges = (makeEdge[#1, p, normal] & ) /@ Flatten[Position[p, 0]]}, 
    Union[halfEdges, (2*p - #1 & ) /@ halfEdges]]]
    
makePolygon[p_] := 
  With[{edges = makeEdges[p]}, With[{verts = 
      Union[Flatten[endPoints /@ edges, 1]]},
      cycleSort[ (pToXYZ /@ verts) ]
      ]]
       
cycleSort[verts_] := 
  With[{centroid = Plus @@ verts/Length[verts]}, 
   sortBy[ArcTan[(#1 - centroid) . (verts[[1]] - centroid), 
      (#1 - centroid) . (verts[[2]] - centroid)] & , verts]]
      
sortBy[fn_, list_] := (#1[[2]] & ) /@ Sort[({fn[#1], #1} & ) /@ list]

zedron[poly_] :=
Block[{}, 
  setStar[poly];
  Polygon /@ Map[makePolygon[#]&, allFaces]
  ]
     
faceCenters[poly_] := Map[cent[#]&, Cases[Level[poly,{-3}],{{_,_,_}..}] ]

edgeMidpoints[poly_] := 
  Union[Flatten[((#1 + RotateLeft[#1])/2 & ) /@ Cases[Level[poly,{-3}],{{_,_,_}..}], 1]]
   
vertices[poly_] := Union[Flatten[parse[poly], 1]]

Other

signvector::usage = "signvector[vecs], where vecs is a list
of n 3-vectors, returns a list of n-vectors with
components -1, 0, or 1. The number of zeroes in any one of the
sign vectors equals the number of sides of the corresponding
face of the zonohedron determined by the n vectors.";

signvector[v_]:=
Block[{dex,x,hf},
dex=hemi[Length[v]];
x=Map[mycross[ Sequence@@v[[ # ]] ]&, dex];
hf=Table[
Map[ Sign[Chop[Dot[x[[k]], #] ] ]&, v],
{k,Length[x]}];
Union[hf,-hf]
]

edgefind::usage = "edgefind[poly] returns the edges
of equilateral regular polytopes.";

edgefind[poly_]:=
Block[{a,b},
b=First@Union[ Map[mag[First@poly-#]&, Rest@poly],
SameTest->(Chop[#1-#2,.000001]==0&)];
Flatten[Table[
a=Drop[poly,i-1];
{First@a,#}& /@ Select[Rest@a, Abs[ mag[First@a-#] - b ]<.001&],
{i,Length[poly]}],1]
]

(*gives the Union of a list {center of zonohedron, zonohedron}*)

tol=.001;
won1[v_]:= Union[v, SameTest->(Max[Abs[ #1[[1]] - #2[[1]] ]] < tol &)] 
hedraunion[v_]:= Module[{r},
r = Split[Sort[Chop[v,.001]], (Abs[#1[[1,1]] - #2[[1,1]]] < tol &)];
Flatten[won1 /@ r, 1]]

getpolycent[poly_]:=
Chop[ cent[ vunion[ Flatten[Cases[poly, Polygon[x_]:>x],1] ] ] ]

(*gives a list of the equally spaced generators of a cone*)

vectors[n_Integer,pitch_,fac_Integer,ratio_]:=
Block[{e},
Table[If[Mod[i,fac]==0,e=ratio,e=1];
	N[{e Cos[Degree pitch] Cos[2Pi i/n],
	e Cos[Degree pitch] Sin[2Pi i/n],
	-e Sin[Degree pitch]}],
	{i,n}]]
	
vectors[n_:12,pitch_:35.2643,fac_:1,ratio_:1]:=
vectors[n, pitch, fac, ratio]

(*restrain vectors to half-plane*)
(*when pitch=0, create regular 2n-gon*)

halfvectors[n_Integer,pitch_,fac_Integer,ratio_]:=
Block[{e},
Table[If[Mod[i,fac]==0,e=ratio,e=1];
	N[{e Cos[Degree pitch] Cos[Pi i/n],
	e Cos[Degree pitch] Sin[Pi i/n],
	-e Sin[Degree pitch]}],
	{i,n}]]
	
halfvectors[n_:12,pitch_:35.2643,fac_:1,ratio_:1]:=
halfvectors[n, pitch, fac, ratio]

Z[v_List,spot_]:=
Block[{n=Length[v],fin,s,a,b,c,gon},
fin=hemi[n];(*make a list of indices for half the faces*)
Polygon /@ Flatten[Table[(*now create all the faces in a Table*)
(*1. create a single face*)
gon=quad[{v[[ fin[[i,1]] ]],v[[ fin[[i,2]] ]]},spot];
(*2. get surface normal*)
s=surf[ gon ];
(*find the vectors not involved in the face*)
a=Complement[v, v[[ fin[[i]] ]] ];
(*sort them according to the dot product with surface normal*)
b=Select[a,Dot[s,#]<0&];
c=Select[a,Dot[s,#]>0&];
(*now add the sets of vectors to the face, creating two faces*)
{Map[(Plus @@ b)+#&,gon],Map[(Plus @@ c)+#&,gon]},
{i,Length[fin]}],1]
]

hemi[n_Integer]:=(*indices, all pairs of n objects*)
	Block[{i,j},Flatten[Table[{i,j},{i,n-1},{j,i+1,n}],1]]
	
demi[n_Integer]:=(*indices, all triples of n objects*)
	Block[{i,j,k},Flatten[Table[{i,j,k},{i,n-2},{j,i+1,n-1},{k,j+1,n}],2]]

quad[x_,ref_]:={ref, ref+x[[1]], ref+x[[1]]+x[[2]], mark=ref+x[[2]]}

mycross[ {ax_, ay_, az_}, {bx_, by_, bz_} ] :=
	{ay bz - az by, az bx - ax bz, ax by - ay bx}

surf[f_List] := mycross[ f[[2]] - f[[1]], f[[3]] - f[[2]] ]

mag[v_]:=Sqrt[v.v](*magnitude of vector 'v'*)
unit[v_]:=v/mag[v](*unit vector in direction of 'v'*)
cent[v_List]:=Plus@@v/Length[v](*centroid of point list*)
tosegs[gon_]:=Transpose[ {gon, RotateLeft[gon]} ](*polygon to edges*)

toPoly[e_]:=With[{e1=(e /. Line[x_]:> x)},
Drop[Join[First@e1,Map[#[[2]]&, Rest@e1]],-1]
]

(*Carl Woll's vector union*)
tol=.000001;
vunion1[v_]:= Union[v, SameTest->(Max[Abs[#1 - #2]] < tol &)] 
vunion[v_]:= Module[{r},
r = Split[Sort[Chop[v,.000001]], (Abs[#1[[1]] - #2[[1]]] < tol &)];
Flatten[vunion1 /@ r, 1]]

(* turn n-gon into list of triangles *)
toTriangles[vlist_List] /; Length[vlist] == 3 := {vlist}

(* treat 4gons specially for efficiency *)
toTriangles[vlist_List] /; Length[vlist] == 4 :=
	{vlist[[{1, 2, 3}]], vlist[[{3, 4, 1}]]}

(* general case: use center of gravity as additional vertex *)
toTriangles[vlist_List] :=
    Module[{bary = (Plus @@ vlist)/Length[vlist], circ},
        circ = Partition[ Append[vlist, First[vlist]], 2, 1 ];
        Apply[ {#1, #2, bary}&, circ, {1} ]
    ]

parse[poly_] := Cases[poly,Polygon[x_]:>x,Infinity]

volumeSort[poly_]:=
Block[{q1,q2,q3,p},
p=Cases[poly,{Polygon[_]..},Infinity];
q1=Table[{p[[i]], hullvol[ p[[i]] ] },{i,Length[p]}];
q2=Split[Sort[q1,#1[[2]]<#2[[2]]&],
TrueQ[(Abs[#1[[2]]-#2[[2]]]<.0001)]&];
Print["You will need at least ",Length[q2],
" colors in the colorlist below."];(*how many colors?*)
Table[Cases[q2[[i]],{{_,_,_}..},Infinity],{i,Length[q2]}]
]

explode[poly_, dist_]:=
Block[{c,d},
c=Map[ #[[1]]&, poly];
d=Map[parse[#]&, poly];
Table[Polygon /@ Map[(dist*c[[j]])+#&, d[[j]], {2}], {j,Length[c]}]
]

(*the volume of a convex polyhedron maybe not
containing the origin within*)
hullvol[poly_]:=
Block[{p,c,tr,a},
p=parse[poly];
c=cent[Flatten[p,1]];
tr=Map[Chop[#-c]&, p, {2}];(*translate to centered on origin*)
(*first divide faces into triangles*)
a=Flatten[Map[toTriangles[#]&, tr ],1];
(*sum one-sixth the volumes of the determined parallelpipeds*)
N[Plus @@ Flatten[Map[(Abs[Minors[#, 3]]/6)&, a]]]
]

(*find indices of edges into list of vertices*)
edgedex[vertlist_, edgelist_] := Block[{e2}, 
    e2 = Partition[Flatten[Table[
 {Position[vertlist, x_ /; Chop[x-edgelist[[i,1]] ]=={0,0},{1}], 
  Position[vertlist, x_ /; Chop[x-edgelist[[i,2]] ]=={0,0},{1}]},
  {i, Length[edgelist]}]], 2]; 
    Union[Sort /@ e2]]

rmoveto[pt1_,pt2_]:= (pt=pt1+pt2)
rlineto[pt1_,pt2_]:= {pt1,(pt=pt1+pt2)}

rlineto[pt1_,pt2_]:= Line[{pt1,(pt=pt1+pt2)}]

(*remove coincident faces; by Xah Lee*)
dePair[poly_]:=
Block[{a},
a=Sort[((Chop[Rationalize[(cent@#),10^-3]])@#&/@ (poly)),
OrderedQ[{Head@#1,Head@#2}]&];
Nest[
	(
	(DeleteCases[#,Null,{1}]&) @ (Flatten[#,1]&) @
	((If[SameQ@@(Head/@#),Null,#]&) /@
		Partition[Join[RotateRight[#],{Null,Null,Null}],2])
	)&,a,2]/.({___}[li_List]:>li)
]

zonotile3D

This awkward function requires two arguments: plns, the list of planes, and d, the number of parallel planes. The argument d is left global from the creation of the planes themselves. Planes are given as unit vectors {x, y, z} with a distance from the origin: {x, y, z, distance}. The set of sign vectors is found and first reduced in length, usually significantly, by Union. Then, if d is greater than 1, a further reduction strips away degenerate polyhedra (in which all the faces lie in a single plane). Finally, duplicate zonohedra are discarded.

There are two global outputs to the function: allpoly and polyflat. In allpoly we have a list of lists of zonohedra; the sublists are formed from zonohedra having the same number of faces. These are usually not zonohedra of the same shape. I find it instructive to display subsets of the zonohedra composing a zonotile; see the examples. The other global output, polyflat, is just a flattened list of all the zonohedra in the zonotile. This form is useful for exploding the zonohedra apart, or sorting them by volume.

Global polyflat actually has individual members in the form {center of zonohedron, zonohedron}. This is useful when displaying subsets or exploding the zonotile.

It is easy to make this function fail; for instance, if the planes are from the icosahedral or octahedral symmetry groups, they must be “balanced” across the origin. This is just one of many problems begging for solutions. It seems to hinge upon the ordering of the planes, but I am not entirely sure of this.

zonotile3D[plns_]:=
Block[{sgn1,s,s1,s2,plns2,r,pt,zaq,
t1,npt,vecs,c,t3,pl,tl,z,xaq,polyunion},

r=Length[plns];

(*****MAKE SIGN VECTORS*****)
Share[sgn1=Flatten[Table[
      Sign[Chop[Det[ plns[[{i,j,k,h}]] ] ] ],
    {i,r-2},
    {j,i+1,r-1},
    {k,j+1,r},
    {h,r}],2];
    Print["Before union, ", Length[sgn1], " sign vectors"];
    Print["After union, ", Length[s=Union[sgn1]], " sign vectors"];
    
    (*to obtain degenerate zonohedra, comment out this part*)
    If[d>1,(*if parallel planes exist*)
    s1=Map[Partition[#,d]&, s];(*break sign vectors into lists of length d*)
s2=DeleteCases[s1, {___,Table[0,{d}],___} ];(*get rid of degenerates*)
s=Map[Flatten[#,1]&, s2];(*regroup list of sign vectors*)
Print["After discard of degenerates, ", Length[s],  " sign vectors"];
]
];

(*****CREATE ZONOTILE*****)
(*create list {center of zonohedron, zonohedron}*)
plns2=Map[#[[{1,2,3}]]&, plns];

Share[
zaq=Table[
pt={0,0,0};
Do[If[ s[[k,j]]==1,
     rmoveto[ pt, plns2[[j]] ],
     rmoveto[ pt, -plns2[[j]] ]
     ],
 {j,r}];

t1=Cases[
Table[
   If[ s[[k,j]]==0,
     rlineto[ pt,2*plns2[[j]] ]
   ],
   {j,r}],{{_,_,_}..},Infinity];
   
   npt=t1[[1,1]];
   vecs=Map[#[[2]]-#[[1]]&, t1];
   c=Chop[ npt+((Plus@@vecs)/2) ];
   t3=
   If[Length[Union[Map[Count[#,0]&, signvector[vecs]]]]==1,
(*if only parallelograms bound zonohedron*)
   Z[vecs, npt],(*make the zonohedron rooted at npt*)
(*else if higher 2n-gons*)
   Map[c+#&,zedron[vecs],{-2}](*make the zonohedron, translate to position*)
   ];

{c,t3},(*center of zonohedron, zonohedron*)
{k,Length[s]}]];


(*****POST PROCESSING AND DISPLAY*****)
Print["There are ",Length[zaq]," zonohedral tiles with duplicates"];
pl=Map[Length[#[[2]] ]&, zaq];(*#[[1]] is centers, #[[2]] is zonohedra*)
Print["There are zonohedra with ",tl=Union[pl]," faces"];
Print["Before removing duplicates there are ",
Table[Count[pl,tl[[k]]],{k,Length[tl]}], " zonohedra in each set"];

(*make lists of tiles by number of faces*)
z=Table[Select[zaq,Length[ #[[2]] ]==tl[[k]]&],
{k,Length[tl]}];
xaq=Drop[z,1];(*we don't need to worry about parallelepipeds*)

(*remove duplicate higher zonohedra, keeping poly centers*)
polyunion=Join[
{z[[1]]},
Table[
hedraunion[ xaq[[k]] ],
{k,Length[xaq]}]
];

(*the unique zonohedra, ordered by number of determining vectors*)
allpoly=Map[#[[2]]&, polyunion, {2} ];(*global*)
Print["There are ", Length[allpoly]," sets altogether"];
Print["There are ",Length /@ allpoly," in each set, after removing duplicates"];

(*flattened list of zonohedra*)
polyflat=Flatten[polyunion,1];(*global*)
Print["There are ", Length[polyflat]," zonohedra altogether"];

](*end Block*)

(*test version*)
newzonotile3D[plns_]:=
Block[{sgn1,s,s1,s2,plns2,r,pt,zaq,
t1,npt,vecs,c,t3,pl,tl,z,xaq,polyunion},

r=Length[plns];

(*****MAKE SIGN VECTORS*****)
Share[sgn1=Flatten[Table[
      Sign[Chop[Det[ plns[[{i,j,k,h}]] ] ] ],
    {i,r-2},
    {j,i+1,r-1},
    {k,j+1,r},
    {h,r}],2];
    Print["Before union, ", Length[sgn1], " sign vectors"];
    Print["After union, ", Length[s=Union[sgn1]], " sign vectors"];
    
    (*to obtain degenerate zonohedra, comment out this part*)
If[d>1,(*if parallel planes exist*)
s=DeleteCases[s, {__,Sequence@@Table[0,{d+1}],__} ];
Print["After discard of degenerates, ", Length[s],  " sign vectors"]
];
];

(*****CREATE ZONOTILE*****)
(*create list {center of zonohedron, zonohedron}*)
plns2=Map[#[[{1,2,3}]]&, plns];

Share[
zaq=Table[
pt={0,0,0};
Do[If[ s[[k,j]]==1,
     rmoveto[ pt, plns2[[j]] ],
     rmoveto[ pt, -plns2[[j]] ]
     ],
 {j,r}];

t1=Cases[
Table[
   If[ s[[k,j]]==0,
     rlineto[ pt,2*plns2[[j]] ]
   ],
   {j,r}],{{_,_,_}..},Infinity];
   
   npt=t1[[1,1]];
   vecs=Map[#[[2]]-#[[1]]&, t1];
     
   vchk1=toTriangles[vecs];
If[ Chop[N[Plus @@ Flatten[Map[(Abs[Minors[#, 3]])&, vchk1]]]]>0,

(*Z[vecs, npt],*)
c=Chop[ npt+((Plus@@vecs)/2) ];
   t3=
   If[Length[Union[Map[Count[#,0]&, signvector[vecs]]]]==1,
(*if only parallelograms bound zonohedron*)
(*Print["using Z"];*)
   Z[vecs, npt],(*make the zonohedron rooted at npt*)
(*else if higher 2n-gons*)
   Map[c+#&,zedron[vecs],{-2}](*make the zonohedron, translate to position*)
   ];

{c,t3},

{}],
{k,Length[s]}]];

zaq=DeleteCases[ zaq, {}];

(*****POST PROCESSING AND DISPLAY*****)
Print["There are ",Length[zaq]," zonohedral tiles with duplicates"];
pl=Map[Length[#[[2]] ]&, zaq];(*#[[1]] is centers, #[[2]] is zonohedra*)
Print["There are zonohedra with ",tl=Union[pl]," faces"];
Print["Before removing duplicates there are ",
Table[Count[pl,tl[[k]]],{k,Length[tl]}], " zonohedra in each set"];

(*make lists of tiles by number of faces*)
z=Table[Select[zaq,Length[ #[[2]] ]==tl[[k]]&],
{k,Length[tl]}];
xaq=Drop[z,1];(*we don't need to worry about parallelepipeds*)

(*remove duplicate higher zonohedra, keeping poly centers*)
polyunion=Join[
{z[[1]]},
Table[
hedraunion[ xaq[[k]] ],
{k,Length[xaq]}]
];

(*the unique zonohedra, ordered by number of determining vectors*)
allpoly=Map[#[[2]]&, polyunion, {2} ];(*global*)
Print["There are ", Length[allpoly]," sets altogether"];
Print["There are ",Length /@ allpoly," in each set, after removing duplicates"];

(*flattened list of zonohedra*)
polyflat=Flatten[polyunion,1];(*global*)
Print["There are ", Length[polyflat]," zonohedra altogether"];

](*end Block*)

Stellation Diagrams

(*Peter Messer's intersection of two planes*)
inter[{{a1_,b1_,c1_},d1_},{{a2_,b2_,c2_},d2_}]:=
Block[{v,sq,DC,DB,AD},
v=Cross[{a1,b1,c1},{a2,b2,c2}];
sq=v.v;
DC= (c1 d2)-(c2 d1);
DB= b1 d2-b2 d1;
AD= a2 d1-a1 d2;
{{(v[[2]]*DC - v[[3]]*DB)/sq, -(v[[1]]*DC + v[[3]]*AD)/sq, (v[[1]]*DB + v[[2]]*AD)/sq},
v}
]

(*intersection of two 2D segments*)
lineinter[{{x1_, y1_},{x2_,y2_}},{{x3_, y3_},{x4_,y4_}}]:=
Block[{num,denom,u,x,y},
num=((x4-x3)*(y1-y3)) - ((y4-y3)*(x1-x3));
denom=Chop[ ((y4-y3)*(x2-x1)) - ((x4-x3)*(y2-y1)) ];
If[ denom==0, {},
u=num/denom;
x=x1+(u*(x2-x1));
y=y1+(u*(y2-y1));
{x,y}]
]

(*find the intersection point, if it exists, of three planes*)
(*planes given as normal vectors {p1, p2, p3},
followed by their magnitudes,{d1,d2,d3}*)
(*I suspect this reverses signs on the points*)
inter3planes[{p1_, p2_, p3_},{d1_, d2_, d3_}]:=
Block[{num, denom,pt},
num=( -d1*(Cross[p2,p3])-d2*(Cross[p3,p1])-d3*(Cross[p1,p2]) );
denom=Chop[( Dot[p1, Cross[p2,p3] ] )];
pt=If[denom==0,Null,Chop[num/denom]]
]

Variants

(*divide-by-zero protected version, Peter Messer's intersection of two planes*)
newinter2planes[{{a1_,b1_,c1_,d1_},{a2_,b2_,c2_,d2_}}]:=
Block[{v,sq,DC,DB,AD},
v=Cross[{a1,b1,c1},{a2,b2,c2}];
sq=v.v;
If[Chop[sq]!=0,
DC= (c1 d2)-(c2 d1);
DB= (b1 d2)-(b2 d1);
AD= (a2 d1)-(a1 d2);
Chop[{{(v[[2]]*DC - v[[3]]*DB)/sq, -(v[[1]]*DC + v[[3]]*AD)/sq, (v[[1]]*DB + v[[2]]*AD)/sq},
v}],Null]
]

(*find the intersection point, if it exists, of three planes*)
newinter3planes[{{a1_,b1_,c1_,d1_},{a2_,b2_,c2_,d2_},{a3_,b3_,c3_,d3_}}]:=
Block[{p1,p2,p3,num,denom,pt},
p1={a1,b1,c1};p2={a2,b2,c2};p3={a3,b3,c3};
num=( d1*(Cross[p2,p3]) + d2*(Cross[p3,p1]) + d3*(Cross[p1,p2]) );
denom=Chop[( Dot[p1, Cross[p2,p3] ] )];
pt=If[denom==0,Null,Chop[num/denom]]
]

Zones

mysurf[poly_]:= 
Block[{f},
f=poly[[1]];
cross[ f[[2]] - f[[1]], f[[3]] - f[[2]] ]
]

mycent[poly_]:= 
Block[{f},
f=poly[[1]];
Plus@@f/Length[f]
]

countsides[shape_] :=
Map[ Length[#]&, parse[shape]]

Zonotile Examples

What is Going On?

My habit in programming is to do things in stages, rather than to wrap a long procedure up into a function and then just let fly at the input, with all due mystery as to how the result is reached. Thus in much of this notebook, several different input cells are evaluated in turn, before a zonotile is ever even displayed; and often, various display variants are offered, separately.

I will not even to attempt to describe the zonotile3D function itself, now. There is a subsection below (Step-by-Step Zonotile Construction) which touches upon how it does what it does. Here, though, it may be useful for those not very familiar with programming in Mathematica to see what is going on, when we create a set of planes, in the form {x, y, z, d}, to use as input to zonotile3D.

Look at the code-block below. Most programming languages provide for ‘comments’, which are not part of the code base itself, but help the programmer document what is happening. In Mathematica we use the syntax (*comment*) for comments. Thus (*****MAKE ARRANGEMENT OF PLANES*****) is a comment, with extra asterisks added for visual impact.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1;(*set Range limits: k=0, 1/2, 1, 3/2,etc.*)
dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];

plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

Another important bit of the Mathematica language is the semicolon. It suppresses non-graphics output. So, to better understand some particular line or lines of Mathematica code, one can delete the semicolons. However, one tiny expression may generate a huge output; for instance, zonotile3D generates lists of lists of polygons, and there may be thousands of these. Hence some care should be taken when removing semicolons. Here is a safe example: evaluate both cells below, and compare the results.

v1=N[Vertices[Icosahedron]];(*with semicolon*)

v1=N[Vertices[Icosahedron]](*without semicolon*)

Now to the next stage in our example: one of each pair of opposite vectors is removed, the result named v, and the length of the list of vectors is named n. We need this n for use in a Table, below. Semicolons are removed, to allow output.

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)]
n=Length[v]

These are the six primitive plane normals, upon which the zonotile will be based. Now we will set our distances. I use the built-in Mathematica function Range to make the distances, as follows; semicolons removed to allow output.

k=1(*set Range limits: k=0, 1/2, 1, 3/2,etc.*)
dist=Table[Range[-k,k,1],{i,n}]
d=Length[dist[[1]]]

Range allows a syntax Range[  lowerlimit, upperlimit, step ] in which one goes from lowerlimit to upperlimit in steps step. Note that we could have obtained the same result by explicitly defining our distances, as follows:

dist=Table[{-1, 0, 1},{i,n}]
d=Length[dist[[1]]]

When trying to understand someone else’s programming, it is always well to consider that some kind of mistake, or error, or oddball construct, might inhabit the code. The above provides two examples: since there are but three distances, {-1, 0, 1}, why create six copies of the three? There is indeed a way to use just one copy. But in the course of developing the code to create sets of planes {x, y, z, d} for input to the zonotile3D function, I hacked together various methods which were, yes, clunky, but in the end, they worked.

The other example is the inscrutable assignment, d=Length[dist[[1]]]. This makes d a global variable, and it is actually used within the zonotile3D procedure, which expects it, and will have some trouble (I seem to think) without it. Here d = 3. And yet here, where a comment in the code would be a good thing, there is no comment. Not only that, since I am always writing about planes in the form {x, y, z, d}, is this d that d? No, it is not. Ah, well, how good it is to be human, and to err.

The final stage in creating the planes, each of the form {x, y, z, d}, is to join the vectors v, each in the form {x, y, z}, to the various distances. This is done using Mathematica’s Table function. The planes are named plns. The raw vectors are converted to unit vectors using my unit function. Semicolons removed, we see the output.

plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j,d}],1]
Print["There are ",r=Length[plns]," planes altogether."]

Print gives output with or without a semicolon. Note that we have to use Flatten on the Table, so that we end up with one list of eighteen planes, rather than three lists of six planes. Below, what happens without Flatten:

Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j,d}]

Note that for each one of the six “primitive” plane normals, the distances -1, 0, and 1 are attached.

When you wish to experiment, and make your own sets of planes to pass to zonotile3D, there are several considerations to keep in mind:

1. The planes must be a simple list of 4-vectors {x, y, z, d}. The components can be integers or real numbers, but avoid exact numbers like Sqrt[2] by wrapping within N, as in N[ Sqrt[2] ].

2. How many planes are created? Any more than thirty, and it will take a while for the zonotile to emerge, and, on my computer, with 100 MB of RAM available to the Mathematica Kernel, more than sixty planes is likely to run the Kernel right out of memory, and take a very long time, perhaps as much as an hour, or more.

3. If any problems arise, such as a zonotile with missing zonohedral cells, or cells that interpenetrate, make sure the planes balance across the origin. That is, for any plane normal {x, y, z}, see that distances a mirrored upon the origin, as in {x, y, z, -1}, and {x, y, z, 1}.

Platonic Zonotiles I

Before delving into more complicated examples, here we construct the simplest zonotiles zonotile3D permits for each of the Platonic solids, in turn. For reasons unknown, zonotile3D must have planes “balanced” across the origin. Although the single distance, zero, is allowed, it will give one undivided zonohedron—a “degenerate” zonotile of edge=1. Except in the case of the cube itself, this zonohedron of edge=1 can be dissected into smaller zonohedra of edge=1; but, not using zonotile3D.

Hence we use the two distances {-1, 1} below, and create zonotiles of edge=2.

Tetrahedron

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Tetrahedron]];(*choose Platonic*)
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[{-1, 1},{i,n}];
d=Length[dist[[1]]];

plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];


(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];


(*Show the first of each of the allpoly*)
If[Length[allpoly]>1,

Do[Show[Graphics3D[ allpoly[[k,1]] ],
Boxed->False],
{k, Length[allpoly]}],

(*else*)
Show[Graphics3D[ allpoly ],
Boxed->False]
];

Octahedron

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Octahedron]];(*choose Platonic*)
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[{-1, 1},{i,n}];
d=Length[dist[[1]]];

plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];


(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];


(*Show the first of each of the allpoly*)
If[Length[allpoly]>1,

Do[Show[Graphics3D[ allpoly[[k,1]] ],
Boxed->False],
{k, Length[allpoly]}],

(*else*)
Show[Graphics3D[ allpoly ],
Boxed->False]
];

Cube

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Cube]];(*choose Platonic*)
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[{-1, 1},{i,n}];
d=Length[dist[[1]]];

plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];


(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];


(*Show the first of each of the allpoly*)
If[Length[allpoly]>1,

Do[Show[Graphics3D[ allpoly[[k,1]] ],
Boxed->False],
{k, Length[allpoly]}],

(*else*)
Show[Graphics3D[ allpoly ],
Boxed->False]
];

Icosahedron

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[{-1, 1},{i,n}];
d=Length[dist[[1]]];

plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];


(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];


(*Show the first of each of the allpoly*)
If[Length[allpoly]>1,

Do[Show[Graphics3D[ allpoly[[k,1]] ],
Boxed->False],
{k, Length[allpoly]}],

(*else*)
Show[Graphics3D[ allpoly ],
Boxed->False]
];

Pentagonal Dodecahedron

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Dodecahedron]];(*choose Platonic*)
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[{-1, 1},{i,n}];
d=Length[dist[[1]]];

plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];


(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];


(*Show the first of each of the allpoly*)
If[Length[allpoly]>1,

Do[Show[Graphics3D[ allpoly[[k,1]] ],
Boxed->False],
{k, Length[allpoly]}],

(*else*)
Show[Graphics3D[ allpoly ],
Boxed->False]
];

Platonic Zonotiles II

Discussion

The Platonic solids give interesting zonotiles. Here we generalize upon their construction, and by choosing, first, which Platonic solid we wish to operate upon, and, second, at what distance(s) we will locate our planes, we obtain zonotiles in edge={1, 2, 3, ...}.

We use the vertices of the various Platonics, named in Mathematica Tetrahedron, Octahedron, Cube (or Hexahedron), Icosahedron, and Dodecahedron. In all cases except the Tetrahedron, vertices occur in diametrically opposed pairs, and we discard one of each pair. Practically speaking, we are creating those zonohedra (the zonotile-at-large, variously tiled by smaller zonohedra) determined by the vectors from the centers to the vertices of the Platonic solids, which are as follows:

Tetrahedron (four vectors) -> Kepler’s Rhombic Dodecahedron.
Octahedron (three vectors) -> Cube.
Cube (or Hexahedron) (four vectors) -> Kepler’s Rhombic Dodecahedron.
Icosahedron (six vectors) -> Kepler’s Rhombic Triacontahedron.
Dodecahedron (ten vectors) -> Rhombic Enneacontahedron.

v1=N[Vertices[Tetrahedron]];
vertstet=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
Print["There are ", Length[v1], " vertices."];
Print["Tetrahedron->Kepler's Rhombic Dodecahedron"];
Show[Graphics3D[zedron[vertstet]]];

v1=N[Vertices[Octahedron]];
vertsoct=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
Print["There are ", Length[v1], " vertices."];
Print["Octahedron->Cube"];
Show[Graphics3D[zedron[vertsoct]]];

v1=N[Vertices[Cube]];
vertscube=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
Print["There are ", Length[v1], " vertices."];
Print["Cube->Kepler's Rhombic Dodecahedron"];
Show[Graphics3D[zedron[vertscube]]];

v1=N[Vertices[Icosahedron]];
vertsicosa=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
Print["There are ", Length[v1], " vertices."];
Print["Icosahedron->Kepler's Rhombic Triacontahedron"];
Show[Graphics3D[zedron[vertsicosa]]];

v1=N[Vertices[Dodecahedron]];
vertsdod=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
Print["There are ", Length[v1], " vertices."];
Print["Dodecahedron->Rhombic Enneacontahedron"];
Show[Graphics3D[zedron[vertsdod]]];

Since here 3D zonotiles arise from arrangements of planes, in the list above we should substitute the word ‘planes’ for the word ‘vectors’, and recognize that it is the face planes of the dual which actually determine the zonotile. The rhombic enneacontahedron is, under this aegis, better related to the stellation of the icosahedron, than to the vertices of the pentagonal dodecahedron—but clearly both are involved.

Having selected a Platonic solid, it remains to choose the distances from the origin at which the planes are set. This is effected using Range.

First, consider that, due to the mysterious (to me) behavior of the Platonic symmetry groups, planes must be equally “balanced” across the origin. For instance, distances of {-1, 1} are permissible, and would give a zonotile of edge=2. The single distance of {0} is also permissible, and gives a zonotile of edge=1, although “degenerate” in that it is not dissected into smaller zonohedra. However, “unbalanced” distances, such as {-2, 1, 0, 1}, will lead to a failure of some kind, and an improper zonotile will arise.

Suppose that, as below, the “step” specification in Range is fixed at 1. Then Range[-0, 0, 1] gives the one distance {0}, Range[-.5, .5, 1] gives the two distances {-.5, .5}, Range[-1, 1, 1] gives the three distances {-1, 0, 1}, Range[-1.5, 1.5, 1] gives the four distances {-1.5, -.5, .5, 1.5}; and so on. Note that whenever there is an even number of distances, there is no distance of 0 among them; but whenever there is an odd number of distances, there is a distance of 0 among them.

Therefore, under the rules of construction followed here, whenever there are an odd number of distances, so that one set of planes intersects at that one point, the origin, {0,0,0}, there is a smaller version of the zonotile-at-large in the very center of the zonotile (or, if there is but the one distance, {0}, the zonotile is but that one zonohedron determined by the Platonic vectors, centered upon the origin).

Thus, in the “Nested Triacontahedra” example below, the Platonic solid is the Icosahedron, and Range variable k is set to equal 1; hence the distances are {-1, 0, 1}, one set of planes passes through the origin, the zonotile-at-large is Kepler’s Rhombic Triacontahedron in edge=3, and at the very center of the zonotile is a smaller triacontahedron, of edge=1.

The “Randomized Triacontahedra” example demonstrates use of multiple random distances, still balanced across the origin, but forbidding any intersections of more than three planes; hence only hexahedra arise within the zonotile.

Perhaps the most interesting of all of these is the rhombic enneacontahedron, derived from the Platonic Dodecahedron. With its ten planes, and k=1, so that each plane lies at three distances, thirty planes altogether, the creation of the zonotile is slow, but worth every second (or minute, in my case, with a slow computer). It is especially nice to look at the subset of zonohedra close to the z=0 plane, so that the central enneacontahedron itself is visible.

Recalling my conjectures about “meta-geometry,” it is fruitful to think of the vertices of the pentagonal dodecahedron as the face-planes of the icosahedron, and vice versa, and so on with the other Platonics. Hence in the example below (if k is set to 1/2) we should find bijection between the three stellations of the dodecahedron, and the three “layers” of zonohedra, from the center to the periphery, of that zonotile.

Example: Nested Triacontahedra

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1;(*set Range limits: k=0, 1/2, 1, 3/2,etc.*)
dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];

plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{3}]] ],
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Boxed->False];

(*select all zonohedra within a given radius from the origin*)
radius=6;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[ shell ],
ViewPoint->{0,-1.5,1},
Boxed->False];

Example: Randomized Triacontahedra

I did not realize until recently that, in balancing distances across the origin for Platonic zonotiles, any random distances can be used, so long as for any one primitive plane normal {x, y, z}, we use both d and -d; that is, in this example, we may create, for instance, twelve different distances, two for each primitive, with negatives, four for each primitive, and twenty-four altogether. I have uncommented the distance-list dist so that these twelve random distances can be seen.

This randomization ensures that in the resulting zonotile (a triacontahedron of edge=4), only rhombic hexahedra arise, 1280 of them.

However, the combinatorics of the plane arrangement can still be influenced. For instance, let the inner pair of random distances be within the range {-1, 1}, and let the outer pair be greater than, oh, {-5, 5}. Then a triacontahedron of edge=2 arises at the center of the triacontahedron of edge=4, as we might expect. Also, holding the inner distances to a narrow range (in the example, the random numbers are constrained between .5 and .6) seems to ensure that twenty equal hexahedra make a star at the very center of the zonotile. This must mean that those twenty points of intersection near the center of the plane arrangement are close enough to those generated by all equal distances, that the same group of zonohedra develops.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[
k=Random[Real,{.5,.6}];
i=Random[Real,{1.3,1.4}];
 {-i,-k, k,i}, {n}]
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Boxed->False];

Example: Randomized Truncated Icosidodecahedron

Here we use the edge-centers of the Icosahedron (or Dodecahedron) to obtain the fifteen unique, non-collinear vectors which determine the Archimedean Truncated Icosidodecahedron. By using different random real numbers for the two balanced distances for each of the fifteen primitives, and adding the third balanced distance 0, we construct a zonotile determined by forty-five planes. Fifteen intersect at the origin and determine an inner Truncated Icosidodecahedron of edge=1, while the zonotile-at-large is also a Truncated Icosidodecahedron, of edge=3. We might expect, from the randomization of distances, that the entire region outside the inner Truncated Icosidodecahedron would be tiled with hexahedra. However, that would be an improper tiling, and instead we obtain a proper tiling, in which arrays of prisms sit atop each face of the inner Truncated Icosidodecahedron, and the gaps between these are sealed with hexahedra. Of course, the thirty square faces of the inner Truncated Icosidodecahedron support various square prisms, which are hexahedra themselves.

This zonotile, generated by 45 planes, takes a long while to construct. Warning: trying to display the entire sublist of hexahedra, allpoly[[ 4 ]], may cause problems. There are many many hexahedra in face-to-face contact. This places great demands upon the Mathematica rendering engine.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=Map[cent[#]&, edgefind[N[Vertices[Icosahedron]]]];

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[
k=Random[Real,{.5,.6}];
{-k, 0, k}, {n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the inner zonohedron*)
Show[Graphics3D[ allpoly[[{4}]] ],
Boxed->False];

(*Show the inner zonohedron and the decagonal prisms*)
Show[Graphics3D[ allpoly[[{3,4}]] ],
Boxed->False];

(*Show the inner zonohedron and the hexagonal prisms*)
Show[Graphics3D[ allpoly[[{2,4}]] ],
Boxed->False];

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2.1;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Boxed->False];

Example: Randomized Truncated Cuboctahedron

Here, an example quite analogous to that above, except we are in the octahedral symmetry group, and set planes at five “types” of distances. There are nine primitive plane normals, formed from the vertices and edge-centers of the octahedron, and one set passes through the origin, giving a truncated cuboctahedron of edge=1. The other four planes are randomized within narrow limits, and the outer planes are far enough away from the inner, that a distinct medial truncated cuboctahedron of edge=3 occurs, while the zonotile-at-large is a truncated cuboctahedron of edge=5.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=Join[N[Vertices[Octahedron]],Map[cent[#]&, edgefind[N[Vertices[Octahedron]]]]];
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[
i=Random[Real,{15,16}];(*outer plane distances*)
k=Random[Real,{.5,.6}];(*inner plane distances*)
{-i, -k, 0, k, i},
{n}];
d=Length[dist[[1]]];

plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the inner zonohedron*)
Show[Graphics3D[ allpoly[[{4}]] ],
Boxed->False];

(*Show the inner zonohedron and the octagonal prisms*)
Show[Graphics3D[ allpoly[[{3,4}]] ],
Boxed->False];

(*Show the inner zonohedron and the hexagonal prisms*)
Show[Graphics3D[ allpoly[[{2,4}]] ],
Boxed->False];

(*Show the inner zonohedron and the hexagonal prisms*)
Show[Graphics3D[ allpoly[[ 1 ]] ],
Boxed->False];

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2.1;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Boxed->False];

Platonic Zonotiles III

Discussion

Here we use a different set of Platonic symmetry vectors for our planes: the centers of the edges. Now, the Platonic solids occur in “dual” pairs: Tetrahedron<->inverted Tetrahedron; Octahedron<->Cube; and Icosahedron<->Dodecahedron. When in their dual positions, with respect to one another, and of equal mid-radii, their edges actually intersect at right angles. Hence for our purposes there are only three sets of edges: Tetrahedron, Octahedron/Cube, and Icosahedron/Dodecahedron. The same considerations apply here as in the section above, with regard to distances—they must be balanced across the origin (for instance, {0}, or {-1,1} or {-1,0,1},  are balanced across the origin). From the three Platonic edge-sets we obtain three types of zonotile-at-large:

Tetrahedron (six edges, three vectors) -> Cube.
Octahedron (twelve edges, six vectors) -> Truncated Octahedron.
Icosahedron (thirty edges, fifteen vectors) -> Truncated Icosidodecahedron.

The three zonohedra may be constructed quickly, bypassing the zonotile algorithm (which would give the same result, using the single plane distance {0}).

v1=Map[cent[#]&, edgefind[N[Vertices[Tetrahedron]]]];
edgetet=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
Print["There are ", Length[v1], " edge-centers."];
Print["Tetrahedron->Cube"];
Show[Graphics3D[zedron[edgetet]]];

v1=Map[cent[#]&, edgefind[N[Vertices[Octahedron]]]];
edgeoct=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
Print["There are ", Length[v1], " edge-centers."];
Print["Octahedron->Truncated Octahedron"];
Show[Graphics3D[zedron[edgeoct]]];

v1=Map[cent[#]&, edgefind[N[Vertices[Icosahedron]]]];
edgeicosa=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
Print["There are ", Length[v1], " edge-centers."];
Print["Icosahedron->Truncated Icosidodecahedron"];
Show[Graphics3D[zedron[edgeicosa]]];

It so happens that a rhombic polyhedron may be associated to each of the three dual pairs; for, let the pairs intersect in dual position, so that the edges cross at right angles: then each edge-cross forms the diagonal of a rhomb (and the centers of the edges coincide with the centers of the rhombs). The three rhombic polyhedra are as follows:

Tetrahedron -> Cube (the rhombs are squares).
Octahedron (Cube) -> Kepler’s Rhombic Dodecahedron.
Icosahedron (Dodecahedron) -> Kepler’s Rhombic Triacontahedron.

Hence, the zonotiles determined by the vectors to the edge-centers of the Platonic solids are (equivalently) determined by the face-planes of the associated rhombic polyhedra. Furthermore, if there are but two plane distances, and these balanced across the origin (say, {-1,1}), the arrangement of planes is exactly the face-planes of the associated rhombic polyhedron; it is this same arrangement from which we derive the stellations of that rhombic polyhedron. The stellations of a convex polyhedron, and the zonotile determined by the arrangement of planes given by its faces, are closely related.

The zonotiles derived from the thirty edges of the Icosahedron, giving fifteen planes, will take much time and kernel memory, for three or more distances.

The example below gives a Truncated Octahedron of edge=3, since each plane is at the three distances {-1, 0, 1}.

Example: Nested Truncated Octahedra

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=Map[cent[#]&, edgefind[N[Vertices[Octahedron]]]];
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1;(*set Range limits: k=0, 1/2, 1, 3/2, etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{3}]] ],
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Boxed->False];

(*select all zonohedra within a given radius from the origin*)
radius=7;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-1.5,1},
Boxed->False];

Platonic Zonotiles IV

Here we combine the vertices of the Icosahedron and Dodecahedron. These number thirty-two in all, but we select just one from each pair of opposites, giving sixteen plane-vectors. Note that the whole set of thirty-two, when the two polyhedra are in reciprocal position and have a common mid-radius, form the vertices of Kepler’s Rhombic Triacontahedron. The Platonics as constructed in Mathematica are in reciprocal position, and with equal mid-radii; hence here we operate upon the face-planes of the dual of the Triacontahedron, the Archimedean Icosidodecahedron. This polyhedron can be regarded as the solid of intersection of the Icosahedron and Dodecahedron. When we use but the two distances, {-1/2, 1/2}, as in the example below, we are creating a zonotile closely related to the stellations of the Icosidodecahedron; we are setting our sixteen “types” of plane in pairs, one to either side of the origin, restoring the original number of thirty-two.

One subtlety here: we do not convert the vectors to unit vectors—if we did so, they would not be the proper vertices of the Triacontahedron (which has two different circum-radii), and the face-planes would not be those of the Icosidodecahedron.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[ Join[ Vertices[Icosahedron], Vertices[Dodecahedron] ] ];(*choose Platonics*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1/2;(*set Range limits: k=0, 1/2, 1, 3/2, etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ v[[i]],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{3,4}]] ],
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Boxed->False];

(*select all zonohedra within a given radius from the origin*)
radius=3;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[ shell ],
ViewPoint->{0,-1.5,1},
Boxed->False];

Archimedean Region: Cubes, Truncated Octahedra, and Truncated Cuboctahedra

The thirteen Archimedean solids are bounded by mixtures of different regular polygons, and “Archimedean” space-fillings are built from mixtures of different regular or quasi-regular polyhedra. For instance, the Platonic tetrahedron and octahedron fit together to fill space in one Archimedean space-filling. Some Archimedean packings are made of zonohedra.

Here, the central region of the zonotile forms the Archimedean space-filling of cubes, truncated octahedra, and truncated cuboctahedra. A larger central array may be constructed (takes time!) by changing the line

f*N[Range[-1, 1, 1]], {i, n}];
to
f*N[Range[-2, 2, 1]], {i, n}];.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=Join[N[Vertices[Octahedron]],
Map[cent[#]&, edgefind[N[Vertices[Hexahedron]]]] ];
d1=Union[Map[mag[#]&, v1],SameTest->(Abs[#1-#2]<.001&)];

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];(*discard opposites*)
Print["There are planes of ", n=Length[v], " types"];

dist=Table[ If[ Abs[ mag[ v[[i]] ] - d1[[1]] ]<.0001, f=d1[[1]], f=d1[[2]]];
f*N[Range[-1,1,1]],{i,n}];

Print["The planes are at ", d=Length[dist[[1]]], " distances"];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ", r=Length[plns], " planes altogether"];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*select all zonohedra within a given radius from the origin*)
radius=8;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-.5,2},
Boxed->False];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{5}]] ],
PlotRange->All,
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

Archimedean Region: Octagonal Prisms and Truncated Cuboctahedra

The thirteen Archimedean solids are bounded by mixtures of different regular polygons, and “Archimedean” space-fillings are built from mixtures of different regular or quasi-regular polyhedra, including “uniform” polyhedra. For instance, the Platonic tetrahedron and octahedron fit together to fill space in one Archimedean space-filling. Some Archimedean packings are made of zonohedra.

Here, the central region of the zonotile forms the Archimedean space-filling of octagonal prisms and truncated cuboctahedra. A larger central array may be constructed (takes time!) by changing the line

f*N[Range[-1, 1, 1]], {i, n}];
to
f*N[Range[-2, 2, 1]], {i, n}];.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=Join[N[Vertices[Octahedron]],
Map[cent[#]&, edgefind[N[Vertices[Hexahedron]]]] ];
d1=Union[Map[mag[#]&, v1],SameTest->(Abs[#1-#2]<.001&)];

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];(*discard opposites*)
Print["There are planes of ", n=Length[v], " types"];

dist=Table[ If[ Abs[ mag[ v[[i]] ] - d1[[1]] ]<.0001, f=d1[[2]], f=d1[[1]]];
f*N[Range[-1,1,1]],{i,n}];

Print["The planes are at ", d=Length[dist[[1]]], " distances"];
plns=Flatten[Table[
f=If[EvenQ[i],1,-1];
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d }],1];
Print["There are ", r=Length[plns], " planes altogether"];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*select all zonohedra within a given radius from the origin*)
radius=8;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-.5,2},
Boxed->False];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{5,6,7}]] ],
PlotRange->All,
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

3D “Penrose” Zonotile

The famous Penrose tiling of the plane, by rhombs of two types, has its analogues in three dimensions. Here an n=5 polar zonohedron forms the zonotile, bounded by twenty rhombic regions, each divided into a 5X5 grid of 25 equal rhombs; and it is dissected into 1250 rhombic hexahedra, of two types. Three different viewing options are provided here; and in the Color Zonohedra by Volume section of this notebook, particular colors may be assigned to the different types of zonohedra, to a nice effect. Be sure to experiment with different n. However, use caution; a large (>30) number of planes can arise, even if n is set to 8, and it may be good to make the Range of distances smaller; the total number of planes is the product of n with the number of distances. Hence if n=8 and we use Range[ -2, 2, 1 ], we will have forty planes. This will take both time and kernel memory.

(*****MAKE ARRANGEMENT OF PLANES*****)
n=5;
p=N[(1/Degree)ArcTan[.5]];
v=vectors[n,p];
dist=Table[ N[Random[Real,{.1,.2}]+Range[-2,2,1]], {i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
f=If[EvenQ[i],1,-1];
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ", Length[plns], " planes altogether"];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=1;
shell=Map[#[[2]]&, Select[polyflat, Abs[ #[[1,3]] ]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->{0,0,2},
Boxed->False];

(*explode the selected region slightly*)
q1=explode[Select[polyflat, Abs[#[[1,3]]]<distance&], .2];
Show[Graphics3D[q1],
ViewPoint->{0,0,2},
Boxed->False];(*exploded*)

Non-equilateral Zonotiles

Here we create the same zonotile as in “Archimedean Region: Cubes, Truncated Octahedra, and Truncated Cuboctahedra,” above, but multiply every third vector by four. With some care, a more symmetrical choice of vectors might have been made, for scaling. Howsoever, we do obtain a proper zonotile, of two edge lengths, one four times the other.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=Join[N[Vertices[Octahedron]],
Map[cent[#]&, edgefind[N[Vertices[Hexahedron]]]] ];
d1=Union[Map[mag[#]&, v1],SameTest->(Abs[#1-#2]<.001&)];

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];(*discard opposites*)
Print["There are planes of ", n=Length[v], " types"];

dist=Table[ If[ Abs[ mag[ v[[i]] ] - d1[[1]] ]<.0001, f=d1[[1]], f=d1[[2]]];
f*N[Range[-1,1,1]],{i,n}];

Print["The planes are at ", d=Length[dist[[1]]], " distances"];
plns=Flatten[Table[
sca=If[ Mod[i-1,3]==0, 4, 1];
Chop[Join[ sca*unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ", r=Length[plns], " planes altogether"];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*select all zonohedra within a given radius from the origin*)
radius=18;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-.5,2},
Boxed->False];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{9}]] ],
PlotRange->All,
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

Polar Zonohedron Tiled by Hexahedra

Here we construct the most ordinary of zonotiles, a polar zonohedron tiled by rhombic hexahedra. By making the distances of the planes non-zero random real numbers, we ensure that only rhombic hexahedra arise within the zonotile. In effect, no more than three planes can intersect at any one point.

(*****MAKE ARRANGEMENT OF PLANES*****)
n=12;(*an n=12 polar zonohedron*)
p=N[(1/Degree)*ArcTan[.5^.5]];
hex=Binomial[n,3];
v=vectors[n,p];
Print["There are ",n," planes."];
Print["There are ",hex," hexahedra."];

plns=Table[(*add random real numbers to unit distances*)
Chop[Join[ unit[ v[[i]] ],{Random[Real,{.1,.2}]+1}]],
{i,n}];
d=1;
r=Length[plns];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

(*select all zonohedra within a given radius from the origin*)
radius=3;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->{0,0,2},
Boxed->False];

Polar Zonohedron Tiled by Hexahedra, Dodecahedra, and k-zonohedra

A mixture of hexahedra, dodecahedra, and triacontahedra form this zonotile. The two triacontahedra lie on the symmetry axis of the n=12 polar zonohedron. The device used to obtain this mix is merely to place the planes, alternately, to either side of the origin.

(*****MAKE ARRANGEMENT OF PLANES*****)
n=12;(*n=12 polar zonohedron*)
p=N[(1/Degree)*ArcTan[.5^.5]];
v=vectors[n,p];
Print["There are ",n," planes."];

plns=Table[
f=If[Mod[i,2]==0,1,-1];
Chop[Join[ unit[ v[[i]] ],{f}]],
{i,n}];
d=1;
r=Length[plns];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{3}]] ],
PlotRange->All,
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

Polar Zonohedron Tiled by Hexahedra, Dodecahedra, and k-zonohedra in edge=2

Here, an n=12 polar zonohedron of edge 2 (each face divided into four rhombs) is created, such that two n=12 polar zonohedra of edge 1 are strung along the symmetry axis, like beads. Rhombic dodecahedra and hexahedra make up the rest of the zonotile.

(*****MAKE ARRANGEMENT OF PLANES*****)
n=12;
p=N[(1/Degree)*ArcTan[.5^.5]];
hex=Binomial[n,3];
v=vectors[n,p];
Print["There are ",n," plane types."];

plns=Flatten[Table[
f=If[Mod[i,2]==0,1,-1];
Chop[Join[ unit[ v[[i]] ],{f*j}]],
{i,n},{j,-1,1,2}],1];
d=2;
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{3}]] ],
PlotRange->All,
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

(*select all zonohedra within a given radius from the origin*)
radius=8;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-1.5,1},
Boxed->False];

Truncated Icosidodecahedron in edge=3

This is my favorite 3D zonotile thus far. The bounding zonohedron is the Archimedean Truncated Icosidodecahedron of edge 3. A truncated icosidodecahedron of edge 1 is at the very center of the zonotile. This zonotile takes a long time to develop, and a still more elaborated zonotile can be made using even more time, by setting the plane distances dist to Range[ -2, 2, 1].

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=Map[cent[#]&, edgefind[N[Vertices[Dodecahedron]]]];
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];
dist=Table[N[Range[-1,1,1]],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{5,6}]] ],
PlotRange->All,
Boxed->False];

(*select all zonohedra with centers less than some distance from the y=0 plane*)
distance=3;
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,2]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-2,0},
Boxed->False];

(*select all zonohedra within a given radius from the origin*)
radius=12;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-1.5,1},
Boxed->False];

Rhombic Enneacontahedron in edge=3

The Platonic solids give interesting zonotiles. Here we use ten of the twenty vertices of the Pentagonal Dodecahedron to give a special zonotile which is (taken as a convex hull, not a tiling) an isometric orthogonal shadow, cast into three dimensions, of a ten-dimensional hypercube. By letting the fourth “distance” components of the 4-vectors run through {-1,0,1}, we have planes at three distances, and realize an enneacontahedron of edge=3. Since one set of ten planes passes all through the origin, at the very center of the zonotile is an enneacontahedron of edge=1.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Dodecahedron]];
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];
dist=Table[N[Range[-1,1,1]],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{5}]] ],
PlotRange->All,
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

(*select all zonohedra with centers less than some distance from the y=0 plane*)
distance=3;
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,2]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-2,0},
Boxed->False];

(*select all zonohedra within a given radius from the origin*)
radius=12;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-1.5,1},
Boxed->False];

n=12 polar zonohedron with symmetrical equatorial neighbors

This zonotile has an equatorial region, centered upon the z=0 plane, in which a central n=12 polar zonohedron is surrounded, alternately, by n=6 and n=4 polar zonohedra, after the manner of the plane tiling of dodecagons, hexagons, and squares. It is likely that if more planes were added, by using Range[ -2, -2, 1], this Archimedean area of polar zonohedra would expand.

(*****MAKE ARRANGEMENT OF PLANES*****)
p=N[(1/Degree)*ArcTan[.5^.5]];
v=vectors[12,p];(*this works, one dodeca, six hexa*)
n=Length[v];
dist=Table[If[EvenQ[i],f=-1,f=Sqrt[3]];
N[f*Range[-1,1,1]],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
f=If[EvenQ[i],1,-1];
Chop[Join[ unit[ v[[i]] ],{f*dist[[i,j]]}]],
{i,Length[v]},{j, d }],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
Timing[newzonotile3D[plns]][[1]]

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{3,4,5}]] ],
PlotRange->All,
Boxed->False];

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=1;
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-.5,2},
Boxed->False];

2n-gonal prisms

Here we test a new zonotile3D function, named, appropriately, newzonotile3D. It differs from the old function in that it uses a different method to weed out degenerate sign vectors before ever the sign vectors are passed to the loops which create the individual zonohedra. Then, since experience has shown that degenerates creep in regardless, each zonohedral tile is tested for non-zero volume before being admitted to the final list. This new function succeeds in properly constructing some zonotiles for which the old function failed. It still fails in some cases for unknown reasons; the principal failure seems to be creation of tiles (zonohedra) which interpenetrate, or, rarely, the outright failure to create a tile where one is needed.

Various types of “prismatic” zonotiles are created, closely related to 2D zonotiles.

Star Prisms

(*Recipe for polar projection of an n-polar zonohedron*)
(*****MAKE ARRANGEMENT OF PLANES*****)
n1=7;(*set n*)
n=If[EvenQ[n1],n1/2,n1];
dist=If[EvenQ[n1],{-1,1},{1}];

v=Join[halfvectors[n,0],{{0,0,1}}];
Print["There are ",Length[v]," types of planes."];

rawplns=Flatten[
	Table[
		sgn=If[EvenQ[i],1,-1];
		Chop[Join[ unit[ v[[i]] ], {sgn*dist[[j]]} ]],
	{i,Length[v]},{j, Length[dist]}],1];

plns=If[EvenQ[n1], Drop[rawplns, -1], rawplns];

d=Length[dist];

(*****CREATE ZONOTILE*****)
newzonotile3D[plns]

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
ViewPoint->2.5*unit[{0,-.5,1}],
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

Tall Prism I

(*Recipe for prismatic 2D zonotile of height=5*)
(*****MAKE ARRANGEMENT OF PLANES*****)
n=8;
dist=Table[If[EvenQ[i], f=1, f=Cos[Pi/n] ];
N[f*Range[-1,0,1]],{i,n}];
d=Length[dist[[1]]];

plns1=Flatten[Table[
Chop[N[
{ Cos[i*Pi/n], Sin[i*Pi/n], 0, dist[[i,j]] }
]],
{i,n}, {j, Length[dist[[1]]]}],1];
plns=Join[plns1, {{0,0,5,1}}];(*this choice of z determines height*)
Print["There are ", r=Length[plns]," planes altogether"];

(*****CREATE ZONOTILE*****)
newzonotile3D[plns];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ],
Axes->True,
ViewPoint->2.5*unit[{0,-2,3}]];

Tall Prism II

(*many of this type fail, this succeeds*)
n=5;
dist=Table[If[i==1 || i==2, f=1, f=2 ];
N[(f)*Range[-2,1,1]],{i,n}];
d=Length[dist];
plns1=Flatten[Table[
Chop[N[
{ Cos[i*Pi/n], Sin[i*Pi/n], 0, dist[[i,j]] }
]],
{i,n}, {j, Length[dist[[1]]]}],1];
plns=Join[plns1, {{0,0,5,1}}];(*this choice of z determines height*)
Print["There are ", r=Length[plns]," planes altogether"];

(*****CREATE ZONOTILE*****)
newzonotile3D[plns];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ],
Axes->True,
ViewPoint->2.5*unit[{0,-2,3}]];

Polar Zonohedra with equatorial hexagons

Here we add one vector {0,0,-1} to the ordinary n vectors which determine a polar zonohedron. When n is even, this stretches the equatorial rhombs into hexagons. It is a curiosity that, in this construction, two ordinary n-polar zonohedra lie on the symmetry axis, and nestle into dodecahedra of mixed rhombs and hexagons at the very center of the zonotile-at-large. This central array may be viewed by itself, in the last example, below.

(*****MAKE ARRANGEMENT OF PLANES*****)
n=8;
p=N[(1/Degree)*ArcTan[.5]];
v=vectors[n,p];(*this works, one dodeca, six hexa*)
dist=Table[If[EvenQ[i],f=-1,f=1];
N[Range[-1,1,2]],{i,n}];
d=Length[dist[[1]]];
plns1=Flatten[Table[
f=If[EvenQ[i],1,-1];
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
plns=Join[plns1, {{0,0,-1,0}}];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
newzonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

Show[Graphics3D[ allpoly[[2]] ],
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

(*select all zonohedra within a given radius from the origin*)
radius=4;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->{0,-2,1},
Boxed->False];

Nested and Beaded Zonotiles

By cunningly combining planes at certain distances, one obtains “nested” zonohedra. Here we use the six vectors which determine an n=6 polar zonohedron—a triacontahedron bounded by thirty rhombs—as our plane unit normals, and let distances run over {-22, -1, 0, 1, 22}. Thus the bounding zonohedron—the zonotile-at-large—is a triacontahedron. Since there is an odd number of distances, the central distance, 0, means that there will be an undivided triacontahedron of edge=1 at the center of the zonotile. Since there are five distances, the outermost triacontahedron has edge=5. Between the innermost triacontahedron of edge=1, and the all-bounding triacontahedron of edge=5, is a middle triacontahedron of edge=3. The distances are set somewhat arbitrarily; here, having one distance=0 forces a central zonohedron; but the other distances, 1 and 22, were chosen merely to ensure (at a guess) that the outermost planes at distance=22 would not intersect any of the planes at distance=0 or distance=1 in a point. If we had used distances {-1000, -.1, 0, .1, 1000}, we would have obtained the identical nested zonotile.

By the way, it is merely fortuitous that the “middle” zonotiles are selected and shown so easily. In other examples, one would either drop component zonohedra, or admit zonohedra which did not belong.

In the second example, “beaded” zonotiles are created.

Nested Polar Zonohedra

(*****MAKE ARRANGEMENT OF PLANES*****)
n=6;
p=N[(1/Degree)*ArcTan[.5^.5]];
v=vectors[n,p];
Print["There are ",n," plane types."];

dist=Table[{-22, -1, 0, 1, 22},{i,n}];(*set distances here*)

d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

Equatorial Region

With a little hard staring, one can see all three triacontahedra in this equatorial quasi-section of the zonotile.

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=1.25;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Axes->True,
Boxed->True];

Innermost Triacontahedron

(*select all zonohedra within a given radius from the origin*)
radius=1;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Middle Triacontahedron

(*select all zonohedra within a given radius from the origin*)
radius=9;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Outermost Triacontahedron

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Nested Truncated Octahedra

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=Map[cent[#]&, edgefind[N[Vertices[Octahedron]]]];
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];
Print["There are ",n," plane types."];

dist=Table[{-22, -1, 0, 1, 22},{i,n}];(*set distances here*)

d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

Equatorial Region

It is difficult to see all three truncated octahedra in this equatorial quasi-section of the zonotile.

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=1.25;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Axes->True,
Boxed->True];

Innermost Truncated Octahedron

(*select all zonohedra within a given radius from the origin*)
radius=1;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Middle Truncated Octahedron

(*select all zonohedra within a given radius from the origin*)
radius=9;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Outermost Truncated Octahedron

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Beaded

Here again we turn to the vectors which determine an n=6 polar zonohedron. Were these planes derived from one of the fiendish symmetry groups of the Platonic solids, we could not admit the distances {0, 1}, but must have used {-1, 1} to achieve the same result.

Please compare the three distances {0, 1, 100}, with {0, 1, 2}. You will find that {0, 1, 100} gives three triacontahedra beaded along the axis of 6-fold symmetry, while {0, 1, 2} gives the three beads, and six others framing in the central bead.

(*****MAKE ARRANGEMENT OF PLANES*****)
n=6;
p=N[(1/Degree)*ArcTan[.5^.5]];
v=vectors[n,p];
Print["There are ",n," plane types."];

dist=Table[{0, 1},{i,n}];(*set distances here*)

d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

Bead Triacontahedra

Two equal triacontahedra are on the axis of 6-fold symmetry.

Show[Graphics3D[ allpoly[[{3}]] ],
ViewPoint->2.5*unit[{0,-2,0}],
Axes->True,
Boxed->True];

Framing Dodecahedra

Twenty-four dodecahedra frame the triacontahedra (not shown).

Show[Graphics3D[ allpoly[[{2}]] ],
ViewPoint->2.5*unit[{0,-1,2}],
Axes->True,
Boxed->True];

Equatorial Region

Twelve of the twenty-four dodecahedra are visible in this equatorial quasi-section of the zonotile. At the center, the “cup” into which one of the two triacontahedra fits.

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=1.25;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Axes->True,
Boxed->True];

Central Array of Dodecahedra

(*select all zonohedra within a given radius from the origin*)
radius=3;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Outermost Triacontahedron

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Random Zonotiles

Here we make a random assortment of planes, and use the same device as with the Platonics to insulate the zonotile3D procedure from improper, intersecting zonohedral tiles: we balance the planes across the origin. We hold our randomized z coordinates within narrower bounds than x or y. Every time the cell below is evaluated, a different zonotile arises.

(*****MAKE ARRANGEMENT OF PLANES*****)
n=5;(*number of plane-types*)
v=Table[{Random[Real,{-1,1}],Random[Real,{-1,1}], Random[Real,{.5,.7}]},{n}];
Print["There are ",n," plane types."];

dist=Table[Range[-1,1,2],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

(*select all zonohedra within a given radius from the origin*)
radius=3;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Step-by-Step Zonotile Construction

Here we examine output from the various steps in an especially simple zonotile creation: an n=4 polar zonohedron, composed by four rhombic hexahedra. We create four planes, as four perpendicular vectors of the form {x, y, z, distance}, where distance is the distance of the plane from the origin. Knowing beforehand that only hexahedra will arise, we use Binomial[ n, 3] to find their number

(*****MAKE ARRANGEMENT OF PLANES*****)
n=4;(*an n=4 polar zonohedron*)
p=N[(1/Degree)*ArcTan[.5^.5]];
hex=Binomial[n,3];
v=vectors[n,p];
Print["There are ",n," planes."];
Print["There are ",hex," hexahedra."];

plns=Table[
Chop[Join[ unit[ v[[i]] ],{i}]],
{i,n}]
d=1;
r=Length[plns];

The sign vectors are created.

(*****MAKE SIGN VECTORS*****)
Share[sgn1=Flatten[Table[
      Sign[Chop[Det[ plns[[{i,j,k,h}]] ] ] ],
    {i,r-2},
    {j,i+1,r-1},
    {k,j+1,r},
    {h,r}],2];
    Print["Before union, ", Length[sgn1], " sign vectors"];
    Print["After union, ", Length[s=Union[sgn1]], " sign vectors"];
    
    (*to obtain degenerate zonohedra, comment out this part*)
    If[d>1,(*if parallel planes exist*)
    s1=Map[Partition[#,d]&, s];(*break sign vectors into lists of length d*)
s2=DeleteCases[s1, {___,Table[0,{d}],___} ];(*get rid of degenerates*)
s=Map[Flatten[#,1]&, s2];(*regroup list of sign vectors*)
Print["After discard of degenerates, ", Length[s],  " sign vectors"];
]
];
s

We form the zonotile, using just {x, y, z} of our original planes {x,  y, z, distance}. After postprocessing, we obtain two lists as output: allpoly and polyflat. The former is a list of lists of polygons—in this case, each cell of the zonotile is a rhombic hexahedron, there are four such cells, and therefore four lists of six rhombic faces each. The latter list, polyflat, is in the former {center of zonohedron, zonohedron}; pre-calculating the centers facilitates exploding a zonotile apart.

(*****CREATE ZONOTILE*****)
(*create list {center of zonohedron, zonohedron}*)
plns2=Map[#[[{1,2,3}]]&, plns];
r=Length[plns2];
Timing[
Share[
zaq=Table[
pt={0,0,0};
Do[If[ s[[k,j]]==1,
     rmoveto[ pt, plns2[[j]] ],
     rmoveto[ pt, -plns2[[j]] ]
     ],
 {j,r}];

t1=Cases[
Table[
   If[ s[[k,j]]==0,
     rlineto[ pt,2*plns2[[j]] ]
   ],
   {j,r}],{{_,_,_}..},Infinity];
   
   npt=t1[[1,1]];
   vecs=Map[#[[2]]-#[[1]]&, t1];
   c=Chop[ npt+((Plus@@vecs)/2) ];
   t3=
   If[Length[Union[Map[Count[#,0]&, signvector[vecs]]]]==1,
(*if only parallelograms bound zonohedron*)
   Z[vecs, npt],(*make the zonohedron rooted at npt*)
(*else if higher 2n-gons*)
   Map[c+#&,zedron[vecs],{-2}](*make the zonohedron, translate to position*)
   ];

{c,t3},(*center of zonohedron, zonohedron*)
{k,Length[s]}]] ][[1]]


(*****POST PROCESSING AND DISPLAY*****)
Print["There are ",Length[zaq]," zonohedral tiles with duplicates"];
pl=Map[Length[#[[2]] ]&, zaq];(*#[[1]] is centers, #[[2]] is zonohedra*)
Print["There are zonohedra with ",tl=Union[pl]," faces"];
Print["Before removing duplicates there are ",
Table[Count[pl,tl[[k]]],{k,Length[tl]}], " zonohedra in each set"];

(*make lists of tiles by number of faces*)
z=Table[Select[zaq,Length[ #[[2]] ]==tl[[k]]&],
{k,Length[tl]}];
xaq=Drop[z,1];(*we don't need to worry about parallelepipeds*)

(*remove duplicate higher zonohedra, keeping poly centers*)
polyunion=Join[
{z[[1]]},
Table[
hedraunion[ xaq[[k]] ],
{k,Length[xaq]}]
];

(*the unique zonohedra, ordered by number of determining vectors*)
allpoly=Map[#[[2]]&, polyunion, {2} ];(*global*)
Print["There are ", Length[allpoly]," sets altogether"];
Print["There are ",Length /@ allpoly," in each set, after removing duplicates"];

(*flattened list of zonohedra*)
polyflat=Flatten[polyunion,1];(*global*)
allpoly (*output the lists of zonohedra*)
polyflat (*output the flattened list of centers and zonohedra*)

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ],Axes->True];

(*animate explosion of the whole zonotile*)
r=4;(*maximum plotrange*)
Do[
Show[Graphics3D[ explode[ polyflat, i] ],
PlotRange->{{-r,r},{-r,r},{-r,r}}],
{i, 0, 2, .1}];

One Failure of  Zonotile Construction

Here we seem to be on the right track to creating a rhombic triacontahedron in edge=1. We allow output of the six planes, each of the form {x, y, z, d}, where d is distance from the origin. All the {x, y, z} are unit vectors. All are in the same half-space, z<0. Yet interpenetrating rhombic icosahedra arise, along with ten rhombic hexahedra, where I myself would expect the output to consist of only twenty rhombic hexahedra, of two types.

(*****MAKE ARRANGEMENT OF PLANES*****)
n=5;(*an n=12 polar zonohedron*)
p=N[(1/Degree)*ArcTan[.5]];
hex=Binomial[n,3];
v=vectors[n,p];

plns1=Table[
Chop[Join[ unit[ v[[i]] ],If[EvenQ[i], {-1}, {1}] ]],
{i,n}];
d=1;
plns=Join[plns1, {{0,0,-1,-1}}]

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*a subset of the whole zonotile*)
Show[Graphics3D[ allpoly[[2]] ],
Boxed->False];

Bitten Zonotiles

This is an experiment to see whether “bitten” 3D zonotiles can be constructed after the fashion of bitten 2D zonotiles. My conjecture is that a bitten polar zonohedron ought to be a rhombic spirallohedron. So far I am unable to actually remove the bites, or even be sure if they are properly formed. The variable bitenumber, below, must be a divisor of n (but greater than 1 and less than n). For instance, when n=12, we can have bitenumbers 2, 3, 4, and 6. When the bitenumber is 2, and were the bites removed, we would (in theory) have nothing left. For the moment, unable to actually take the two bites, the packing of zonohedra within the zonotile does exhibit a kind of 2-fold symmetry.

In two dimensions, perfect bitten zonogons are obtained by setting the distances of the planes to Mod[i-1,kmod]-((kmod-1)/2) within the loop. In three dimensions, better results seem to arise with Mod[i-1,kmod]^2.

(*****MAKE ARRANGEMENT OF PLANES*****)
n=12;(*bounded by an n=12 polar zonohedron*)
bitenumber=4;
kmod=n/bitenumber;

d=1;
plns=Table[
Chop[N[{ Cos[(2*i*Pi/n)], Sin[(2*i*Pi/n)], N[Sqrt[.5]], Mod[i-1,kmod]^2 } ]],
{i,n}];

Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

Show[Graphics3D[ allpoly ],
ViewPoint->2.5*unit[{0,-1,0}],
Axes->True,
Boxed->True];

Show[Graphics3D[ allpoly[[2]] ],
ViewPoint->2.5*unit[{0,0,1}],
Axes->True,
Boxed->True];

Equatorial Region

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=1.0;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Axes->True,
Boxed->True];

Bitten Intersection Diagrams

Here we create a “stellation diagram” for the same set of 4-vectors, plns, which is used as input to the zonotile3D function. The example below uses the planes for a 3-bite, n=9 zonotile. The zonotile is created, and an equatorial quasi-section displayed. Then we proceed through the same principle steps followed in the “Stellation Diagrams” section of the notebook, but use somewhat modified versions of the functions which find the lines of intersection of two planes, or the points of intersection of three planes. The modifications were made to allow easier use of the 4-vectors {x, y, z, d}, but as a result some awkwardness arises elsewhere, for instance, ViewPoint->2.5*unit[ rawplanes[[ 1, {1,2,3} ]] ].

To see all the lines and points of intersection of the planes helps to understand the combinatorics of these “bitten” zonotiles. Some other integer pairs of {n, bitenumber} to try include {6, 3}, {8, 4}, {10, 2}, {10, 5}, {12, 2}, {12, 3}, {12, 4}, and {12, 6}. However, when n is any greater than 10, the lines and points are many, probably too many.

Since every zonotile example in this notebook uses a plane arrangement named plns, one can bring the planes from some other zonotile here, and see how they intersect.

Here, the default ViewPoint is on the positive z axis, since this is the axis of n-fold symmetry of the n vectors. This perspective allows the “bite-fold” symmetry to emerge into view, as well. It may help to crop and enlarge the central part of the array of all lines and all points of intersection, to see this bite-fold symmetry.

Create Planes, Zonotile

(*****MAKE ARRANGEMENT OF PLANES*****)
n=9;(*bounded by an n=9 polar zonohedron*)
bitenumber=3;(*three bites*)
kmod=n/bitenumber;

d=1;
plns=Table[
Chop[N[{ Cos[(2*i*Pi/n)], Sin[(2*i*Pi/n)], N[Sqrt[.5]], Mod[i-1,kmod]+((kmod-1)/2) } ]],
{i,n}];

Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

Show Equatorial Region

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=1.0;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Axes->True,
Boxed->True];

Find Some Lines and All Points of Intersection

rawplanes=plns;
(*get rid of any parallel to the first*)
planes=Select[rawplanes, Chop[ rawplanes[[1,{1,2,3}]]+#[[{1,2,3}]] ]!={0,0,0}&];
Print["We find the intersections of ", Length[planes]-1, " planes with the first plane."];

(*find lines of intersection with first plane in list*)
stel=DeleteCases[Table[(*intersect all legal planes with the first*)
newinter2planes[ {planes[[k]], planes[[1]]} ],
{k,1,Length[planes]}],Null];
(*make lines of intersection*)
Print[Length[stel]," lines of intersection were found."];
lns=Map[{ #[[1]]+(100*(#[[2]])),#[[1]]-(100*(#[[2]])) }&, stel];


(*now for the points of intersection between all planes*)
dex=demi[ Length[rawplanes] ];(*indices of all plane triples*)
rawpts=DeleteCases[
Table[
newinter3planes[ plns[[ dex[[k]] ]] ],
{k, Length[dex]}],
Null];
Print[Length[rawpts], " raw points of intersection were found."];
Print[Length[pts=vunion[rawpts]], " distinct points of intersection were found."];

(*and now, knowing the points, we know the PlotRange*)
plndist=Union[ Map[ (Abs[#.plns[[1,{1,2,3}]]])&, pts], SameTest->(#1 - #2 < .001 &)];
Print["The points lie on ", Length[spheres=Union[ Map[ mag[#]&, pts], SameTest->(#1 - #2 < .001 &)]], " concentric spheres."];
maxrad=Max[spheres];

Display Lines of Intersection With First Plane, and All Points

r=maxrad;(*PlotRange*)
Show[Graphics3D[
{Line /@ lns,
AbsolutePointSize[5], Point /@ pts}
],
PlotRange->{{-r,r},{-r,r},{-r,r}},
ViewPoint->2.5*unit[ {0,0,1} ],(*rawplanes[[1,{1,2,3}]]*)
Axes->True,
Boxed->False];

Find All Lines of Intersection

(*now, all lines of intersection between planes*)
dex=hemi[ Length[rawplanes] ];(*indices of all plane pairs*)
allsegs=Chop[
DeleteCases[
Table[
newinter2planes[ plns[[ dex[[k]] ]] ],
{k, Length[dex]}],
Null]
];
(*make lines of intersection*)
lnsall=Map[{ #[[1]]+(2*r*(#[[2]])),#[[1]]-(2*r*(#[[2]])) }&, allsegs];
Print[Length[lnsall], " lines of intersection were found altogether."];

Display All Lines and All Points

r=maxrad+2;(*PlotRange*)
Show[Graphics3D[
{Line /@ lnsall,
AbsolutePointSize[5], Point/@ pts
}
],
PlotRange->{{-r,r},{-r,r},{-r,r}},
ViewPoint->2.5*unit[ {0,0,1} ],(*rawplanes[[1,{1,2,3}]]*)
Axes->True,
Boxed->False];

Selecting Regions of Zonotiles

We operate upon a flattened list of all the zonohedra in the zonotile, polyflat, in the form {center of zonohedron, zonohedron). At first we select subsets of the zonohedra composing the zonotile somewhat crudely; in the subsections farther below, more exact methods are used.

Note that in one of the “crude” examples involving a coordinate plane, the z=h plane is used, by the device Select[polyflat, Abs[#[[1,3]]-h]<distance&], where we obtain all zonohedra within distance of the specified plane. For x=h we'd write #[[ 1, 1 ]] and for y=h, #[[ 1, 2 ]]; for we are selecting with reference to the centers, {x, y, z}, of the zonohedra.

The other subsections contain methods to select zonohedra by exact distance from the origin, or by exact distances from some arbitrary plane, possibly perpendicular to a symmetry vector of the zonotile. For instance, we obtain “quasi-sections” along the axes of two-fold, three-fold, and five-fold symmetry of the icosahedral group, or, if a zonotile created elsewhere has octahedral symmetry, using the axes of two-fold, three-fold, and four-fold symmetry. Cylindrical regions around symmetry vectors are found, and finally a combination of cylindrical and plane regions.

Further possibilities which occur to me are selecting hollow cylinders of zonohedra, or cutting a zonotile by a sphere outside the zonotile, or making a kind of polyhedral lattice of zonohedra—for instance, given a zonotile of icosahedral symmetry, selecting all zonohedra within some radius of the edges of the pentagonal dodecahedron (or their equivalents, within the structure of the zonotile). These types of selections would be most practical and pleasing using zonotiles composed by thousands of zonohedra, the more the merrier. To make such a zonotile requires at least fifty or sixty planes, and thus very much time and Kernel RAM.

Crude Methods

Create Sample Zonotile

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1;(*set Range limits: k=0, 1/2, 1, 3/2, etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

Central Region

(*select all zonohedra within a given radius from the origin*)
radius=6;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
If[shell!={},
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True],
Print["No zonohedra selected!"]
];

Zonohedra Near a Coordinate Plane

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2;
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
If[shell!={},
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True],
Print["No zonohedra selected!"]
];

Zonohedra Near a Coordinate Plane; A Toroidal Region

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2;
shell=Map[#[[2]]&, Select[polyflat, (Abs[#[[1,3]]]<distance && mag[ #[[1]] ] > 5)&] ];
If[shell!={},
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True],
Print["No zonohedra selected!"]
];

Zonohedra to One Side of a Coordinate Plane

(*select all zonohedra with centers below the z=distance plane*)
distance=0;
shell=Map[#[[2]]&, Select[polyflat, #[[1,3]]<distance&] ];
If[shell!={},
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True],
Print["No zonohedra selected!"]
];

Zonohedra Near a Plane Parallel to a Coordinate Plane

(*select all zonohedra with centers less than some distance from the z=h plane*)
distance=2;(*set distance*)
h=4;(*set z coordinate, 'height' at which to make quasi-section*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]-h]<distance&] ];
If[shell!={},
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True],
Print["No zonohedra selected!"]
];

Sorting by Exact Distance from Origin

This requires some attention: having constructed a zonotile, we have the global list polyflat, of the form {center of zonohedron, zonohedron}. Evaluating the first cell below, we find distances from the origin, for all centers of zonohedra. Then, knowing these distances, and how many there are, we select only those zonohedra greater than some one distance, and less than some other distance.

Of course, one could simplify things by selecting all zonohedra at less than some distance, or all at more than some distance.

It is interesting to make animations of a zonotile being gradually built up by the additions of zonohedra, from the center, out.

The zonotile upon which the cells below were based is that triacontahedron made using the vertices of the Icosahedron (i.e., the face-planes of the Pentagonal Dodecahedron), and distances { -1, 0, 1}, giving edge=3, with a triacontahedron of edge=1 at the center.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1;(*set Range spaces: k=0, 1/2,1,1.5,etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*distances of zonohedral centers from origin*)
rawrad=Map[mag[ #[[1]] ]&, polyflat];
tol=.001;(*tolerance for Union*)
rad=Union[rawrad, SameTest->(#1 - #2 < tol &)]
Print["Zonohedra lie at ",Length[rad], " distances from the center of the zonotile."];

(*knowing the distances, choose lower and upper bounds*)
shell=Map[#[[2]]&, Select[polyflat,
(Chop[ mag[#[[1]]]-rad[[4]] ]>=0 && Chop[ mag[#[[1]]]-rad[[6]] ]<=0)&] ];
Print["There are ", Length[shell], " zonohedra selected."];
Show[Graphics3D[shell],
ViewPoint->{0,-.25,2},
Axes->True,
Boxed->True];

(*knowing the distances, choose one distance*)
shell=Map[#[[2]]&, Select[polyflat, (Chop[ mag[#[[1]]]-rad[[30]] ]==0)&] ];
Print["There are ", Length[shell], " zonohedra selected."];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

(*knowing the distances, choose all within one distance*)
shell=Map[#[[2]]&, Select[polyflat, (Chop[ mag[#[[1]]]-rad[[2]] ]<=0)&] ];
Print["There are ", Length[shell], " zonohedra selected."];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Sort Zonohedra by Exact Distance

(*Make a list of zonohedra at each distance*)
shells=Table[
Map[#[[2]]&, Select[polyflat, (Chop[ mag[#[[1]]]-rad[[k]] ]==0)&] ],
{k, Length[rad]}];
Print["Zonohedra lie at ",Length[rad], " distances from the center of the zonotile."];
Print["There are ",Length /@ shells, " zonohedra at the ", Length[rad], " distances from the center of the zonotile."];

(*Show one of the "shells"*)
Show[Graphics3D[ shells[[5]] ],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

(*Show a Range of "shells"*)
(*Range[3,5,1] == {3,4,5}*)
Show[Graphics3D[ shells[[ Range[3,5,1] ]] ],
ViewPoint->2.5*unit[{0,0,2}],
Axes->True,
Boxed->True];

Quasi-sections Along Symmetry Axes

Here we select zonohedra according to the distances of their centers from some arbitrary plane. In the example, a nested triacontahedral zonotile is created, with one triacontahedron at the exact center. This central one is allpoly[[ 4 ]]. Three sets of symmetry vectors are established: the axes of five-fold, three-fold, and two-fold symmetry of the icosahedral symmetry group. A fourth random vector is also provided, along with symmetry vectors for the octahedral group.

We choose one of these seven types, by a simple assignment: vec  =  icosa3 is the default example.

Then we loop through the zonohedral centers, taking their dot products with the chosen axis, vec. These are the signed distances from a plane through the origin, perpendicular to vec. We take their Union and name that list rad. 

We also use vec as the ViewPoint when we Show a subset of the zonohedra. That is, we look at the zonohedra lying in or near a plane perpendicular to vec from a multiple of vec.

After testing the efficacy of the method by selecting all zonohedra near one of the planes-containing-zonohedra, we make a list of zonohedra with centers lying in those planes. This allows animation of the building of the zonotile, layer by layer, parallel to the plane perpendicular to vec.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1;(*set Range limits: k=0, 1/2, 1, 3/2,etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the entire zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];

(*Show a subset*)
Show[Graphics3D[ allpoly[[4]] ],
Axes->True,
Boxed->True];

Icosahedral & Octahedral Symmetry Vectors

(*axes of 5-fold symmetry*)
icosa5=First@Map[Chop[unit[#]]&, Vertices[Icosahedron] ];
(*axes of 3-fold symmetry*)
icosa3=First@Map[Chop[unit[#]]&, Vertices[Dodecahedron] ];
(*axes of 2-fold symmetry*)
icosa2=First@Map[unit[cent[#]]&, edgefind[N[Vertices[Icosahedron]]]];

(*arbitrary vector*)
arb=unit[ {Random[Real,{0,1}], Random[Real,{0,1}], Random[Real,{0,1}]} ];

(*for zonotiles of octahedral symmetry*)
oct2=unit[{1,0,1}];
oct3=unit[{1,1,1}];
oct4=unit[{0,0,1}];

Choose Vector, Find Face Centers and Distances

vec=icosa3;(*choose which symmetry axis*)

dist=Map[ Chop[ #[[1]].vec ]&, polyflat ];
tol=.001;(*tolerance for Union*)
rad=Union[dist, SameTest->(#1 - #2 < tol &)]
Print["Zonohedra lie at ",Length[rad], " distances"];
Print["The greatest distance is ",maxrad=Max[rad]];

Select Zonohedra Near One Distance

(*select all zonohedra with centers less than some distance from the chosen plane*)
breadth=1.5;(*set breadth of selected region*)
h=0;(*set distance at which to make quasi-section; could be h=rad[[k]]*)
shell=Map[#[[2]]&, Select[polyflat, (Abs[(#[[1]].vec)-h]<breadth)&] ];
If[shell!={},
Show[Graphics3D[shell],
ViewPoint->2.5*unit[ vec ],
Axes->True,
Boxed->True],
Print["No zonohedra selected!"]
];

Find Sets of Zonohedra at Each Distance

(*create all layers*)
layers=Table[
Map[#[[2]]&, Select[polyflat, Chop[ Abs[(#[[1]].vec)-rad[[k]] ] ]==0&] ],
{k, Length[rad]}];
Print["There are ", Length[layers], " sets of zonohedra"];

Show[Graphics3D[ layers[[15]] ],
ViewPoint->2.5*unit[ vec ],
Axes->True,
Boxed->True];

Build Zonotile Layer by Layer: Animation

r=maxrad+2;(*set PlotRange*)
Table[
stage=Polygon /@ dePair[Chop[parse[ Flatten[layers[[ Range[i] ]],1] ] ] ];
Show[Graphics3D[ stage ],
PlotRange->{{-r,r},{-r,r},{-r,r}},
Boxed->False,
Background->GrayLevel[1],
ViewPoint->2.25*unit[ vec ]],(*from perpendicular to planes-of-selection*)
{i, Length[layers]}];

Cylindrical Regions Along Symmetry Axes

Create Zonotile

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=2;(*set Range limits: k=0, 1/2, 1, 3/2,etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

Icosahedral & Octahedral Symmetry Vectors

(*axes of 5-fold symmetry*)
icosa5=First@Map[Chop[unit[#]]&, Vertices[Icosahedron] ];
(*axes of 3-fold symmetry*)
icosa3=First@Map[Chop[unit[#]]&, Vertices[Dodecahedron] ];
(*axes of 2-fold symmetry*)
icosa2=First@Map[unit[cent[#]]&, edgefind[N[Vertices[Icosahedron]]]];

(*arbitrary vector*)
arb=unit[ {Random[Real,{0,1}], Random[Real,{0,1}], Random[Real,{0,1}]} ];

(*for zonotiles of octahedral symmetry*)
oct2=unit[{1,0,1}];
oct3=unit[{1,1,1}];
oct4=unit[{0,0,1}];

Choose Vector, Select Cylindrical Region, Display

vec=icosa3;(*choose which symmetry axis*)
side=First@NullSpace[ {vec} ];(*for ViewPoint*)
radius=4;(*set radius of cylindrical region*)

sceptre=Map[#[[2]]&, Select[polyflat,(mag[(vec*(#[[1]].vec))-#[[1]]]<radius)&] ];
Print[Length[sceptre], " zonohedra were selected within the cylinder of radius ", radius];

If[sceptre!={},
Show[Graphics3D[sceptre],
ViewPoint->2.5*unit[ side ],
Axes->True,
Boxed->True],
Print["No zonohedra selected!"]
];

Cylindrical Region & Perpendicular Plane Region

Here we combine a disk-like planar region with a cylindrical region, the one perpendicular to an axis vec, the other surrounding that axis. To evaluate the example, first create a zonotile and declare symmetry vectors as in  Cylindrical Region subsection above.

ViewPoint is set to the intermediate direction vec - side, where vec is the chosen axis, and side is a vector perpendicular to vec. These types of zonohedral arrays are somewhat difficult to see properly from any direction; real-time rotation helps, which can be had by evaluating <<RealTime3D`. To return to regular display evaluate <<Default3D`.

vec=icosa5;(*choose axis*)
side=First@NullSpace[ {vec} ];
(*select all zonohedra with centers less than some distance from the chosen plane*)
breadth=1.5;(*set breadth of selected region*)
h=0;(*set distance at which to make quasi-section*)
radius=3;(*and all zonohedra within a cylinder of this radius around axis vec*)

shell=Map[#[[2]]&, Select[polyflat,
(Abs[(#[[1]].vec)-h]<breadth)
||
(mag[(vec*(#[[1]].vec))-#[[1]]]<3)&]
];
If[shell!={},
Show[Graphics3D[shell],
ViewPoint->2.5*unit[ vec-side ],
PlotRange->All,
Axes->True,
Boxed->True],
Print["No zonohedra selected!"]
];

Building Zonotiles: Animation

I have been making animations of zonotile construction for many years. The examples below are simple. It is possible to have each set of zonohedra come spinning in from space and settle into their own little nests, gradually changing color as they do so; and I have done just this! Or, every zonohedron in the zonotile may be exported to POV-Ray, where degrees of transparency may be assigned, and the zonohedra may spin into view out of nothingness, and then assume their proper colors, and gravitate to their positions.

Animation I: Build Zonotile From Center to Periphery

We animate the building of a zonotile from the central zonohedra, outward, and use the dePair function to remove coincident faces. This drastically reduces render times, without affecting the outward appearance. By way of example we use the Triacontahedron of edge=3, derived from the vertices of the Platonic Icosahedron (hence, face-planes of the Platonic Dodecahedron). This zonotile is composed by 233 zonohedra.

1. Create Zonotile

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1;(*set Range spaces: k=0, 1/2,1,1.5,etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

2. Find distances of zonohedra

(*distances of zonohedral centers from origin*)
rawrad=Map[mag[ #[[1]] ]&, polyflat];
tol=.001;(*tolerance for Union*)
rad=Union[rawrad, SameTest->(#1 - #2 < tol &)]
Print["Zonohedra lie at ",Length[rad], " distances from the center of the zonotile."];
maxrad=Max[rad];(*center of farthest zonohedron*)

3. Animate

r=maxrad+3;(*set plotrange, adjust as necessary*)

Do[

rawshell=Map[#[[2]]&, Select[polyflat, (Chop[ mag[#[[1]]]-rad[[k]] ]<=0)&] ];
shell=Polygon /@ dePair[Chop[parse[ Flatten[rawshell,1] ] ] ];
(*Print["Number of zonohedra is ", Length[rawshell]];*)

	Show[
		Graphics3D[shell],
		PlotRange->{{-r,r},{-r,r},{-r,r}},
		Boxed->False,
		ViewPoint->2.25*unit[{0,0,2}]
	],
	
{k, 1, Length[rad], 1}];

Animation II: Build Zonotile From Lowest to Highest

We animate the building of a zonotile from the bottom, upward. First, the Union of all z coordinates of zonohedral centers is found. Then, a Do loop effects the animation. To reduce polygon count, the dePair function strips away any polygons which have fallen into coincidence, so that all we see is the exterior of the zonohedral array—which is all we would see, anyway.

1. Create Zonotile

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1;(*set Range spaces: k=0, 1/2,1,1.5,etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

2. Find z coordinates of zonohedral centers

(*z coords of zonohedral centers*)
zrawrad=Map[#[[1,3]]&, polyflat];
tol=.001;(*tolerance for Union*)
zrad=Union[zrawrad, SameTest->(#1 - #2 < tol &)];
Print["Zonohedral centers have ",Length[zrad], " unique z coordinates."];
maxrad=Max[rad];(*center of highest/lowest zonohedron*)

3. Animate

(*knowing how many z coords, and what they are, choose all below a given one*)
r=maxrad+3;(*plotrange, adjust as necessary*)

Do[
rawshell=Map[#[[2]]&, Select[polyflat, (Chop[ #[[1,3]]-zrad[[k]] ]<=0)&] ];
shell=Polygon /@ dePair[Chop[parse[ Flatten[rawshell,1] ] ] ];

Show[Graphics3D[shell],
PlotRange->{{-r,r},{-r,r},{-r,r}},
ViewPoint->2.5*unit[{0,-1,.5}],
Boxed->False],
{k, 1, Length[zrad], 1}];

Animation III: Build Zonotile by Expanding Cylindrical Regions

1. Create Zonotile

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1;(*set Range spaces: k=0, 1/2,1,1.5,etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

Icosahedral & Octahedral Symmetry Vectors

(*axes of 5-fold symmetry*)
icosa5=First@Map[Chop[unit[#]]&, Vertices[Icosahedron] ];
(*axes of 3-fold symmetry*)
icosa3=First@Map[Chop[unit[#]]&, Vertices[Dodecahedron] ];
(*axes of 2-fold symmetry*)
icosa2=First@Map[unit[cent[#]]&, edgefind[N[Vertices[Icosahedron]]]];

(*arbitrary vector*)
arb=unit[ {Random[Real,{0,1}], Random[Real,{0,1}], Random[Real,{0,1}]} ];

(*for zonotiles of octahedral symmetry*)
oct2=unit[{1,0,1}];
oct3=unit[{1,1,1}];
oct4=unit[{0,0,1}];

3. Choose Axis, Find Radii

vec=icosa5;(*choose which symmetry axis*)
side=First@NullSpace[ {vec} ];(*for ViewPoint*)

dist=Map[ ( mag[(vec*(#[[1]].vec))-#[[1]]] )&, polyflat ];
tol=.001;(*tolerance for Union*)
Print["The radii are ", rad=Union[dist, SameTest->(#1 - #2 < tol &)]];
Print["Zonohedra lie at ",Length[rad], " distances from the chosen axis"];
Print["The greatest distance is ",maxrad=Max[rad]];

radius=rad[[3]];(*set radius of cylindrical region*)

sceptre=Map[#[[2]]&, Select[polyflat,(mag[(vec*(#[[1]].vec))-#[[1]]]<radius)&] ];
Print[Length[sceptre], " zonohedra were selected within the cylinder of radius ", radius];

If[sceptre!={},
Show[Graphics3D[sceptre],
ViewPoint->2.5*unit[ side ],
Axes->True,
Boxed->True],
Print["No zonohedra selected!"]
];

4. Animate

r=maxrad+2;(*set PlotRange*)
Do[
sceptre=Map[#[[2]]&, Select[polyflat,
(Chop[ mag[(vec*(#[[1]].vec))-#[[1]]] -rad[[k]] ]<=0)&] ];

shell=Polygon /@ dePair[Chop[parse[ Flatten[sceptre,1] ] ] ];

Show[Graphics3D[shell],
PlotRange->{{-r,r},{-r,r},{-r,r}},
ViewPoint->2.5*unit[ side ],
Boxed->False],
{k,Length[rad]}];

Exploding Zonotiles

The explode function accepts a list in the form {{center, zonohedron}, ...}, and translates each zonohedron along the vector to its center by a multiple of that vector; as in explode[ zonolist, multiple].

Explode the entire zonotile

Show[Graphics3D[ explode[ polyflat, .2] ]];

Explode a selected region

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2;
zonolist=Select[polyflat, Abs[#[[1,3]]]<distance&];
Show[Graphics3D[ Map[#[[2]]&, zonolist] ],
ViewPoint->{0,0,2.5}];(*unexploded*)

Show[Graphics3D[explode[zonolist, .2]],
ViewPoint->{0,0,2.5}];(*exploded*)

Animation

(*how big is the slightly exploded zonotile?*)
(*use this info to set the plotrange 'r' in the next cell, to animate*)
Show[Graphics3D[ explode[ polyflat, .2] ], Axes->True];

(*animate explosion of the whole zonotile*)
r=26;(*maximum plotrange; if i goes from 0 to 1, set to twice initial i=0 radius*)
Do[
Show[Graphics3D[ explode[ polyflat, i] ],
PlotRange->{{-r,r},{-r,r},{-r,r}}],
{i, 0, 1, .1}];

Color Zonohedra by Volume

The volumeSort function accepts a list of zonohedra, and sorts it into sublists according to volume. Usually, but not always, this sorts the zonohedra into their different types; that is, it rarely does happen that two different zonohedra have the same volume. The resulting list of lists of zonohedra can then be assigned colors and displayed.

A zonotile example is provided, so that all examples below will work properly. Choices must be made at various stages. However, any zonotile created elsewhere in the notebook will work here, by skipping to Stage 2, Choose What to Sort. In fact, if some subset of a zonotile has already been created elsewhere, perhaps, who knows, named shell, one can substitute shell for polyflat in Stage 2. The volumeSort function operates upon a list of polyhedra, whether it be in the form found in polyflat, {{center of zonohedron, zonohedron}, ...}, or simply {{zonohedron},{ zonohedron}, ...}.

1. Construct Zonotile

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Dodecahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1/2;(*set Range limits: k=0, 1/2, 1, 3/2,etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

2. Choose What to Sort ...

(*the whole zonotile, unexploded*)
s=volumeSort[ polyflat ];

(*the whole zonotile, exploded*)
s=volumeSort[ explode[ polyflat, .2] ];

(*a region of the zonotile, exploded*)
distance=2;
s=volumeSort[ explode[Select[polyflat, Abs[#[[1,3]]]<distance&], .2] ];

3. Choose Color List

(*random colors, if more than 10 colors are needed*)
colorlist=Table[RGBColor[Random[Real,{0,1}],
		Random[Real,{0,1}],
		Random[Real,{0,1}]],{Length[s]}];

(*if the number of colors needed is 10 or less*)
colorlist={
RGBColor[1,.25,.25],
RGBColor[.25,.25,1],
RGBColor[1,1,0],
RGBColor[.5,.5,.5],
RGBColor[1,0,1],
RGBColor[.5,0,0],
RGBColor[0,0,.5],
RGBColor[.5,.5,0],
RGBColor[.5,0,.5],
RGBColor[0,.5,0]
};

4a. Display All Colored Zonohedra

(*display colored zonotile*)
qaq=Graphics3D[
{EdgeForm[GrayLevel[0]],
	Table[ {SurfaceColor[colorlist[[i]],GrayLevel[.5], 5],
		Polygon /@ s[[ i ]]}, {i,Length[s]}]
}];
Show[qaq,
Boxed->False,
Background->GrayLevel[0],
LightSources -> 
   {
   {{0,0,2}, RGBColor[.5, .5, .5]},
   {{-2,-.5,1}, RGBColor[1, 1, 1]}
   },
   PlotRange->All,
ViewPoint->2.25*unit[{0,-.5,1}]];

4b. Display One Type of Colored Zonohedra

(*display one volume-group of colored zonotile*)
(*which one?*)
k=6;
qaq=Graphics3D[
{EdgeForm[GrayLevel[0]],
	{SurfaceColor[colorlist[[k]],GrayLevel[.5], 5],
		Polygon /@ s[[ k ]]}}];
Show[qaq,
Boxed->False,
Background->GrayLevel[0],
LightSources -> 
   {
   {{0,0,2}, RGBColor[.5, .5, .5]},
   {{-2,-.5,1}, RGBColor[1, 1, 1]}
   },
   PlotRange->All,
ViewPoint->2.25*unit[{0,-1,1}]];

4c. Display Two or More Types of Colored Zonohedra

Length[s](*check how many sets?*)

(*display some subset of the colored zonotile*)
(*declare indices*)
dex={1,7};
qaq=Graphics3D[
	Table[
		{EdgeForm[GrayLevel[0]],
		{SurfaceColor[colorlist[[ dex[[k]] ]],GrayLevel[.5], 5],
		Polygon /@ s[[ dex[[k]] ]]}},
	{k, Length[dex]}]
	];
Show[qaq,
Boxed->False,
Background->GrayLevel[0],
LightSources -> 
   {
   {{0,0,2}, RGBColor[.5, .5, .5]},
   {{-2,-.5,1}, RGBColor[1, 1, 1]}
   },
   PlotRange->All,
ViewPoint->2.25*unit[{0,-1,1}]];

Color Zonohedra by Distance From Origin

We sort the zonohedra by distance from the center of the zonotile (the origin), and then assign a random color to each set.

When the zonotile is determined from the vertices of a Platonic—hence, the face-planes of its dual—and the plane distances are set to {-1/2, 1/2}, as in the example below, we obtain “something like” the reciprocal of the set of all mainline stellations of that polyhedron. In the example below, we use the vertices of the Dodecahedron, hence the face-planes of the Icosahedron, which has eight mainline stellations. We find that the zonohedra are divided into eight sets, according to their distances from the origin. In Step 6 below, Build Entire Zonotile, we add each shell or layer of zonohedra in turn, that is, first shells[[1]], then shells[[{1,2}]], then shells[[{1,2,3}]], etc. Perhaps these eight stages of construction correspond most closely to the eight mainline stellations. Also note that in the stellation diagram for the Icosahedron, points of intersection lie at eight distances from the origin.

1. Create Zonotile

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Dodecahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

k=1/2;(*set Range limits: k=0, 1/2, 1, 3/2,etc.*)

dist=Table[Range[-k,k,1],{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

2. Find Distances

(*distances of zonohedral centers from origin*)
rawrad=Map[mag[ #[[1]] ]&, polyflat];
tol=.001;(*tolerance for Union*)
rad=Union[rawrad, SameTest->(#1 - #2 < tol &)]
Print["Zonohedra lie at ",Length[rad], " distances from the center of the zonotile."];
Print["The greatest distance is ",maxrad=Max[rad], "."];

3. Sort Zonohedra by Distance

(*Make a list of zonohedra at each distance*)
shells=Table[
Map[#[[2]]&, Select[polyflat, (Chop[ mag[#[[1]]]-rad[[k]] ]==0)&] ],
{k, Length[rad]}];
Print["Zonohedra lie at ",Length[rad], " distances from the center of the zonotile."];
Print["There are ",Length /@ shells, " zonohedra at the ", Length[rad], " distances."];

4. Color Zonohedral Sets

(*random colors*)
colorlist=Table[RGBColor[Random[Real,{0,1}],
		Random[Real,{0,1}],
		Random[Real,{0,1}]],{Length[shells]}];

coloredshells=Table[
Graphics3D[
{EdgeForm[GrayLevel[0]],
{SurfaceColor[colorlist[[i]],GrayLevel[.5], 5],
		shells[[ i ]]}}],
		
{i,Length[shells]}];

5. Display One Set

Show[coloredshells[[ 6 ]],
Boxed->False,
Background->GrayLevel[0],
LightSources -> 
   {
   {{0,0,2}, RGBColor[.5, .5, .5]},
   {{-2,-.5,1}, RGBColor[1, 1, 1]}
   },
ViewPoint->2.25*unit[{0,-.5,1}]];

6. Build Entire Zonotile

r=maxrad+2;(*set PlotRange*)
Table[
Show[coloredshells[[ Range[i] ]],
PlotRange->{{-r,r},{-r,r},{-r,r}},
Boxed->False,
Background->GrayLevel[0],
LightSources -> 
   {
   {{0,0,2}, RGBColor[.5, .5, .5]},
   {{-2,-.5,1}, RGBColor[1, 1, 1]}
   },
ViewPoint->2.25*unit[{0,-.5,1}]],
{i, Length[coloredshells]}];

7. Display Each Stage

r=maxrad+2;(*set PlotRange*)
Table[
Show[coloredshells[[ i ]],
PlotRange->{{-r,r},{-r,r},{-r,r}},
Boxed->False,
Background->GrayLevel[0],
LightSources -> 
   {
   {{0,0,2}, RGBColor[.5, .5, .5]},
   {{-2,-.5,1}, RGBColor[1, 1, 1]}
   },
ViewPoint->2.25*unit[{0,-.5,1}]],
{i, Length[coloredshells]}];

Stellation Diagrams

Since some of the arrangements of planes treated above are just exactly those which give stellations of convex polyhedra, I provide some examples of stellation diagrams, first in 3D and then projected into a plane. Then intersection points are found in the plane, and the radii of those points, and circles are drawn containing those points.

The method used here will work on any polyhedron such that all faces are perpendicular to the lines joining the center of the polyhedron to the centers of the faces. For instance, the Platonic solids. With modifications, it could produce stellation diagrams for at least some of the Archimedean solids; for in those cases, there are two or three types of polygon, and thus two or three different diagrams.

Unmodified, it will work on polyhedra such as Kepler’s Rhombic Dodecahedron and Triacontahedron. To do so, one must change the code slightly at the beginning of Stage 1 of the Single Stellation subsection, Find Lines and Points of Intersection. The Pentagonal Dodecahedron is the default example.

In the last stage, Mark Intersections and Find Radii, we have reduced the stellation diagram to a plane (as one usually sees it; the 3D diagram created here is unusual). The points of intersection in the plane are found and marked by red dots. The distances of these points from the center of the diagram are found, and marked by blue circles. In the case of the Pentagonal Dodecahedron, there are three distances (radii), and three circles comprehend all points of intersection. Note that these three correspond exactly to the three concentric shells of zonohedra in the corresponding zonotile. To make this zonotile, go to Zonotile Examples, subsection “Platonic Zonotiles I” and set k to 1/2.

Kepler’s Rhombic Triacontahedron

The KRT is the dual of the Archimedean Icosidodecahedron, and is bounded by thirty equal rhombs. Its stellations number in the hundreds of millions. We create a KRT here and scale it properly for use below, in the Single Stellation subsection, where the first lines must then be commented out.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[{0},{i,n}];(*all planes intersect at the origin*)
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the zonohedron*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];

rawfaces=parse[ allpoly ];
sca=1/Chop[mag[cent[ rawfaces[[1]] ] ] ];(*scale factor*)
polyfaces1=sca*rawfaces;(*scale the KRT to unit inradius*)

Single Stellation of the Pentagonal Dodecahedron

1. Find Lines and Points of Intersection

(*comment out the three lines below if using the KRT*)
rawfaces=parse[ Polyhedron[Dodecahedron] ];
sca=1/Chop[mag[cent[ rawfaces[[1]] ] ] ];(*scale factor*)
polyfaces1=sca*rawfaces;

rawplanes=Map[Chop[cent[#]]&, polyfaces1];(*all face-centers*)
(*get rid of any parallel to the first*)
planes=Select[rawplanes, Chop[rawplanes[[1]]+#]!={0,0,0}&];
Print["We find the intersections of ", Length[planes]-1, " planes with the first plane."];

(*find lines of intersection with first plane in list*)
stel=Table[(*intersect all legal planes with the first*)
inter[{unit[planes[[k]]], mag[ planes[[k]] ]}, {unit[planes[[1]]], mag[ planes[[1]] ]} ],
{k,2,Length[planes]}];

(*make lines of intersection*)
lns=Map[{ #[[1]]+(30*(#[[2]])),#[[1]]-(30*(#[[2]])) }&, stel];

(*now for the points of intersection between all planes*)
plns=Map[unit[#]&, rawplanes];(*convert normals to unit vectors*)
dex=demi[ Length[rawplanes] ];(*indices of all plane triples*)
rawpts=DeleteCases[
Table[(*use unit vectors plus actual distances*)
inter3planes[ plns[[ dex[[k]] ]], mag /@ rawplanes[[ dex[[k]] ]] ],
{k, Length[dex]}],
Null];
Print[Length[pts=vunion[rawpts]], " distinct points of intersection were found."];

(*and now, knowing the points, we find plane-distances and PlotRange*)
(*we must divide by actual inradius*)
inrad=mag[ rawplanes[[1]] ];(*we assume all planes equal*)
plndist=Union[ Map[ (Abs[#.plns[[1]]]/inrad)&, pts], SameTest->(#1 - #2 < .001 &)];
Print["The points lie on ", Length[spheres=Union[ Map[ mag[#]&, pts], SameTest->(#1 - #2 < .001 &)]], " concentric spheres."];
maxrad=Max[spheres];(*minimum PlotRange*)

2. 3D Stellation Diagram

r=maxrad+1;(*PlotRange*)
Show[Graphics3D[
{Polygon /@ polyfaces1,
Line /@ lns,
AbsolutePointSize[5], Point/@ pts}
],
PlotRange->{{-r,r},{-r,r},{-r,r}},
ViewPoint->2.5*unit[ rawplanes[[1]] ],
Axes->True,
Boxed->True];

3. Plane Projection of Diagram

gbas=NullSpace[ {N[ unit[planes[[1]]] ]} ];(*basis vectors*)

lns1=Map[ Dot[ gbas, # ]&, lns, {2}];(*project lines*)

diagram=Graphics[{AbsoluteThickness[1],Line /@ lns1}];(*project lines and points*)

r=maxrad;(*PlotRange*)
Show[diagram,
AspectRatio->Automatic,
PlotRange->{{-r,r},{-r,r}},
Frame->True];

4. Mark Intersections and Find Radii

n=Length[lns1];
dex=Flatten[Table[{i,j},{i,n-1},{j,i+1,n}],1];

points=DeleteCases[
Map[ Chop[ lineinter[ lns1[[ #[[1]] ]], lns1[[ #[[2]] ]] ] ]&, dex],
{}];
rad=Union[ Map[mag[#]&, points], SameTest->(#1 - #2 < tol &)];

redpts=Graphics[{AbsolutePointSize[8],RGBColor[1,0,0], Point /@ points}];
bluecircles=Graphics[{RGBColor[0,0,1], Map[Circle[{0,0}, #]&, rad]}];

r=maxrad+1;(*PlotRange*)
Show[{diagram,redpts,bluecircles},
AspectRatio->Automatic,
PlotRange->{{-r,r},{-r,r}},
Frame->True];

Multiple Stellations

Here is an arcane subject—doubly arcane, for we explore “double stellations” of polyhedra. Of course, we do not create the stellations themselves, only the arrangements of planes, and some stellation diagrams. We nestle one polyhedron inside another larger one—the default example uses two concentric Platonic dodecahedra. Then we create the zonotiles(s) determined by such planar arrangements.

This routine suffers from confusion about the inradius, and how to find the distances at which the face-planes of the second dodecahedron will intersect points of intersection from the first. Eventually various assumptions are made (faces all equal, all equal inradii) which make the routine of little general use. In the example, the Platonic dodecahedron is doubly stellated—two concentric dodecahedra are created such that the outer just kisses some k-ic points of intersection of the first. Now k+1 planes meet there, and a (k+1)-ic zonohedron is determined at its proper relative position.

Kepler’s Rhombic Triacontahedron

The KRT is the dual of the Archimedean Icosidodecahedron, and is bounded by thirty equal rhombs. Its stellations number in the hundreds of millions. We create a KRT here and scale it properly for use below, in the Double Stellation subsection, where the first lines must then be commented out.

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[{0},{i,n}];(*all planes intersect at the origin*)
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the zonohedron*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];

rawfaces=parse[ allpoly ];
sca=1/Chop[mag[cent[ rawfaces[[1]] ] ] ];(*scale factor*)
polyfaces1=sca*rawfaces;(*scale the KRT to unit inradius*)

Double Stellation of the Platonic Pentagonal Dodecahedron

1. Find Lines and Points of Intersection

(*comment out the three lines below if using the KRT*)
rawfaces=parse[ Polyhedron[Dodecahedron] ];
sca=1/Chop[mag[cent[ rawfaces[[1]] ] ] ];(*scale factor*)
polyfaces1=sca*rawfaces;

rawplanes=Map[Chop[cent[#]]&, polyfaces1];(*all face-centers*)
(*get rid of any parallel to the first*)
planes=Select[rawplanes, Chop[rawplanes[[1]]+#]!={0,0,0}&];
Print["We find the intersections of ", Length[planes]-1, " planes with the first plane."];

(*find lines of intersection with first plane in list*)
stel=Table[(*intersect all legal planes with the first*)
inter[{unit[planes[[k]]], mag[ planes[[k]] ]}, {unit[planes[[1]]], mag[ planes[[1]] ]} ],
{k,2,Length[planes]}];

(*make lines of intersection*)
lns=Map[{ #[[1]]+(30*(#[[2]])),#[[1]]-(30*(#[[2]])) }&, stel];

(*now for the points of intersection between all planes*)
plns=Map[unit[#]&, rawplanes];(*convert normals to unit vectors*)
dex=demi[ Length[rawplanes] ];(*indices of all plane triples*)
rawpts=DeleteCases[
Table[(*use unit vectors plus actual distances*)
inter3planes[ plns[[ dex[[k]] ]], mag /@ rawplanes[[ dex[[k]] ]] ],
{k, Length[dex]}],
Null];
Print[Length[pts=vunion[rawpts]], " distinct points of intersection were found."];

(*and now, knowing the points, we find plane-distances and PlotRange*)
(*we must divide by actual inradius*)
inrad=mag[ rawplanes[[1]] ];(*we assume all planes equal*)
plndist=Union[ Map[ (Abs[#.plns[[1]]]/inrad)&, pts], SameTest->(#1 - #2 < .001 &)];
Print["The points lie on ", Length[spheres=Union[ Map[ mag[#]&, pts], SameTest->(#1 - #2 < .001 &)]], " concentric spheres."];
maxrad=Max[spheres];(*minimum PlotRange*)

2. 3D Stellation Diagram

r=maxrad+1;(*PlotRange*)
Show[Graphics3D[
{Polygon /@ polyfaces1,
Line /@ lns,
AbsolutePointSize[5], Point/@ pts}
],
PlotRange->{{-r,r},{-r,r},{-r,r}},
ViewPoint->2.5*unit[ rawplanes[[1]] ],
Axes->True,
Boxed->True];

Create and Display All Lines of Intersection

allsegs=Table[
rplns=RotateRight[rawplanes,j-1];(*cycle the list for outer loop*)
planes=Select[rplns, Chop[rplns[[1]]+#]!={0,0,0}&];


stel=Table[(*intersect all legal planes with the first*)
inter[{unit[planes[[k]]], mag[ planes[[k]] ]}, {unit[planes[[1]]], mag[ planes[[1]] ]} ],
{k,2,Length[planes]}];(*inner loop*)

(*make lines of intersection*)
Map[{ #[[1]]+(100*(#[[2]])),#[[1]]-(100*(#[[2]])) }&, stel],
{j, Length[rawplanes]}];(*outer loop*)
lnsall=Flatten[allsegs, 1];

r=maxrad+2;(*PlotRange*)
Show[Graphics3D[
{Polygon /@ polyfaces1,
Line /@ lnsall,
AbsolutePointSize[5], Point/@ pts}
],
PlotRange->{{-r,r},{-r,r},{-r,r}},
ViewPoint->2.5*unit[ rawplanes[[1]] ],
Axes->True,
Boxed->True];

3. Create Outer Polyhedron

Print["The ", Length[pts], " intersection points are at the ", Length[plndist], " plane distances ", Chop[plndist]];

The ratios between one plane-distance and another are the most significant factor in the determined zonotiles. In our default example, using the Pentagonal Dodecahedron, the distance-pair {1,2} is equivalent to {2,4} because the two ratios are equal. In fact, when doubly stellating the dodecahedron, there are only three possible outcomes, using distance-pairs {2, 3}, {2, 4}, and, say, {2, any random number greater than plane-distance[[ 2 ]]}. For instance: {{1., 2.23607}, {1., 4.23607}, {1, 574.531}}.

Similarly, in the other provided example, double stellation of Kepler’s Rhombic Triacontahedron, we may ignore all plane-distances less than 1, except 0. The unique distance 0, when applied to a zonotile, means that one set of fifteen planes intersects at the origin, and thirty other planes intersect outside the origin. Instead of the four plane-distances of the “usual” double stellation, {-d2, -d1, d1, d2}, we have {-d1, 0, d1}, and d1 can be set to any non-zero number whatsoever without affecting the zonotile. A zonotile of edge=3 is created.

Setting aside this special type of double stellation, we ignore all plane-distances less than 1, and consider only all distinct ratios between pairs; for each distinct ratio gives a distinct zonotile. To do this, suppose there are just k distances, k>=1: then the number of k things taken two at a time equals how many pairs can be assembled. But some pairs will give equal ratios to others. By weeding these out, we are left with some number of pairs of distances embodying unique ratios. These will give unique zonotiles. Among these we will find the unique ratios given by {1, any other distance>1}.

polyfaces2=polyfaces1*plndist[[3]];(*this works on the default example, change as desired*)

r=maxrad+2;(*PlotRange*)
Show[Graphics3D[
{Polygon /@ polyfaces2,
Line /@ lns,
AbsolutePointSize[5], Point/@ pts
}
],
PlotRange->{{-r,r},{-r,r},{-r,r}},
ViewPoint->2.5*unit[ rawplanes[[1]] ],
Boxed->False];

4. Create Zonotile

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=rawplanes;

v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];

dist=Table[{-plndist[[3]], -1, 1, plndist[[3]]},{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];

(*Show the whole zonotile*)
Show[Graphics3D[ allpoly ],
Axes->True,
Boxed->True];

(*Show just some of the subsets*)
Show[Graphics3D[ allpoly[[{4}]] ],
Boxed->False];

(*explode the whole zonotile*)
Show[Graphics3D[ explode[ polyflat, .2] ]];

(*select all zonohedra with centers less than some distance from the z=0 plane*)
distance=2.5;(*set distance*)
shell=Map[#[[2]]&, Select[polyflat, Abs[#[[1,3]]]<distance&] ];
Show[Graphics3D[shell],
ViewPoint->2.5*unit[{0,0,1}],
Boxed->False];

(*select all zonohedra within a given radius from the origin*)
radius=6;
shell=Map[#[[2]]&, Select[ polyflat, mag[ #[[1]] ]<radius&] ];
Show[Graphics3D[ shell ],
ViewPoint->{0,-1.5,1},
Boxed->False];

Zones of Zonohedra

Zonotiles are tiled by zonohedra, zonohedra are bounded by zonogons, and for any one such zonogon there is an equal and opposite zonogon on the other side of the zonohedron. When zonohedra close-pack, they share common zonogons. The idea here is to select a “zone” of zonohedra, such that the same zonogon arises again and again, connecting the zonohedra. These equal zonogons are called the “interfaces” of the solid zone.

1. Construct Zonotile, Preview Zonohedra

The method is somewhat as follows. Having created a zonotile, we have the list polyflat, {center of zonohedron, zonohedron}, where the second member “zonohedron” is itself a list of zonogons. We also have the crudely-sorted list allpoly, in which the sublists contain zonohedra of an equal number of faces. Those zonohedra within any one list may yet be, and most likely are, of more than one type. The first sublist is always the hexahedra. First, we display the first member of each sublist in turn (if the zonotile is composed entirely by hexahedra, there is but the one sublist).

(*****MAKE ARRANGEMENT OF PLANES*****)
v1=N[Vertices[Icosahedron]];(*choose Platonic*)
v=Union[v1,SameTest->(Chop[#1+#2]=={0,0,0}&)];
n=Length[v];
dist=Table[{-1, 0, 1},{i,n}];
d=Length[dist[[1]]];
plns=Flatten[Table[
Chop[Join[ unit[ v[[i]] ],{dist[[i,j]]}]],
{i,Length[v]},{j, d}],1];
Print["There are ",r=Length[plns]," planes altogether."];

(*****CREATE ZONOTILE*****)
zonotile3D[plns];


(*Show the first of each of the allpoly*)
If[Length[allpoly]>1,

Do[Show[Graphics3D[ allpoly[[k,1]] ],
Boxed->False],
{k, Length[allpoly]}],

(*else*)
Show[Graphics3D[ allpoly ],
Boxed->False]
]

2. Choose k-gonal Interface

Set k to 4, 6, 8, etc., to find rhombs, hexagons, octagons, etc.

k=4;(*what type of face do we want?*)
(*find k-gonal faces are*)
poly=Map[parse[#]&, polyflat];
types=Map[ Length[#]&, poly, {2}];
kgonal=Union[First /@ Position[types, k] ];

3. Find Zonohedron Containing a k-gonal face near the center of the zonotile.

It may be that one of our k-gons bounds a zonohedron at the center of the zonotile. Or perhaps not. We search within a radius of 3. If none are found, increase the radius. Or, look for a zonohedron near the periphery of the zonotile: set radius higher, and change ‘< radius’ to ‘>radius’. No matter which zonohedron is chosen, a solid zone connects it with others.

Using the first of the list, we obtain the index into polyflat of a zonohedron within our solid zone, and then display. However, we may miss the desired zonohedron in this way.

radius=5;
hedra=Select[ polyflat[[ kgonal ]], mag[#[[1]]]<radius&];
If[Length[hedra]<1, Print["No zonohedra found, change radius"],

Print[Length[hedra], " zonohedra found"];
cendex=kgonal[[ First@Flatten[Position[ polyflat[[ kgonal ]], hedra[[1,1]] ]] ]];

Show[Graphics3D[ polyflat[[cendex,2]] ],
ViewPoint->{0,1,2},
Boxed->False];]

We may not have found the desired zonohedron upon which to root our solid zone. We can examine each one in turn, and note which member of the list makes our target. Caution: if many zonohedra compose the list hedra, we may overstrain available FrontEnd RAM.

Do[
Show[Graphics3D[ hedra[[k,2]] ],
PlotLabel->ToString[k]],
{k, Length[hedra]}];

Suppose the 5th zonohedron looks good.

cendex=kgonal[[ First@Flatten[Position[ polyflat[[ kgonal ]], hedra[[5,1]] ]] ]];

Show[Graphics3D[ polyflat[[cendex,2]] ],
ViewPoint->{0,1,2},
Boxed->False];

4. Find and Display k-gon; Find Opposite k-gon

Having found a central zonohedron, we now find the k-gon and its opposite.

selected=Flatten[Position[countsides[ polyflat[[cendex]] ], k]]

Show[Graphics3D[ polyflat[[ cendex,2,selected[[1]] ]] ],
Boxed->False];

chos=mysurf[ polyflat[[ cendex,2,selected[[1]] ]] ];
opps=Select[polyflat[[cendex,2]], (Chop[ mysurf[#] + chos]=={0,0,0} || Chop[ mysurf[#] - chos]=={0,0,0})&]

5. Construct Solid Zone of Zonohedra and Display

Two While loops find the zonohedra within the zone, first to one side, then to the other side, of the central zonohedron. A list of indices into polyflat is returned from each While, and then the two are combined into the complete list, bothdex.

(*initialize*)
refpt=mycent[ opps[[1]] ];
newdex=cendex;
alldex={};
seldex={0,0};

While[Length[seldex]>1,

seldex=DeleteCases[
Table[
If[
Select[polyflat[[k,2]], ( Chop[mycent[#]-refpt]=={0,0,0} )&] != {},
k],
{k, Length[polyflat]}], Null];

If[Length[seldex]>1,
newdex=First@Complement[seldex,{newdex}];
alldex=Union[alldex,seldex];
refpt=Chop[ (polyflat[[newdex,1]]-refpt)+polyflat[[newdex,1]] ],

(*else*)
Print["Empty set reached after ", Length[alldex], " zonohedra found"]
]

];
fdex=alldex

(*initialize*)
refpt=mycent[ opps[[2]] ];
newdex=cendex;
alldex={};
seldex={0,0};

While[Length[seldex]>1,

seldex=DeleteCases[
Table[
If[
Select[polyflat[[k,2]], ( Chop[mycent[#]-refpt]=={0,0,0} )&] != {},
k],
{k, Length[polyflat]}], Null];

If[Length[seldex]>1,
newdex=First@Complement[seldex,{newdex}];
alldex=Union[alldex,seldex];
refpt=Chop[ (polyflat[[newdex,1]]-refpt)+polyflat[[newdex,1]] ],

(*else*)
Print["Empty set reached after ", Length[alldex], " zonohedra found"]
]

]
sdex=alldex

bothdex=Union[fdex,sdex]
Length[bothdex]

Show[Graphics3D[ polyflat[[bothdex,2]] ],
Axes->True,
Boxed->True,
ViewPoint->2.5*unit[{0,0,3}]];

References

The latest versions of both this notebook (3D_Zonotiles.nb) and my 2D zonotiles notebook (2D_Zonotiles.nb) can be found using the link below.

http://library.wolfram.com/infocenter/MathSource/1197/

The Regular Polytopes notebook can be found at the link below. Create hidden-detail-removed projections of convex four-dimensional polytopes which are reminiscent of 3D zonotiles.

http://library.wolfram.com/infocenter/MathSource/593/

My Rhombic Spirallohedra notebook can be found at the link below. These non-convex polyhedra, when considered to arise from the stripping-away of rhombic hexahedra from polar zonohedra, exemplify what I call “bitten” zonotiles, which seem to arise in any n-dimensional space; and some, at least, are space-fillers. An attempt is made in the Zonotile Examples section, above, to create bitten zonotiles from (of course) an arrangement of planes. Bitten 2D zonotiles are created from line arrangements in my 2D zonotiles notebook (link above).

http://library.wolfram.com/infocenter/MathSource/613/

Some other references bearing upon zonotopal tilings:

J. RICHTER-GEBERT AND G. ZIEGLER, Zonotopal Tilings and the Bohne-Dress Theorem,
in Jerusalem Combinatorics 1993 (H. Barcelo and G. Kalai, eds.), pp. 211–232,
Contemporary Mathematics 178, Amer. Math. Soc. 1994.

A. Bjørner, M. Las Vergnas, B. Sturmfels, N. White and G.M. Ziegler, Oriented Matroids,
Encyclopedia of Mathematics and Its Applications 46, Cambridge University Press, Cambridge, 1993.

J. RAMBAU AND G. M. ZIEGLER, Projections of Polytopes and the Generalized Baues Conjecture, Discrete Comput. Geom. 16 (1996), 215–237.

N.G. de Bruijn. Algebraic Theory of Penrose's Non-periodic Tilings of the Plane.
Kon. Nederl. Akad. Wetensch. Proc. Ser. A, 43(84), 1981.

N.G. de Bruijn. Dualization of Multigrids. J. Phys. France, C(47):3{9, 1981.