Minimum-strain symmetrization of Bravais lattices

Bravais lattices are the most fundamental building blocks of crystallography. They are classiﬁed into groups according to their translational, rotational, and inversion symmetries. In computational analysis of Bravais lattices, fulﬁllment of symmetry conditions is usually determined by analysis of the metric tensor, using either a numerical tolerance to produce a binary (i.e., yes or no) classiﬁcation or a distance function which quantiﬁes the deviation from an ideal lattice type. The metric tensor, though, is not scale invariant, which complicates the choice of threshold and the interpretation of the distance function. Here, we quantify the distance of a lattice from a target Bravais class using strain. For an arbitrary lattice, we ﬁnd the minimum-strain transformation needed to fulﬁll the symmetry conditions of a desired Bravais lattice type; the norm of the strain tensor is used to quantify the degree of symmetry breaking. The resulting distance is invariant to scale and rotation, and is a physically intuitive quantity. By symmetrizing to all Bravais classes, each lattice can be placed in a 14-dimensional space, which we use to create a map of the space of Bravais lattices and the transformation paths between them. A software implementation is available online under a permissive license.


I. INTRODUCTION
Many properties of a crystalline material are governed by the geometry of its Bravais lattice and changes in the Bravais lattice, both static and dynamic.For example, modification of the Bravais geometry by the imposition of hydrostatic and/or shear strains forms the basis of elastic strain engineering [1], a field the successes of which include improved photoluminescence and electronic spectra in semiconductors [2][3][4], vibrational properties in micromechanical oscillators [5], and magnetic properties in multiferroics [6].Alternatively, the design of materials which can sustain repeated cyclic changes in the Bravais lattice geometry is the fundamental aim of shape memory alloy research [7][8][9].
In addition to the practical consequences for materials engineering, there is a fundamental theoretical interest in the classification of Bravais lattices.This has motivated the development of many algorithms for lattice analysis.A traditional lattice analysis method is straightforward to describe.
(1) Find a canonical description of the lattice.
(2) Classify the lattice according to the symmetries of its cell parameters.
The most commonly used canonical description is the Niggli-reduced [10,11] form of the lattice, which is defined using the metric tensor.For a unit cell with lattice vectors a, * peter.mahler.larsen@gmail.comPublished by the American Physical Society under the terms of the Creative Commons Attribution 4.0 International license.Further distribution of this work must maintain attribution to the author(s) and the published article's title, journal citation, and DOI.b, and c, the metric tensor is given by G = ⎡ ⎣ g 11 g 12 g 13 g 12 g 22 g 23 g 13 g 23 g 33 Niggli reduction specifies a set of constraints which the metric tensor must satisfy.After reduction, the Bravais class of a lattice can be determined by inspecting the elements of its metric tensor.For example, the Niggli-reduced form of the primitive cubic lattice is the same as the conventional setting, where the metric tensor satisfies g 11 = g 22 = g 33 , (2a) Similar conditions can be specified for each of the Bravais lattice types.While constraints of this form provide a clean theoretical description of Bravais classification, practical classification is complicated by the presence of noise.Lattice parameters which have been determined by experiment are subject to measurement errors, and even lattice parameters obtained from computer simulations are subject to numerical errors resulting from the use of floating-point arithmetic.Both of these sources of error mean that symmetry conditions are rarely fulfilled exactly.
Noise is treated in one of two ways, using either a tolerance parameter or a distance calculation.A tolerance parameter, , is used to permit approximate, rather than exact, fulfillment of symmetry conditions.In this case, the cubic symmetry conditions in Eq. ( 2) are relaxed: The tolerance parameter effectively defines a boundary between symmetry classes.Exactly where this boundary should lie is decided by the user, but involves a degree of arbitrariness: given a sufficiently large number of structures, there will be lattices on either side of any boundary, the classification of which is dependent on the exact choice of threshold.This use of a tolerance parameter is the approach used by Grosse-Kunstleve et al. [12] and in the lattice determination component of the Spglib [13] space-group classification library.
To overcome the need for an a priori selection of a tolerance parameter, a distance can be calculated, which effectively measures the degree of constraint violation: By measuring the distances from each Bravais class, the choice of distance threshold can be made a posteriori, and informed by the available options.A variety of distance functions on the metric tensor have been used in the literature [14][15][16][17].
The commonality of existing methods is the use of the metric tensor.While the metric tensor is a rotationally invariant description of a lattice basis, it is not invariant to scale, which complicates the choice of threshold and the interpretation of the distance function.
In this paper we determine the distance of a lattice from a chosen Bravais type using strain.The problem we aim to solve is the following: given an observed lattice which we wish to symmetrize, find the closest lattice in the target Bravais class.The distance is determined by the minimum strain needed to elastically deform a lattice such that the symmetry conditions of the target Bravais class are satisfied.Strain has the advantage of being a rotationally invariant and physically intuitive quantity, and can be easily made scale invariant.Furthermore, since the transformation between any two lattices is a linear map, strain also represents the most natural distance measure on lattices.In addition to lattice classification, the other immediate application of this paper is simply to symmetrize a lattice in a well-defined manner.
The rest of this paper is structured as follows: we define a similarity measure for fixed Bravais lattices using the strain tensor in Sec.II.We extend this to variable lattices in Sec.III, which enables us to find minimum-strain symmetrizations.The symmetrization method is illustrated in Sec.IV, and we use it to present a map of the Bravais lattices in Sec.V. Concluding remarks are given in Sec.VI.

II. QUANTIFICATION OF LATTICE SIMILARITY
The symmetrization procedure we develop in this paper is based on quantification of the deformation from one lattice to another.In this section we describe the necessary deformation theory, which is a cornerstone of continuum mechanics.

A. Lattice basis comparison
A Bravais lattice, A , is described by a lattice basis, A ∈ R 3×3 , which consists of three lattice vectors, described by the columns of A. The lattice points lie at all integer combinations of the lattice basis vectors, that is, a lattice point p satisfies p = A h where h ∈ Z 3 .Here we wish to compare two lattices geometrically.
Let A and B be Bravais lattices, with lattice bases A ∈ R 3×3 and B ∈ R 3×3 .We will compare the bases by quantifying the minimum deformation necessary to map one basis onto the other.There is a linear map, F ∈ R 3×3 , also called the deformation gradient, which exactly transforms B into A: A = FB. ( We can remove the rotational dependence of F with the use of a polar decomposition: Here U is a symmetric matrix containing pure stretches only, known as the right stretch tensor.The Doyle-Ericksen strain tensors [18] (or sometimes, the Seth-Hill strain tensors [19,20]) are a generalized class of strain tensors of the form where m typically takes an integer value.In this paper we choose m = 1, in which case E is the Biot strain tensor.Combining Eqs. ( 5)-( 7), we define a distance function which quantifies the strain tensor norm: where X F = √ i, j X 2 i j is the Frobenius norm.The distance function is invariant to orthogonal transformations of A and B (trivially for A, and by the spectral theorem for B).To provide some intuition for the distance function, we can express it in terms of the principal stretches, ν i , of F: The principal stretches are equal to the eigenvalues of U, and the distance is equal to the norm of the principal strains.Alternatively, we can use the relationship between the polar decomposition and the singular value decomposition to express the distance function in terms of the singular values, σ i , of F: The principal stretches and singular values are in fact identical and differ only in their origin and interpretation; the former formulation is from continuum mechanics and has a physical meaning, while the latter is a more general formulation from linear algebra which we use in this paper.The distance function shown here is described in further detail by Koumatos and Muehlemann [21], who also discuss different choices of m.

B. Lattice correspondences
The similarity function, d, compares two lattice bases.However, every lattice is generated by an infinite number of bases.Any two bases, A and A , which generate the same FIG. 1.The target lattice, A (shown on the left), can be described by an infinite number of bases, which differ by a unimodular matrix, L. Three such lattice bases, AL 1 , AL 2 , and AL 3 , are shown in orange, green, and red.A lattice, B , is shown in the center, the distance of which from A we wish to determine.A single basis, B, is shown.The correspondence AL 1 = FB has the lowest associated strain.Shown on the right is the best rigid match (prior to the stretch operation).
lattice A are related by A = AL, where L ∈ SL ± 3 (Z) is a unimodular matrix, that is, L ∈ Z 3×3 where det (L) = ±1.We call L the correspondence matrix.
In order to compare two lattices, rather than specific lattice bases, we optimize the function: which compares all possible bases of the first lattice against a single basis of the second, and in doing so finds the optimal lattice correspondence.This process is illustrated in Fig. 1.For any two lattices A and B , it is easily shown that the minimum distance shown in Eq. ( 11) is invariant to the choice of lattice bases.An efficient algorithm for determining the optimal lattice correspondence is described by Chen et al. [22].Their algorithm proceeds by bounding the maximum norm of the optimal correspondence matrix, and testing all correspondences that lie within that sphere.We refer to their paper for further details, as it applies to correspondences between fixed lattices only.

III. QUANTIFICATION OF SYMMETRY BREAKING
In this section we extend the concept of the distance function for fixed lattices to variable lattices, and demonstrate its application to the symmetrization of Bravais lattices.After describing a function for quantification of symmetry breaking, we analyze some properties of the function, and provide some visual examples.

A. Variable lattice bases
The distance function in Eq. ( 8) compares two lattices with fixed parameters.In order to quantify symmetry breaking, we introduce lattices with variable lattice parameters.As stated above, the problem we aim to solve is the following: given an observed lattice which we wish to symmetrize, find the closest lattice in the target Bravais class.We will first consider this problem for a fixed correspondence.
Let B ∈ R 3×3 be a lattice basis which we wish to symmetrize, and let Z ∈ R 3×3 be a variable lattice basis.The TABLE I. Semidefinite program to calculate the minimum-strain symmetrization of a lattice basis, B. Here, Z is a lattice basis matrix the geometry of which respects the symmetry conditions of the target Bravais type.The symmetry conditions are enforced using additional linear constraints on the Gramian matrix, G.
Subject to: symmetrized basis is found by solving the problem Here, d B denotes the distance function in Eq. ( 8) subject to geometric constraints on Z which enforce the symmetry conditions of a chosen Bravais lattice type (B) (described below).We can rearrange the distance function in terms of the Frobenius and nuclear norms of Let G = Z T Z be the Gramian [23] matrix of Z.Then, Using the Gramian form, the minimization problem in Eq. ( 12) can be expressed as a semidefinite program (SDP).Semidefinite programming is a form of convex optimization which generalizes linear programming [24]; in addition to a linear objective and linear constraints, a SDP permits positive semidefiniteness constraints on matrix variables.
Table I shows a SDP for solving Eq. ( 12).The input parameter is the lattice basis we wish to symmetrize (P1).The matrix variables are the variable lattice basis (P2) and the Gramian matrix (P3).The objective is to minimize the square of the distance function (P4).By the Schur complement condition for positive semidefiniteness [25,26], the linear matrix inequality in (P5) enforces the condition G − Z T Z 0 (where the relation X 0 means X is positive semidefinite).At the optimum value of Z this inequality is tight, i.e., G = Z T Z.In addition to the constraint (P5), we impose a set of linear constraints on G in order to control the geometry of Z and thereby enforce the symmetry conditions of a chosen Bravais lattice type.The constraints for each Bravais type are shown in Fig. 2.
The SDP shown above finds a symmetrized basis, Z * , by simultaneously finding the lattice parameters and the lattice rotation.This solution is the closest lattice basis to B under the defined distance function.The basis respects the specified symmetry constraints, and represents the minimumstrain symmetrization of the basis B. The resulting distance quantifies the symmetry breaking of the target Bravais type.

B. Lattice correspondence search
By formulating the distance function as a SDP, we have shown that the distance function is convex for a fixed lattice correspondence.In order to find the minimum-strain solution, however, we need to optimize over all lattice correspondences, of which there are an infinite number.For fixed lattices, the optimal correspondence can be found efficiently, owing to the existence of an upper bound on the norm of the optimal correspondence matrix [22].For variable lattices no such bound is known.Instead, we use a steepest-descent search over lattice correspondences.
Starting from a good candidate correspondence, we search the surrounding correspondences, iteratively, until no better solution can be found.For the initial "good candidate" correspondence, we use a right-handed Minkowski-reduced [27] (MR) basis of the B .A MR basis is a lattice basis the lattice vectors of which are the shortest possible.It can be computed efficiently [28].
We define the neighborhood of a correspondence as all unimodular matrices with positive determinant and the elements of which have magnitude at most 1: The steepest-descent approach is shown in the algorithm in Table II.The inputs to the algorithm are a MR basis (B) of the lattice we wish to symmetrize, the target Bravais type (B), and the correspondence neighborhood (L).The algorithm keeps track of the minimum distance (line 2) and the optimal correspondence (line 3), and maintains a set of visited correspondences (line 4), i.e., correspondences which have been tested.In each iteration (line 6) the search tests all correspondences in the neighborhood of the current best solution (line 9).Each time a correspondence is tested (line 12), the best solution (lines 13 and 14) is updated as needed.
Here d B denotes the distance function with the appropriate linear constraints on G for a chosen Bravais type, B. The algorithm iteratively tests the neighborhood of the current best solution until no better solution can be found.
The neighborhood we define in Eq. ( 14) is motivated by the need for a distance function which is continuous in the presence of a continuous lattice deformation.The distance function is continuous if we consider, at a minimum, all lattice bases consisting of Voronoi-relevant vectors [29].Voronoirelevant vectors are illustrated in Fig. 3.A theorem due to Minkowski [30] states that the Voronoi-relevant vectors are a nonstrict subset of all integer combinations of {−1, 0, 1} of the MR basis vectors.The neighborhood defined in Eq. ( 14) contains all right-handed unimodular matrices of this form.
With the lattice correspondence search in place, we can perform minimum-strain symmetrization.This is illustrated for a selection of lattice types in Fig. 4. Of the lattices shown, only the monoclinic symmetrization has a small strain tensor norm.Nonetheless, even for very dissimilar lattice types with a correspondingly large strain tensor, we can still compute a minimum-strain symmetrization.

C. Scale invariance and range of the symmetrization distance
We have described minimum-strain symmetrization for the different Bravais lattice types.Here, we show that the distance function is scale invariant, and that the scale invariance sets an upper bound on the range of the distance function.
When comparing two fixed lattices, d (A, B) is dependent on the scales of A and B. However, when using variable lattices to quantify symmetry breaking, the resulting distance is scale invariant.For a scaling factor, k ∈ R, we have argmin Due to the lattice parameters of Z being optimized, any change in scale of B is accommodated by a corresponding change in Z.This is a particularly useful property when comparing symmetry breaking in cells of different sizes.
The scale invariance also sets an upper bound on the maximum strain.For a chosen Bravais type and a correspondence matrix L, the range of the distance function ).The zero-strain solution is trivial.To show that the interval is bounded from above by √ 2, we observe that even the most constrained Bravais types (the cubic lattices) accommodate a change in scale of Z. Representing this by a parameter s ∈ R, the distance function (minimized over all lattice correspondences) is given by min Let σ = [σ 1 , σ 2 , σ 3 ] be the singular values of ALB −1 .It can be shown that the optimal scaling parameter is given by Expressed in terms of singular values, the distance function is which attains a maximum (at √ 2) when σ = [σ 1 , 0, 0].The rank of a matrix is equal to the number of nonzero singular values.Let B −1 be a rank-1 matrix.Then ALB −1 is a rank-1 matrix for every L ∈ SL 3 (Z) and any matrix A, and σ = [σ 1 , 0, 0].Since a rank-1 matrix is not invertible, this vector of singular values can only be attained asymptotically.
The analysis shown here establishes the upper bound of the distance function.A similar analysis for Bravais types with more degrees of freedom is significantly more complicated.The additional degrees of freedom, however, can only serve to reduce the range of the distance function.

D. Implementation
We have presented the symmetrization procedure as a semidefinite program.In addition to being the most compact description, formulation as a semidefinite program is the most convenient for analysis of the scale invariance and range of the distance function.For a practical implementation, however, it is less than ideal to rely upon semidefinite programming, as SDP solvers are slow.Fortunately, any optimization approach which achieves the same result as the SDP is equally valid.
By employing a quaternion parametrization, symmetrization can be done by solution of a constrained multivariate polynomial, which can be solved efficiently using sequential quadratic programming [31].Further details are given in the Appendix.Compared to the existing Bravais classification methods this procedure has significantly higher computational requirements, taking on the order of 0.3 s per cell.Despite the increased computation time, the relevant comparison is not between classification methods, but to the time required to perform an experiment or an ab initio structure calculation.In this regard, the increased computation time is of little importance.

IV. SYMMETRIZATION DISTANCE VECTORS
The procedure we have described above finds a minimumstrain symmetrization to a target Bravais type, and an associated symmetrization distance.By calculating the symmetrization distance from every Bravais type we obtain a vector of length 14, which we will call a "distance vector."We can develop some intuition for the distance vectors by observing how they change in the two-parameter Bravais lattices.Figure 5 shows the distance vectors for the tetragonal Bravais lattices, along a path containing all possible permitted geometries permitted by their respective symmetry conditions.
The geometry of the primitive tetragonal lattice is the simpler of the two.At φ = 0, the lattice is degenerate, with lattice parameters c = 1 and a = 0; we denote this lattice the "tall square," since, in the limit φ → 0, the base is square.As mentioned above, symmetrization cannot be performed for degenerate lattices, so this degenerate state is achieved only asymptotically.At φ = π/4, c = a and the lattice is identical to the primitive cubic lattice, which is reflected by the associated zero symmetrization distance.At φ = π/2, the lattice is again degenerate, with lattice parameters c = 0 and a = 1; we denote this lattice the "flat square."Here, the symmetrization distance from the cubic Bravais types approaches the maximum value of √ 2, as discussed in Sec.III C. The path of the body-centered tetragonal lattice is similar to that of the primitive tetragonal lattice, but rather than passing through the primitive cubic lattice it contains the Bain transformation [32,33], from the fcc to the bcc lattice The distances from the primitive triclinic lattice type are zero everywhere, since all lattices are trivially triclinic.Similarly, it can be seen that the symmetrization distance for rhombohedral symmetrization is never greater than the fcc distance.Both lattice types have lattice vectors of equal length, but the rhombohedral lattice has a variable angle.The extra degree of freedom serves to reduce the strain.This illustrates the important concept that Bravais lattices exist in a hierarchy, whereby the symmetries of some Bravais types are subsets of other types.The subset relationships are illustrated in Fig. 6.In general we can state the following: for two lattice types, B and B , if the symmetries of B are a subset of those FIG. 5. Symmetrization distances for all geometries of the primitive and body-centered tetragonal lattices (modulo scale).Since the symmetrization distance is scale invariant, these lattices can be parametrized using a single parameter: using the lattice parameters of the conventional setting a = sin φ and c = cos φ.Zero-distance symmetrizations are not shown.In the primitive tetragonal figure, the body-centered orthorhombic symmetrization distances are identical to the body-centered tetragonal distances and are not shown.Similarly, in the body-centered tetragonal figure, the primitive monoclinic symmetrization distances are identical to the base-centered monoclinic distances, and the primitive orthorhombic symmetrization distances are identical to the primitive tetragonal distances. of B , then min This fact should be considered when classifying a lattice.After imposing a threshold on the symmetrization distance, multiple potential Bravais types can be assigned.In this case, the Bravais type with the highest number of symmetries should be selected.
FIG. 6. Hierarchy of symmetry relationships between Bravais types.The Bravais types are arranged vertically according to the number of degrees of freedom in the lattice.Presence of an upwards path indicates a subset relationship.For clarity, Bravais types are denoted by their Pearson symbols: a, triclinic; m, monoclinic; o, orthorhombic; t, tetragonal; h, hexagonal; c, cubic; P, primitive; C, base-centered; R, rhombohedral; I, body-centered; F, face-centered.

V. MAP OF THE BRAVAIS LATTICES
As shown above, by calculating the symmetrization distance from each Bravais type, we can assign every lattice a 14-dimensional distance vector.It is instructive to consider the space of possible distance vectors.To do so we use truncated principal component analysis (PCA). Figure 7 shows a two-dimensional PCA projection of distance vectors.Two projections have been used to illustrate the distance vector space.The PCA projection in Fig. 7(a) is calculated using the randomly sampled Bravais lattices, including many degenerate lattices.While this projection is useful for studying the complete space of Bravais geometries, it wastes a lot of space on physically uninteresting lattices with degenerate geometries.The PCA projection in Figs.7(b)-7(h) is calculated using structures from the Crystallography Open Database [34] (COD).Since only experimentally observed structures are used (as filtered in earlier work by some of the authors [35]), the region of interest for practical applications occupies a significantly larger area.
The extreme vertices in Fig. 7(a) consist of the three cubic Bravais types and six degenerate Bravais lattices.Within each of the cubic Bravais types, lattices differ by a scale factor only.Since the symmetrization distance is invariant to scale, all lattices within a cubic Bravais type have the same distance vector.Four of the other six extreme vertices are degenerate states of the primitive tetragonal and primitive hexagonal Bravais lattices: the tall square and flat square lattices (described above), and the "tall hexagonal" and "flat hexagonal" lattices, which are similar to the square variants but with hexagonal bases.The two remaining extreme vertices, the "flat rectangular" and "tall rectangular" lattices, are degenerate states of the primitive orthorhombic lattice, the cell vectors of which have sufficiently different lengths that two degrees of freedom are insufficient to accommodate a low-strain symmetrization.
The tetragonal and hexagonal Bravais types have two degrees of freedom.As shown above, the scale invariance of the symmetrization distance means that the distance vectors of these Bravais types exist on a line.The Bain transformation (from fcc to bcc), which can be achieved by tracing the body-FIG.8. Fraction of variance explained by each PCA component (dark blue), and the cumulative sum of the explained variance (light blue).The projection matrix is calculated using experimental structures taken from the Crystallography Open Database.centered tetragonal transition line, is widely studied due to its energetic accessibility in iron-based materials.It is noteworthy that the small strain required for this transformation is reflected by the proximity of the fcc and bcc distance vectors.
The orthorhombic lattices have three degrees of freedom.Distance vectors of these lattice types exist on a twodimensional manifold.The monoclinic lattices have four degrees of freedom and exist on a three-dimensional manifold.The two-dimensional PCA projection captures the symmetry subset relationships shown in Fig. 6, but also introduces a degree of spurious overlap.The overlap is worse for lattices with more degrees of freedom.For example, the path traced by the body-centered tetragonal lattice appears to self-intersect at multiple points.These crossings are resolved in a threedimensional projection.Similarly, the spaces occupied by the primitive and base-centered monoclinic lattices appear to overlap, despite being disjoint.
To explain the overlap, we observe that five parameters are sufficient to describe a Bravais lattice if scale is ignored, and that the space of distance vectors is therefore a nonlinear five-dimensional manifold embedded in R 14 .It is therefore inevitable that a two-dimensional projection contains spurious overlap.The loss of information resulting from PCA projection is shown in Fig. 8. PCA is a linear projection and therefore requires more than five dimensions to capture the full distance vector space.Nonetheless, even a two-dimensional projection explains most of the variance and is useful for visualizing the space of lattices, and particularly for studying phase transitions.
Although the hcp crystal structure is not a Bravais lattice, we have marked it on the projection for reference.In structure maps based on energetic descriptors [36], the hcp and fcc structures are typically close together, due to their small energetic differences.Geometrically, however, the primitive hexagonal and fcc lattices are far apart, since a large strain is required to transform one lattice into the other.The difference between an energetic map and our map can be described in terms of transformation paths: an energetic map encodes differences in the energetic end points of the transition, whereas our map encodes the lattice deformation of the transition itself.

VI. CONCLUSION
We have described a method for symmetrization of Bravais lattices using strain, which is a physically intuitive quantity that also has the attractive features of rotation and scale invariance.A distance from each Bravais type is determined by quantifying the strain necessary for symmetrization.This allows classification to be performed without the need for an a priori selection of a tolerance parameter.
By projecting the distance vectors using PCA, we obtain insight into the positions of lattices in an abstract Bravais space, where distances are determined by displacive deformations only.The symmetrization procedure finds minimumstrain solutions, rather than minimum-energy solutions, and considers the Bravais lattice only.Indeed, a correct analysis of a material's lattice type requires that a primitive unit cell can be identified.Factors which complicate this, such as perturbations of the atomic basis coordinates within a unit cell, are not accounted for when considering the Bravais lattice only.Nonetheless, we envisage that the method will be useful for a range of experimental and computational applications.Ongoing work aims to combine the symmetrization distance presented here with an atomic basis to create a distance on space groups.
A software implementation (C++ with PYTHON wrappers) is available online [37].

ACKNOWLEDGMENT
This work was supported by Grant No. 7026-00126B from the Danish Council for Independent Research.

APPENDIX: NUMERICAL SOLUTION
Here we rework the distance expression in Eq. (P4) into a multivariate polynomial.We perform a substitution where Q ∈ SO( 3) is a right-handed orthogonal matrix, and {T i | i ∈ 1 . . .n} is a set of n matrices (a "template") which, in conjunction with a vector x ∈ R n , parametrizes a Bravais lattice type.The template maintains the desired symmetry by construction.For example, a primitive orthorhombic basis has a template In order to express the distance function as a multivariate polynomial we parametrize the orthogonal matrix using maxima of this equation are symmetrically equivalent and have equal solution values (see Fig. 9).
Equation (A10) can be solved using stepwise iteration [41].This is guaranteed to converge to a local maximum (rather than a saddle point) but only at a linear convergence rate.On the other hand, by solving the Lagrangian form of Eq. (A10) we can use Newton-Raphson iteration to achieve quadratic convergence.We have found that a good compromise is to perform stepwise iteration for 10-20 iterations (to get close to the maximum), and then switch to Newton-Raphson (to get fast convergence).
The total number of iterations required is dependent on the number of template variables, but even for monoclinic templates convergence is typically achieved within 20 iterations.For cubic templates the problem reduces to a polar decomposition, and convergence is achieved in a single stepwise iteration.A larger number of correspondences must be investigated for highly skewed cells, but for cells of typical crystal structures the running time does not exhibit large deviations.Calculation of the symmetrization vectors for 22 000 structures in the Crystallography Open Database takes ≈2 h running on a single core of a 2014 Macbook Pro.On a modern HPC server node this is reduced to a few minutes.

FIG. 2 .
FIG.2.Linear constraints on the Gramian matrices for each of the 14 three-dimensional Bravais lattices.The constraints are constructed such that the symmetry conditions of the Bravais lattice are respected.The triclinic lattice type does not have any constraints on the form of the Gramian; all lattices are trivially triclinic, since there are no symmetry conditions to satisfy.

FIG. 3 .
FIG.3.Voronoi cells (Wigner-Seitz cells) in a 2D oblique lattice.Voronoi-relevant vectors are those which pass through a Voronoi face to an adjacent lattice point; they constitute a subset of the red lattice points, which are combinations of {−1, 0, 1} of the MR basis vectors.

FIG. 4 .
FIG. 4. Successively aggressive lattice symmetrizations.A triclinic cell (blue) is symmetrized to each of the six primitive Bravais lattice types.The symmetrized cells (orange) can be mapped onto the triclinic cell with a pure stretch.The lattice parameters and rotations of the symmetrized cells represent minimum-strain solutions.

FIG. 7 .
FIG.7.Two-dimensional maps of the space of Bravais lattices, obtained using PCA projections of the 14-dimensional distance vector space.Theoretical, randomly sampled structures are shown in gray.Experimentally observed lattices from the Crystallography Open Database (COD) are marked with blue points.Two projections are used to illustrate the data: in (a) the projection matrix is calculated using the randomly sampled structures; in (b) the projection matrix is calculated using the COD structures.Regions occupied by each of the monoclinic and orthorhombic lattice types in the COD projection are shown in (c)-(h).
templates can be constructed for the other Bravais lattice types.

TABLE II .
A steepest-descent algorithm for finding the optimal lattice correspondence.