Reverse Monte Carlo with geometric analysis – RMC+GA

Wells, S.; Dove, M.; Tucker, M.

doi:10.1107/S0021889804008957

research papers

JOURNAL OF
APPLIED
CRYSTALLOGRAPHY

ISSN: 1600-5767

Volume 37| Part 4| August 2004| Pages 536-544

doi:10.1107/S0021889804008957

Reverse Monte Carlo with geometric analysis – RMC+GA

Stephen Wells,^a ^* Martin Dove ^b and Matthew Tucker ^b

^aRoyal Institution, Albemarle Street, London, UK, and ^bDepartment of Earth Sciences, Downing Street, Cambridge, UK
^*Correspondence e-mail: stephen@ri.ac.uk

(Received 11 March 2004; accepted 14 April 2004)

An analysis of simulated framework structures based on the rigid-unit model quantifies the distortion of polyhedra from their ideal geometric forms and decomposes the motions of the structure into components of rigid-unit displacement, rigid-unit rotation and distortion. Case studies analysing the behaviour of quartz and of other silicates demonstrate that the method provides a novel way of extracting information from reverse Monte Carlo simulations. A program called GASP has been developed to perform this analysis and is freely available to researchers. Rotations are handled using the rotor method of geometric algebra.

Keywords: simulated framework structures; flexible framework structures; rigid-unit model; reverse Monte Carlo simulations; geometric analysis; computer programs.

1. Introduction

While crystal structures might be expected to be highly rigid, framework structures such as aluminosilicates in fact display considerable internal flexibility via cooperative motions of polyhedral groups of atoms. This fact has important consequences for the modelling of framework structures. We present a methodology, namely reverse Monte Carlo with geometric analysis, that quantifies the significance of such flexibility in the dynamics of framework structures.

The question of the rigidity or flexibility of framework structures in three dimensions remains poorly understood. Constraint counting as used by Maxwell fails to capture the flexibility of structures in which symmetry makes certain constraints redundant. This shortcoming is particularly relevant to crystal framework structures such as silicates, which possess many degrees of freedom and in which the SiO₄ tetrahedra rotate and translate without distortion (Dove et al., 2000 ). These rigid-unit modes (RUMs) can manifest themselves as soft modes in displacive phase transitions (Giddy et al., 1993 ; Hammonds et al., 1996 ), as sources of negative thermal expansion (Heine et al., 1999 ; Welche et al., 1998 ) and as mechanisms for large-amplitude thermal motion of vertex atoms at low energy cost (Gambhir et al., 1999 ).

RUMs can be identified in reciprocal space using the `split-atom' model (Hammonds et al., 1994 ). While this approach can successfully identify candidate soft modes for phase transitions, and quantify to a degree the rigid-unit flexibility of a structure, the model does not capture the real-space motion of the structure. Firstly, this approach operates in the limit of infinitesimal amplitude; secondly, rotations are inherently anharmonic, as the motion of atoms is not linear in the rotation angle; and thirdly, rotations are not commutative. A motion of atoms due to a large number of active RUMs will therefore differ from a superposition of RUM eigenvectors.

Reverse Monte Carlo (RMC) modelling allows us to simulate configurations of atoms based directly on experimental data (McGreevy & Pusztai, 1988 ; McGreevy, 1995 ). The RMC technique used to generate the configurations in this paper was modified to account explicitly for both long-range and short-range order (Tucker et al., 2000 ; Tucker, Dove & Keen, 2001 ; Tucker, Keen & Dove, 2001 ). Our configurations are consistent with the total scattering intensity, the pair distribution factors and the intensities of the Bragg peaks.

We have been working for some time on a method to extract information on RUM motion in real space from such techniques (Wells, Dove & Tucker, 2002 ; Wells, Dove, Tucker & Trachenko, 2002 ; Wells, 2003 ) in order to complement the reciprocal-space treatment of rigidity described by Hammonds et al. (1994). The principal problem is to model the behaviour of the polyhedra formed by real atoms in such a way as to relate it to the behaviour of the ideal polyhedra that are the objects of the rigid-unit model. Our analysis now not only quantifies the total atomic motion but also allows for both rotational and displacive motion of the polyhedra, whereas the earlier approach of Wells, Dove & Tucker (2002) and Wells, Dove, Tucker & Trachenko (2002) was restricted to rotational motion only.

We use the language of geometric algebra (GA; Lasenby et al., 2000 ), a form of Clifford algebra, in which the basis vectors of three-dimensional space and Hamilton's quaternions exist within the same algebra. It is thus particularly well suited to the description of rotations. None of our results depends in detail on the form of the rotation operator used (it is simply computationally convenient for us to use the GA formalism) and so we will only quote those expressions necessary for the interpretation of our results.

The combination of RMC modelling with geometric analysis represents a new methodology, RMC+GA, for obtaining information on the behaviour of framework structures. This geometric analysis can also be applied to configurations generated by other methods, such as molecular dynamics (MD).

2. Finding best-fit polyhedral motions

Suppose that we have two slightly different forms of the same structure, in which each polyhedron has undergone some displacement, rotation and distortion due to the motion of the atoms but there has been no change in topology, i.e. no bonds have been broken or have formed. Specifically, let us take two independently generated RMC fits to the same data; since RMC includes a stochastic element, these fits will differ in detail as to the positions of the atoms. As argued by Wells, Dove, Tucker & Trachenko (2002), we take each fit as an instantaneous snapshot of the system in dynamic disorder. Two independent (uncorrelated) fits can be taken as two snapshots separated by an infinite time interval.

We could also compare two different frames from a molecular dynamics simulation, either successive frames separated by a single time step or two frames with a wider temporal separation.

By comparison of the two structures, polyhedron by polyhedron, we can identify the displacements and rotations of the polyhedra, and the remaining distortion, namely the difference in shape between the two forms of the same polyhedron. We proceed by defining `mismatch scores', representing the differences between the two structures, and then minimizing these scores as a function of the rigid-unit motion of the polyhedra. We can also compare the polyhedra in a simulated structure with those in an idealized structure with perfectly regular polyhedra, so as to quantify the distortion of the polyhedra or to identify an RUM acting as a soft mode for a phase transition.

The steps in the comparison process are as follows: comparison of atomic positions between the two configurations; displacement of polyhedra so that their centres coincide; and rotation of the polyhedra so that their orientations match as closely as possible. The process is illustrated in Fig. 1. Effectively, we decompose the motions of all the atoms of the polyhedron into a bodily displacement (given by the displacement of the central atom), a bodily rotation and a residual distortion.

Figure 1
The process of fitting; the comparison of structures proceeds via decomposition into polyhedra, displacement for coincidence of centres, and rotation, leaving some residual distortion. Comparison of the magnitude of the residual distortions with the magnitudes of the atomic motions reveals the significance of rigid-unit motion.

2.1. Definition of mismatch scores

We define an initial mismatch score to be the sum of the squares of the displacements of all the atoms from one configuration to the other. Thus, if the atoms in one configuration have positions r_i and each atom's position in the second configuration is $[r_i + \Delta r_i \,]$ , the initial mismatch score before any fitting takes place is now defined to be

$[M_{\rm before} = \textstyle \sum\limits_i \Delta r_i^2. \eqno (1)]$

Since this score is simply a sum over atomic contributions, we can either consider the mismatch score for the structure as a whole or decompose the score by atomic species; for each species, we obtain a mean-square atomic motion

$[\left| \Delta r_{\rm{spec}}^2 \right| = {n_{\rm{spec}}^{-1} \textstyle\sum\limits_{i_{\rm{spec}}} \Delta r_i^2 \eqno (2)]$

and an r.m.s. atomic displacement

$[\overline{\Delta r_{\rm{spec}} } = { \left| \Delta r_{\rm{spec}}^2 \right|^{1/2} }, \eqno (3)]$

while, of course, the r.m.s. atomic displacement taken over all atoms is simply

$[\overline{\Delta r } = M_{\rm before} ^{1/2}. \eqno (4)]$

Since M_before scales linearly with the size of the system, comparison between different systems is easier if we divide M by the size of the system – that is, the number of polyhedra – to obtain the value of M per polyhedron. The mismatch values in the case studies that follow are given in this form.

Having obtained these measures of the atomic motion, we now consider the structure as an assemblage of polyhedra, decomposing it as illustrated in Fig. 1. We make the approximation that the displacement of each polyhedron is equal to the displacement of its central atom. We now wish to identify the rotational motions of structural polyhedra. Let us consider the structure bond by bond. For each bond, there is a vector from the central atom to an atom q at a vertex of the polyhedron p; for one form of the structure, we call this vector pq, for the other form, $[{\bf pq'}]$ . Our assumption, from the rigid-unit picture, is that there exists a rotation of the polyhedron p that takes the set of vectors $[{\bf pq}]$ very close to the set $[{\bf pq'}]$ . We can write down a vector $[{\bf pq''}({\bf pq},B)]$ , representing the vector pq after a rotation B; this allows us to write a mismatch vector, m, given by $[{\bf m}_q = {\bf pq''} - {\bf pq'}]$ . Minimizing the total magnitude of all the mismatches for the polyhedron, $[ M_p = \textstyle\sum_q {\bf m}_q^2]$ , with respect to the parameters of the rotation B will give the rotation that best fits pq onto $[{\bf pq'}]$ . The value of M_p before minimization with respect to B represents the degree of mismatch between the polyhedra when they have been displaced so that their centres coincide but they have not yet been allowed to rotate. We therefore designate this value as M_{p,fit centres}. The final value of M_p, representing the residual distortion of the polyhedra, we write as M_p,after.

M_after, the average value of M_p,after for all the polyhedra in the structure, can be compared with the original measure of the degree of atomic motion, M_before per polyhedron, in order to obtain a measure of the importance of rigid-unit motion. For concision, we frequently refer to M_before as the disorder score and to M_after as the distortion score.

Full mathematical details of the process of finding the best-fit rotation are given by Wells, Dove & Tucker (2002). However, note that the quantity defined in that paper as `M_before' or as `total disorder' is, in fact, the sum over all polyhedra of our M_{p,fit centres} and is not related to M_before or the disorder score as defined in the present paper.

2.2. Interpretation of rotors

We describe a rotation using three independent variables, B_x,y,z, which define a `bivector' object describing a `rotor' operator; these are concepts of geometric algebra. To interpret the results in the case studies that follow, it suffices to interpret B in terms of an axis $[\hat{b}]$ and an angle $[\varphi]$ . The size of the rotation angle $[\varphi]$ , in radians, is given by

$[\varphi = 2 \arcsin ({{\left| B \right|}/{2}}). \eqno (5)]$

The axis of the rotation is given by a unit vector, $[\hat{b}]$ , whose Cartesian components, b_i, are given by

$[\textstyle\sum\limits_i b_i^2 \, = \,1, \eqno (6)]$

$[b_x: b_y: b_z \, = \, B_x: B_y: B_z. \eqno (7)]$

The case B = 0 is the identity operation, that is, R(0) = 1. To first order, therefore, the values of B_x,y,z are the rotations about the x,y,z axes, in radians.

2.3. How important is rigid-unit motion? Comparison of mismatch scores

We wish to define a measure of the importance of rigid-unit motion in the dynamics of a structure. Perfect accommodation (zero distortion score M_after) cannot be this criterion, since even when RUMs are dominant, higher-frequency distortive modes will also be active. Some distortion of the polyhedra is therefore inevitable. In the disorder score, each atom contributes the square of its displacement, $[\Delta r^2]$ . For each atomic species we find an r.m.s. atomic displacement

$[\overline{\Delta r_{\rm{species}} } = { \left| \Delta r_{\rm{species}}^2 \right| }^{1/2}. \eqno (8)]$

The distortion score is the sum of the mean-square displacements of the vertex atoms from their geometrically ideal positions; when different types of polyhedra are present, this value will, of course, be different for each. Thus, from the distortion score, we can obtain the r.m.s. distance of vertex atoms from ideal vertices, d_{atom vertex}, and these distances can be compared with the r.m.s. atomic displacements in order to gauge the importance of rigid-unit motion. Our identification of the displacement of central atoms with the displacement of polyhedra means that central atom species make no direct contribution to the distortion score. The significant comparison is therefore between the motion of vertex atom species and the distortion of polyhedra. In the case of quartz discussed below, for example, the r.m.s. motion of O atoms is twice as large as the r.m.s. distortion, d, at room temperature, and four times as large at high temperatures in the $[\beta]$ phase.

The value of M_{fit centres} is not particularly significant in itself. We do, however, note that, when rotational amplitudes are large, M_{fit centres} lies at or above the value of M_before. If, on the other hand, displacements are more significant than rotations, M_{fit centres} will lie below M_before. Thus far, we find rotational-dominant behaviour in the crystalline framework silicates such as quartz, while silica glass also shows significant displacive motion.

Evidently, the fitting process generates a rotor for every polyhedron in the structure. The r.m.s. magnitude of the rotor, or of a component of the rotor, provides an indication of the magnitude of rigid-unit motion in the structure, and it may be of interest to observe this magnitude as a function of temperature. Furthermore, it is possible to plot the distribution of rotor components for all the polyhedra, rather than finding an r.m.s. value; this process may yield information on the nature of the rigid-unit motion. We use both forms of analysis in our discussion of quartz below.

3. Applying a force model via geometric modelling

The geometric model on which this analysis is based, in which each polyhedron exists both as a set of real atoms and as an ideal geometric form, can itself be used to generate configurations of atoms by establishing the ideal forms first and then constraining the real atoms to fit them. This process provides a simple means of studying the geometric response of frameworks to substitutional defects or compression, and of rapidly relaxing idealized structures (for example, silica networks with 180° Si—O—Si angles) to obtain realistic bond lengths and angles.

The model requires several iterations over two distinct steps. In the first step, we fit a geometrically ideal polyhedron over each real polyhedron using the GA rotor-fitting approach, so that the vertices of the ideal polyhedra provide a close match to the positions of the vertex atoms. The initial fitting of the ideal polyhedron involves only three variables (the parameters describing the orientation of the polyhedron), and this step is easily performed using the method given by Wells, Dove & Tucker (2002). In the second step, each atom moves within a simple-harmonic force model based on the positions of the vertices. The displacement of the atom relative to the position of an ideal vertex can be divided into the components parallel and perpendicular to the direction of the ideal vertex–centre (e.g. Si—O) bond. The parallel component represents the extension of the bond, while the perpendicular component represents distortion of the internal bond angles of the polyhedron. To each component we associate a spring force, tending to pull the atom towards the ideal position. The three-body terms constraining the O—Si—O bond angles are thus replaced by constraints on the displacement of the vertex atoms from their ideal positions. In keeping with the principle of using two-body rather than three-body potentials, the centre–vertex–centre (e.g. Si—O—Si) bond angle constraints are represented as constraints on the centre–centre distances.

Having relaxed the atoms within this simple force model, we then refit the ideal polyhedra over the real polyhedra using rotor fitting and iterate until the structure is fully geometrically relaxed. The geometrically ideal polyhedra therefore represent a form of multibody potential in which the evaluation of bond angles is unnecessary, as angular information is implicit in the shape of the polyhedron. It is, of course, possible to include constraints on the relaxation by, for example, pinning certain atoms or polyhedra in place. This allows for a direct connection between the real-space and reciprocal-space pictures of polyhedral motion, since the ideal polyhedra of the model are the objects of RUM theory. Thus we can both obtain the RUM component of a given set of atomic motions from the motions of the ideal polyhedra, and apply one or more RUM eigenvectors to the ideal polyhedra and observe the response of the real polyhedra.

This geometric modelling process has been used successfully to model the high-P deformation mechanism of a zeolite framework with EDI topology (Gatta & Wells, 2004 ).

4. Program availability and capability

We hope that the analysis will be of interest to the crystallography community. We have therefore implemented it in a self-contained program, called GASP (geometric analysis of structural polyhedra). The source code is written entirely in Fortran 90 and is available with documentation from http://www.esc.cam.ac.uk/minsci/downloads/gaspfiles.tar.gz . Questions, comments, bug reports and fixes should be reported to the corresponding author. The program in its current version (1.5) performs three principal functions. The first function, known as preconditioning, analyses a structure and identifies the polyhedra within it. During this process a connectivity file is generated, recording which atoms are bonded together, and an ideal polyhedron is fitted over each real polyhedron, quantifying the distortion of the polyhedra from their geometrically ideal forms. If the connectivity is already known, it can be provided as input. Tetrahedra, octahedra and planar triangles can all be present.

The comparison function takes as its input two structures that have been preconditioned, that is, decomposed into their component polyhedra. This function performs the analysis described in §2.1, quantifying the difference between the two structures and analysing it into components of rigid-unit displacement, rigid-unit rotation and distortion. The two structures could be two independent RMC fits, an ideal (average) structure and a fitted (disordered) structure, two frames from an MD simulation, or representations of a structure before and after the introduction of a defect.

The third function makes use of the modelling technique described in §3 to relax a structure towards geometric ideality. The program and analysis are designed for the comparison of structures containing very large numbers (of order 10²–10⁵) of atoms. The program is driven by a text file containing geometric information and program commands. Structures can be given either in a native .dat format or in the common .xtl format. Full information on file formats and commands is available with the code.

5. Case study

Our principal case study is the analysis of a series of quartz configurations generated by reverse Monte Carlo fitting to total neutron scattering data. Frameworks of corner-linked SiO₄ tetrahedra such as these are the classic objects of the RUM model, and such frameworks are therefore a suitable first application for our real-space analysis. Our analysis provides information on the RUM contribution to the dynamic disorder and the effects of phase transitions on the RUM behaviour. We quantify the significance of RUM motions and show that they are very important in the dynamic disorder, particularly at high temperatures. By analysis of the $[\alpha]$ and $[\beta]$ phases, we show that the disorder in the high phase can be explained on the basis of the excitation of RUMs. We show the effect of the phase transition on the RUM behaviour in the $[\alpha]$ phase, demonstrating the link between the symmetry change due to the transition and the amplitudes of rotations about different axes.

While quartz shows considerable geometric flexibility, it is not the most flexible crystalline silicate. We compare our results for quartz with those for cristobalite and zeolite Y, demonstrating that larger amplitudes of atomic motion are obtained in the more RUM-flexible frameworks.

5.1. Quartz

In a recent total neutron scattering study of quartz (Tucker et al., 2000; Tucker, Keen & Dove, 2001), measurements were performed over a wide range of temperatures below the phase transition temperature and for some temperatures above the phase transition. The RMC configurations therefore provide information about changes in both short-range and long-range order as a function of temperature on heating through the phase transition. The configurations contained 6000 SiO₄ polyhedra. Sections of some of the configurations so generated are shown in Fig. 2. It is clear that, as temperature increases, the variation about the `average' structure (the structure found from the Bragg peaks) increases markedly. Since RUMs provide a mechanism for large reorientations of the polyhedra at small energy cost, we would expect RUMs to play an important role in this variation. By visual inspection, it would seem that the polyhedra are undergoing rigid-unit motion. The rotational motion and polyhedron distortion were considered by Wells, Dove, Tucker & Trachenko (2002). We can now complete that analysis by considering also the total atomic motion. Fortunately, data were available for a large number of points covering both the low and the high phase, and Fig. 3 shows the degree of disorder in quartz as a function of temperature. This graph shows the mismatch scores defined in §2, in Å² per polyhedron, between two independent RMC fits to the scattering data at each temperature.

Figure 2
Sections through RMC configurations of quartz. With increasing temperature, the variation of atomic positions about the positions in the average structure (insets) increases dramatically. Note the large atomic displacements and the apparently rigid-unit character of the polyhedral motion.

Figure 3
Average mismatch scores, M, in Å² per polyhedron, for quartz. At each temperature, we are comparing one RMC configuration with another, where both configurations are fitted to the same data. The difference between the disorder M_before and the tetrahedral distortion M_after is due to rigid-unit motion. Since the disorder is much greater than the distortion, the dynamic disorder in quartz is largely accounted for by rigid-unit motions of the tetrahedra.

Note that the tetrahedral distortion is much less than the total disorder. The large difference between the two indicates the contribution of RUMs to the dynamic disorder. The value of M_{fit centres} lies slightly above the value of the distortion score M_before, indicating that the RUM motions are predominantly rotational rather then displacive, as discussed in §2.3. Note particularly that the distortion score is relatively insensitive to temperature and to the phase transition, while the total disorder increases sharply with temperature and shows the phase transition quite clearly.

The tetrahedral distortion can be further divided into the terms due to variation in bond lengths and due to bending of the O—Si—O bond angle. It is clear that the distortion is primarily due to bond bending, with a much smaller bond-stretching component. The nature of the distortions does not appear to change over the temperature range of this study. Wells, Dove, Tucker & Trachenko (2002) showed that there is no correlation between the distortions in the two different configurations.

Fig. 4 shows the breakdown of the mismatch scores into r.m.s. motion values for Si and O atoms, plus distortion values for vertex O atoms in SiO₄ polyhedra. The motion of O atoms is twice the distortion value at room temperature and in the $[\beta]$ phase reaches four times the distortion value. This comparison is very revealing of the importance of rigid-unit motion. Since the displacement of the Si atoms is identified with the displacement of the polyhedra, they do not contribute directly to the distortion value.

Figure 4
The mismatch values for quartz, interpreted as r.m.s. motions of O and Si atoms, and r.m.s. distortion of SiO₄ polyhedra. It is clear that the vertex O atoms are far more mobile than the Si atoms, as we would expect from rigid-unit motion. The r.m.s. distortion (distance of O atoms from geometrically ideal vertices) lies well below the r.m.s. motion of O atoms, indicating that RUMs make a dominant contribution to the amplitude of the dynamic disorder.

5.2. Order parameter in quartz

In terms of the average structure, the $[\alpha]$ phase of quartz is obtained from the $[\beta]$ phase by the condensation of an RUM as a soft mode (Vallade et al., 1992 ). This RUM involves rotations of the polyhedra about the $[\langle 100 \rangle]$ axes, with the rotation angle reaching ∼17° at room temperature. It should be possible to detect this rotation using geometric analysis. We therefore compare the fitted configurations with an ideal (zero-order parameter) high-temperature form and observe the distribution of rotations of the fitted polyhedra. We expect to observe rotations of the polyhedra about the $[\langle 100 \rangle]$ axes, which in terms of GA means the following; if we plot the X and Y components of the rotors for the polyhedra, relative to the high phase, against each other, the results should cluster about three points. One point will be at $[B_y = 0, B_x \simeq \varphi]$ (where $[\varphi \simeq 0.3]$ radians or 17° at room temperature), while the other two points will lie at 120° to the first point in the XY plane. If more than one domain were present in the sample, there would be six equally spaced clusters instead of three.

This prediction is amply born out by the results given in Figs. 5 –8.

Figure 5
X and Y rotor components for polyhedra in a fitted configuration at 150 K, relative to their orientations in the ideal $[\beta]$ phase. The distribution of the rotations of the polyhedra, relative to the average positions in the $[\beta]$ phase, clearly shows the rotations about the $[\langle 100 \rangle]$ axes by which the $[\alpha]$ phase is derived from the $[\beta]$ phase.

Figure 6
X and Y rotor components for polyhedra in a fitted configuration at 473 K, relative to their orientations in the ideal $[\beta]$ phase. As the temperature rises, the distribution of the rotations broadens. The polyhedra display more variation about the mean positions as RUM amplitude increases in the dynamic disorder.

Figure 7
X and Y rotor components for polyhedra in a fitted configuration at 793 K, relative to their orientations in the ideal $[\beta]$ phase. At high temperatures near the phase transition, the distribution is so broad that the average rotation can barely be discerned. Note, however, that the distribution retains its threefold symmetry.

Figure 8
X and Y rotor components for polyhedra in a fitted configuration at 1073 K, relative to their orientations in the ideal $[\beta]$ phase. The distribution of orientations has become isotropic about zero and is very broad, indicating the large amplitude of rigid-unit motion of the polyhedra. There is no evidence for the formation of domains of the $[\alpha]$ phase.

With increasing temperature, the distribution of rotations relative to the high phase becomes broader as the rigid-unit component of the dynamic disorder increases in amplitude. In the $[\beta]$ phase, the distribution of rotations becomes isotropic at about zero and is very broad. It appears that the disorder in $[\beta]$ -quartz arises largely from the excitation of the RUM spectrum.

The value of the tetrahedral rotation angle, $[\varphi]$ , as obtained by fitting Gaussians over the rotor distributions and extracting the peak position, is plotted in Fig. 9. At room temperature, we find an average rotation of 17° relative to the $[\beta]$ phase, as expected.

Figure 9
The average rotation of the polyhedra in an RMC configuration of quartz, relative to their positions in the ideal $[\beta]$ phase. The rotation displays a dependence on the order parameter. Unfortunately the distribution of angles near the transition is so broad that the parameter cannot be obtained with any accuracy in this region.

5.2.1. Amplitude of rotational motion

The $[\alpha]$ – $[\beta]$ phase transition is clearly visible in the dynamic disorder, which drops dramatically as the temperature falls below the transition temperature. There exist phonon modes that are RUMs in the high-symmetry $[\beta]$ phase but develop a component of tetrahedral distortion in the low-symmetry $[\alpha]$ phase. This change is accompanied by a dramatic decrease in amplitude as the energy cost of the distortion increases the frequency of the mode. We attribute the changes in dynamic disorder with temperature to the effects of the symmetry change on the RUM phonons. We investigate this effect by considering the distribution of polyhedral rotations in the dynamic disorder, as obtained from the analysis in §5.1. For each polyhedron, there is a rotation describing its reorientation from one configuration to the other. The distribution of rotor components B_x,y,z describing these rotations is a Gaussian centred on zero. The r.m.s. value of the rotor component is a measure of the amplitude of rotational motion, as shown in Fig. 10. It appears from the behaviour of the rotor components that, in the $[\beta]$ phase, the polyhedra are freer to rotate about axes lying in the xy plane than they are about the z axis, a fact that reflects the hexagonal symmetry of the structure. The effect of the phase transition, however, with its imposition of an average rotation of the polyhedra about the $[\langle 100 \rangle]$ axes, is to constrain the x and y rotations and reduce them to the same amplitude as rotations about the z axis.

Figure 10
The r.m.s. rotor components (amplitudes of rotation of polyhedra) in the dynamic disorder of quartz. The r.m.s. values of the rotor components show different behaviour for the XY and the Z components with temperature. Rotations about the Z axis appear more constrained than those about the X and Y axes in the high-symmetry phase.

5.2.2. RUM frequencies and amplitudes and their relation to the order parameter

We now look more closely at this effect, whereby the symmetry-breaking in the phase transition affects the rotational behaviour, by considering the link between the order parameter, the frequencies of the RUMs and their amplitudes. A simple Einstein model for the rotational oscillations of the polyhedra suggests that

$[\langle \theta^2 \rangle \propto \textstyle\sum\limits_{\rm modes} {{T}/{ \omega_{\rm mode}^2 }}, \eqno (9)]$

and we therefore plot $[{{\langle \theta^2 \rangle}/{T}}]$ , a measure of $[\textstyle\sum_{\rm modes} {{1}/{ \omega_{\rm mode}^2 }}]$ , against T to observe the behaviour of the RUM frequencies. The split-atom model (Giddy et al., 1993, Hammonds et al., 1994) would suggest that, as the average tilt, $[\left| \varphi \right|]$ , of the polyhedra relative to their positions in the $[\beta]$ phase increases from zero, modes that are RUMs in the $[\beta]$ phase gain frequency in proportion to $[\sin \left| \varphi \right|]$ . Since the phase transition in quartz is of almost tricritical character (Carpenter et al., 1998 ), and the data do not approach the phase transition closely enough to allow us to observe any first-order discontinuity, we take it that

$[\left| \varphi \right| \propto (1 - {{T}/{T_c}})^{{{1}/{4}}}. \eqno (10)]$

As a simple model, we suppose that the value of $[{{\langle \theta^2 \rangle}/{T}}]$ for the high phase is a result of a population of RUMs having a typical frequency $[\omega_{\rm RUM}]$ , and that $[{{\langle \theta^2 \rangle}/{T}} \propto {{1}/{\omega_{\rm RUM}^2}}]$ . We extrapolate linearly from the values in the high phase to obtain the quantity $[{{\langle \theta^2 \rangle}/{T}}_{\rm RUM}]$ , that is, the predicted value of $[{{\langle \theta^2 \rangle}/{T}}]$ if all the RUMs from the high phase remained RUMs in the low phase.

We now suppose that RUMs in the low phase develop a component of tetrahedral distortion and gain frequency $[\Delta \omega]$ in proportion to the growth of the order parameter: $[{{\Delta \omega}/{\omega}} = K (1 - {{T}/{T_c}})^{{{1}/{4}}}]$ . These modes will therefore have their amplitude reduced and our prediction for $[{{\langle \theta^2 \rangle}/{T}}]$ then becomes

$[{{\langle \theta^2 \rangle}/{T}} = ({{\langle \theta^2 \rangle}/{T}}_{\rm{RUM}}) \left [{{1 + K (1 - {{T}/{T_c}})^{{{1}/{4}}}}} \right]^{-2}. \eqno (11)]$

We perform two separate fits, one to the average of the X and Y rotor components, and another to the Z rotors. The results are plotted in Fig. 11. For the X and Y rotors, the model fits the data well, with a value of K = 0.75. For the Z rotors, however, the fit is best with K = 0.19. This outcome is understandable in terms of our model in which rotations about the z axis are more constrained than those about the x and y axes in the high phase.

Figure 11
Mean-square rotor amplitude divided by temperature ( $[{{\langle \Theta^2 \rangle}/{T}}]$ ) versus temperature for quartz. The fit is to an Einstein model incorporating the dependence of RUM frequencies on the order parameter according to the split-atom model.

This result confirms the impression obtained from the rotor components in Fig. 10; rotations about the z axis are constrained throughout, while rotations about axes in the xy plane are freer in the high phase but become constrained because of the symmetry change in the phase transition.

5.3. Comparison with other crystalline silicates

As we have seen, quartz displays a considerable degree of geometrical flexibility; however, it is not the most flexible of the crystalline silicates.

Tucker, Squires et al. (2001 ) demonstrated that the structural disorder in the high phase of cristobalite could be accounted for on the basis of rigid-unit motion, without the need for a domain model. Wells, Dove, Tucker & Trachenko (2002) considered the flexibility of the cristobalite structure from the point of view of RUM theory, showing that it is greater than that of quartz.

Zeolites are another class of framework structures for which an RUM analysis is fruitful. The high symmetry and large unit cells of some zeolite structures give them a large number of RUMs, and in some cases, the number of RUMs can even be a finite fraction of the total number of modes (Hammonds et al., 1998 ).

We have performed geometric analyses on cristobalite configurations fitted to total neutron scattering data at five temperature points and on siliceous zeolite Y configurations fitted to total neutron scattering data at two temperature points, using the same RMC technique as for the quartz. A geometric analysis of the cristobalite configurations, considering only the rotational motions, was performed by Wells, Dove, Tucker & Trachenko (2002). We now present the results of our full analysis in the form of the r.m.s. motions of the Si and O atoms, and the r.m.s. distortion of the SiO₄ polyhedra (Fig. 12). The corresponding values for quartz, from Fig. 4, are also given for comparison.

Figure 12
Comparative graph of the r.m.s. motions of Si and O atoms, and the r.m.s. distortion distance in SiO₄ polyhedra, for quartz (open symbols), cristobalite (black symbols) and zeolite Y (grey). The magnitude of the motions in cristobalite and the zeolite are considerably greater than in quartz, reflecting the greater geometric flexibility of the structures. For all three silicates, the ratio of the r.m.s. oxygen motion to the r.m.s. distortion distance is large, indicating that rigid-unit motion is a dominant contribution to the amplitude of dynamic disorder.

The zeolite and the cristobalite display very similar levels of atomic motion and of tetrahedral distortion; these values may therefore be typical of highly flexible crystalline silica. The degree of atomic motion is considerably higher than that seen in quartz, reflecting the greater geometric flexibility of these two structures. For the zeolite, cristobalite and quartz in the $[\beta]$ phase, the ratios of the r.m.s. oxygen motion to the r.m.s. tetrahedral distortion are similar and indicate the dominance of rigid-unit motion in the dynamic disorder, while in $[\alpha]$ quartz at lower temperatures, this ratio decreases, as the symmetry change from the phase transition decreases the geometric flexibility.

6. Conclusions

We have developed an analysis connecting the motion of atoms in real space to the concepts of rigid-unit theory. Using this approach, we can quantify the importance of rigid-unit motion in dynamic disorder. We have shown that the method of `reverse Monte Carlo with geometric analysis' (RMC+GA) is extremely useful for extracting information from RMC simulations of framework structures. The method is implemented in a Fortran 90 code, GASP, which is freely available.

Acknowledgements

The authors thank Kostya Trachenko, Alexandra Simperler and Diego Gatta for helpful discussions and suggestions. Funding for the thesis upon which this work is largely based was provided by the EPSRC.

References

Carpenter, M. A., Salje, E. K. H., Graeme-Barber, A., Wruck, B., Dove, M. T. & Knight, K. S. (1998). Am. Mineral. 83, 2–22. Google Scholar
Dove, M. T., Trachenko, K. O., Tucker, M. G. & Keen, D. A. (2000). Rev. Mineral. Geochem. 39, 1–33. CrossRef Google Scholar
Gambhir, M., Dove, M. T. & Heine, V. (1999). Phys. Chem. Miner. 26, 484–495. Web of Science CrossRef CAS Google Scholar
Gatta, G. D. & Wells, S. A. (2004). Phys. Chem. Miner. Submitted. Google Scholar
Giddy, A. P., Dove, M. T., Pawley, G. S. & Heine, V. (1993). Acta Cryst. A49, 697–703. CrossRef CAS Web of Science IUCr Journals Google Scholar
Hammonds, K. D., Dove, M. T., Giddy, A. P. & Heine, V. (1994). Am. Mineral. 79, 1207–1209. CAS Google Scholar
Hammonds, K. D., Dove, M. T., Giddy, A. P, Heine, V. & Winkler, B. (1996). Am. Mineral. 81, 1057–1079. CAS Google Scholar
Hammonds, K. D., Heine, V. & Dove, M. T. (1998). J. Phys. Chem. B, 102, 1759–1767. Web of Science CrossRef CAS Google Scholar
Heine, V., Welche, P. R. L. & Dove, M. T. (1999). J. Am. Ceram. Soc. 82, 1793–1802. CrossRef CAS Google Scholar
Lasenby, J., Lasenby, A. N. & Doran, C. J. L. (2000). Philos. Trans. R. Soc. London Ser. A, 358, 21–39. CrossRef Google Scholar
McGreevy, R. L. (1995). Nucl. Instrum. Methods Phys. Res. Sect. A, 354, 1–16. CrossRef CAS Web of Science Google Scholar
McGreevy, R. L. & Pusztai, L. (1988). Mol. Simul. 1, 359–367. Web of Science CrossRef Google Scholar
Tucker, M. G., Dove, M. T. & Keen, D. A. (2000). J. Phys. Condens. Matter, 12, 723–730. Web of Science CrossRef Google Scholar
Tucker, M. G., Dove, M. T. & Keen, D. A. (2001). J. Appl. Cryst. 34, 630–638. Web of Science CrossRef CAS IUCr Journals Google Scholar
Tucker, M. G., Keen, D. A. & Dove, M. T. (2001). Mineral. Mag. 65, 489–507. Web of Science CrossRef CAS Google Scholar
Tucker, M. G., Squires, M. D., Dove, M. T. & Keen, D. A. (2001). J. Phys. Condens. Matter, 13, 403–423. Web of Science CrossRef CAS Google Scholar
Vallade, M., Berge, B. & Dolino, G. (1992). J. Phys. I, 2, 1481–1495. CrossRef CAS Google Scholar
Welche, P. R. L., Heine, V. & Dove, M. T. (1998). Phys. Chem. Miner. 26, 63–77. Web of Science CrossRef CAS Google Scholar
Wells, S. A. (2003). PhD thesis, University of Cambridge, UK. Google Scholar
Wells, S. A., Dove, M. T. & Tucker, M. G. (2002). J. Phys. Condens. Matter, 14, 4567–4584. Web of Science CrossRef CAS Google Scholar
Wells, S. A., Dove, M. T., Tucker, M. G. & Trachenko, K. (2002). J. Phys. Condens. Matter, 14, 4645–4657. Web of Science CrossRef CAS Google Scholar

© International Union of Crystallography. Prior permission is not required to reproduce short quotations, tables and figures from this article, provided the original authors and source are cited. For more information, click here.

JOURNAL OF
APPLIED
CRYSTALLOGRAPHY

ISSN: 1600-5767

Volume 37| Part 4| August 2004| Pages 536-544

doi:10.1107/S0021889804008957