Zonohedrification

George W. Hart,                               August, 1996
http://www.georgehart.com

Reference

"Zonohedrification" submitted to Mathematica Journal 
for publication, by George W. Hart, August 1996.

Discussion

Runs under Mathematica Kernel 2.1, Front end 2.2
See above reference for explanation of algorithm.

Implementation

SetDirectory["c:/mathemat/packages/polyhedr"]

<< "Graphics`Polyhedra`"

builtInPoly[polyhedron_] := 
      Map[Vertices[polyhedron][[#1
            ]]&, Faces[polyhedron], {2}]

vertices[poly_] := Union[Flatten[poly, 1]]

view[poly_] := Show[Graphics3D[Polygon/@poly], Boxed -> False]

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

tolerance = 0.000001;

collinear[v1_, v2_] := 
      And@@(Abs[#1] < tolerance&)/@cross[v1, v2]

setStar[vlist_] := 
      Module[{selected = {}}, 
        Scan[Function[v, 
            If[v ≠ {0, 0, 0} && 
                Select[selected, collinear[v, #1]&] == {}, 
              AppendTo[selected, v]]], vlist]; 
        Print[Length[selected], " zonal directions."]; 
        gStar = (selected/2)]

pToXYZ[p_] := Plus@@MapThread[Times, {gStar, p}]

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}

allPairs[n_] := 
      With[{indices = Table[i, {i, n}]}, 
        Select[Flatten[Outer[List, indices, indices], 
            1], #1[[1]] < #1
                [[2]]&]]

makeFace[{i_, j_}] := 
      With[{normal = 
            cross[gStar[[i]], 
              gStar[[j]]]}, 
        (approxSign[normal . #1]&)/@gStar]

approxSign[x_] := If[Abs[x] < tolerance, 0, Sign[x]]

allFaces := 
      With[{halfTheFaces = makeFace/@allPairs[Length[gStar]]}, 
        halfTheFaces ⋃ -halfTheFaces]

makeEdge[i_, p_, faceNormal_] := 
      With[{normal = 
            cross[gStar[[i]], 
              faceNormal]}, 
        MapThread[
          If[#1 ≠ 0, #1, approxSign[#2 . normal]]&, {p, gStar}]]

makeEdges[p_] := 
      With[{normal = faceNormal[p]}, 
        With[{halfEdges = 
              (makeEdge[#1, p, normal]&)/@Flatten[Position[p, 0]]}, 
          halfEdges ⋃ (2 p - #1&)/@halfEdges]]

makePolygon[p_] := 
      (pToXYZ[#1]&)/@
        With[{edges = makeEdges[p]}, 
          With[{vertices = Union[Flatten[endPoints/@edges, 1]]}, 
            Select[cycleSort[Join[edges, vertices]], Count[#1, 0] == 0&]]]

allEdges := Union[Flatten[makeEdges/@allFaces, 1]]

lineDrawing := (Line[pToXYZ/@endPoints[#1]]&)/@allEdges;

cycleSort[list_] := 
      Module[{result = {list[[1]]}, 
          remaining = Rest[list]}, 
        While[remaining ≠ {}, 
          choices = 
            Select[remaining, 
              incident[#1, 
                  result[[1]]]&]; 
          If[choices ≠ {}, 
            next = choices[[1]], 
            Print["zono-error!:", list]; Break[]]; 
          PrependTo[result, next]; 
          remaining = Complement[remaining, {next}]]; result]

incident[p1_, p2_] := Plus@@Abs[p1 - p2] == 1

Null

zonohedrify[poly_] := 
      Module[{plist}, setStar[vertices[poly]]; plist = allFaces; 
        printSummary[plist]; makePolygon/@plist]

printSummary[p2list_] := 
      With[{maxSize = Max[(Count[#1, 0]&)/@p2list]}, 
        Do[With[{howMany = Length[Select[p2list, Count[#1, 0] == i&]]}, 
            If[howMany > 0, 
              Print["     ", 2 i, "-gons:  ", howMany]]], {i, 
            maxSize, 2, -1}]; Print["Total faces:  ", Length[p2list]]]

polarStar[n_] := 
      N[Table[{Sin[((2 π i)/n)], Cos[((2 π i)/n)], 1}, {i, 1, 
            n}]];

faceCenters[poly_] := (Plus@@#(1/Length)[#1]&)/@poly

edgeMidpoints[poly_] := 
      Union[Flatten[((1/2) (#1 + RotateLeft[#1])&)/@poly, 1]]

unit[v_] := (v/(Sqrt[v . v]))

viewStar[star_] := 
      Show[Graphics3D[(Line[{{0, 0, 0}, #1}]&)/@star, 
          Boxed -> False]];

Use the functions

view[builtInPoly[Icosahedron]]

view[zonohedrify[builtInPoly[Octahedron]]]

view[zonohedrify[{edgeMidpoints[builtInPoly[Cube]]}]]

truncatedCuboctahedron = 
      zonohedrify[{
          unit/@Join[faceCenters[builtInPoly[Cube]], 
              edgeMidpoints[builtInPoly[Cube]]]}];

view[truncatedCuboctahedron]

view[zonohedrify[{edgeMidpoints[builtInPoly[Icosahedron]]}]]

view[zonohedrify[builtInPoly[Cube]]]

star = 
      N[{{1, 0, 0}, {0, 1, 0}, {((Sqrt[2])/2), ((Sqrt[2])/2), 0}, {((Sqrt[2])/2), -((Sqrt[2])/2), 
            0}, {0, 0, 1}}]; view[zonohedrify[{star}]]

view[Nest[zonohedrify, builtInPoly[Octahedron], 3]]

star = 
      N[{{1, 0, 1}, {-1, 0, 1}, {0, 1, 1}, {0, -1, 1}, {0, 0, 2}}]; 
    view[zonohedrify[{star}]]

rhombicEnneacontahedron = zonohedrify[builtInPoly[Dodecahedron]]; 
    view[rhombicEnneacontahedron]

rhombicTriacontahedron = zonohedrify[builtInPoly[Icosahedron]]; 
    view[rhombicTriacontahedron]

RTstar = gStar; view[zonohedrify[{Drop[RTstar, -1]}]]

view[zonohedrify[{Drop[RTstar, -2]}]]

view[zonohedrify[{Drop[RTstar, -3]}]]

starFig3a = {{1, 0, 0}, {0, 1, 0}, {0, 0, 1}, {-(1/2), (1/2), 
          (1/2)}}; view[zonohedrify[{starFig3a}]]

starFig3b = {{1, 0, 0}, {0, 1, 0}, {0, 0, 1}, {1, 1, 0}}; view[zonohedrify[{starFig3b}]]

setStar[starFig3a]; Show[Graphics3D[lineDrawing, Boxed -> False]]

setStar[starFig3b]; Show[Graphics3D[lineDrawing, Boxed -> False]]

viewStar[starFig3a]

viewStar[starFig3b]

view[zonohedrify[rhombicTriacontahedron]]

view[
      zonohedrify[{
          unit/@Join[vertices[builtInPoly[Icosahedron]], 
              faceCenters[builtInPoly[Icosahedron]], 
              edgeMidpoints[builtInPoly[Icosahedron]]]}]]

view[zonohedrify[truncatedCuboctahedron]]

view[zonohedrify[{polarStar[10]}]]

view[
      zonohedrify[{
          unit/@Join[vertices[builtInPoly[Icosahedron]], 
              edgeMidpoints[builtInPoly[Icosahedron]]]}]]

view[zonohedrify[rhombicEnneacontahedron]]