Representing and Learning
Functions Invariant Under Crystallographic Groups

This section sketches our results in a purely heuristic way; proper definitions follow in Section 2.
Crystallographic symmetry. Crystallographic groups are groups that tile a Euclidean space $\mathbb{R}^{n}$ with a convex shape. Suppose we place a convex polytope $\Pi$ in the space $\mathbb{R}^{n}$, say a square or a rectangle in the plane. Now make a copy of $\Pi$, and use a transformation ${\phi:\mathbb{R}^{2}\rightarrow\mathbb{R}^{2}}$ to move this copy to another location. We require that $\phi$ is an isometry, which means it may shift, rotate or flip $\Pi$, but does not change its shape or size. Here are some examples, where the original shape $\Pi$ is marked in red:

The descriptors p1, p2, and p2mm follow the naming standard for groups developed by crystallographers [31], and the symbol “F” is inscribed to clarify which transformations are used. The transformations in these examples are horizontal and vertical shifts (in p1), rotations around the corners of the rectangle in (p2), and reflections about its edges (p2mm). Suppose we repeat one of these processes indefinitely so that the copies of $\Pi$ cover the entire plane and overlap only on their boundaries. That requires a countably infinite number of transformations, one per copy. Collect these into a set $\mathbb{G}$. If this set forms a group, this group is called crystallographic. Such groups describe all possible symmetries of crystals, and have been thoroughly studied in crystallography. For each dimension $n$, there is—up to a natural notion of isomorphy that we explain in Section 2—only a finite number of crystallographic groups: Two on $\mathbb{R}$, 17 on $\mathbb{R}^{2}$, 230 on $\mathbb{R}^{3}$, and so forth. Those on $\mathbb{R}^{2}$ are also known as wallpaper groups, and those on $\mathbb{R}^{3}$ as space groups.
The objects of interest. A function $f$ is invariant under $\mathbb{G}$ if it satisfies

A simple way to construct such a function is to start with a tiling as above, define a function on $\Pi$, and then replicate it on every copy of $\Pi$. Here are two examples on $\mathbb{R}^{2}$, corresponding to (ii) and (iii) above, and an example on $\mathbb{R}^{3}$:

However, as the examples illustrate, functions obtained this way are typically not continuous. Our goal is to construct smooth invariant functions, such as these:

We identify two representations of such functions, one linear and one nonlinear. Working with either representation algorithmically requires a data structure representing the invariance constraint. We construct such a structure, which we call an orbit graph, in Section 4. This graph is constructed from a description of the group (which can be encoded as a finite set of matrices) and of $\Pi$ (a finite set of vectors).
Linear representations: Invariant Fourier transforms. We are primarily interested in two and three dimensions, but a one-dimensional example is a good place to start: In one dimension, a convex polytope is always an interval, say ${\Pi=[0,1]}$. If we choose $\mathbb{G}$ as the group $\mathbb{Z}$ of all shifts of integer length, it tiles the line $\mathbb{R}$ with $\Pi$. In this case, an invariant function is simply a periodic function with period $1$. Smooth periodic functions can be represented as a Fourier series,

for sequences of scalar coefficients $a_{i}$ and $b_{i}$. Note each sine and cosine on the right is $\mathbb{G}$-invariant and infinitely often differentiable. Now suppose we abstract from the specific form of these sines and cosines, and only regard them as $\mathbb{G}$-invariant functions that are very smooth. The series representation above then has the general form

where the $e_{i}$ are smooth, $\mathbb{G}$-invariant functions that depend only on $\mathbb{G}$ and $\Pi$, and the $c_{i}$ are scalar coefficients that depend on $f$. (In the Fourier series, $e_{i}$ is a cosine for odd and a sine for even indices.) In Section 5, we obtain generalizations of this representation to crystallographically invariant functions. To do so, we observe that the Fourier basis can be derived as the set of eigenfunctions of the Laplace operator: The sine and cosine functions above are precisely those functions ${e:\mathbb{R}\rightarrow\mathbb{R}}$ that solve

for some ${\lambda\geq 0}$. (The negative sign is chosen to make the eigenvalues $\lambda$ non-negative.) The periodicity constraint is equivalent to saying that $e$ is invariant under the shift group ${\mathbb{G}=\mathbb{Z}}$. The corresponding problem for a general crystallographic group $\mathbb{G}$ on $\mathbb{R}^{n}$ is hence

Theorem 7 shows that this problem has solutions for any dimension $n$, convex polytope ${\Pi\subset\mathbb{R}^{n}}$, and crystallographic group $\mathbb{G}$ that tiles $\mathbb{R}^{n}$ with $\Pi$. As in the Fourier case, the solution functions ${e_{1},e_{2},\ldots}$ are very smooth.

If we choose ${\Pi\subset\mathbb{R}^{2}}$ as the square ${[0,1]^{2}}$ and $\mathbb{G}$ as the group $\mathbb{Z}^{2}$ of discrete horizontal and vertical shifts—that is, the two-dimensional analogue of the example above—we recover the two-dimensional Fourier transform. The function $e_{0}$ is constant; the functions ${e_{1},\ldots,e_{5}}$ are shown in Figure 4. If the group also contains other transformations, the basis looks less familiar. These are the basis functions ${e_{1},\ldots,e_{5}}$ for a group (p3) containing shifts and rotations of order three:

The same idea applies in any finite dimension $n$. For ${n=3}$, the $e_{i}$ can be visualized as contour plots. For instance, the first five non-constant basis elements for a specific three-dimensional group, designated I4₁ by crystallographers, look like this:

Our results show that any continuous invariant function can be represented by a series expansion in functions $e_{i}$. As for the Fourier transform, the functions form a orthonormal basis of the relevant $\mathbf{L}_{2}$ space. The functions $e_{i}$ can hence be seen as a generalization of the Fourier transform from pure shift groups to crystallographic groups. All of this is made precise in Section 5.
Nonlinear representations: Factoring through an orbifold. The second representation, in Section 6, generalizes an idea of David MacKay [45], who constructs periodic functions on the line as follows: Start with a continuous function ${h:\mathbb{R}^{2}\rightarrow\mathbb{R}}$. Choose a circle of circumference $1$ in $\mathbb{R}^{2}$, and restrict $h$ to the circle. The restriction is still continuous. Now “cut and unfold the circle with $h$ on it” to obtain a function on the unit interval. Since this function takes the same value at both interval boundaries, replicating it by shifts of integer length defines a function on $\mathbb{R}$ that is periodic and continuous:

More formally, MacKay’s approach constructs a function ${\rho:\mathbb{R}\rightarrow\text{circle}\subset\mathbb{R}^{2}}$ such that

We show how to generalize this construction to any finite dimension $n$, any crystallographic group $\mathbb{G}$ on $\mathbb{R}^{n}$, and any convex polytope with which $\mathbb{G}$ tiles the space: For each $\mathbb{G}$ and $\Pi$, there is a continuous, surjective map

such that

This is Theorem 15. Section 6.1 shows how to compute a representation of $\rho$ using multidimensional scaling.

Consider a group $\mathbb{G}$ of isometries of $\mathbb{R}^{n}$. (See Appendix A for a brief review of definitions.) Every isometry $\phi$ of $\mathbb{R}^{n}$ is of the form

Let ${M\subset\mathbb{R}^{n}}$ be a set. We say that $\mathbb{G}$ tiles the space $\mathbb{R}^{n}$ with $M$ if the image sets ${\phi M}$ completely cover the space so that only their boundaries overlap:

Each set $\phi M$ is a tile, and the collection ${\mathbb{G}M:={\{\phi M|\phi\in\mathbb{G}\}}}$ is a tiling of $\mathbb{R}^{n}$.

By a convex polytope, we mean the convex hull of a finite set of points [70]. Let ${\Pi\subset\mathbb{R}^{n}}$ be an $n$-dimensional convex polytope. The boundary $\partial\Pi$ consists of a finite number of ${(n-1)}$-dimensional convex polytopes, called the facets of $\Pi$. Thus, if $\mathbb{G}$ tiles $\mathbb{R}^{n}$ with $\Pi$, only points on facets are contained in more than one tile.

Some properties of $\mathbb{G}$ can be read right off the definition: Since $\mathbb{G}$ tiles the entire space with a set $\Pi$ of finite diameter, we must have ${|\mathbb{G}|=\infty}$. Since $\Pi$ is $n$-dimensional and convex, it contains an open metric ball of positive radius. Each tile contains a copy of this ball, and these copies do not overlap. It follows that

A group of isometries that satisfies (4) for some ${\varepsilon>0}$ is called discrete, in contrast to groups which contain, e.g., continuous rotations. Discreteness implies $\mathbb{G}$ is countable, but not all countable groups of isometries are discrete (the group $\mathbb{Q}^{n}$ of rational-valued shifts is a non-example). In summary, every crystallographic group is an infinite, discrete (and hence countable) subgroup of the Euclidean group on $\mathbb{R}^{n}$.

Suppose we choose one of the tilings in Figure 1, and rotate or shift the entire plane with the tiling on it. Informally speaking, that changes the tiling, but not the tiling mechanism, and it is natural to consider the two tilings isomorphic. More formally, two crystallographic groups $\mathbb{G}$ and $\mathbb{G}^{\prime}$ are isomorphic if there is an orientation-preserving, invertible, and affine (but not necessarily isometric) map ${\gamma:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n}}$ such that ${\mathbb{G}^{\prime}=\gamma\mathbb{G}}$, where ${\gamma\mathbb{G}:={\{\gamma\phi\gamma^{-1}\,|\,\phi\in\mathbb{G}\}}}$.

To construct a $\mathbb{G}$-invariant function, we may start with a function $h$ on $\Pi$ and “replicate it by tiling”. For that to be possible, $h$ must in turn be the restriction of a $\mathbb{G}$-invariant function to $\Pi$. It must then satisfy ${h(\phi x)=h(x)}$ if both $\phi x$ and $x$ are in $\Pi$. We hence define the relation

We note immediately that ${x\sim y}$ implies each point is also contained in an adjacent tile, so both must be on the boundary $\partial\Pi$ of $\Pi$. The requirement

is therefore a periodic boundary condition. If it holds, the function

is well-defined on $\mathbb{R}^{n}$, and is $\mathbb{G}$-invariant. Conversely, every $\mathbb{G}$-invariant function $f$ can be obtained this way (by choosing $h$ as the restriction $f|_{\Pi}$). Informally, (7) says that we stitch together function segments on tiles that are all copies of $h$, and these segments overlap on the tile boundaries. The boundary condition ensures that wherever such overlaps occur, the segments have the same value, so that (7) produces no ambiguities. The special case of (6) for pure shift groups—where $A_{\phi}$ is the identity matrix for all ${\phi\in\mathbb{G}}$—is known as a Born-von Karman boundary condition (e.g., Ashcroft and Mermin [7]).

An alternative way to express invariance is as follows: A function is $\mathbb{G}$-invariant if and only if it is constant on each set of the form

The set $\mathbb{G}(x)$ is called the orbit of $x$. We see immediately that each orbit of a crystallographic group is countably infinite, but locally finite: The definition of discreteness in (4) implies that every bounded subset of $\mathbb{R}^{n}$ contains only finitely many points of each orbit. We also see that each point ${x\in\mathbb{R}^{n}}$ is in one and only one orbit, which means the orbits form a partition of $\mathbb{R}^{n}$. The assignment ${x\mapsto\mathbb{G}(x)}$ is hence a well-defined map

The orbit set $\mathbb{R}^{n}/\mathbb{G}$ is also called the quotient set or just the quotient of $\mathbb{G}$, and $q$ is called the quotient map (e.g., Bonahon [15]). Since the orbits are mutually disjoint, we can informally think of $q$ as collapsing each orbit into a single point, and $\mathbb{R}^{n}/\mathbb{G}$ is the set of such points.

Quotient spaces are abstract but useful tools for expressing invariance properties: For any function ${f:\mathbb{R}^{n}\rightarrow\mathbb{R}}$, we have

since each point of ${\mathbb{R}^{n}/\mathbb{G}}$ represents an orbit and $f$ is invariant iff it is constant on orbits. We can also use the quotient to express continuity, by equipping it with a topology that satisfies

There is exactly one such topology, called the quotient topology in the literature. Its definition can be made more concrete by metrizing it:

Since $\mathbb{G}$ is discrete, the infimum in $d_{\mathbb{G}}$ is a minimum. The distance of two orbits (considered as points in $\mathbb{R}^{n}/\mathbb{G}$) is hence the shortest Euclidean distance between points in these orbits (considered as sets in $\mathbb{R}^{n}$), see Figure 2 (right). If $x$ and $y$ are points in the polytope $\Pi$, we have

Informally speaking, $d_{\mathbb{G}}$ implements the periodic boundary condition (6). The metric space $(\mathbb{R}^{n}/\mathbb{G},d_{\mathbb{G}})$ is also called the quotient space or orbit space of $\mathbb{G}$. A very important property of crystallographic groups is that they have compact quotient spaces:

Since orbit spaces are abstract objects, we can only work with them implicitly. One way to do so is by representing each orbit by one of its points in $\mathbb{R}^{n}$. A subset of $\mathbb{R}^{n}$ that contains exactly one point of each orbit is called a transversal. In general, transversals can be exceedingly complex sets [9], but crystallographic groups always have simple transversals. Algorithm 2 in the next section constructs a transversal explicitly. In the following, we will always write $\tilde{\Pi}$ to mean

Given such a transversal, we can define the projector ${p:\mathbb{R}^{n}\rightarrow\Pi}$ as

If we think of each point in $\tilde{\Pi}$ as a concrete representative of an element of ${\mathbb{R}^{n}/\mathbb{G}}$, then $p$ is similarly a concrete representation of the quotient map $q$, and we can translate the identities above accordingly: The projector is by definition $\mathbb{G}$-invariant, since we can write $f$ in (7) as ${f=h\circ p}$. That shows

Although $p$ is not continuous as a function ${\mathbb{R}^{n}\rightarrow\Pi}$, continuity only fails at the boundary, and $p$ behaves like a continuous function when composed with $h$:

Algorithmically, an orbit graph can be constructed as follows: Constructing an $\varepsilon$-net is a standard problem in computational geometry and can be solved efficiently (e.g., Haussler and Welzl [33]). Having done so, the problem we have to solve is:

Since $\Pi$ is a polytope, its diameter

can also be evaluated computationally. By definition of $d_{\mathbb{G}}$, we have

That shows the minimum is always attained for a point $\phi y$ on a tile $\phi\Pi$ that lies within distance $\text{diam}(\Pi)$ of $x$. The set of transformations that specify these tiles is

This set is always finite, since $\mathbb{G}$ is discrete and the ball of radius $\text{diam}(\Pi)$ is compact. We can hence evaluate the quotient metric as

which reduces the construction of $E$ to a finite search problem. In summary:

The construction is illustrated in Figure 3.

Recall that the faces of a polytope are its vertices, edges, and so forth; the facets are the ${(n-1)}$-dimensional faces. The polytope itself is also a face, of dimension $n$. See [70] for a precise definition. Given $\Pi$ and $\mathbb{G}$, we will call two faces $S$ and $S^{\prime}$ $\mathbb{G}$-equivalent if ${S^{\prime}=\phi S}$ for some ${\phi\in\mathbb{G}}$. Thus, if ${S=\Pi}$, its equivalence class is ${\{\Pi\}}$. If $S$ is a facet, it is equivalent to at most one distinct facet, so its equivalence class has one or two elements. The equivalence classes of lower-dimensional faces may be larger—if $\mathbb{G}$ is p1 and $\Pi$ a square, for example, all four vertices of $\Pi$ are $\mathbb{G}$-equivalent.

Since $\mathbb{G}$ is crystallographic, it contains shifts in $n$ linearly independent directions, and these shifts hence specify a coordinate system of $\mathbb{R}^{n}$. More precisely: There are $n$ elements ${\phi_{1},\ldots,\phi_{n}}$ of $\mathbb{G}$ that (1) are pure shifts (satisfy ${A_{\phi_{i}}=\mathbb{I}}$), (2) are linearly independent, and (3) are the shortest such elements (in terms of the Euclidean norm of $b_{\phi_{i}}$). Up to a sign, each of these elements is uniquely determined. We refer to the vectors ${\phi_{1},\ldots,\phi_{n}}$ as the shift coordinate system of $\mathbb{G}$.

For any open set ${M\subseteq\mathbb{R}^{n}}$, we define the Laplace operator on twice differentiable functions ${h:M\rightarrow\mathbb{R}}$ as

Now consider specifically functions ${e:\mathbb{R}^{n}\rightarrow\mathbb{R}}$. Fix some ${\lambda\in\mathbb{R}}$, and consider the constrained partial differential equation

Clearly, there is always a trivial solution, namely the constant function ${e=0}$. If (13) has a non-trivial solution $e$, we call this $e$ a $\mathbb{G}$-eigenfunction and $\lambda$ a $\mathbb{G}$-eigenvalue of the linear operator $-\Delta$. Denote the set of solutions by

Since $0$ is a solution, and any linear combination of solutions is again a solution, ${\mathcal{V}(\lambda)}$ is a vector space, called the eigenspace of $\lambda$. Its dimension

is the multiplicity of $\lambda$.

The standard Fourier bases for periodic functions on $\mathbb{R}^{n}$ can be obtained as the special cases of Theorem 7 for shift groups: Fix some edge width ${c>0}$, and choose $\Pi$ and $\mathbb{G}$ as

For these groups, all eigenvalue multiplicities are ${k(\lambda_{i})=2n}$ for each ${i\in\mathbb{N}}$. For ${n=2}$, the group $\mathbb{G}$ is p1 (see Figure 1). Its eigenfunctions are shown in Figure 4.

To clarify the relationship in more detail, consider the case ${n=1}$: Since $\Delta$ is a second derivative, the functions ${e(x)=\cos(\nu x)}$ and ${e(x)=\sin(\nu x)}$ satisfy

and are hence eigenfunctions of ${-\Delta}$ with eigenvalue ${\lambda=\nu^{2}}$. For this choice of $\Pi$ and $\mathbb{G}$, the invariance constraint in (13) holds iff ${e(x)=e(x+c)}$ for every ${x\in\mathbb{R}}$. That is true iff

The eigenspaces are therefore the two-dimensional vector spaces

Any continuous function $f$ that is $\mathbb{G}$-invariant (or, equivalently, $c$-periodic) can be expanded as

In the notation of Theorem 7, the coefficients are ${c_{2i}=a_{i}}$ and ${c_{2i+1}=b_{i}}$, and

Note that the unconstrained equation has solutions for all $\lambda$ in the uncountable set ${[0,\infty)}$. The invariance constraint limits possible values to the countable set ${\lambda_{1},\lambda_{2},\ldots}$. If $f$ was continuous but not invariant, the expansion (14) would hence require an integral on the right. Since $f$ is invariant, a series suffices.

The eigenfunctions in Theorem 7 can be approximated by eigenvectors of a suitable graph Laplacian of the orbit graph as follows. We first compute an orbit graph ${\mathcal{G}=(\Gamma,E)}$ as described in Section 4. We weight each edge $(x,y)$ of the graph by

The normalized Laplacian of the weighted graph is

and ${D}$ is the diagonal matrix containing the sum of each row of $W$. See e.g., Chung [20] for more on the matrix $L$. Our estimates of the eigenvalues and -functions of $\Delta$ are the eigenvalues and eigenvectors of $L$,

These approximate the spectrum of $\Delta$ in the sense that

see Singer [56]. Once an eigenvector $\hat{e}_{i}$ is computed, values of $e_{i}$ at points ${x\not\in\Gamma}$ can be estimated using standard interpolation methods.

Alternatively, the basis can be computed using a Galerkin approach, which is described in Section 9.5. The functions in Figures 4, 5 and 6 are computed using the Galerkin method.