Neutron scattering study of the orientational disorder in potassium cyanide

We report the results of a combined neutron diffraction and total scattering study of the orientational order-disorder phase transitions in potassium cyanide, KCN. The diffraction data are analysed in terms of the spontaneous strains that accompany the phase transitions. The total scattering data are analysed using the Reverse Monte Carlo method, which gives direct access to the distribution of atomic positions and hence molecular orientations in each phase. Incorporating information from diffuse scattering in this way provides a means to measure the coefficients of the orientation distribution function up to almost arbitrarily high orders, and furthermore has the advantage that this function is naturally positive everywhere. The results for the cubic phase show that the distribution of orientations never exceeds 25% difference from an isotropic distribution.


Orientationally disordered crystals
Orientational disorder is a common phenomenon in molecular crystals or inorganic crystals. 6 Often, at least in the first phases found on cooling below the melting point, the disorder is dynamic, with the molecules or molecular ions undergoing rotational motion across a continuous distribution of molecular orientations, or jump rotations between a well-defined set of distinct orientations. Examples range from molecular crystals [1] such as fullerenes [2] to polyatomic ions in metal-organic frameworks [3] and even hemes in proteins [4]. Orientational disorder may lead to anomalous materials properties, such as reversible amorphisation under pressure [5] and giant barocaloric effects [6]. For this reason, understanding orientational disorder in materials plays an important role in understanding and tuning [7] such properties of materials.
One simple and intensively studied case is that of the cyanide ion, CN − , in metal cyanides [8]. CN − exhibits two types of orientational disorder in metal cyanides. In the cristobalite-like group-12 metal cyanides Zn(CN) 2 [9] and Cd(CN) 2 [10], the cyanide ions can be treated as electric dipoles with static head-to-tail orientational disorder. Apparently this disorder persists to low temperature in Zn(CN) 2 , but recently reorientational motion and short-range order have been identified in Cd(CN) 2 [10]. On the other hand, in the alkali metal cyanides NaCN, KCN and RbCN with the rocksalt structure, the cyanide molecular ions have a much broader distribution of orientations and show rotational diffusion [11] rather than jump motions [12]. In fact, as we will Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. 6 In the past orientationally-disordered materials were called 'plastic crystals' as a result of their mechanical properties, but the scientific community later preferred the term we use here. However, as orientationally disordered crystals have come back into fashion, sometimes people are rediscovering the older name! For the good reasons why the name 'plastic crystals' was replaced by the name 'orientionally disordered crystals' we choose to work with the latter name. see below, the first phase transition on cooling after freezing involves ordering of the elastic quadrupole moments rather than of the molecular electrical dipole moments.

Methods to quantify the spatial part of orientational disorder from experiment
The atomic structure of crystals with continuous orientational disorder is not easy to study using conventional methods such as Bragg diffraction, because it is not reasonable to identify discrete sites for the positions of atoms. Instead, the atomic structure must be described in terms of a continuous distribution of molecular orientations. A common approach is to expand these distribution functions in terms of an appropriate basis of 'rotator functions' such as symmetry-adapted spherical harmonics [13,14]. This process reduces the problem of describing the structure to that of determining the appropriate coefficients for the basis functions, which might be done in analysis of Bragg diffraction data using a standard least-squares refinement of the Bragg diffraction data. However, this method has two related problems. First is that it is only practical to relatively low order, and second, as a consequence the resulting distribution function is not necessarily positive for all possible orientations, although there have been some practical suggestions how to solve this by fitting instead an effective rotational potential [15,16].
An alternative approach is to analyse a set of large molecular configurations. In the past this typically meant analysis of configurations generated by the molecular dynamics simulation method, but from an experimental approach this is now possible using neutron total scattering data together with the Reverse Monte Carlo (RMC) method [17][18][19]. One recent example of this technique for studying orientational disorder was an investigation of crystalline adamantane [20], but in this case the molecules were found to be confined to be in either of two preferred orientations rather than with a continuous distribution. If the orientational probability distribution function is expanded in terms of orthonormal functions, it follows that the coefficients c ℓ can be calculated directly from the statistical averages = á W ñ c K ( ) ℓ ℓ evaluated from a set of independent configurations. This is not subject to the problem of only allowing low-order terms and thus can be developed to arbitrary accuracy. To date this approach has never been applied to analysis of orientational disorder from total scattering and RMC data. This analysis will be discussed in more detail below.
In fact it is now relevant to develop this experimental approach, because some of the new optoelectronic materials that are attracting a lot of interest contain molecular ions such as methylammonium or formamidinium that are free to rotate. One example is the hybrid perovskite solar cell material methylammonium lead iodide (MAPbI 3 ), which has been investigated in a huge number of studies, yet still facing the problem of characterising the orientational disorder from experimental data. MAPbI 3 undergoes two phase transitions at ambient pressure involving changes in the degree of orientational order, and it has been shown that these correlate with its temperature-dependent dielectric properties [21]. The methylammonium includes a C-N bond, a direct analogue of the cyanide ion in KCN, for example, and indeed a recent study by molecular dynamics simulation of MAPbI 3 compared the ordering to that of KCN from single-crystal diffraction, suggesting that the methylammonium ions are somewhat more isotropic, although again the spherical harmonic expansion was terminated at very low order [22].
For these reasons it is timely to return to the example of KCN. Here we take this material as a case study in the use of the neutron total scattering method to provide information about the distribution of orientations-and hence of electric dipole and elastic quadrupole moments-in an orientationally-disordered or partially-ordered crystal.

The phase transitions in KCN
Below its melting point (907 K), potassium cyanide exists in three distinct crystalline phases [23]. From 168 K up to its melting point, KCN has the well-known rock salt structure with space group Fm m 3 [23,24]. The cyanide anions occupy sites of the full cubic point symmetry (m m 3 ) and therefore there is three-dimensional orientational disorder. At 168 K is a phase transition involving considerable but not complete ordering of the cyanide orientations [23]. The crystal structure of this phase is orthorhombic, space group Immm. The cyanide anions have their long axes parallel to the orthorhombic [010] axis (equivalent to the á ñ 110 axes in the cubic phase) but with head-to-tail disorder of the electric dipoles [24,25]. The ordering can be described as involving the molecular elastic quadrupole moments. It is well established that the ordering couples linearly to the softening of the C 44 shear elastic constant as a result of coupling between the translations and rotations [26][27][28][29][30][31], meaning that the phase transition is a proper ferroelastic phase transition. At 83 K there is a second phase transition to a fully-ordered structure of orthorhombic space group Pmnm [23,32]. This has ordering of the electric dipole moments in a proper antiferroelectric arrangement [23,33].
The atomic structures of the three crystalline phases of KCN are illustrated in figure 1. This figure shows the relationship between the phases, and illustrates the change in orientational ordering.

Outline of this paper
Although the main thrust of this paper is towards a new quantification of the orientational distribution in KCN, we take a wider look too at the sequence of phase transitions. At the time when KCN was being studied before, there was little opportunity to measure the structure across a wide range of temperatures, and in addition to our total scattering measurements presented here we also consider in more detail than before the nature of the structural changes as a function of temperature. In particular, we can analyse spontaneous strains which give high-precision information about the changes in order, and it seems that until now there are only sketchy reports of the crystal structure of the low-temperature ordered phase.
In the next section we consider the methods we used, both experiment and analysis. We then consider the crystal structure analysis in more detail in section 3. This is followed in section 4 by our analysis of the orientation distribution functions from neutron total scattering and RMC analysis.

Sample
The samples of KCN were obtained commercially from (as then) the Aldrich company. Quoted purity by standard assay was more than 98%, with the largest cited contaminant being Na as a replacement for K at the level of 0.5%. Samples were ground within a glove box in a dry atmosphere, and then immediately loaded into the sample can (see below) whilst still within the glove box.

Neutron scattering measurements
Two sets of neutron scattering experiments were performed in this study. A set of neutron powder diffraction measurements were performed using a 10-detector array on the E3 spectrometer/diffractometer on the NRU nuclear reactor at the Chalk River Laboratories [34]. The second set were a combination of neutron diffraction and total scattering measurements performed on the GEM diffractometer at the ISIS pulsed spallation neutron source in the UK [35]. The data from ISIS are publicly available 7 . The data reduction for the NRU data was performed using in-house software, and for the ISIS data using the MANTID software [36].
In both cases the sample was contained within thin-walled vanadium cans. Temperatures were controlled by standard cryostats.

Analysis of diffraction data
The crystal structures of KCN were refined using the the GSAS/EXPGUI programs [37][38][39]. Starting coordinates for the intermediate temperature phase were taken from literature values with disordered occupancy of the C/N sites. For the low-temperature phase the structure was taken as an antiferroelectic ordering of the electric dipole moments of the intermediate-temperature phase. The disordered hightemperature phase has no well-defined positions for the C/N atoms, but nevertheless such positions are and (c) the high-temperature phase (above 168 K) with black, cyan and purple balls representing carbon, nitrogen and potassium atoms respectively. The grey lines represent a pseudo-rhombohedral or a cubic cell containing the same set of atoms. In the low-temperature phase, the cyanide ions align along the y axis with the same direction within one layer (shown as red dashed lines) and adjacent layers point in opposite directions. This is known as head-to-tail order. In the intermediatetemperature phase, the cyanide ions lose head-to-tail order while still aligning along the b axis. In the high-temperature phase, the cyanide anions are no longer restrained to align along b and orient randomly over the full unit sphere. required for standard Rietveld analysis and so were selected to lie parallel to the cubic unit cell axes. Backgrounds were fitted using Chebychev polynomials.

Total scattering data and the reverse Monte Carlo analysis
The total scattering data were transformed to the appropriate scattering functions i(Q) and pair distribution functions D(r) using standard procedures; we have described these functions in several previous publications, including one where they are compared with other functions to be found in the literature [40,47]. The function D(r) is the sine Fourier transform of the scattering function Qi(Q). To extract i(Q) from the scattering data requires account to be taken of the background scattering within the instrument, and the scattering and beam attenuation by the sample container and sample environment. These are obtained in separate measurements, together with a normalisation measurement of a vanadium rod, and the final normalised scattering function i (Q) was formed using the GUDRUN package [41,48].
The Reverse Monte Carlo simulation was performed using the RMCprofile code [42]. Configurations of the low-temperature phase used a 12×10×8 supercell based on the refined crystal structure, set up using the data2config/RMCcreate code [43]. The configuration for the intermediate-temperature phase was set up in a similar manner, but using an option to randomly rotate the cyanide molecular ions by 180°about the normal to the molecular axis. The configuration for the cubic phase was a 8×8×8 supercell of the cubic phase, using an option to give each cyanide molecular ion a random orientation. The RMC simulations were carried out for 2 315 263 steps in total until convergence (see figure 2(d), giving 70.1 accepted moves per atom and an average move acceptance rate of 18.6%. For each temperature we generated 100 independent configurations for analysis.
3. Rietveld analysis of KCN: lattice parameters and spontaneous strain 3.1. Crystal structure refinements The atomic structure of the ordered Pmnm phase has atomic coordinates K u , , ( ) and N y , , The atomic structure of the partially-ordered Immm phase has atomic coordinates K (0, 0, 0) and  ( ) with fractional occupancy of 1/12 for both atom types, as an over-simplified model in order to be able extract lattice parameters. Results for the crystal structure of all three phases refined from the GEM data are given in tables 1-2. We believe that this is the first time that the crystal structure of the Pmnm phase has been reported in detail.

Spontaneous strains from the lattice parameters
The variations of the lattice parameters of the cubic and orthorhombic phases with temperature are shown in figure 3(a), using the scaled parameters ¢ = a a 2 and ¢ = b b 2 for the orthorhombic phase in order to give values that will tend towards the value of the lattice parameter of the cubic phase. To interpret the lattice Table 1. Crystal structure data for the low-temperature ordered phase of KCN, space group symmetry Pmnm (number 59). Atomic coordinates are of the form y , , 1 4 1 4 ( ) for atom type K, and of the form y , , 3 4 1 4 ( ) for atom types C and N. These results were obtained by Rietveld refinement of the diffraction data from the GEM diffractometer at ISIS.  (2) 1.90 (7) 1.14 (7) 1.13 (7) (2) 1.94 (7) 1.00 (7) 1.36 (7) (2) 2.17 (7) 1.28 (7) 1.38 (7) (7) 1.14 (7) 1.49 (7)   ( ) , both with occupancy 0.5. These results were obtained by Rietveld refinement of the diffraction data from the GEM diffractometer at ISIS.   (4) parameters, it is instructive to look at the crystal structure viewed down the orthorhombic c axis, figure 3(b). The cyanide anions orient along one of the cubic á ñ 110 axes, which corresponds to the orthorhombic b axis. This leads to an expansion of the b lattice parameter on cooling compared with its value in the cubic phase, with a corresponding reduction in the size of the a and c axes.
It is striking that the volume of the unit cell appears to show no change in passing through the cubicorthorhombic phase transition. We can therefore define the parameter a 0 =(a′ b′ c) 1/3 , such that a 0 acts as the extrapolation of the cubic a lattice parameter to low temperature. We define the three axial spontaneous strains as The axial strains can be combined to create three symmetry-adapted strains [44]. The first is the totallysymmetric volume strain: The second is the tetragonal strain, namely an expansion along b with a shrinkage in the two orthogonal directions: The third is the shear strain, which in the axial setting of the orthorhombic phase is defined as = -   5 6 2 1

( )
These three strains are shown in figure 4(b). As expected, ò a is small across the whole range of temperatures. And also, as expected, ò 6 is the largest, but ò t is not insignificant in size.
In figure 5 we show that over the range of the intermediate-temperature Immm phase ò t ∝ò 6 . This linearity is slightly broken in the low-temperature phase. As will be discussed below, this linear relation is unexpected and unexplained.

Landau free energy function
In order to understand the strain behaviour associated with the phase transition, we obtained the form of the Landau expansion of the free energy associated with the Fm m 3 -Immm phase transition from the ISOTROPY software [45]: where Q is the order parameter, which at this point does not need to be physically defined but will be associated with the orientational order of the cyanide molecular anions. The C ij 0 parameters are the bare elastic constants, noting that the actual values of the elastic constants themselves will vary with temperature, including, significantly, a softening of C 44 . At equilibrium, the three strains can be obtained by minimisation of the free energy:    Substitution back into equation (6) gives From these equations we see that the linear coupling between the shear strain ò 6 and order parameter causes an increase in the phase transition temperature, and that quadratic coupling between the other two strains and order parameter lead to a decreased value of the fourth-order coefficient. The significant size of this coupling leads to the phase transition becoming first-order (discontinuous). The analysis suggests that we expect -µ   t 6 2 , consistent with the fact that ò 6 captures the full symmetry change from Fm m 3 to Immm as a proper ferrroelastic phase transition, whereas the strain ò t would preserve a 4-fold rotation axis in addition to other symmetry operations. This analysis is not consistent with the linear relationship between ò 6 and −ò t seen in the experimental data, figure 5; it is not at all clear to us how such a clear linear relationship can exist.

4.
Orientational order and local atomic structure from neutron total scattering and RMC simulation 4

.1. Atomic distributions
The spatial distribution of individual atoms in the unit cell was constructed by collapsing each configuration onto one unit cell, and merging the results from many independent RMC simulations. Figures 6(a)-(c) shows the density profile for each atom type at the cubic phase. Figure 6(a) reflects the large amount of thermal motion of the K atom as seen in the size of the atomic displacement parameter from the Rietveld analysis, table 2. The distribution of C and N atoms is wider in extent. What is interesting from figures 6(b) and (c) are that the distribution of C and N atoms are centred in the middle of the unit cell, rather than the distribution of C and N atoms forming a spherical shell. Thus the orientational disorder is accompanied by a high degree of translational disorder.
The corresponding atomic distributions for the intermediate-temperature (100 K) partially-ordered Immm phase are shown in figures 7(a)-(c), and for the low-temperature (20 K) ordered Pmnm phase in figures 7(d)-(f). The reduction in translation disorder is clear, and now we see also the distinct sites for C and N atoms. In the intermediate phase there is head-to-tail disorder of the cyanide molecular ions, which is clearly seen in the distribution of C and N atoms.

Information from pair distribution functions
The partial PDFs g(r) for K-K, K-(C/N) and (C/N)-(C/N) are shown in figures 8-10 respectively. For the K-K PDF (figure 8), the high-temperature data show broad peaks at positions consistent with the face-centred cubic arrangement of ions, and in the lower-temperature phases the transition to the orthorhombic phases results in the splitting of these peaks, with the positions being very similar in the two orthorhombic phases but with the expected sharpening of the peaks on cooling through the phase transition from 100 K to 10 K. The peaks up to a distance of around 7 Å in the PDF for the two lower-temperature phases correspond to peaks in the PDF of the high-temperature phase, with the expected broadening, but after 7 Å the peaks positions of the cubic phase are clearly out of register with those in the orthorhombic phases, reflecting the change in symmetry. The first peak in the PDF of the orthorhombic phase is actually comprised of three peaks, at mid-point distances of around 4.2, 4.4 and 4.6 Å. The first peak in the PDF of the cubic phase is centred on position a 2 , and it is broad to the extent that it encompasses the distribution of peaks in the low-temperature phases. The second distinct peak in the PDF for the data at 10 K becomes a shoulder in the data for 100 K, and is subsumed within the tail of the first peak in the cubic phase. The fact that the first peak in the PDF of the cubic phase encompasses the distribution of peaks within the low-temperature orthorhombic peaks suggests that there are local distortions of the unit cell that reflect the local orientations of the CN molecular anions. The same can be said of the second peak too, The K-C and K-N g(r) results (figure 9) for 10 K show two peaks at around 3.1 and 3.8Å, which correspond to two cases in which the K atom lies equidistant from both atoms in the CN molecular anion and where the K.. CN and K...NC connections are almost straight lines. It is interesting that we see the same basic form of this g(r) for the cubic phase too for both the same first pair of peaks and beyond, albeit with expected broadening. Clearly  the first peak should be the same because it is the nearest-neighbour contact between a K atom and either a C or N atom, but it is interesting that the second peak remains into the high-temperature disordered phase, but broadened and with a midpoint shifted to lower distance. The remaining peaks for the cubic phase are consistent with the positions expected for the cation-anion distance in the rocksalt structure with broadening due to the orientational disorder of the cyanide anions.
The C-C, C-N and N-N g(r) results are shown in figure 10. For the cubic phase, there is no difference between C and N, reflecting the orientational disorder, and the g(r) function for both is similar to that for the K-K g(r) (figure 8) with broadening of peaks reflecting the orientational disorder. In the low-temperature phase, the C-C and N-N distributions are almost identical, and due to the ordering these are not the same as the C-N distribution.

Molecular orientational distribution function
As the simplest type of molecule possible, the cyanide molecular anions consist of only two atoms and therefore containing only one chemical bond. The orientational distribution can be described quantitatively by the bond orientational distribution function P(Ω). Here Ω represents the polar angles (θ, f), where θ is the zenith angle (0θπ) and f is the azimuthal angle (0f2π). P(Ω) describes the probability of a C-N bond lying within a given element of solid angle dΩ=sin θ dθ df.  We expand P(Ω) in an orthonormal basis of 'rotator functions' S i (Ω): where the S i are themselves linear combinations of spherical harmonics W Y m ( ) ℓ [31]. It will always be convenient to take real combinations such as where the 'c' refers to the fact that this function is proportional to f m cos( ) 8 [46]. Further combinations may be required, and particular values of m and ℓexcluded altogether, depending on the molecular site symmetry. Here we give only a brief explanation of the mathematics relevant to this particular system; full details for all crystallographic point groups are tabulated in the reference by Bradley and Cracknell [46]. With appropriate normalisation applied, we have c 0 =1.
In  [14,47], and are given for convenience in the appendix.
As noted earlier, the coefficients in the bond orientational distribution function can be obtained directly by forming the average value of the relevant harmonic from the RMC configurations, in contrast to fitting to data as in diffraction experiments [14]. Here we can compare the calculated distribution function directly with that formed from histograms of bond orientations in the configurations. The form of P(Ω) for the disordered cubic phase for two arcs of Ω are shown in figures 11(a) and (b) for temperatures 200 K and 250 K respectively. These show the function calculated from the histogram and from the form of P(Ω) of equation (1). Values of the coefficients c ℓ were evaluated from averaging over around 100 independent configurations, and are given in table 3, and compared to values obtained from molecular dynamics simulations [25]. All values are small compared with c 0 =1, indicative of a high degree of disorder as seen in figures 11(a) and (b), where the fluctuations from the uniform value P(Ω)=1/4π do not exceed 10%. There is a very slight preference for ordering along the cubic á ñ 100 directions, and a very slight preference away from the á ñ 111 directions, with no obvious preference for the á ñ 110 directions, the directions in which the cyanide molecular anions order in the two lower-temperature orthorhombic phases. The results appear to show a slightly less uniform distribution at 250 K than at 200 K, but the difference is very small and may not be significant. The RMC results for the coefficients c ℓ are similar to those from the simulations [25], particularly in so far as the coefficients for ℓ>0are remarkably small. Although the coefficients determined by a single-crystal neutron diffraction study [48] presented in table 3 appear slightly different-albeit with a relatively high error-the bond orientational distribution function has some similarities to the results presented here but with larger variance from uniform distribution.
The coefficients for the low-temperature Immm phase at 20 K and 50 K, and for the intermediatetemperature Pmnm phase at 100 K and 140 K, are given in tables 4 and 5 respectively. The bond orientation functions are shown in figure 12, showing a much higher degree of orientational order as expected. The degree of orientational order/disorder described by the orientation distribution function of each phases is highlighted by calculating the three-dimensional representations as shown in figure 13, together with views of representative RMC configurations.
The picture that emerges from this analysis is that there is almost uniform distribution of orientations of the cyanide anions in the high-temperature cubic phase. The small preference for or against particular orientations is only of order of 10% of the uniform distribution. On the other hand, we see a much greater degree of quadrupolar ordering (that is orientational ordering of the bond ignoring the direction of the dipole moment) in the lower-symmetry phases, and complete ordering on the low-temperature phase. Furthermore, in the lowtemperature phase the widths of the distribution function reflect and increase in librational motion on heating from 20 K to 100 K. 8 We follow the results and notation of Bradley and Cracknell [46], but note that since they do not include the Condon-Shortley phase factor in their definition of W Y m ( ) ℓ , neither do they need or include the factor of (−1) m . The equation is given here in the form that we believe will be most useful to a modern reader; in particular this form is compatible with the spherical harmonic functions built in to most popular scientific computing systems, including Mathematica and the SciPy Python library, and with the International tables for Crystallography.

Conclusions
In this paper we have used neutron scattering methods-Bragg scattering analysed by the Rietveld method, and total scattering analysed by the Reverse Monte Carlo method-to study the orientational-ordering phase transition in KCN in more detail than was previously possible. Bragg scattering has primarily given information about the lattice parameters, which in turn give information on the spontaneous strains that accompany the phase transitions. The strains are coupled to each other and to the degree of order as determined by the symmetry of the free energy function, enabling us to quantify the change in the degree of order within the intermediate and low-temperature phases. That said, there is a linear relationship between the shear and the tetragonal strain that is very difficult to understand.
The orientational order in each phase of KCN has been evaluated from the atomic configurations generated by the RMC method. We have shown that the distribution of orientations of the CN bonds does not vary significantly from random, with the largest deviation being around 25% in favour of orientations in the cubic á ñ 100 directions and showing a reduction of similar size in the á ñ 111 directions. We were able to obtain good Figure 11. Bond orientation distribution function for the CN molecular anion in KCN at two temperatures in the disordered cubic phase at two temperatures, a) 200 K and b) 250 K. In both cases the points represent averages from many RMC configurations, and the red lines are calculated from the Kubic harmonic expansion with coefficients calculated from the configurations. The left panels show orientations in the directions from á ñ 001 to á ñ 110 , and the right panels show the directions from á ñ 110 to á ñ 100 . The broken black line at f (Ω)=1/4π shows the value for a uniform distribution. agreement between histograms of the bond orientation distribution function formed from the atomic configurations with an expansion of the bond orientation distribution function in terms of symmetry-adapted orientational functions. We are able to take the expansion of the distribution function to high order and therefore maintain positive-definiteness, a problem that is encountered when fitting to diffraction data. For future studies of orientational order, we believe that the experience described here should prompt the use of similar methods. Table 3. Calculated values of the coefficients c ℓ of the Kubic harmonic expansion of the bond orientation distribution function P(Ω) defined in equation (1). RMC results are compared with values calculated in a molecular dynamics simulation [25] and single crystal neutron diffraction study (SCND) [48].    Figure 12. A slice view of spherical harmonic expansion and cyanides orientation distribution for KCN at 20 K, 50 K, 100 K and 140 K. Each subplot demonstrates the direction from á ñ 010 to á ñ 001 . The black dots represent the orientation distribution directly counted from the RMC output configuration and the red curve represents the orientation distribution function calculated from cubic harmonics. There is a broken black line at f (Ω)=1/4π indicating the uniformly distributed states.