Planar Network Systems

Name: Maxim Piskunov
Instructor: Todd Rowland
Wolfram Science Summer School 2014

Homework Solution

OuterCA[rule_, init_, steps_] := CellularAutomaton[{rule, {2, {{2, 2, 2}, {2, 1, 2}, {2, 2, 2}}}, {1, 1}}, init, {{{steps}}, All, All}]
InitCA[] := ReplacePart[Array[0&, {1000, 1000}], {{500, 500}  1, {501, 500}  1, {500, 501}  1, {501, 502}  1}];
PlotRule28503[] := ArrayPlot[OuterCA[28503, InitCA[], 4000], PixelConstrained  1]
PlotRule28503[]

Project Description

To figure out the fundamental theory of physics, we need to explore network systems (see chapter 9 of NKS).
In this project I studied a few rules for network systems which preserve planarity:







starting from a hexagonal grid:
Assuming that networks come into an equilibrium after a certain amount of steps, the goal was to compute some properties of this equilibrium, such as the distribution of face sizes, and local dimensionalities.

Code

SetDirectory["/Users/maxitg/Projects/A New Kind of Science/Summer School"];
Needs["GeneralUtilities`"];

Rules

◼

Hexagonal Grid

◼
HexagonalGrid[width_?EvenQ, height_?EvenQ] := OrientedGraph​​ AssociationTableIfEvenQ[i + Ceiling[i / width]],​​ Int[i]  ​​ If[i - width ≥ 1, {Int[i - width], 3}, {Ext[i], 0}],​​ If[i + width ≤ width height, {Int[i + width], 1}, {Ext[Mod[i - 1, width] + 1 + width], 0}],​​ IfMod[i, width]  0, Ext
Quotient[i, width] + 4 width + height
2
, 0, {Int[i + 1], 2}​​ ,​​ Int[i]  ​​ If[i - width ≥ 1, {Int[i - width], 2}, {Ext[i], 0}],​​ IfMod[i, width]  1, Ext
Quotient[i, width] + 4 width + 1
2
, 0, {Int[i - 1], 3},​​ If[i + width ≤ width height, {Int[i + width], 1}, {Ext[Mod[i - 1, width] + 1 + width], 0}]​​ , {i, 1, width height}​​  /. Missing[]  {Ext[RandomInteger[{0,
64
2
-1}]], 0},​​ Association @ Table​​ Ifi ≤ width, Ext[i]  {Int[i], 1},​​ Ifi ≤ 2 width, Ext[i]  {Int[width (height - 2) + i], Mod[i, 2] + 2},​​ Ifi ≤ 2 width +
height
2
, Ext[i]  {Int[(i - 2 width) 2 width - width + 1], 2},​​ Ext[i]  Inti - 2 width -
height
2
2 width, 3​​  ​​ ​​ ,​​ {i, 1, 2 width + height}​​

Balls

◼
Neighborhood[OrientedGraph[intEdges_, extEdges_], v_, radius_] := Block{queue = <|v  0|>, position = 1},​​ IfMissingQ @ Catch @ Whilequeueposition < radius,​​ Do​​ If!KeyExistsQ[queue, intEdges[Keys[queue]〚position〛]〚i, 1〛], ​​ IfintEdges[Keys[queue]〚position〛]i, 1, 0 === Ext,​​ Throw[Missing[]],​​ queue[intEdges[Keys[queue]〚position〛]〚i, 1〛] = queueposition + 1​​ ​​ ;​​ , {i, 1, 3};​​ position++;​​ , Missing["Not available. Out of range."], queue​​
NeighborhoodSize[graph_, v_, radius_] := With[{nhood = Neighborhood[graph, v, radius]}, If[MissingQ[nhood], nhood, Length @ nhood]]

Enumeration

◼
EdgeLess[Int[_], {Ext[_], _}, {Int[_], _}] := False;​​EdgeLess[Int[_], {Int[_], _}, {Ext[_], _}] := True;​​​​EdgeLess[Int[_], {Ext[v1_], 0}, {Ext[v2_], 0}] := v1 < v2;​​​​EdgeLess[Int[v_], {Int[v_], _}, {Int[v2_], _}] /; v ≠ v2 := True;​​EdgeLess[Int[v_], {Int[v1_], _}, {Int[v_], _}] /; v ≠ v1 := False;​​​​EdgeLess[Int[v_], {Int[v_], i1_}, {Int[v_], i2_}] := i1  2 && i2  1;​​​​EdgeLess[Int[v_], {Int[v1_], i1_}, {Int[v2_], i2_}] /; v ≠ v1 && v ≠ v2 := If[v1 < v2, True,​​ If[v2 < v1, False,​​ If[i1 < i2, True, False​​]]];
HexagonalGridSubgraph[n_, i_] := Block[{v, randomVertex, j, index = i, selected = {}, grid = HexagonalGrid[4n, 4n], it = 1, internal, external},​​ AppendTo[selected, Int[2n(2n + 1)]];​​ For[v = 1, v ≤ n-1, v++,​​​​ randomVertex = Mod[index, Length @ selected] + 1; index = Quotient[index, Length @ selected];​​ j = Mod[index, 3] + 1; index = Quotient[index, 3];​​ If[FreeQ[selected, grid〚1〛[selected〚randomVertex〛]〚j, 1〛],​​ AppendTo[selected, grid〚1〛[selected〚randomVertex〛]〚j, 1〛],​​ Return[Missing["Impossible."]];​​ ]​​ ​​ ];​​ internal = Association[Normal[KeyTake[grid〚1〛, selected]] /. {x_Int?(Not[KeyExistsQ[KeyTake[grid〚1〛, selected], #]]&), ind_}  {Ext[it++], 0}];​​ external = Association[MapIndexed[Ext[#2〚1〛]  {#1〚1, 1〛, #1〚2〛}&, Position[internal, Ext[_]]]];​​ OrientedGraph[internal, external]​​]
Get["HexagonalGridSubgraphs.m"]
HexagonalGridSubgraphs[n_] := HexagonalGridSubgraphs[n] = ​​ Union ​​ @ Select[Not @* MissingQ]​​ @ Table[HexagonalGridSubgraph[n, itt], {itt, 0,
n - 1
3
(n - 1)! - 1}]​​;
HexagonalGridSubgraphsSave[] := Save["HexagonalGridSubgraphs.m", HexagonalGridSubgraphs]

Visualization

◼
ExtendColorSequence[{color1_, color2_}, n_] := TableRGBColorList@@color2
i - 1
n - 1
+ List@@color1
n - i
n - 1
, {i, 1, n}
OrientedGraph /: Graph[OrientedGraph[edges_, extEdges_], layout_:"SpringElectricalEmbedding"] := Graph[​​ Union[Keys[edges], Keys[extEdges]],​​ ​​ Select[EdgeLess[Int[RandomInteger[2^64]], Sequence @@ #] &] @​​ Flatten @​​ Table[Thread[Table[{v, i}, {i, 1, 3}]  edges[v]], {v, Keys[edges]}]​​  /. {v_, i_Integer}  v,​​ GraphLayout  layout​​];
Options[OrientedGraphPlot] = ​​ EdgePlotPoints  2,​​ VertexLabels  False,​​ GraphLayout  "SpringElectricalEmbedding",​​ Coloring  "None",​​ DimensionalityRange  10,​​ HighlightList  {},​​ ImageSize  Automatic​​;​​OrientedGraphPlot[OrientedGraph[edges_, extEdges_], OptionsPattern[]] := Block{vertexCoordinates, fullNameList, vertexIndexes, edgeLine},​​ fullNameList = Union[Keys[edges], Keys[extEdges]];​​ vertexCoordinates = VertexCoordinates /. AbsoluteOptions[Graph[OrientedGraph[edges, extEdges], OptionValue[GraphLayout]], VertexCoordinates];​​ vertexIndexes = Association[Table[fullNameList〚i〛  i, {i, 1, Length @ fullNameList}]];​​ edgeLine[{v1_, i1_}, {v2_, i2_}] := LineIfOptionValue[EdgePlotPoints]  2,​​ {vertexCoordinates〚vertexIndexes[v1]〛, vertexCoordinates〚vertexIndexes[v2]〛},​​ BezierFunction​​ vertexCoordinatesvertexIndexes[v1],​​ vertexCoordinatesvertexIndexes[v1] + IfHead@v1 === Int, .5 FlattenRotationMatrix
2 π
3
(i1 - 1) . {{0}, {1}}, 0,​​ vertexCoordinatesvertexIndexes[v2] + IfHead@v2 === Int, .5 FlattenRotationMatrix
2 π
3
(i2 - 1) . {{0}, {1}}, 0,​​ vertexCoordinatesvertexIndexes[v2]​​  /@ Range0, 1,
1.
OptionValue[EdgePlotPoints] - 1
,​​ SwitchOptionValue[Coloring],​​ "Orientation", VertexColors  ExtendColorSequence[​​ {i1, i2} /. {{0, 0}  {RGBColor@@{.5, .5, .5}, RGBColor@@{.5, .5, .5}}, 1  Red, 2  Green, 3  Blue, 0  RGBColor@@{1., 1., 1.}},​​ OptionValue[EdgePlotPoints]​​ ],​​ "Dimensionality", VertexColors  WithdimValues = Table​​ Mean @ ​​
Log[N @ NeighborhoodSize[OrientedGraph[edges, extEdges], If[Head[v1] === Int, v1, v2], OptionValue[DimensionalityRange]]]
Log[OptionValue[DimensionalityRange] + 1.]
,​​
Log[N @ NeighborhoodSize[OrientedGraph[edges, extEdges], If[Head[v2] === Int, v2, v1], OptionValue[DimensionalityRange]]]
Log[OptionValue[DimensionalityRange] + 1.]
​​ ,​​ {OptionValue[EdgePlotPoints]}​​ , dimValues /. d_Real  If1  3., Hue[.5], Hue
d - 1.
3. - 1.
,​​ "Highlight", VertexColors  Table[​​ If[FreeQ[OptionValue[HighlightList], v1] || FreeQ[OptionValue[HighlightList], v2], Black, Red]​​ , {OptionValue[EdgePlotPoints]}​​ ],​​ "BallHighlight", VertexColors  Table​​ RGBColor​​ Min
1
Log[1. Lookup[OptionValue[HighlightList], v1, ∞] + ]
,
1
Log[1. Lookup[OptionValue[HighlightList], v2, ∞] + ]
,​​ 1 - Min
1
Log[1. Lookup[OptionValue[HighlightList], v1, 0] + ]
,
1
Log[1. Lookup[OptionValue[HighlightList], v2, 0] + ]
,​​ 0​​ ,​​ {OptionValue[EdgePlotPoints]}​​ ,​​ "None", VertexColors  Table[Black, {OptionValue[EdgePlotPoints]}] ​​ ​​ ;​​ Graphics[​​ Flatten @ Table[​​ Table[edgeLine[{v, i}, edges[v]〚i〛], {i, 1, Length @ edges[v]}]​​ , {v, Keys[edges]}​​ ],​​ Table[edgeLine[{v, 0}, extEdges[v]], {v, Keys[extEdges]}],​​ If[OptionValue[VertexLabels], Table[Text[fullNameList〚i〛, vertexCoordinates〚i〛, {Left, Top}], {i, 1, Length @ fullNameList}], {}],​​ , ImageSize  OptionValue[ImageSize]]​​;

Evolution

◼
GraphExtrapolate[​​ graph_,​​ subgraph_,​​ {subV_, subI_}  {superV_, superI_}​​] := Block{stack = {{subV, subI}}, isomorphism = <|{subV, subI}  {superV, superI}|>},​​ WhileLength[stack] > 0, Block{v = First @ Last @ stack, i = Last @ Last @ stack},​​ stack = Most @ stack;​​ If[Head[v] === Int && !KeyExistsQ[isomorphism, {v, Mod[i, 3] + 1}],​​ If[Head[isomorphism[{v, i}] // First] === Ext, Return[Missing["Impossible. Out of range."]]];​​ isomorphism[{v, Mod[i, 3] + 1}] = {First @ isomorphism[{v, i}], Mod[Last @ isomorphism[{v, i}], 3] + 1};​​ AppendTo[stack, {v, Mod[i, 3] + 1}];​​ ];​​ IfHead[v] === Int && !KeyExistsQ[isomorphism, (First @ subgraph)[v]〚i〛],​​ If[Head[isomorphism[{v, i}] // First] === Ext, Return[Missing["Impossible"]]];​​ isomorphism[(First @ subgraph)[v]〚i〛] = (First @ graph)[isomorphism[{v, i}] // First]isomorphism[{v, i}] // Last;​​ AppendTo[stack, (First @ subgraph)[v]〚i〛];​​ ;​​ ;​​ If[DuplicateFreeQ[Values[isomorphism]], isomorphism, Missing["Impossible. Different topology."]]​​
GraphEvolve[graph_, inGraph_  outGraph_] := Block{isomorphism, newIntEdges = First @ graph, newExtEdges = Last @ graph},​​ evolutionCounter++;​​ If[Mod[evolutionCounter, 100]  0, slowEvolutionCounter = evolutionCounter];​​ While[MissingQ[isomorphism = GraphExtrapolate[​​ graph,​​ inGraph,​​ {First @ Keys[inGraph〚1〛], 1}  {(Keys @ graph〚1, {RandomInteger[{1, Length @ graph〚1〛}]}〛)〚1〛, RandomChoice @ {1, 2, 3}}​​ ]]];​​​​ beforeLastEvolution = First /@ Values[isomorphism]; ​​​​ (* removing old stuff *)​​ Do​​ Do​​ IfHead @ First @ isomorphism @ (First @ inGraph)[vertex]i === Int,​​ newIntEdges[(*EchoHeld @ *)First @ isomorphism @ (First @ inGraph)[vertex]〚i〛] = ReplacePart​​ newIntEdges[First @ isomorphism @ (First @ inGraph)[vertex]〚i〛],​​ Last @ isomorphism @ (First @ inGraph)[vertex]i  {}​​ ,​​ newExtEdges[First @ isomorphism @ (First @ inGraph)[vertex]〚i〛] = {}​​ ;​​ (*Print["--------------------- New Step ---------------------"];​​ Print[newIntEdges];​​ Print[newExtEdges];*)​​ , {i, 1, 3};​​ , {vertex, Keys[inGraph〚1〛]};​​ ​​ Do[​​ newIntEdges[First @ isomorphism[{vertex, 1}]] =.;​​ (*Print["--------------------- Deleting ---------------------"];​​ Print[newIntEdges];​​ Print[newExtEdges];*)​​ , {vertex, Keys[inGraph〚1〛]}];​​​​ (* adding new stuff *)​​ Do[Block[{newName = Int @ RandomInteger[{0,
64
2
-1}]},​​ Do[​​ isomorphism[{vertex, i}] = ReplacePart[isomorphism[{vertex, i}], 1  newName];​​ , {i, 1, 3}];​​ ], {vertex, Keys[outGraph〚1〛]}];​​ ​​ afterLastEvolution = First /@ Values[isomorphism]; ​​​​ Do[​​ newIntEdges[First @ isomorphism[{vertex, 1}]] = {0, 0, 0};​​ Do[​​ newIntEdges[First @ isomorphism[{vertex, i}]] = ReplacePart[​​ newIntEdges[First @ isomorphism[{vertex, i}]],​​ Last @ isomorphism[{vertex, i}]  isomorphism[(First @ outGraph)[vertex]〚i〛]​​ ];​​ , {i, 1, 3}];​​ (*Print["--------------------- Inserted ---------------------"];​​ Print[newIntEdges];​​ Print[newExtEdges];*)​​ , {vertex, Keys[outGraph〚1〛]}];​​​​ Do​​ Do​​ IfHead @ First @ isomorphism @ (First @ outGraph)[vertex]i === Int,​​ newIntEdges[First @ isomorphism @ (First @ outGraph)[vertex]〚i〛] = (*EchoHeld@*)ReplacePart​​ newIntEdges[First @ isomorphism @ (First @ outGraph)[vertex]〚i〛],​​ Last @ isomorphism @ (First @ outGraph)[vertex]i  isomorphism[{vertex, i}]​​ ,​​ newExtEdges[First @ isomorphism @ (First @ outGraph)[vertex]〚i〛] = isomorphism[{vertex, i}];​​ ;​​ (*Print["--------------------- New Step ---------------------"];​​ Print[newIntEdges];​​ Print[newExtEdges];*)​​ , {i, 1, 3};​​ , {vertex, Keys[outGraph〚1〛]};​​ ​​ OrientedGraph[newIntEdges, newExtEdges]​​ 
Get["GraphEvolve.m"];
GraphEvolve[graph_, inGraph_  outGraph_, steps_] := GraphEvolve[graph, inGraph  outGraph, steps] = ​​ evolutionCounter = 0;​​ started = AbsoluteTime[];​​ slowEvolutionCounter = 1;​​ PrintTemporary​​ steps,​​ Dynamic[evolutionCounter],​​ DynamicDateString
AbsoluteTime[] - started
slowEvolutionCounter / steps
+ started,​​ DynamicDateDifferenceDateList[], DateString
AbsoluteTime[] - started
slowEvolutionCounter / steps
+ started, {"Year", "Day", "Hour", "Minute", "Second"}​​ ;​​ Nest[GraphEvolve[#, inGraph  outGraph] &, graph, steps]​​
GraphEvolveSave[] := Save["GraphEvolve.m", GraphEvolve]

Faces

◼
FaceWalk[OrientedGraph[intEdges_, extEdges_], {v_?(Head[#] === Ext&), i_}] := Missing["Not available. Out of range."];
FaceWalk[OrientedGraph[intEdges_, extEdges_], {v_?(Head[#] === Int&), i_}] := ​​ First @ intEdges[v]i,​​ Mod[Last @ intEdges[v]〚i〛, 3] + 1​​;
GraphFace[graph_, {v_, i_}] := With[​​ {face = RotateRight @ NestWhileList[FaceWalk[graph, #]&, FaceWalk[graph, {v, i}], # =!= {v, i} && !MissingQ[#]&]},​​ If[MissingQ[First @ face], First @ face, face]​​]
GraphFacesList[graph_] := Union @ RotateLeft[#, Ordering[#, 1, EdgeLess[Int[0], #1, #2]&] - 1]& /@ Select[!MissingQ[#]&] @ Flatten[Table[​​ GraphFace[graph, {v, i}],​​ {v, Keys[graph〚1〛]},​​ {i, 1, 3}​​], 1];
GraphFacesLengths[graph_] := Sort[Length /@ GraphFacesList[graph]]

Dimentionality

◼
BallVolumes[graph_, n_] := Select[Not @* MissingQ] @ (NeighborhoodSize[graph, #, n]& /@ Keys[graph〚1〛])
BallDimensionalities[graph_, n_] := Log[BallVolumes[graph, n]] / Log[n + 1]

Visualizations

Evolution after the first few steps (for the first rule):
Equilibrium states (coloring stands for local dimensionality:
):