Rietveld refinement

For general application of the Rietveld method, irrespective of the software used, the observed Bragg peaks in a powder diffraction pattern are best described by the so-called peak shape function (PSF). The PSF is a convolution of three functions: the instrumental broadening ${\displaystyle \Omega (\theta )}$, wavelength dispersion ${\displaystyle \Lambda (\theta )}$, and the specimen function ${\displaystyle \Psi (\theta )}$, with the addition of a background function, ${\displaystyle b(\theta )}$. It is represented as follows:

${\displaystyle PSF(\theta )=\Omega (\theta )\otimes \Lambda (\theta )\otimes \Psi (\theta )+b(\theta )}$,

where ${\displaystyle \otimes }$ denotes convolution, which is defined for two functions ${\displaystyle f}$ and ${\displaystyle g}$ as an integral:

${\displaystyle f(t)\otimes g(t)=\int _{-\infty }^{\infty }f(\tau )g(t-\tau )d\tau =\int _{-\infty }^{\infty }g(\tau )f(t-\tau )d\tau }$

The instrumental function depends on the location and geometry of the source, monochromator, and sample. Wavelength function accounts for the distribution of the wavelengths in the source, and varies with the nature of the source and monochromatizing technique. The specimen function depends on several things. First is dynamic scattering, and secondly the physical properties of the sample such as crystallite size, and microstrain.

A short aside: unlike the other contributions, those from the specimen function can be interesting in materials characterization. As such, average crystallite size, ${\displaystyle \tau }$, and microstrain, ${\displaystyle \varepsilon }$, effects on Bragg peak broadening, ${\displaystyle \beta }$ (in radians), can be described as follows, where ${\displaystyle k}$ is a constant:

${\displaystyle \beta ={\frac {\lambda }{\tau \cdot \cos \theta }}}$ and ${\displaystyle \beta =\kappa \cdot \varepsilon \cdot \tan \theta }$.

Returning to the peak shape function, the goal is to correctly model the Bragg peaks which exist in the observed powder diffraction data. In the most general form, the intensity, ${\displaystyle Y(i)}$, of the ${\displaystyle i^{\text{th}}}$ point (${\displaystyle 1\leq i\leq n}$, where ${\displaystyle n}$ is the number of measured points) is the sum of the contributions ${\displaystyle y_{k}}$ from the ${\displaystyle m}$ overlapped Bragg peaks (${\displaystyle 1\leq k\leq m}$), and the background, ${\displaystyle b(i)}$, and can be described as follows:

${\displaystyle Y(i)=b(i)+\sum _{k=1}^{m}{I_{k}[y_{k}(x_{k})]}}$

where ${\displaystyle I_{k}}$ is the intensity of the ${\displaystyle k^{\text{th}}}$ Bragg peak, and ${\displaystyle x_{i}=2\theta _{i}-2\theta _{k}}$. Since ${\displaystyle I_{k}}$ is a multiplier, it is possible to analyze the behaviour of different normalized peak functions ${\displaystyle y(x)}$ independently of peak intensity, under the condition that the integral over infinity of the PSF is unity. There are various functions that can be chosen to do this with varying degrees of complexity. The most basic functions used in this way to represent Bragg reflections are the Gauss, and Lorentzian functions. Most commonly though, is the pseudo-Voigt function, a weighted sum of the former two (the full Voigt profile is a convolution of the two, but is computationally more demanding). The pseudo-Voigt profile is the most common and is the basis for most other PSF's. The pseudo-Voigt function can be represented as:

${\displaystyle y(x)=V_{p}(x)=n*G(x)+(1-n)*L(x)}$,

where

${\displaystyle G(x)={\frac {C_{G}^{\frac {1}{2}}}{{\sqrt {\pi }}H}}e^{-C_{G}x^{2}}}$

and

${\displaystyle L(x)={\frac {C_{L}^{\frac {1}{2}}}{{\sqrt {\pi }}H'}}\left(1+C_{L}x^{2}\right)^{-1}}$

are the Gaussian and Lorentzian contributions, respectively.

Thus,

${\displaystyle V_{p}(x)=\eta {\frac {C_{G}^{\frac {1}{2}}}{{\sqrt {\pi }}H}}e^{-C_{G}x^{2}}+(1-\eta ){\frac {C_{L}^{\frac {1}{2}}}{{\sqrt {\pi }}H'}}(1+C_{L}x^{2})^{-1}.}$

where:

• ${\displaystyle H}$ and ${\displaystyle H'}$ are the full widths at half maximum (FWHM)
• ${\displaystyle x={\frac {2\theta _{i}-2\theta _{k}}{H_{k}}}}$ is essentially the Bragg angle of the ${\displaystyle i^{\text{th}}}$point in the powder pattern with its origin in the position of the ${\displaystyle k^{\text{th}}}$peak divided by the peak's FWHM.
• ${\displaystyle C_{G}=4\ln 2}$, ${\displaystyle C_{L}=4}$ and ${\textstyle {\frac {C_{G}^{\frac {1}{2}}}{{\sqrt {\pi }}H}}}$ and ${\textstyle {\frac {C_{L}^{\frac {1}{2}}}{{\sqrt {\pi }}H'}}}$ are normalization factors such that ${\textstyle \int _{-\infty }^{\infty }G(x)dx=1}$ and ${\textstyle \int _{-\infty }^{\infty }L(x)dx=1}$ respectively.
• ${\displaystyle H^{2}=U\tan ^{2}\theta +V\tan \theta +W}$, known as the Caglioti formula, is the FWHM as a function of ${\displaystyle \theta }$ for Gauss, and pseudo-Voigt profiles. ${\displaystyle U}$, ${\displaystyle V}$, and ${\displaystyle W}$ are free parameters.
• ${\textstyle H'={\frac {X}{\cos \theta }}+Y\tan \theta }$ is the FWHM vs. ${\displaystyle 2\theta }$ for the Lorentz function. ${\displaystyle X}$ and ${\displaystyle Y}$ are free variables
• ${\displaystyle \eta =\eta _{0}+\eta _{1}2\theta +\eta _{2}\theta ^{2}}$, where ${\displaystyle 0\leq \eta \leq 1}$ is the pseudo-Voigt mixing parameter, and ${\displaystyle \eta _{0,1,2}}$ are free variables.

The pseudo-Voigt function, like the Gaussian and Lorentz functions, is a centrosymmetric function, and as such does not model asymmetry. This can be problematic for non-ideal powder XRD data, such as those collected at synchrotron radiation sources, which generally exhibit asymmetry due to the use of multiple focusing optics.

The Finger–Cox–Jephcoat function is similar to the pseudo-Voigt, but has better handling of asymmetry, which is treated in terms of axial divergence. The function is a convolution of pseudo-Voigt with the intersection of the diffraction cone and a finite receiving slit length using two geometrical parameters, ${\displaystyle S/L}$, and ${\displaystyle H/L}$, where ${\displaystyle S}$ and ${\displaystyle H}$ are the sample and the detector slit dimensions in the direction parallel to the goniometer axis, and ${\displaystyle L}$ is the goniometer radius.

The shape of a powder diffraction reflection is influenced by the characteristics of the beam, the experimental arrangement, and the sample size and shape. In the case of monochromatic neutron sources the convolution of the various effects has been found to result in a reflex almost exactly Gaussian in shape. If this distribution is assumed then the contribution of a given reflection to the profile ${\displaystyle y_{i}}$ at position ${\displaystyle 2\theta _{i}}$ is:

${\displaystyle y_{i}=I_{k}\exp \left[{\frac {-4\ln \left(2\right)}{H_{k}^{2}}}\left(2\theta _{i}-2\theta _{k}\right)^{2}\right]}$

where ${\displaystyle H_{k}}$ is the full width at half peak height (full-width half-maximum), ${\displaystyle 2\theta _{k}}$ is the center of the reflex, and ${\displaystyle I_{k}}$ is the calculated intensity of the reflex (determined from the structure factor, the Lorentz factor, and multiplicity of the reflection).

At very low diffraction angles the reflections may acquire an asymmetry due to the vertical divergence of the beam. Rietveld used a semi-empirical correction factor, ${\displaystyle A_{s}}$ to account for this asymmetry:

${\displaystyle A_{s}=1-\left[{\frac {P\left(2\theta _{i}-2\theta _{k}\right)^{2}}{\tan \theta _{k}}}\right]}$

where ${\displaystyle P}$ is the asymmetry factor and ${\displaystyle s}$ is ${\displaystyle +1}$, ${\displaystyle 0}$, or ${\displaystyle -1}$ depending on the difference ${\displaystyle 2\theta _{i}-2\theta _{k}}$ being positive, zero, or negative respectively.

At a given position more than one diffraction peak may contribute to the profile. The intensity is simply the sum of all reflections contributing at the point ${\displaystyle 2\theta _{i}}$.

For a Bragg peak ${\displaystyle (hkl)}$, the observed integrated intensity, ${\displaystyle I_{hkl}}$, as determined from numerical integration is

${\displaystyle I_{hkl}=\sum _{i=1}^{j}(Y_{i}^{obs}-b_{i})}$,

where ${\displaystyle j}$ is the total number of data points in the range of the Bragg peak. The integrated intensity depends on multiple factors, and can be expressed as the following product:

${\displaystyle I_{hkl}=K\times p_{hkl}\times L_{\theta }\times P_{\theta }\times A_{\theta }\times T_{hkl}\times E_{hkl}\times |F_{hkl}|^{2}}$

where:

• ${\displaystyle K}$: scale factor
• ${\displaystyle p_{hkl}}$: multiplicity factor, which accounts for symmetrically equivalent points in the reciprocal lattice
• ${\displaystyle L_{\theta }}$: Lorentz multiplier, defined by diffraction geometry
• ${\displaystyle P_{\theta }}$: polarization factor
• ${\displaystyle A_{\theta }}$: absorption multiplier
• ${\displaystyle T_{hkl}}$: preferred orientation factor
• ${\displaystyle E_{hkl}}$: extinction factor (often neglected as it is usually insignificant in powders)
• ${\displaystyle F_{hkl}}$: structure factor as determined by the crystal structure of the material

The width of the diffraction peaks are found to broaden at higher Bragg angles. This angular dependency was originally represented by

${\displaystyle H_{k}^{2}=U\tan ^{2}\theta _{k}+V\tan \theta _{k}+W}$

where ${\displaystyle U}$, ${\displaystyle V}$, and ${\displaystyle W}$ are the half-width parameters and may be refined during the fit.

The fitting method used in Rietveld refinement is the non-linear least squares approach. A detailed derivation of non-linear least squares fitting will not be given here. Further detail can be found in Chapter 6 of Pecharsky and Zavalij's text 12 . There are a few things to note however. First, non-linear least squares fitting has an iterative nature for which convergence may be difficult to achieve if the initial approximation is too far from correct, or when the minimized function is poorly defined. The latter occurs when correlated parameters are being refined at the same time, which may result in divergence and instability of the minimization. This iterative nature also means that convergence to a solution does not occur immediately for the method is not exact. Each iteration depends on the results of the last which dictate the new set of parameters used for refinement. Thus, multiple refinement iterations are required to eventually converge to a possible solution.

Using non-linear least squares minimization, the following system is solved:

${\displaystyle {\begin{pmatrix}Y_{i}^{\text{calc}}=kY_{i}^{\text{obs}}\\\vdots \\Y_{n}^{\text{calc}}=kY_{n}^{\text{obs}}\end{pmatrix}}}$

where ${\displaystyle Y_{i}^{\text{calc}}}$is the calculated intensity and ${\displaystyle Y_{i}^{\text{obs}}}$is the observed intensity of a point ${\displaystyle i}$ in the powder pattern, ${\displaystyle k}$, is a scale factor, and ${\displaystyle n}$ is the number of measured data points. The minimized function is given by:

${\displaystyle \Phi =\sum _{i=1}^{n}w_{i}(Y_{i}^{\text{obs}}-Y_{i}^{\text{calc}})^{2}}$

where ${\displaystyle w_{i}}$is the weight, and ${\displaystyle k}$from the previous equation is unity (since ${\displaystyle k}$ is usually absorbed in the phase scale factor). The summation extends to all ${\displaystyle n}$ data points. Considering the peak shape functions and accounting for the overlapping of Bragg peaks because of the one-dimensionality of XRD data, the expanded form of the above equation for the case of a single phase measured with a single wavelength becomes:

${\displaystyle \Phi =\sum _{i=1}^{n}w_{i}{\biggl (}Y_{i}^{\text{obs}}-{\Bigl (}b_{i}+K\sum _{j=1}^{m}I_{j}y_{j}(x_{j}){\Bigr )}{\biggr )}^{2}}$

where:

• ${\displaystyle b_{i}}$ is the background at the ${\displaystyle i^{\text{th}}}$ data point.
• ${\displaystyle K}$ is the phase scale factor.
• ${\displaystyle m}$ is the number of Bragg reflections contributing to the intensity of the ${\displaystyle i^{\text{th}}}$ reflection.
• ${\displaystyle I_{j}}$ is the integrated intensity of the ${\displaystyle j^{\text{th}}}$ Bragg peak.
• ${\displaystyle y_{i}(x_{i})}$ is the peak shape function.

For a material that contains several phases (${\displaystyle p}$), the contribution from each is accounted for by modifying the above equation as follows:

${\displaystyle \Phi =\sum _{i=1}^{n}w_{i}{\biggl (}Y_{i}^{\text{obs}}-{\Bigl (}b_{i}+\sum _{l=1}^{p}K_{l}\sum _{j=1}^{m}I_{l,j}y_{l,j}(x_{l,j}){\Bigr )}{\biggr )}^{2}}$

It can easily be seen from the above equations that experimentally minimizing the background, which holds no useful structural information, is paramount for a successful profile fitting. For a low background, the functions are defined by contributions from the integrated intensities and peak shape parameters. But with a high background, the function being minimized depends on the adequacy of the background and not integrated intensities or peak shapes. Thus, a structure refinement cannot adequately yield structural information in the presence of a large background.

It is also worth noting the increased complexity brought forth by the presence of multiple phases. Each additional phase adds to the fitting, more Bragg peaks, and another scale factor tied to corresponding structural parameters, and peak shape. Mathematically they are easily accounted for, but practically, due to the finite accuracy and limited resolution of experimental data, each new phase can lower the quality and stability of the refinement. It is advantageous to use single phase materials when interested in finding precise structural parameters of a material. However, since the scale factors of each phase are determined independently, Rietveld refinement of multi phase materials can quantitatively examine the mixing ratio of each phase in the material.

Generally, the background is calculated as a Chebyshev polynomial. In GSAS and GSAS-II they appear as follows. Again, background is treated as a Chebyshev polynomial of the first kind ("Handbook of Mathematical Functions", M. Abramowitz and IA. Stegun, Ch. 22), with intensity given by:

${\displaystyle I_{i}=\sum _{j=1}^{N}P_{j}T'_{j-1}}$

where ${\displaystyle T'_{j-1}}$are the coefficients of the Chebyshev polynomial taken from Table 22.3, pg. 795 of the Handbook. The coefficients have the form:

${\displaystyle T'_{n}=\sum _{m=0}^{i-1}C_{m}X^{m}}$

and the values for ${\displaystyle C_{m}}$are found in the Handbook. The angular range (${\displaystyle 2\theta }$) is converted to ${\displaystyle X}$ to make the Chebyshev polynomial orthogonal by

${\displaystyle X={\frac {2}{\tau }}-1}$

And, the orthogonal range for this function is –1 to +1.

Now, given the considerations of background, peak shape functions, integrated intensity, and non-linear least squares minimization, the parameters used in the Rietveld refinement which put these things together can be introduced. Below are the groups of independent least squares parameters generally refined in a Rietveld refinement.

• Background parameters: usually 1 to 12 parameters.
• Sample displacement: sample transparency, and zero shift corrections. (move peak position)
• Multiple peak shape parameters.
  • FWHM parameters: i.e. Caglioti parameters (see section 3.1.2)
  • Asymmetry parameters (FCJ parameters)
• Unit cell dimensions
  • one to six parameters (a, b, c, α, β, γ), depending on the crystal family/system, for each present phase.
• Preferred orientation, and sometimes absorption, porosity, and extinction coefficients, which can be independent for each phase.
• Scale factors (for each phase)
• Positional parameters of all independent atoms in the crystal model (generally 0 to 3 per atom).
• Population parameters
  • Occupation of site positions by atoms.
• Atomic displacement parameters
  • Isotropic and anisotropic (temperature) parameters.

Each Rietveld refinement is unique and there is no prescribed sequence of parameters to include in a refinement. It is up to the user to determine and find the best sequence of parameters for refinement. It is worth noting that it is rarely possible to refine all relevant variables simultaneously from the beginning of a refinement, nor near the end since the least squares fitting will be destabilized or lead to a false minimum. It is important for the user to determine a stopping point for a given refinement. Given the complexity of Rietveld refinement it is important to have a clear grasp of the system being studied (sample, and instrumentation) to ensure that results are accurate, realistic, and meaningful. High data quality, a large enough ${\displaystyle 2\theta }$ range, and a good model – to serve as the initial approximation in the least squares fitting – are necessary for a successful, reliable, and meaningful Rietveld refinement.