Kenyon’s squarespirals

The other day by chance I happened to look at Richard Kenyon’s web page, and was struck by a very beautiful animated image there. The image is of a region tiled by colored squares, which are slowly rotating. As the squares rotate, they change size in such a way that the new (skewed, resized) squares still tile the same region. I thought it might be fun to try to guess how the image was constructed, and to produce my own version of his image.

I already know a little bit about square tilings. This is a subject with a history, going back at least to the work of Tutte and his colleagues. The basic problem is just to tile a rectangular region by squares. Easy enough, you say.

Well, yes; if the rectangle has sides which are rationally related, in can be filled up by squares with commensurable side lengths pretty easily. Here a 4 by 8 rectangle is filled with 7 squares of edge length 1, 1 square of edge length 3, and 1 square of edge length 4. It’s more amusing to look for a tiling in which all the squares have different lengths. One well-known tiling, found by Tutte’s colleague Stone, is as follows:

Obviously the problem becomes more interesting and challenging if one starts in advance with the combinatorics of a square tiling, and then tries to assign edge lengths to squares in such a way that they fit together nicely. One elegant method, developed by Brooks, Smith, Stone and Tutte, assigns a directed graph to the tiling, with one vertex for each vertical edge (say) and one directed edge for each square. Here’s the graph associated to Stone’s tiling:

The condition that the sum of square lengths on either side of a vertical edge sum to the length of that edge implies that the incoming edge weights and the outgoing edge weights at each vertex sum to the same value (except for at the leftmost and rightmost vertices). On the other hand, the fact that the squares are all square implies that each edge weight (as above) is equal to the length of its projection to a horizontal line; this means that the sum of edge weights around each loop in the graph (with sign changed when the orientation disagrees with the orientation on the loop) is equal to zero. These two conditions are precisely Kirchoff’s two laws for the current flowing through an electrical network where every edge has resistance 1, and the voltage difference between left and right vertices is the width of the rectangle. There is a unique solution; it might have some weights negative, in which case we can reverse the orientation of the edge so that the weights are all positive, and determine a square tiling with slightly different combinatorics. By the way, the uniqueness of the solution has an interesting (and well-known) consequence: since Kirchoff’s laws both impose linear conditions on the edge weights, the space of solutions is a rational affine space (in units for which the width is equal to 1). Since this space of solutions consists of a single point, this point has rational coordinates; this implies in particular that the height of the rectangle is a rational multiple of the width, and so are the widths of the squares.

In more homological language, the assignment of weights to edges is a (simplicial) 1-chain. The condition that the incoming and outgoing edge weights at each vertex have equal sum says that this 1-chain is actually a (relative) 1-cycle; i.e. that it is closed. The condition that the sum around every loop is zero says that if we think of this 1-chain as a 1-cochain it is actually a 1-cocycle; i.e. it is co-closed. A (co)-chain which is both closed and co-closed is said to be harmonic, and the uniqueness of a solution corresponds to the uniqueness of a harmonic representative of a (relative co-) homology class.

Incidentally, if we form the graph with one vertex for each vertical horizontal line and one edge for each square, this will be the (planar) dual to the graph above. Edges in one graph correspond to edges in the other, and the closed condition for one set of edge weights becomes the co-closed condition for the other, and vice versa.

Now instead of considering a square tiling of a rectangle, let’s consider a square tilings of a Euclidean torus. A combinatorial tiling gives us a graph, and a harmonic 1-cycle gives us a square tiling with the desired combinatorics. Changing the 1-cycle by rescaling it just rescales the torus and all the squares by the same factor, which is not very interesting. However, there is something interesting we can do. The homology of a torus is 2-dimensional, so we can consider a 1-parameter family of homology classes whose projective classes are changing, and a 1-parameter family of harmonic 1-cycles and of square tilings.

Let’s start with the simplest possible example. We fix a graph G embedded in the torus. Since we want G to be able to carry every homology class, we need at least two edges. So let’s take as G the graph with one vertex and two edges, one of which wraps horizontally once around the torus, and one of which wraps vertically around. Any assignment of weights to the edges will be both closed and co-closed, so a 1-parameter family is given by taking weights cos(t), sin(t) for t in the unit circle. The resulting square tilings of the torus have two squares, one of side length cos(t) and one of side length sin(t). The total area of the torus is thus normalized to be 1. The pattern of tilings “rotates” with t as follows:

(click on the image to see it rotate)

OK, how about a more complicated example? Let’s let G be some complicated embedded graph on the torus (so that it can carry any homology class). For the sake of concreteness, let’s let G be the following graph:

G has 10 edges (corresponding to 10 squares in the tiling), 5 vertices and 5 complementary faces. There are 5 vertex conditions and 5 face conditions; however, this system of 10 equations is redundant, and has a 2 dimensional space of solutions.

Weights on the edges of G form a vector space, and there is an inner product on this space which is just the ordinary Euclidean inner product with co-ordinates the weights on each each edge. We want to normalize our weights to have length (i.e. square root of their inner product with themselves) equal to 1, so that the resulting torus will have area 1. All we need to do is find two orthogonal weights M and L which are closed and co-closed, orthogonal to each other (i.e. the inner product of M and L is zero) and of length 1, and then we can form the family cos(t)M + sin(t)L of weights, and the associated square tilings.

The resulting rotating family of tilings is as follows:

(click on image to see it rotate)

Something else is needed to get the “spiraling” evident in Kenyon’s picture. For our square tilings of a torus above, the result of laying down a sequence of squares that winds once around a loop in the torus is to displace the tiling by a translation of the plane; this translation is called the holonomy around the loop, and only depends on its homotopy class (actually: on its homology class). Essentially, this is the result of integrating the (dual) 1-form associated to the weight. An educated guess is that in Kenyon’s picture, the holonomy is not a translation, but rather a dilation of the plane, centered at some point. At the level of homology, one can think of the dilation factor around a loop as a representation of the fundamental group, and we need to consider (harmonic) 1-cycles with coefficients twisted by this representation.

How to translate this into the language of square tilings and weights? Instead of thinking of a weight on the graph G, let’s let G~ denote the lift of G to the universal cover of the torus; i.e. G~ is a periodic graph in the plane. A twisted weight on G with coefficients in a representation is the same thing as a weight on G~ that transforms according to the given representation. For the sake of simplicity, let’s work with the graph G with one vertex and two edges as in the first example above, so that G~ has one vertex, one horizontal edge, and one vertical edge for each pair of integers. Pick a pair of edges H,V of G~, going to the right and up incoming to the vertex (0,0) respectively and let h,v be the weights on these edges.

If we let A denote the multiplication factor for horizontal translation, and B the multiplication factor for vertical translation, the vertex equation at (0,0) is

The vertex equations at every other vertex are obtained from this one by scaling by power of A and B, so they are satisfied if this one is. The face equation for the face with vertices (-1,-1), (0,-1), (0,0), (-1,0) is

Eliminating h from this pair of equations and dividing out by v gives

In order to enforce spiraling, we would like moving “horizontally” some fixed number of steps to be the same as moving “vertically” some (other) fixed number of steps; this can be imposed by setting for some coprime integers p,q. With these constraints, there is a unique solution h,v in complex numbers, up to scale. The real  part of any such solution gives a “spiral” tiling, and the 1-parameter family obtained by multiplying by before taking the real part gives a rotating spiral.

Let’s try an example. Taking p=2,q=1 gives and . There is a totally real solution, giving rise to the following “degenerate” spiral:

Since this solution is totally real, it can’t be “rotated”. Hmm, I wasn’t expecting that. OK, taking p=3,q=1 gives and .

(click on image to see it rotate)

Success!

Getting more squares in the picture is a matter of spiraling slower, which can be achieved by taking p and q bigger. Let’s try p=7,q=1.

(click on image to see it rotate)

If you want to have a play with this yourself, the source of the .eps file that generated these figures is below. To change the amount of spiraling, change the values of A and B, subject to the constraint that . The resulting .eps file can be transformed to a layered .pdf (eg using Preview on a Mac) then to a .gif (eg in gimp). The case q=1 is pretty easy, since then A is the root of with smallest (nonzero) argument, and . Wolframalpha will cough up the values of A and B if you coax it long enough.

(Update January 16): Just for fun, here’s the tiling with p=101, q=1 (warning: the .gif file is quite large!)

(click on image to see it rotate)

 

(Update May 8 2015): There is a lot of literature out there on what might be called “discrete complex analysis”, which includes such topics as circle packing, square packing,  discrete Dirichlet boundary value problems on planar graphs, etc. One recent example is this paper by Sa’ar Hersonsky.


%!PS-Adobe-2.0 EPSF-2.0
%%BoundingBox: 0 0 400 400
gsave
400 400 scale
1 20 div setlinewidth
1 setlinejoin
0.5 0.5 translate
/square{4 dict begin
/z exch def
/y exch def
/x exch def
gsave
newpath
rand 10 mod 10 div rand 10 mod 10 div rand 10 mod 10 div setrgbcolor
x y moveto
x z add y lineto
x z add y z add lineto
x y z add lineto
closepath
fill
stroke
grestore
end } def
/simple_edge_squares{4 dict begin
/v exch def
/n v length def
0 1 n 1 sub{
/i exch def
0 0 v i get square
0 v i get translate
} for
end} def
/rcmul{2 dict begin
/t exch def
/z exch def
[ z 0 get t mul z 1 get t mul]
end} def
/ccmul{2 dict begin
/w exch def
/z exch def
[
z 0 get w 0 get mul z 1 get w 1 get mul sub
z 0 get w 1 get mul z 1 get w 0 get mul add
]
end} def
/cconj{1 dict begin
/z exch def
[
z 0 get 0 z 1 get sub
]
end} def
/cnorm{1 dict begin % |z|^2
/z exch def
z 0 get dup mul z 1 get dup mul add
end} def
/ccdiv{2 dict begin % w/z = w*zbar/|z|^2
/z exch def
/w exch def
w z cconj ccmul 1 z cnorm div rcmul
end} def
/ccadd{2 dict begin
/w exch def
/z exch def
[ z 0 get w 0 get add z 1 get w 1 get add ]
end} def
/creal{1 dict begin
/z exch def
z 0 get
end} def
/cimag{1 dict begin
/z exch def
z 1 get
end} def
0 5 355 {
/t exch def
0 srand
gsave
/A [0.71469 -0.870643] def % A is root of x^14+x^8-4x^7+x^6+1=0
/B [0.432505 -0.0429583] def % B = A^-7
% check: A^3B=1
/h [t cos t sin] def
/v h A [-1 0] ccadd ccmul [1 0] B -1 rcmul ccadd ccdiv def %
% h*(A-1)/(1-B)
/Ainv [1 0] A ccdiv def
/Aser [1 0] Ainv ccadd Ainv Ainv ccmul ccadd Ainv Ainv ccmul Ainv ccmul ccadd Ainv Ainv ccmul Ainv ccmul Ainv ccmul ccadd Ainv Ainv ccmul Ainv ccmul Ainv ccmul Ainv ccmul ccadd def
/Acom [1 0] Ainv -1 rcmul ccadd def
/htran h Acom ccdiv def
/vtran v Acom ccdiv def
t rotate
htran creal vtran creal -1 mul translate
[h A ccmul creal v B ccmul creal h [-1 0] ccmul creal v [-1 0] ccmul creal] simple_edge_squares
1 1 50{
h -1 rcmul creal v creal translate
/h h A ccdiv def
/v v A ccdiv def
[h A ccmul creal v B ccmul creal h [-1 0] ccmul creal v [-1 0] ccmul creal] simple_edge_squares
} for
showpage
grestore
} for
grestore
%eof