ChemWiki - User contributions [en]

Rep:Mod:MgOben14

2017-03-23T09:31:17Z

Ben14: /* Experimental */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon which is discussed and presented here through the example of magnesium oxide, MgO. One such example is the the construction of any large structure where changes in dimension can occur due to temperature variations, engineers must account for this thermal expansion in its design and operation, though due to the often very large size of these structures, only a crude approximation is required. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is fundamental in enabling the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Helmholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell had been selected as an optimised function of accuracy and computational time. As the Helmholtz Free Energy is a property dependent on the temperature and volume, this then enabled the prediction of how the material expands when heated, through the computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic (or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available<ref name=":4">C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>.

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models used and for the running of simulations are listed below. The Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the molecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively known as the Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based on the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> A purely harmonic model (as modelled by a parabola) is unable to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur). Therefore, the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds, hence termed ''quasi-harmonic''.<ref name=":1" /> It should be noted that the model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed in the output of the calculations. For the purpose of these calculations, this is disregarded though will contribute to the error in the calculations, when compared to literature.

The Helmholtz Free Energy in the Quasi-Harmonic Model is given by Equation 2 below, (''where F is Helmholtz free energy, E is the zero potential energy, T is temperature (K) '').
[[File:BN Helmholtz free energy equation.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy ]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the ''Results and Discussion'' section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties. Its accuracy and applicability for the system under investigation is evaluated in the Results and Discussion section below.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed in the Figures below.

The primitive cell is the smallest cell that can be used to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref name=":4" /> It should be noted that the densities of these cells do not differ, the volume of the convensional cell is simply four times that of the primitive cell.
[[File:BN Convensional and primitive cell inverted colour.jpg|none|thumb|718x718px|Figure 1 to show:
Left - primitive cell Right- conventional cell
]]
==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in Figure 2 below. The aim was to investigate how vibration affected expansion. The figure shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter, and a* is the magnitude of the Brillouin Zone (the unit cell of the reciprocal lattice). A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|437x437px|Figure 2 to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, Figure 3 below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|434x434px|Figure 3 to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the subsequent section, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken arbitrarily by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) showing only a translation but leaving the model unaffected) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels (DOS) was studied to determine the required ''basis'' for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out. This was then evaluated against the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions in the x, y and z plane (1x1x1..2x2x2...64x64x64) yielded systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. Thus ultimately the computational time was not a significant factor in this optimisation process, each submitted calculation ran well within an acceptable period for the timescale of the lab (maximum 5 minutes for the longest). From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure 4 to show DOS Graphs for MgO using various Grid Sizes ]]To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in Table 1 below. With increasing grid dimension, there is an incremental decrease in the Helmholtz free energy, smaller at each step until the value is optimised.
{| class="wikitable"
|+Table 1 to show Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting (for a less symmetric system, for example).

It is interesting to examine also how this optimisation process would have to be altered for other materials.Given three contrasting examples: an oxide similar to MgO (CaO), a zeolite (faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin zone, a* as shown in Equation 3 below.
[[File:BN Brillouin zone.PNG|none|thumb|Equation 3 to show conversion to reciprocal space]]
The literature lattice constants for the materials discussed above are presented in Table 2 below, where an approximate relative optimised grid size is given qualitatively, followed by a rationale in the discussion below.<ref>H. Zhang and M. S. T. Bukowinski, Phys. Rev. B, 1991, 44, 2495–2503.</ref><ref>O. Madelung, U. Rössler and M. Schulz, Eds., in II-VI and I-VII Compounds; Semimagnetic Compounds, Springer Berlin Heidelberg, Berlin, Heidelberg, 1999, pp. 1–3.</ref><ref>Clays and Clay Minerals, 197,3, Vol. 21, pp. 379-389</ref><ref>S. H. Kellington, D. Loveridge, J. M. Titman, H. F. H. and A. B. L, K. G. and M. A, N. S. K, T. J. M. and K. S. H and V. H, J. Phys. D. Appl. Phys., 1969.</ref>
{| class="wikitable"
|+Table 2: Lattice Constants, with Comparison to MgO
!Material
!Lattice Constant, Å
!Ideal Grid Size (relative to MgO)
|-
|CaO
|4.800
|similar to MgO
|-
|Faujasite
|24.638-24.650
|smaller than MgO
|-
|Lithium
|3.5107
|smaller than MgO
|}
As the CaO is of a similar sized lattice to MgO, a similar number of k points (and therefore grid dimensions) would be needed to provide the optimum output. The zeolite, faujasite, has a lattice constant much greater than MgO, and due to the inverse relationship discussed above, a much narrower Brillouin zone. This leads to a smaller grid size being sufficient for the sampling of the necessary number of k points for an accurate model.

Initially, it may seem that due to the smaller lattice constant of Li, and the reciprocal relationship, the grid size should be larger than for MgO (which would lead to long computational times). However, this is in fact not the case due to the difference in bonding present in the lithium metal (ionic) compared to magnesium oxide (covalent). The ionic nature of the Li bonds results in the presence of delocalized electrons, which in turn leads to less repulsion between the metal cations.<ref name=":2">C. Kittel, Solid State Phys., 2005, 703.</ref> Thus a narrower energy band (a less disperse Density of States function) would be observed, and accordingly a smaller grid size would be adequate.

=== The Thermal Energy of MgO ===

==== The Quasi-Harmonic Approximation ====
The results of the optimisation with respect to the free energy computed within the quasi-harmonic approximation at a number of temperatures is presented and discussed in this section. This is quantified through the use of the Helmholtz free energy, as shown in Equation 4 below, with an entropic term that sums over all k-values so is a useful indication of accuracy.
[[File:BN Helmholtz free energy equation.JPG|none|thumb|689x689px|Equation 4 to show the calculation of the Helmholtz Free energy in the Quasi-Harmonic Model]]

The charts below show the plot of the Helmholtz free energy (Figure 5 a)) as well as the lattice constant (Figure 5 b)) and cell volume against temperature (Figure 5 c)) for the temperature range 0 -1000 K in 100 K increments. [[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1091x1091px|Figure 5 (From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume of MgO against Temperature]]The trend in the variation of lattice constant is as expected: the greater the the temperature, the greater the thermal energy and therefore there is a larger average distance between the atoms. This is directly responsible for the same trend observed in the cell volume.

Within the temperature range 300-1000 K there is a linear relationship with temperature. A linear regression can be applied to this region, with the gradient of the curve showing the change in volume with respect to the change in temperature (at constant pressure, it is a solid and thus this can be assumed). Employing Equation 1 (see ''Introduction)'' for the thermal expansion coefficient, a value of 2.65442E-05 K-1 is obtained. This can be compared as to the literature value of 3.16E-05 K-1 (for the temperature range 303 - 1273 K) to produce an error of 16%. This error can be accounted to the approximations made in the Quasi-Harmonic model, however this is within an acceptable range. .

The linear trend stated above however does not accurately model the low temperature regime, where the thermal expansion shows deviation from linearity due to a degree of temperature dependence. The Helmholtz free energy computed at 0 K is denoted as the ''Zero Point Energy, ''which is effectively a residual energy arising from quantum mechanics''. ''

This can be explained through examination of the underlying principles in the quasi-harmonic approximation model. As introduced in the ''Experimental'' section of the report, this is a quantum-mechanical based model that builds on the harmonic oscillator, however introduces anharmonicity to allow thermal expansion. This cannot happen in a purely harmonic system as the bond length cannot vary with the population of thermally excited states. The quasi-harmonic approach allows for the variation of bond length, and thus thermal expansion with the population of higher energy states as a result of greater thermal energy. Figure 6, showing the harmonic and anharmonic oscillator, illustrates this further.
[[File:BN Harmonic anharmonic.png|none|thumb|439x439px|Figure 6 to show the Harmonic and Anharmonic Oscillator ]]In the upper limit of 1000 K the Quasi-harmonic approximation continues to be valid, due to the melting point of MgO being much greater than this, at 3125 K<ref name=":2" />. At that temperature, a change in the free energy trend would be expected, since it would no longer be a crystal structure, with vibrations but rather a much more disordered system. The Molecular Dynamics section below discusses this high temperature regime further.

The Helmholtz free energy decreases with increasing temperature and this can simply be accounted for using the thermodynamic definition of Helmholtz Free Energy, provided in the equation below (Equation 5).
[[File:Hemholtz thermodynamic.JPG|none|thumb|231x231px|Equation 5 to show the Thermodynamic Free Energy Equation]]

As can be seen from the equation, with an increase in temperature, the -TS term becomes increasingly negative, in addition to the increase in entropy due to the greater thermal disordering. The non-linearity once again arises from the quasi-harmonic approximation.

==== Molecular Dynamics Calculations ====
Following the quasi-harmonic calculations, an alternative approach was carried out for comparison. The molecular dynamics calculations model vibration as random motion, based on classical mechanics. The graph presented below in Figure 7 shows this, alongside a comparison to the LD approach. As can be seen from the plots, the MD approach is not accurate in the low temperature regime, though appears to converge to the quasi-harmonic approximation in the high temperature limit of this 300- 1000 K region. The discussion below provides rationale for the observed results and compares the thermal expansion coefficient, α obtained with literature. [[File:BN Graph to compare LD with MD.JPG|none|thumb|754x754px|Figure 7 to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]The MD approach yields a thermal expansion coefficient of 3.1853E-05 K-1 (over the same 300-1000 K temperature range as above), which is a variation of less than 1% from the literature value stated above, thus in excellent accordance. The improved accuracy of the molecular dynamics model in this 300- 1000 K temperature range can be attributed to the fact that at these higher temperatures random motion due to thermal excitation is increasingly significant.

At temperatures below this, the model becomes increasingly inaccurate because there is less thermal motion. At 0 K, it is unable to provide an approximation - the thermal motion would be zero, rather than the zero-point energy that the LD model correctly predicts, showing the limits of the approach. Harmonicity and a quantum mechanical treatment is required here: due to intrinsic limitations imposed on the material, as expressed by the Heisenberg Uncertainty Principle, the molecular dynamics approach becomes invalid.

In the high temperature limit, the model would continue to work well, because in the melted state the system is very well described by random thermal motion. The model would predict a persistent thermal expansion of the volume, as the material changes from solid to liquid to gas (provided the thermal expansion coefficient were positive, which is not always the case, such as for pure water below 3.98 °C).<ref name=":3">K. Röttger, A. Endriss, J. Ihringer, S. Doyle and W. F. Kuhs, Acta Crystallogr. Sect. B Struct. Sci., 2012, 68, 91–91.</ref>

The MD calculations made use of 32 MgO molecules (i.e. 64 atoms) which is much less than that for the Quasi-Harmonic calculations, but the greater increase in computational time required to run these calculations justifies this.

== Conclusion ==
It can be seen from the discussion presented in the sections above that the two approaches (Quasi-Harmonic and Molecular Dynamics) implemented for the modeling of the thermal expansion of MgO are each more suited to the low and high temperature regions, respectively. To compare, a linear fit was made on the 300 - 1000 K region of both plots, to yield thermal expansion coefficients of 2.65442E-05 K-1 for the Quasi-Harmonic Calculations (16% error) and 3.1853E-05 K-1 (0.8% error) for the Molecular Dynamics calculations. The methodology utilised in each step was rationalised through the an explanation of the relevant theory and by considering the computational calculation process to identify where shortcuts (such as reduction in calculations due to symmetry) were possible. The use of other materials for such calculations were also explored, reaching the conclusion that both metals and large structures such as zeolites would suffice with smaller grid dimensions (and hence shorter computational times) than MgO, although for distinct reasons. To build further on the computational work, it would be of interest to model other structures, such as ice perhaps, which exhibit ''negative thermal expansion (NTE)'' in both very low (-200 °C) temperatures and also in its liquid form, when below 3.98 °C.<ref name=":3" />The cell size required for molecular dynamics could also be optimised for further work, as was done for the quasi-harmonic model.

== References ==
<references />

Rep:Mod:MgOben14

2017-03-23T09:05:39Z

Ben14:

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon which is discussed and presented here through the example of magnesium oxide, MgO. One such example is the the construction of any large structure where changes in dimension can occur due to temperature variations, engineers must account for this thermal expansion in its design and operation, though due to the often very large size of these structures, only a crude approximation is required. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is fundamental in enabling the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Helmholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell had been selected as an optimised function of accuracy and computational time. As the Helmholtz Free Energy is a property dependent on the temperature and volume, this then enabled the prediction of how the material expands when heated, through the computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic (or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available<ref name=":4">C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>.

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models used and for the running of simulations are listed below. The Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the molecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively known as the Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based on the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> A purely harmonic model (as modelled by a parabola) is unable to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur). Therefore, the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds, hence ''quasi-harmonic''.<ref name=":1" /> It should be noted that the model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed in the output of the calculations.

The Helmholtz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholtz free energy, E is the zero potential energy, T is temperature (K) '').
[[File:BN Helmholtz free energy equation.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy ]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the ''Results and Discussion'' section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed in the Figures below.

The primitive cell is the smallest cell that can be used to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref name=":4" />
[[File:BN Convensional and primitive cell inverted colour.jpg|none|thumb|718x718px|Figure 1 to show:
Left - primitive cell Right- conventional cell
]]
==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter, and a* is the magnitude of the Brillouin Zone. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|437x437px|Figure 2 to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|434x434px|Figure 3 to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the subsequent section, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken arbitrarily by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) showing only a translation but leaving the model unaffected) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels (DOS) was studied to determine the required ''basis'' for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out. This was then evaluated against the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1..2x2x2...64x64x64) yielded systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. Thus the computational time was not a significant factor in this optimisation process, each submitted calculation ran well within an acceptable period for the timescale of the lab (maximum 5 minutes for the longest). From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure 4 to show DOS Graphs for MgO using various Grid Sizes ]]To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below. With increasing grid dimension, there is an incremental decrease in the Helmholtz free energy, smaller at each step until the value is optimised.
{| class="wikitable"
|+Table 1 to show Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

It is interesting to examine how other materials would alter this optimisation process.Given three contrasting examples: an oxide similar to MgO (CaO), a zeolite (faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin zone, a* as shown in the equation below.
[[File:BN Brillouin zone.PNG|none|thumb|Equation 3 to show conversion to reciprocal space]]
The literature lattice constants for the materials discussed above are presented in the table below, where an approximate relative optimised grid size is given qualitatively, followed by a chemical rationale in the discussion below.<ref>H. Zhang and M. S. T. Bukowinski, Phys. Rev. B, 1991, 44, 2495–2503.</ref><ref>O. Madelung, U. Rössler and M. Schulz, Eds., in II-VI and I-VII Compounds; Semimagnetic Compounds, Springer Berlin Heidelberg, Berlin, Heidelberg, 1999, pp. 1–3.</ref><ref>Clays and Clay Minerals, 197,3, Vol. 21, pp. 379-389</ref><ref>S. H. Kellington, D. Loveridge, J. M. Titman, H. F. H. and A. B. L, K. G. and M. A, N. S. K, T. J. M. and K. S. H and V. H, J. Phys. D. Appl. Phys., 1969.</ref>
{| class="wikitable"
|+Table 2: Lattice Constants, with Comparison to MgO
!Material
!Lattice Constant, Å
!Ideal Grid Size (relative to MgO)
|-
|CaO
|4.800
|similar to MgO
|-
|Faujasite
|24.638-24.650
|smaller than MgO
|-
|Lithium
|3.5107
|smaller than MgO
|}
As the CaO is of a similar sized lattice to MgO, a similar number of k points (and therefore grid dimensions) would be needed to provide the optimum output. The zeolite, faujasite, has a lattice constant much greater than MgO, and due to the inverse relationship discussed above, a much narrower Brillouin zone. This leads to a smaller grid size being sufficient for the sampling of the necessary number of k points for an accurate model.

Initially, it may seem that due to the smaller lattice constant of Li, and the reciprocal relationship, the grid size should be larger than for MgO (which would lead to long computational times). However, this is in fact not the case due to the difference in bonding present in the lithium metal (ionic) compared to magnesium oxide (covalent). The ionic nature of the Li bonds results in the presence of delocalized electrons, which in turn leads to less repulsion between the metal cations.<ref name=":2">C. Kittel, Solid State Phys., 2005, 703.</ref> Thus a narrower energy band (a less disperse Density of States function) would be observed, and accordingly a smaller grid size would be adequate.

=== The Thermal Energy of MgO ===

==== The Quasi-Harmonic Approximation ====
The results of the optimisation with respect to the free energy computed within the quasi-harmonic approximation at a number of temperatures is presented and discussed in this section. This is quantified through the use of the Helmholtz free energy, as shown in the equation below, with an entropic term that sums over all k-values so is a useful indication of accuracy.
[[File:BN Helmholtz free energy equation.JPG|none|thumb|689x689px|Equation 4 to show the calculation of the Helmholtz Free energy in the Quasi-Harmonic Model]]

The charts below show the plot of the Helmholtz free energy (Figure a)) as well as the lattice constant (Figure b)) and cell volume against temperature (Figure c)) for the temperature range 0 -1000 K in 100 K increments. [[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1091x1091px|Figure 5 (From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]The trend in the variation of lattice constant is as expected: the greater the the temperature, the greater the thermal energy and therefore there is a larger average distance between the atoms. This is directly responsible for the same trend observed in the cell volume.

Within the temperature range 300-1000 K there is a linear relationship with temperature. A linear regression can be applied to this region, with the gradient of the curve showing the change in volume with respect to the change in temperature (at constant pressure, it is a solid and thus this can be assumed). Employing Equation 1 (see ''Introduction)'' for the thermal expansion coefficient, a value of 2.65442E-05 K-1 is obtained. This can be compared as to the literature value of 3.16E-05 K-1 (for the temperature range 303 - 1273 K) to produce an error of 16%. This error can be accounted to the approximations made in the Quasi-Harmonic model, however this is within an acceptable range. .

The linear trend stated above however does not accurately model the low temperature regime, where the thermal expansion shows deviation from linearity due to a degree of temperature dependence. The Helmholtz free energy computed at 0 K is denoted as the ''Zero Point Energy, ''which is effectively a residual energy arising from quantum mechanics''. ''

This can be explained through examination of the underlying principles in the Quasi-Harmonic Approximation model. As introduced in the ''Experimental'' section of the report, this is a quantum-mechanical based model that builds on the harmonic oscillator, however introduces anharmonicity to allow thermal expansion. This cannot happen in a purely harmonic system as the bond length cannot vary with the population of thermally excited states. The quasi-harmonic approach allows for the variation of bond length, and thus thermal expansion with the population of higher energy states as a result of greater thermal energy. The figure below, showing the harmonic and anharmonic oscillator illustrates this further.
[[File:BN Harmonic anharmonic.png|none|thumb|439x439px|Figure 6 to show the Harmonic and Anharmonic Oscillator ]]In the upper limit of 1000 K the Quasi-harmonic approximation continues to be valid, due to the melting point of MgO being much greater than this, at 3125 K<ref name=":2" />. At that temperature, a change in the free energy trend would be expected, since it would no longer be a crystal structure, with vibrations but rather a much more disordered system. The Molecular Dynamics section below discusses this high temperature regime further.

The Helmholtz free energy decreases with increasing temperature and this can simply be accounted for using the thermodynamic definition of Helmholtz Free Energy, provided in the equation below.
[[File:Hemholtz thermodynamic.JPG|none|thumb|231x231px|Equation 5 to show the Thermodynamic Free Energy Equation]]

As can be seen from the equation, with an increase in temperature, the -TS term becomes increasingly negative, in addition to the increase in entropy due to the greater thermal disordering. The non-linearity once again arises from the quasi-harmonic approximation.

==== Molecular Dynamics Calculations ====
Following the Quasi-Harmonic calculations, an alternative approach was carried out for comparison. The Molecular Dynamics calculations model vibration as random motion, based on classical mechanics. The graph presented below in the Figure shows this, alongside a comparison to the LD approach. As can be seen from the plots, the MD approach is not accurate in the low temperature regime, though appears to converge to the Quasi-Harmonic Approximation in the high temperature limit of this 300- 1000 K region. The discussion below provides rationale for the observed results and compares the thermal expansion coefficient, α obtained with literature. [[File:BN Graph to compare LD with MD.JPG|none|thumb|754x754px|Figure 7 to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]The MD approach yields a thermal expansion coefficient of 3.1853E-05 K-1 (over the same 300-1000 K temperature range), which is a variation of less than 1% from the literature value stated above, thus in excellent accordance. The improved accuracy of the molecular dynamics model in this 300- 1000 K temperature range can be attributed to the fact that at these higher temperatures random motion due to thermal excitation is increasingly significant.

At temperatures below this, the model becomes increasingly inaccurate because there is less thermal motion. At 0 K, this model is unable to provide an approximation - the thermal motion would be zero, rather than the zero-point energy that the LD model correctly predicts, showing the limits of the approach. Harmonicity and a quantum mechanical treatment is required here: due to intrinsic limitations imposed on the material, as expressed by the Heisenberg Uncertainty Principle, the molecular dynamics approach becomes invalid.

In the high temperature limit, the model would continue to work well, in the melted state, the system is very well described by random thermal motion. The model would predict a persistent thermal expansion of the volume, as the material changes from solid to liquid to gas (provided the thermal expansion coefficient were positive, which is not always the case, such as for pure water below 3.98 °C).<ref name=":3">K. Röttger, A. Endriss, J. Ihringer, S. Doyle and W. F. Kuhs, Acta Crystallogr. Sect. B Struct. Sci., 2012, 68, 91–91.</ref>

The MD calculations made use of 32 MgO molecules (i.e. 64 atoms) which is much less than that for the Quasi-Harmonic calculations, but the greater increase in computational time required to run these calculations justifies this.

== Conclusion ==
It can be seen from the discussion presented in the sections above that the two approaches (Quasi-Harmonic and Molecular Dynamics) implemented for the modeling of the thermal expansion of MgO are each more suited to the low and high temperature regions, respectively. To compare, a linear fit was made on the 300 - 1000 K region of both plots, to yield thermal expansion coefficients of 2.65442E-05 K-1 for the Quasi-Harmonic Calculations (16% error) and 3.1853E-05 K-1 (0.8% error) for the Molecular Dynamics calculations. The methodology utilised in each step was rationalised through the an explanation of the relevant theory and by considering the computational calculation process to identify where shortcuts (such as because of symmetry) were possible. The use of other materials for such calculations were also explored, reaching the conclusion that both metals and large structures such as zeolites would suffice with smaller grid dimensions (and hence shorter computational times) than MgO, although for distinct reasons. To build further on the computational work, it would be of interest to model other structures, such as ice perhaps, which exhibits ''negative thermal expansion (NTE)'' in both very low (-200 °C) temperatures and also in its liquid form, below 3.98 °C.<ref name=":3" />

== References ==
<references />

2017-03-22T12:02:35Z

Ben14: /* Molecular Dynamics Calculations */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon which is discussed and presented here through the example of magnesium oxide, MgO. During the the construction of any large structures where changes in dimension can occur due to temperature variations, engineers must account for this thermal expansion in its design and operation, though due to the often very large size of these structures, only a crude approximation is required. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is fundamental in enabling the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Helmholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell had been selected as an optimised function of accuracy and computational time. This then enabled the prediction of how the material expands when heated, through the computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic (or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>.

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models used and for the running of simulations are listed below. Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> The purely harmonic model is unable to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur). Therefore, the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds, hence ''quasi-harmonic''.<ref name=":1" /> It should be noted that the model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed in the output of the calculations.

The Helmholtz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholtz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the ''Results and Discussion'' section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed in the Figures below.

The primitive cell is the smallest cell that can be used to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>
[[File:BN Convensional and primitive cell inverted colour.jpg|none|thumb|718x718px|Figure to show Left: primitive cell Right: conventional cell ]]
==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter, and a* is the magnitude of the Brillouin Zone. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|437x437px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|434x434px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the subsequent section, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken arbitrarily by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) showing only a translation but leaving the model unaffected) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels (DOS) was studied to determine the required ''basis'' for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out. This was then evaluated against the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1..2x2x2...64x64x64) yielded systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. Thus the computational time was not a significant factor in this optimisation process, each submitted calculation ran well within an acceptable period for the timescale of the lab (maximum 5 minutes for the longest). From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below. With increasing grid dimension, there is an incremental decrease in the Helmholtz free energy, smaller at each step until the value is optimised.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

It is interesting to examine how other materials would alter this optimisation process.Given three contrasting examples: an oxide similar to MgO (CaO), a zeolite (faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin zone, a* as shown in the equation below.
[[File:BN Brillouin zone.PNG|none|thumb|Equation to show conversion to reciprocal space]]
The literature lattice constants for the materials discussed above are presented in the table below, where an approximate relative optimised grid size is given qualitatively, followed by a chemical rationale in the discussion below.<ref>Zhang, H., & Bukowinski, M. S. T. (1991). Modified potential-induced-breathing model of potentials between close-shell ions. Phys. Rev. B, 44(6), 2495–2503. doi:10.1103/PhysRevB.44.2495</ref><ref>Madelung, O., Rössler, U., & Schulz, M. (Eds.). (1999). Calcium oxide (CaO) crystal structure, lattice parameters, thermal expansion. In II-VI and I-VII Compounds; Semimagnetic Compounds (pp. 1–3). Berlin, Heidelberg: Springer Berlin Heidelberg. doi:10.1007/10681719_224</ref><ref>Clays and Clay Minerals, 197,3, Vol. 21, pp. 379-389</ref><ref>Kellington, S. H., Loveridge, D., Titman, J. M., L, H. F. H. and A. B., A, K. G. and M., K, N. S., … H, V. (1969). The lattice parameters of some alloys of lithium. Journal of Physics D: Applied Physics. doi:10.1088/0022-3727/2/8/415</ref>
{| class="wikitable"
!Material
!Lattice Constant, Å
!Ideal Grid Size (relative to MgO)
|-
|CaO
|4.800
|similar to MgO
|-
|Faujasite
|24.638-24.650
|smaller than MgO
|-
|Lithium
|3.5107
|smaller than MgO
|}
As the CaO is of a similar sized lattice to MgO, a similar number of k points (and therefore grid dimensions) would be needed to provide the optimum output. The zeolite, faujasite, has a lattice constant much greater than MgO, and due to the inverse relationship discussed above, a much narrower Brillouin zone. This leads to a smaller grid size being sufficient for the sampling of the necessary number of k points for an accurate model.

Initially, it may seem that due to the smaller lattice constant of Li, and the reciprocal relationship, the grid size should be larger than for MgO (which would lead to long computational times). However, this is in fact not the case due to the difference in bonding present in the lithium metal (ionic) compared to magnesium oxide (covalent). The ionic nature of the Li bonds results in the presence of delocalized electrons, which in turn leads to less repulsion between the metal cations.<ref name=":2">C. Kittel, Solid State Phys., 2005, 703.</ref> Thus a narrower energy band (a less disperse Density of States function) would be observed, and accordingly a smaller grid size would be adequate.

=== The Thermal Energy of MgO ===

==== The Quasi-Harmonic Approximation ====
The results of the optimisation with respect to the free energy computed within the quasi-harmonic approximation at a number of temperatures is presented and discussed in this section. This is quantified through the use of the Helmholtz free energy, as shown in the equation below, with an entropic term that sums over all k-values so is a useful indication of accuracy.
[[File:BN Helmholtz free energy equation.JPG|none|thumb|663x663px|Equation to show the calculation of the Helmholtz Free energy in the Quasi-Harmonic Model]]

The charts below show the plot of the Helmholtz free energy (Figure a)) as well as the lattice constant (Figure b)) and cell volume against temperature (Figure c)) for the temperature range 0 -1000 K in 100 K increments. [[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1091x1091px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]The trend in the variation of lattice constant is as expected: the greater the the temperature, the greater the thermal energy and therefore there is a larger average distance between the atoms. This is directly responsible for the same trend observed in the cell volume.

Within the temperature range 300-1000 K there is a linear relationship with temperature. A linear regression can be applied to this region, with the gradient of the curve showing the change in volume with respect to the change in temperature (at constant pressure, it is a solid and thus this can be assumed). Employing Equation 1 (see ''Introduction)'' for the thermal expansion coefficient, a value of 2.65442E-05 is obtained. This can be compared as to the literature value of 3.16E-05 (temperature range 303 - 1273 K) to produce an error of 16%. This error can be accounted to the approximations made in the Quasi-Harmonic model, however this is within an acceptable range. .

The linear trend stated above however does not accurately model the low temperature regime, where the thermal expansion shows deviation from linearity due to a degree of temperature dependence. The Helmholtz free energy computed at 0 K is denoted as the ''Zero Point Energy, ''which is effectively a residual energy arising from quantum mechanics''. ''

This can be explained through examination of the underlying principles in the Quasi-Harmonic Approximation model. As introduced in the ''Experimental'' section of the report, this is a quantum-mechanical based model that builds on the harmonic oscillator, however introduces anharmonicity to allow thermal expansion. This cannot happen in a purely harmonic system as the bond length cannot vary with the population of thermally excited states. The quasi-harmonic approach allows for the variation of bond length, and thus thermal expansion with the population of higher energy states as a result of greater thermal energy. The figure below, showing the harmonic and anharmonic oscillator illustrates this further.
[[File:BN Harmonic anharmonic.png|none|thumb|439x439px|Figure to show the Harmonic and Anharmonic Oscillator ]]In the upper limit of 1000 K the Quasi-harmonic approximation continues to be valid, due to the melting point of MgO being much greater than this, at 3125 K<ref name=":2" />. At that temperature, a change in the free energy trend would be expected, since it would no longer be a crystal structure, with vibrations but rather a much more disordered system. The Molecular Dynamics section below discusses this high temperature regime further.

The Helmholtz free energy decreases with increasing temperature and this can simply be accounted for using the thermodynamic definition of Helmholtz Free Energy, provided in the equation below.
[[File:Hemholtz thermodynamic.JPG|none|thumb|231x231px|Equation to show the Thermodynamic Free Energy Equation]]

As can be seen from the equation, with an increase in temperature, the -TS term becomes increasingly negative, in addition to the increase in entropy due to the greater thermal disordering. The non-linearity once again arises from the quasi-harmonic approximation.

==== Molecular Dynamics Calculations ====
Following the Quasi-Harmonic calculations, an alternative approach was carried out for comparison. The Molecular Dynamics calculations model vibration as random motion, based on classical mechanics. The graph presented below in the Figure shows this, alongside a comparison to the LD approach. As can be seen from the plots, the MD approach is not accurate in the low temperature regime, though appears to converge to the Quasi-Harmonic Approximation in the high temperature limit. The discussion below provides rationale for the observed results and compares the thermal expansion coefficient, α obtained with literature. [[File:BN Graph to compare LD with MD.JPG|none|thumb|754x754px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]The MD approach yields a thermal expansion coefficient of 3.1853E-05, which is

== Conclusion ==

== References ==
<references />

2017-03-21T16:34:53Z

Ben14:

Rep:Mod:MgOben14

2017-03-21T16:33:43Z

Ben14: /* Molecular Dynamics Calculations */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon which is discussed and presented here through the example of magnesium oxide, MgO. During the the construction of any large structures where changes in dimension can occur due to temperature variations, engineers must account for this thermal expansion in its design and operation, though due to the often very large size of these structures, only a crude approximation is required. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is fundamental in enabling the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Helmholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell had been selected as an optimised function of accuracy and computational time. This then enabled the prediction of how the material expands when heated, through the computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic (or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models used and for the running of simulations are listed below. Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> The purely harmonic model is unable to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur). Therefore, the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds, hence ''quasi-harmonic''.<ref name=":1" /> It should be noted that the model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed in the output of the calculations.

The Helmholtz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholtz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the ''Results and Discussion'' section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed in the Figures below.

The primitive cell is the smallest cell that can be used to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>
[[File:BN Convensional and primitive cell inverted colour.jpg|none|thumb|718x718px|Figure to show Left: primitive cell Right: conventional cell ]]
==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter, and a* is the magnitude of the Brillouin Zone. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the subsequent section, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken arbitrarily by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) showing only a translation but leaving the model unaffected) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels (DOS) was studied to determine the required ''basis'' for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out. This was then evaluated against the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1..2x2x2...64x64x64) yielded systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. Thus the computational time was not a significant factor in this optimisation process, each submitted calculation ran well within an acceptable period for the timescale of the lab (maximum 5 minutes for the longest). From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below. With increasing grid dimension, there is an incremental decrease in the Helmholtz free energy, smaller at each step until the value is optimised.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

It is interesting to examine how other materials would alter this optimisation process.Given three contrasting examples: an oxide similar to MgO (CaO), a zeolite (faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin zone, a* as shown in the equation below.
[[File:BN Brillouin zone.PNG|none|thumb|Equation to show conversion to reciprocal space]]
The literature lattice constants for the materials discussed above are presented in the table below, where an approximate relative optimised grid size is given qualitatively, followed by a chemical rationale in the discussion below.<ref>Zhang, H., & Bukowinski, M. S. T. (1991). Modified potential-induced-breathing model of potentials between close-shell ions. Phys. Rev. B, 44(6), 2495–2503. doi:10.1103/PhysRevB.44.2495</ref><ref>Madelung, O., Rössler, U., & Schulz, M. (Eds.). (1999). Calcium oxide (CaO) crystal structure, lattice parameters, thermal expansion. In II-VI and I-VII Compounds; Semimagnetic Compounds (pp. 1–3). Berlin, Heidelberg: Springer Berlin Heidelberg. doi:10.1007/10681719_224</ref><ref>Clays and Clay Minerals, 197,3, Vol. 21, pp. 379-389</ref><ref>Kellington, S. H., Loveridge, D., Titman, J. M., L, H. F. H. and A. B., A, K. G. and M., K, N. S., … H, V. (1969). The lattice parameters of some alloys of lithium. Journal of Physics D: Applied Physics. doi:10.1088/0022-3727/2/8/415</ref>
{| class="wikitable"
!Material
!Lattice Constant, Å
!Ideal Grid Size (relative to MgO)
|-
|CaO
|4.800
|similar to MgO
|-
|Faujasite
|24.638-24.650
|smaller than MgO
|-
|Lithium
|3.5107
|smaller than MgO
|}
As the CaO is of a similar sized lattice to MgO, a similar number of k points (and therefore grid dimensions) would be needed to provide the optimum output. The zeolite, faujasite, has a lattice constant much greater than MgO, and due to the inverse relationship discussed above, a much narrower Brillouin zone. This leads to a smaller grid size being sufficient for the sampling of the necessary number of k points for an accurate model.

Initially, it may seem that due to the smaller lattice constant of Li, and the reciprocal relationship, the grid size should be larger than for MgO (which would lead to long computational times). However, this is in fact not the case due to the difference in bonding present in the lithium metal (ionic) compared to magnesium oxide (covalent). The ionic nature of the Li bonds results in the presence of delocalized electrons, which in turn leads to less repulsion between the metal cations.<ref>C. Kittel, Solid State Phys., 2005, 703.</ref> Thus a narrower energy band (a less disperse Density of States function) would be observed, and accordingly a smaller grid size would be adequate.

=== The Thermal Energy of MgO ===

==== The Quasi-Harmonic Approximation ====
The results of the optimisation with respect to the free energy computed within the quasi-harmonic approximation at a number of temperatures is presented and discussed in this section. This is quantified through the use of the Helmholtz free energy, as shown in the equation below, with an entropic term that sums over all k-values so is a useful indication of accuracy.
[[File:BN Helmholtz free energy equation.JPG|none|thumb|663x663px|Equation to show the calculation of the Helmholtz Free energy in the Quasi-Harmonic Model]]

The charts below show the plot of the Helmholtz free energy (Figure a)) as well as the lattice constant (Figure b)) and cell volume against temperature (Figure c)) for the temperature range 0 -1000 K in 100 K increments. [[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1091x1091px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]The trend in the variation of lattice constant is as expected: the greater the the temperature, the greater the thermal energy and therefore there is a larger average distance between the atoms. This is directly responsible for the same trend observed in the cell volume.

Within the temperature range 300-1000 K there is a linear relationship with temperature. A linear regression can be applied to this region, with the gradient of the curve showing the change in volume with respect to the change in temperature (at constant pressure, it is a solid and thus this can be assumed). Employing Equation 1 (see ''Introduction)'' for the thermal expansion coefficient, a value of 2.65442E-05 is obtained. There is a ********

The linear trend stated above however does not accurately model the low temperature regime, where the thermal expansion shows deviation from linearity due to a degree of temperature dependence. The Helmholtz free energy computed at 0 K is denoted as the ''Zero Point Energy, ''which is effectively a residual energy arising from quantum mechanics''. ''

This can be explained through examination of the underlying principles in the Quasi-Harmonic Approximation model. As introduced in the ''Experimental'' section of the report, this is a quantum-mechanical based model that builds on the harmonic oscillator, however introduces anharmonicity to allow thermal expansion. This cannot happen in a purely harmonic system as the bond length cannot vary with the population of thermally excited states. The quasi-harmonic approach allows for the variation of bond length, and thus thermal expansion with the population of higher energy states as a result of greater thermal energy. The figure below, showing the harmonic and anharmonic oscillator illustrates this further.
[[File:BN Harmonic anharmonic.png|none|thumb|439x439px|Figure to show the Harmonic and Anharmonic Oscillator ]]In the upper limit of 1000 K the Quasi-harmonic approximation continues to be valid, due to the melting point of MgO being much greater than this, at 3125 K. At that temperature, a change in the free energy trend would be expected, since it would no longer be a crystal structure, with vibrations but rather a much more disordered system. The Molecular Dynamics section below discusses this high temperature regime further.

The Helmoltz free energy decreases with increasing temperature and this can simply be accounted for using the thermodynamic definition of Helmholtz Free Energy, provided in the equation below.

--> CALCULATE THE THERMAL EXPANSION COEFFICIENT AND STATE THE LIT VALUE

-

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

Rep:Mod:MgOben14

2017-03-21T16:13:52Z

Ben14: /* The Quasi-Harmonic Approximation */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon which is discussed and presented here through the example of magnesium oxide, MgO. During the the construction of any large structures where changes in dimension can occur due to temperature variations, engineers must account for this thermal expansion in its design and operation, though due to the often very large size of these structures, only a crude approximation is required. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is fundamental in enabling the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Helmholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell had been selected as an optimised function of accuracy and computational time. This then enabled the prediction of how the material expands when heated, through the computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic (or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models used and for the running of simulations are listed below. Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> The purely harmonic model is unable to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur). Therefore, the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds, hence ''quasi-harmonic''.<ref name=":1" /> It should be noted that the model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed in the output of the calculations.

The Helmholtz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholtz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the ''Results and Discussion'' section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed in the Figures below.

The primitive cell is the smallest cell that can be used to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>
[[File:BN Convensional and primitive cell inverted colour.jpg|none|thumb|718x718px|Figure to show Left: primitive cell Right: conventional cell ]]
==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter, and a* is the magnitude of the Brillouin Zone. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the subsequent section, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken arbitrarily by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) showing only a translation but leaving the model unaffected) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels (DOS) was studied to determine the required ''basis'' for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out. This was then evaluated against the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1..2x2x2...64x64x64) yielded systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. Thus the computational time was not a significant factor in this optimisation process, each submitted calculation ran well within an acceptable period for the timescale of the lab (maximum 5 minutes for the longest). From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below. With increasing grid dimension, there is an incremental decrease in the Helmholtz free energy, smaller at each step until the value is optimised.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

It is interesting to examine how other materials would alter this optimisation process.Given three contrasting examples: an oxide similar to MgO (CaO), a zeolite (faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin zone, a* as shown in the equation below.
[[File:BN Brillouin zone.PNG|none|thumb|Equation to show conversion to reciprocal space]]
The literature lattice constants for the materials discussed above are presented in the table below, where an approximate relative optimised grid size is given qualitatively, followed by a chemical rationale in the discussion below.<ref>Zhang, H., & Bukowinski, M. S. T. (1991). Modified potential-induced-breathing model of potentials between close-shell ions. Phys. Rev. B, 44(6), 2495–2503. doi:10.1103/PhysRevB.44.2495</ref><ref>Madelung, O., Rössler, U., & Schulz, M. (Eds.). (1999). Calcium oxide (CaO) crystal structure, lattice parameters, thermal expansion. In II-VI and I-VII Compounds; Semimagnetic Compounds (pp. 1–3). Berlin, Heidelberg: Springer Berlin Heidelberg. doi:10.1007/10681719_224</ref><ref>Clays and Clay Minerals, 197,3, Vol. 21, pp. 379-389</ref><ref>Kellington, S. H., Loveridge, D., Titman, J. M., L, H. F. H. and A. B., A, K. G. and M., K, N. S., … H, V. (1969). The lattice parameters of some alloys of lithium. Journal of Physics D: Applied Physics. doi:10.1088/0022-3727/2/8/415</ref>
{| class="wikitable"
!Material
!Lattice Constant, Å
!Ideal Grid Size (relative to MgO)
|-
|CaO
|4.800
|similar to MgO
|-
|Faujasite
|24.638-24.650
|smaller than MgO
|-
|Lithium
|3.5107
|smaller than MgO
|}
As the CaO is of a similar sized lattice to MgO, a similar number of k points (and therefore grid dimensions) would be needed to provide the optimum output. The zeolite, faujasite, has a lattice constant much greater than MgO, and due to the inverse relationship discussed above, a much narrower Brillouin zone. This leads to a smaller grid size being sufficient for the sampling of the necessary number of k points for an accurate model.

Initially, it may seem that due to the smaller lattice constant of Li, and the reciprocal relationship, the grid size should be larger than for MgO (which would lead to long computational times). However, this is in fact not the case due to the difference in bonding present in the lithium metal (ionic) compared to magnesium oxide (covalent). The ionic nature of the Li bonds results in the presence of delocalized electrons, which in turn leads to less repulsion between the metal cations.<ref>C. Kittel, Solid State Phys., 2005, 703.</ref> Thus a narrower energy band (a less disperse Density of States function) would be observed, and accordingly a smaller grid size would be adequate.

=== The Thermal Energy of MgO ===

==== The Quasi-Harmonic Approximation ====
The results of the optimisation with respect to the free energy computed within the quasi-harmonic approximation at a number of temperatures is presented and discussed in this section. This is quantified through the use of the Helmholtz free energy, as shown in the equation below, with an entropic term that sums over all k-values so is a useful indication of accuracy.
[[File:BN Helmholtz free energy equation.JPG|none|thumb|663x663px|Equation to show the calculation of the Helmholtz Free energy in the Quasi-Harmonic Model]]

The charts below show the plot of the Helmholtz free energy (Figure a)) as well as the lattice constant (Figure b)) and cell volume against temperature (Figure c)) for the temperature range 0 -1000 K in 100 K increments. [[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1091x1091px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]The trend in the variation of lattice constant is as expected: the greater the the temperature, the greater the thermal energy and therefore there is a larger average distance between the atoms. This is directly responsible for the same trend observed in the cell volume.

Within the temperature range 300-1000 K there is a linear relationship with temperature, however this does not accurately model the low temperature regime, where the thermal expansion shows deviation from linearity due to a degree of temperature dependence. The Helmholtz free energy computed at 0 K is denoted as the ''Zero Point Energy, ''which is effectively a residual energy arising from quantum mechanics''. ''

This can be explained through examination of the underlying principles in the Quasi-Harmonic Approximation model. As introduced in the ''Experimental'' section of the report, this is a quantum-mechanical based model that builds on the harmonic oscillator, however introduces anharmonicity to allow thermal expansion. This cannot happen in a purely harmonic system as the bond length cannot vary with the population of thermally excited states. The quasi-harmonic approach allows for the variation of bond length, and thus thermal expansion with the population of higher energy states as a result of greater thermal energy. The figure below, showing the harmonic and anharmonic oscillator illustrates this further.
[[File:BN Harmonic anharmonic.png|none|thumb|439x439px|Figure to show the Harmonic and Anharmonic Oscillator ]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

File:BN Harmonic anharmonic.png

2017-03-21T16:06:54Z

Ben14:

Rep:Mod:MgOben14

2017-03-21T16:05:57Z

Ben14: /* The Thermal Energy of MgO */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon which is discussed and presented here through the example of magnesium oxide, MgO. During the the construction of any large structures where changes in dimension can occur due to temperature variations, engineers must account for this thermal expansion in its design and operation, though due to the often very large size of these structures, only a crude approximation is required. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is fundamental in enabling the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Helmholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell had been selected as an optimised function of accuracy and computational time. This then enabled the prediction of how the material expands when heated, through the computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic (or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models used and for the running of simulations are listed below. Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> The purely harmonic model is unable to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur). Therefore, the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds, hence ''quasi-harmonic''.<ref name=":1" /> It should be noted that the model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed in the output of the calculations.

The Helmholtz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholtz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the ''Results and Discussion'' section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed in the Figures below.

The primitive cell is the smallest cell that can be used to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>
[[File:BN Convensional and primitive cell inverted colour.jpg|none|thumb|718x718px|Figure to show Left: primitive cell Right: conventional cell ]]
==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter, and a* is the magnitude of the Brillouin Zone. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the subsequent section, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken arbitrarily by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) showing only a translation but leaving the model unaffected) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels (DOS) was studied to determine the required ''basis'' for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out. This was then evaluated against the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1..2x2x2...64x64x64) yielded systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. Thus the computational time was not a significant factor in this optimisation process, each submitted calculation ran well within an acceptable period for the timescale of the lab (maximum 5 minutes for the longest). From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below. With increasing grid dimension, there is an incremental decrease in the Helmholtz free energy, smaller at each step until the value is optimised.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

It is interesting to examine how other materials would alter this optimisation process.Given three contrasting examples: an oxide similar to MgO (CaO), a zeolite (faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin zone, a* as shown in the equation below.
[[File:BN Brillouin zone.PNG|none|thumb|Equation to show conversion to reciprocal space]]
The literature lattice constants for the materials discussed above are presented in the table below, where an approximate relative optimised grid size is given qualitatively, followed by a chemical rationale in the discussion below.<ref>Zhang, H., & Bukowinski, M. S. T. (1991). Modified potential-induced-breathing model of potentials between close-shell ions. Phys. Rev. B, 44(6), 2495–2503. doi:10.1103/PhysRevB.44.2495</ref><ref>Madelung, O., Rössler, U., & Schulz, M. (Eds.). (1999). Calcium oxide (CaO) crystal structure, lattice parameters, thermal expansion. In II-VI and I-VII Compounds; Semimagnetic Compounds (pp. 1–3). Berlin, Heidelberg: Springer Berlin Heidelberg. doi:10.1007/10681719_224</ref><ref>Clays and Clay Minerals, 197,3, Vol. 21, pp. 379-389</ref><ref>Kellington, S. H., Loveridge, D., Titman, J. M., L, H. F. H. and A. B., A, K. G. and M., K, N. S., … H, V. (1969). The lattice parameters of some alloys of lithium. Journal of Physics D: Applied Physics. doi:10.1088/0022-3727/2/8/415</ref>
{| class="wikitable"
!Material
!Lattice Constant, Å
!Ideal Grid Size (relative to MgO)
|-
|CaO
|4.800
|similar to MgO
|-
|Faujasite
|24.638-24.650
|smaller than MgO
|-
|Lithium
|3.5107
|smaller than MgO
|}
As the CaO is of a similar sized lattice to MgO, a similar number of k points (and therefore grid dimensions) would be needed to provide the optimum output. The zeolite, faujasite, has a lattice constant much greater than MgO, and due to the inverse relationship discussed above, a much narrower Brillouin zone. This leads to a smaller grid size being sufficient for the sampling of the necessary number of k points for an accurate model.

Initially, it may seem that due to the smaller lattice constant of Li, and the reciprocal relationship, the grid size should be larger than for MgO (which would lead to long computational times). However, this is in fact not the case due to the difference in bonding present in the lithium metal (ionic) compared to magnesium oxide (covalent). The ionic nature of the Li bonds results in the presence of delocalized electrons, which in turn leads to less repulsion between the metal cations.<ref>C. Kittel, Solid State Phys., 2005, 703.</ref> Thus a narrower energy band (a less disperse Density of States function) would be observed, and accordingly a smaller grid size would be adequate.

=== The Thermal Energy of MgO ===

==== The Quasi-Harmonic Approximation ====
The results of the optimisation with respect to the free energy computed within the quasi-harmonic approximation at a number of temperatures is presented and discussed in this section. This is quantified through the use of the Helmholtz free energy, as shown in the equation below, with an entropic term that sums over all k-values so is a useful indication of accuracy.
[[File:BN Helmholtz free energy equation.JPG|none|thumb|663x663px|Equation to show the calculation of the Helmholtz Free energy in the Quasi-Harmonic Model]]

The charts below show the plot of the Helmholtz free energy (Figure a)) as well as the lattice constant (Figure b)) and cell volume against temperature (Figure c)) for the temperature range 0 -1000 K in 100 K increments. [[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1091x1091px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]The trend in the variation of lattice constant is as expected: the greater the the temperature, the greater the thermal energy and therefore there is a larger average distance between the atoms. This is directly responsible for the same trend observed in the cell volume.

Within the temperature range 300-1000 K there is a linear relationship with temperature, however this does not accurately model the low temperature regime, where the thermal expansion shows deviation from linearity due to a degree of temperature dependence. The Helmholtz free energy computed at 0 K is denoted as the ''Zero Point Energy, ''which is effectively a residual energy arising from quantum mechanics''. ''

This can be explained through examination of the underlying principles in the Quasi-Harmonic Approximation model. As introduced in the ''Experimental'' section of the report, this is a quantum-mechanical based model that builds on the harmonic oscillator, however introduces anharmonicity to allow thermal expansion. This cannot happen in a purely harmonic system as the bond length cannot vary with the population of thermally excited states. The quasi-harmonic approach allows for the variation of bond length, and thus thermal expansion with the population of higher energy states as a result of greater thermal energy. The figure below, showing the harmonic and anharmonic oscillator illustrates this further.

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

Rep:Mod:MgOben14

2017-03-21T16:00:18Z

Ben14:

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon which is discussed and presented here through the example of magnesium oxide, MgO. During the the construction of any large structures where changes in dimension can occur due to temperature variations, engineers must account for this thermal expansion in its design and operation, though due to the often very large size of these structures, only a crude approximation is required. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is fundamental in enabling the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Helmholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell had been selected as an optimised function of accuracy and computational time. This then enabled the prediction of how the material expands when heated, through the computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic (or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models used and for the running of simulations are listed below. Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> The purely harmonic model is unable to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur). Therefore, the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds, hence ''quasi-harmonic''.<ref name=":1" /> It should be noted that the model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed in the output of the calculations.

The Helmholtz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholtz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the ''Results and Discussion'' section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed in the Figures below.

The primitive cell is the smallest cell that can be used to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>
[[File:BN Convensional and primitive cell inverted colour.jpg|none|thumb|718x718px|Figure to show Left: primitive cell Right: conventional cell ]]
==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter, and a* is the magnitude of the Brillouin Zone. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the subsequent section, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken arbitrarily by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) showing only a translation but leaving the model unaffected) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels (DOS) was studied to determine the required ''basis'' for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out. This was then evaluated against the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1..2x2x2...64x64x64) yielded systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. Thus the computational time was not a significant factor in this optimisation process, each submitted calculation ran well within an acceptable period for the timescale of the lab (maximum 5 minutes for the longest). From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below. With increasing grid dimension, there is an incremental decrease in the Helmholtz free energy, smaller at each step until the value is optimised.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

It is interesting to examine how other materials would alter this optimisation process.Given three contrasting examples: an oxide similar to MgO (CaO), a zeolite (faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin zone, a* as shown in the equation below.
[[File:BN Brillouin zone.PNG|none|thumb|Equation to show conversion to reciprocal space]]
The literature lattice constants for the materials discussed above are presented in the table below, where an approximate relative optimised grid size is given qualitatively, followed by a chemical rationale in the discussion below.<ref>Zhang, H., & Bukowinski, M. S. T. (1991). Modified potential-induced-breathing model of potentials between close-shell ions. Phys. Rev. B, 44(6), 2495–2503. doi:10.1103/PhysRevB.44.2495</ref><ref>Madelung, O., Rössler, U., & Schulz, M. (Eds.). (1999). Calcium oxide (CaO) crystal structure, lattice parameters, thermal expansion. In II-VI and I-VII Compounds; Semimagnetic Compounds (pp. 1–3). Berlin, Heidelberg: Springer Berlin Heidelberg. doi:10.1007/10681719_224</ref><ref>Clays and Clay Minerals, 197,3, Vol. 21, pp. 379-389</ref><ref>Kellington, S. H., Loveridge, D., Titman, J. M., L, H. F. H. and A. B., A, K. G. and M., K, N. S., … H, V. (1969). The lattice parameters of some alloys of lithium. Journal of Physics D: Applied Physics. doi:10.1088/0022-3727/2/8/415</ref>
{| class="wikitable"
!Material
!Lattice Constant, Å
!Ideal Grid Size (relative to MgO)
|-
|CaO
|4.800
|similar to MgO
|-
|Faujasite
|24.638-24.650
|smaller than MgO
|-
|Lithium
|3.5107
|smaller than MgO
|}
As the CaO is of a similar sized lattice to MgO, a similar number of k points (and therefore grid dimensions) would be needed to provide the optimum output. The zeolite, faujasite, has a lattice constant much greater than MgO, and due to the inverse relationship discussed above, a much narrower Brillouin zone. This leads to a smaller grid size being sufficient for the sampling of the necessary number of k points for an accurate model.

Initially, it may seem that due to the smaller lattice constant of Li, and the reciprocal relationship, the grid size should be larger than for MgO (which would lead to long computational times). However, this is in fact not the case due to the difference in bonding present in the lithium metal (ionic) compared to magnesium oxide (covalent). The ionic nature of the Li bonds results in the presence of delocalized electrons, which in turn leads to less repulsion between the metal cations.<ref>C. Kittel, Solid State Phys., 2005, 703.</ref> Thus a narrower energy band (a less disperse Density of States function) would be observed, and accordingly a smaller grid size would be adequate.

=== The Thermal Energy of MgO ===

==== The Quasi-Harmonic Approximation ====
The results of the optimisation with respect to the free energy computed within the quasi-harmonic approximation at a number of temperatures is presented and discussed in this section. This is quantified through the use of the Helmholtz free energy, as shown in the equation below, with an entropic term that sums over all k-values so is a useful indication of accuracy.
[[File:BN Helmholtz free energy equation.JPG|none|thumb|663x663px|Equation to show the calculation of the Helmholtz Free energy in the Quasi-Harmonic Model]]

The charts below show the plot of the Helmholtz free energy (Figure a)) as well as the lattice constant (Figure b)) and cell volume against temperature (Figure c)) for the temperature range 0 -1000 K in 100 K increments. [[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1091x1091px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]The trend in the variation of lattice constant is as expected: the greater the the temperature, the greater the thermal energy and therefore there is a larger average distance between the atoms. This is directly responsible for the same trend observed in the cell volume.

Within the temperature range 300-1000 K there is a linear relationship with temperature, however this does not accurately model the low temperature regime, where the thermal expansion shows deviation from linearity due to a degree of temperature dependence. The Helmholtz free energy computed at 0 K is denoted as the ''Zero Point Energy, ''which is effectively a residual energy arising from quantum mechanics''. ''

This can be explained through examination of the underlying principles in the Quasi-Harmonic Approximation model. As introduced in the ''Experimental'' section of the report, this is a quantum-mechanical based model that builds on the harmonic oscillator, however introduces anharmonicity to allow thermal expansion. This cannot happen in a purely harmonic system as the bond length cannot vary with a variation of

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

Rep:Mod:MgOben14

2017-03-21T15:51:54Z

2017-03-21T13:24:19Z

Ben14: /* Grid Size Optimisation */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal Expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon (through the example of magnesium oxide) presented in this report. The construction of any large structures where changes in dimension can occur due to temperature variations must account for this thermal expansion in its engineering and design, though the size of these structures only require a crude approximation. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is required to enable the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Hemholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell to be used, had been selected as an optimised function of accuracy and computational time. This enabled the prediction of how the material expands when heated, through computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic ( or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models and running of simulations are listed below. The use of the Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> As the purely harmonic model isn't able to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur), the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds.<ref name=":1" />The model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed.

The Helmholz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholtz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the Results and Discussion Section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the modelled MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed below.

The primitive cell is the smallest cell to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>

<nowiki> </nowiki>INSERT IMAGES OF PRIMITIVE AND CONVENSIONAL CELL HERE![[File:BN_MgO_snapshot_of_lattice.PNG|none|thumb|251x251px|Figure to show Snapshot of MgO lattice]]

==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the next section of the lab, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) CHECK THIS) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels was studied to determine the required basis for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out, requiring a balance with the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1...64x64x64) yields systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]
<blockquote></blockquote>To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

It is interesting to examine how other materials would alter this optimisation process. Taking the three contrasting examples of a similar oxide (CaO), a zeolite (Faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin Zone, a* as shown in the equation below.
[[File:BN Brillouin zone.PNG|none|thumb|Equation to show conversion to reciprocal space]]
The lattice constants for the materials discussed above are shown in the table below, where the optimised grid size, relative to MgO is proposed, followed by a rationale in the discussion below.<ref>Zhang, H., & Bukowinski, M. S. T. (1991). Modified potential-induced-breathing model of potentials between close-shell ions. Phys. Rev. B, 44(6), 2495–2503. doi:10.1103/PhysRevB.44.2495</ref><ref>Madelung, O., Rössler, U., & Schulz, M. (Eds.). (1999). Calcium oxide (CaO) crystal structure, lattice parameters, thermal expansion. In II-VI and I-VII Compounds; Semimagnetic Compounds (pp. 1–3). Berlin, Heidelberg: Springer Berlin Heidelberg. doi:10.1007/10681719_224</ref><ref>Clays and Clay Minerals, 197,3, Vol. 21, pp. 379-389</ref><ref>Kellington, S. H., Loveridge, D., Titman, J. M., L, H. F. H. and A. B., A, K. G. and M., K, N. S., … H, V. (1969). The lattice parameters of some alloys of lithium. Journal of Physics D: Applied Physics. doi:10.1088/0022-3727/2/8/415</ref>
{| class="wikitable"
!Material
!Lattice Constant, Å
!Ideal Grid Size (relative to MgO)
|-
|CaO
|4.800
|similar; 32x32x32
|-
|Faujasite
|24.638-24.650
|smaller
|-
|Lithium
|3.5107
|smaller
|}
As the CaO is of a similar sized lattice constant to MgO, a similar number of k points (and therefore grid dimensions) would be sufficient to provide the optimum output. The zeolite, Faujasite, has a lattice constant much greater than MgO, and due to the inverse relationship discussed above, a much narrower Brillouin zone. This leads to a smaller grid size being sufficient for the sampling of k points.

Initially, it may seem that due to the smaller lattice constant of Li, the grid size should be larger than for MgO (which would lead to long computational times). However, this is in fact not the case due to the difference in bonding present in the lithium metal (ionic) compared to magnesium oxide (covalent). The ionic nature of these bonds results in the presence of delocalized electrons, which lead to less repulsion between the metal cation and in

==== The Quasi-Harmonic Approximation ====
[[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1105x1105px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

Rep:Mod:MgOben14

2017-03-21T12:05:56Z

Ben14: /* Grid Size Optimisation */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal Expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon (through the example of magnesium oxide) presented in this report. The construction of any large structures where changes in dimension can occur due to temperature variations must account for this thermal expansion in its engineering and design, though the size of these structures only require a crude approximation. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is required to enable the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Hemholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell to be used, had been selected as an optimised function of accuracy and computational time. This enabled the prediction of how the material expands when heated, through computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic ( or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models and running of simulations are listed below. The use of the Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> As the purely harmonic model isn't able to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur), the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds.<ref name=":1" />The model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed.

The Helmholz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholtz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the Results and Discussion Section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the modelled MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed below.

The primitive cell is the smallest cell to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>

<nowiki> </nowiki>INSERT IMAGES OF PRIMITIVE AND CONVENSIONAL CELL HERE![[File:BN_MgO_snapshot_of_lattice.PNG|none|thumb|251x251px|Figure to show Snapshot of MgO lattice]]

==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the next section of the lab, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) CHECK THIS) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels was studied to determine the required basis for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out, requiring a balance with the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1...64x64x64) yields systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]
<blockquote></blockquote>To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

It is interesting to examine how other materials would alter this optimisation process. Taking the three contrasting examples of a similar oxide (CaO), a zeolite (Faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin Zone, a* as shown in the equation below.
[[File:BN Brillouin zone.PNG|none|thumb|Equation to show conversion to reciprocal space]]
The lattice constants for the materials discussed above are shown in the table below.<ref>Zhang, H., & Bukowinski, M. S. T. (1991). Modified potential-induced-breathing model of potentials between close-shell ions. Phys. Rev. B, 44(6), 2495–2503. doi:10.1103/PhysRevB.44.2495</ref><ref>Madelung, O., Rössler, U., & Schulz, M. (Eds.). (1999). Calcium oxide (CaO) crystal structure, lattice parameters, thermal expansion. In II-VI and I-VII Compounds; Semimagnetic Compounds (pp. 1–3). Berlin, Heidelberg: Springer Berlin Heidelberg. doi:10.1007/10681719_224</ref><ref>Clays and Clay Minerals, 197,3, Vol. 21, pp. 379-389</ref><ref>Kellington, S. H., Loveridge, D., Titman, J. M., L, H. F. H. and A. B., A, K. G. and M., K, N. S., … H, V. (1969). The lattice parameters of some alloys of lithium. Journal of Physics D: Applied Physics. doi:10.1088/0022-3727/2/8/415</ref>
{| class="wikitable"
!Material
!Lattice Constant, Å
!Ideal Grid Size (relative to MgO)
|-
|CaO
|4.800
|
|-
|Faujasite
|24.638-24.650
|
|-
|Lithium
|3.5107
|
|}

==== The Quasi-Harmonic Approximation ====
[[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1455x1455px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

File:BN Brillouin zone.PNG

2017-03-21T11:23:44Z

Ben14:

Rep:Mod:MgOben14

2017-03-21T11:23:16Z

Ben14: /* Grid Size Optimisation */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal Expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon (through the example of magnesium oxide) presented in this report. The construction of any large structures where changes in dimension can occur due to temperature variations must account for this thermal expansion in its engineering and design, though the size of these structures only require a crude approximation. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is required to enable the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Hemholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell to be used, had been selected as an optimised function of accuracy and computational time. This enabled the prediction of how the material expands when heated, through computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic ( or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models and running of simulations are listed below. The use of the Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> As the purely harmonic model isn't able to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur), the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds.<ref name=":1" />The model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed.

The Helmholz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholtz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the Results and Discussion Section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the modelled MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed below.

The primitive cell is the smallest cell to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>

<nowiki> </nowiki>INSERT IMAGES OF PRIMITIVE AND CONVENSIONAL CELL HERE![[File:BN_MgO_snapshot_of_lattice.PNG|none|thumb|251x251px|Figure to show Snapshot of MgO lattice]]

==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the next section of the lab, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) CHECK THIS) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels was studied to determine the required basis for the accurate modelling of the system. The Helmholtz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out, requiring a balance with the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1...64x64x64) yields systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]
<blockquote></blockquote>To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholtz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholtz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

It is interesting to examine how other materials would alter this optimisation process. Taking the three contrasting examples of a similar oxide (CaO), a zeolite (Faujasite) and a metal (lithium), the variation in their lattice constants is an excellent starting point, due to the inverse relationship between the length of the unique Brillouin Zone, a* as shown in the equation below.

==== The Quasi-Harmonic Approximation ====
[[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1455x1455px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

Rep:Mod:MgOben14

2017-03-21T11:02:53Z

Ben14: /* Phonon Dispersion Curves */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal Expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon (through the example of magnesium oxide) presented in this report. The construction of any large structures where changes in dimension can occur due to temperature variations must account for this thermal expansion in its engineering and design, though the size of these structures only require a crude approximation. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is required to enable the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Hemholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell to be used, had been selected as an optimised function of accuracy and computational time. This enabled the prediction of how the material expands when heated, through computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic ( or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models and running of simulations are listed below. The use of the Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> As the purely harmonic model isn't able to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur), the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds.<ref name=":1" />The model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed.

The Helmholz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the Results and Discussion Section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the modelled MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed below.

The primitive cell is the smallest cell to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>

<nowiki> </nowiki>INSERT IMAGES OF PRIMITIVE AND CONVENSIONAL CELL HERE![[File:BN_MgO_snapshot_of_lattice.PNG|none|thumb|251x251px|Figure to show Snapshot of MgO lattice]]

==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter. A vibrational band is denoted by the term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the next section of the lab, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) CHECK THIS) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels was studied to determine the required basis for the accurate modelling of the system. The Helmholz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out, requiring a balance with the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1...64x64x64) yields systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]
<blockquote></blockquote>To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

NEXT, EXAMINE OTHER MATERIALS: WOULD THE GRID SIZE NEEDED BE COMPARATIVE

==== The Quasi-Harmonic Approximation ====
[[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1455x1455px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

Rep:Mod:MgOben14

2017-03-18T16:49:46Z

Ben14: /* Grid Size Optimisation */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal Expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon (through the example of magnesium oxide) presented in this report. The construction of any large structures where changes in dimension can occur due to temperature variations must account for this thermal expansion in its engineering and design, though the size of these structures only require a crude approximation. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is required to enable the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Hemholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell to be used, had been selected as an optimised function of accuracy and computational time. This enabled the prediction of how the material expands when heated, through computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic ( or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models and running of simulations are listed below. The use of the Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> As the purely harmonic model isn't able to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur), the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds.<ref name=":1" />The model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed.

The Helmholz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the Results and Discussion Section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the modelled MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed below.

The primitive cell is the smallest cell to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>

<nowiki> </nowiki>INSERT IMAGES OF PRIMITIVE AND CONVENSIONAL CELL HERE![[File:BN_MgO_snapshot_of_lattice.PNG|none|thumb|251x251px|Figure to show Snapshot of MgO lattice]]

==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter. A vibrational band is denoted by term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the next section of the lab, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) CHECK THIS) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels was studied to determine the required basis for the accurate modelling of the system. The Helmholz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out, requiring a balance with the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1...64x64x64) yields systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]
<blockquote></blockquote>To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

NEXT, EXAMINE OTHER MATERIALS: WOULD THE GRID SIZE NEEDED BE COMPARATIVE

==== The Quasi-Harmonic Approximation ====
[[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1455x1455px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

Rep:Mod:MgOben14

2017-03-18T16:48:03Z

Ben14: /* Phonon Dispersion Curves */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal Expansion is the change in shape, area or volume of matter in response to a change in temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this intrinsic property of materials (through normalisation with the V0 term), at constant pressure as given by the equation below, Equation 1<ref name=":0">K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.
[[File:BN Thermal expansion coefficient.JPG|none|thumb|160x160px|Equation 1: Thermal Expansion Coefficient ]]
Isotropic expansion is such that the rate of expansion is equal in all directions, and is more simple to model.<ref name=":0" /> There are many significant applications for the modelling of this phenomenon (through the example of magnesium oxide) presented in this report. The construction of any large structures where changes in dimension can occur due to temperature variations must account for this thermal expansion in its engineering and design, though the size of these structures only require a crude approximation. The choice of metals with minimal expansion coefficients for use in alloys in aerospace engineering is one example where precise understanding of this field is required to enable the safe and correct operation of the aircraft<ref>K. U. Kainer, Metal Matrix Composites: Custom-made Materials for Automotive and Aerospace Engineering, 2006.</ref>.

==== Aims of MgO Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of an ideal, non-defective and periodic crystal of MgO. In particular, phonon dispersion curves showed the crystal's vibrational energy structure.

The vibrational energies were then used to determine computationally the Hemholtz Free Energy of the crystal, once an ideal basis i.e. dimensions of the cell to be used, had been selected as an optimised function of accuracy and computational time. This enabled the prediction of how the material expands when heated, through computation of selected properties of the crystal, such as the lattice parameter and unit cell volume. The Quasi-Harmonic ( or Lattice Dynamics, LD) and Molecular Dynamics (MD) approximations were compared for their accuracy in this modelling of the thermal expansion, alongside a comparison with the literature available (CITATION NEEDED).

== Experimental ==

==== Tools in Use ====
The computational tools employed for the generation of the models and running of simulations are listed below. The use of the Linux operating system provided a suitable environment for numerical simulations. DL Visualize is a graphical user interface which served to run the mollecular modelling code GULP.
# RedHat Linux
# DL Visualize
# GULP (General Utility Lattice Programme)

==== Methodology Used ====

===== The Quasi-Harmonic Approximation Model =====
The quantum mechanical Quasi-Harmonic Approximation model (alternatively, Lattice Dynamics Model) is used for the description of volume dependent thermal effects. It is based in the assumption that the Harmonic Approximation holds for each value of the lattice constant.<ref name=":1">M. T. Dove, Introduction to Lattice Dynamics, Cambridge University Press, 1993.</ref> As the purely harmonic model isn't able to adequately model thermal expansion because equilibrium distance would be independent of temperature (and as such no expansion would occur), the phonon frequencies are allowed to become volume dependent, but such that for each volume, the harmonic approximation holds.<ref name=":1" />The model does not account for any 'phonon-phonon' interactions however, and as such, anharmonicity is observed.

The Helmholz Free Energy in the Quasi-Harmonic Model is given by the Equation below, (''where F is Helmholz free energy, U is internal energy, E is the zero potential energy, V is volume, T is temperature (K) and S is entropy'').[[File:BN Quasi-harmonic free energy.JPG|none|thumb|475x475px|Equation 2 to show the Free Energy as given by the Quasi-Harmonic Approximation]]

===== The Molecular Dynamics Model =====
The Molecular Dynamics model was employed in addition to to the Quasi-Harmonic Approximation discussed above to allow comparisons of the models, as presented in the Results and Discussion Section. This is a classical mechanics model where vibration is represented as random motion. It makes use of Newton's Second Law of Motion (''F = ma'') to determine the force acting on an atom moving along its trajectory over a period of time. The model provides an alternative approach to the computation of thermodynamic properties.

== Results and Discussion ==

==== Representing the MgO Structure ====
The output of the initial calculations provided the parameters of the modelled MgO crystal. The cell can be represented both in terms of the ''primitive cell ''and the ''conventional cell'', displayed below.

The primitive cell is the smallest cell to describe the crystal structure: it is a rhombohedral structure consisting of 1 oxygen at the body centre and 8 x 1/8 Mg (each located on the corners of the shape), yielding a two atom cell. The output of the calculations is provided below.

<code>bond angle </code> <code> 60.0000° </code>

<code>lattice constant</code> <code>2.9783 '''Å ''' </code>

The conventional cell represents the FCC (face centred cubic) nature of the system more clearly. The output of the calculations yielded a lattice constant of 4.212 Å, which is in agreement with the published literature value.<ref>C. Y. Ho and R. E. Taylor, Thermal Expansion of Solids, ASM International, 1998.</ref>

<nowiki> </nowiki>INSERT IMAGES OF PRIMITIVE AND CONVENSIONAL CELL HERE![[File:BN_MgO_snapshot_of_lattice.PNG|none|thumb|251x251px|Figure to show Snapshot of MgO lattice]]

==== Phonon Dispersion Curves ====
The phonon dispersion curve was computed for the modelling of the vibrational energy structure of the crystal, as shown in the figure below. It shows the relationship between the vibrational frequency of the phonon (ωk) and the wavevector k. It is periodic with a*=2π/a where a is the lattice parameter. A vibrational band is denoted by term ''branch ''and these branches are observed in the figure below. The total number of branches is given by the product of the number of atoms in the primitive cell and the number of dimensions. For our primitive cell consisting of 2 atoms (as discussed in the section above), in 3 dimensions, 6 branches are observed.

The merging of branches that can be seen for MgO is due to reduction by symmetry and the resulting degeneracy of the vibrational modes. Optical and acoustic branches are exhibited for this crystal, as it has a primitive cell consisting of more than one atom, which gives rise to both ''coherent movements'' of atoms in the lattice out of their equilibrium positions (acoustic phonons) and ''out of phase movements'' (optical phonons)<ref>L. Saviot, B. Champagnon, E. Duval, I. A. Kudriavtsev and A. I. Ekimov, J. Non. Cryst. Solids, 1996, 197, 238–246.</ref>. The optical modes have a non zero frequency at the centre of the Brillouin zone (denoted by Γ =(0,0,0)).''' ''' [[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|501x501px|Figure to show the Phonon Dispersion Curve Computed for MgO (Grid size = 1x1x1)]]The frequencies of the bands at the single k point (as determined by overlaying with the Phonon Density of States graph, below) used in the 1x1x1 calculation were: 288.49, 288.49, 351.76, 351.76, 676.23, 818.82. For the degenerate vibrational modes, the density of states is greater, to indicate this. [[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|415x415px|Figure to show MgO DOS for 1x1x1]]The density of states (DOS), as shown in the figure above is the energy distribution of the system, i.e. the number of states at each frequency. In the next section of the lab, greater grid sizes are investigated, in order to optimise the calculation in terms of accuracy and computational time. The one k point (taken by the software at 0.5,0.5,0.5 but represented in the phonon dispersion curve as at L (1,1,1) CHECK THIS) is not sufficient for the accurate depiction of the DOS.

==== Grid Size Optimisation ====
The density of states of the vibrational energy levels was studied to determine the required basis for the accurate modelling of the system. The Helmholz free energy was taken at various grid dimensions (as detailed in both the figure and table below) in order to quantify the accuracy to which the calculations should be carried out, requiring a balance with the expense of computational time in which the calculations were to be run. Examining grid sizes of equal dimensions (1x1x1...64x64x64) yields systems with a high degree of symmetry, reducing greatly the calculations needed to model this system, and thereby reducing computational time, which is advantageous when carrying out this work. From the figures presented below, it can be seen than an increase in the number of unit cells used increases the resolution of the DOS because it takes into account more of the possible vibrations of the system.[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]
<blockquote></blockquote>To quantify the level of accuracy to which the above grid dimensions allow the calculations to be run, a physical parameter was extracted from each of the output log files, the Helmholtz free energy, presented in the table below.
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}
As can be seen, the 32x32x32 system is sufficient, as the value of the Helmholz free energy is determined to a convergent value, there is no further benefit in increasing to 64x64x64. The 3x3x3 grid size would be sufficient to provide a calculation accurate to 0.1 meV, should computational time be limiting.

==== The Quasi-Harmonic Approximation ====
[[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1455x1455px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

Rep:Mod:MgOben14

2017-03-18T15:46:18Z

Ben14:

File:Diatomic phonons.png

2017-03-18T15:24:59Z

Ben14: Ben14 uploaded a new version of File:Diatomic phonons.png

By Brews ohare (Own work) [CC BY-SA 3.0 (http://creativecommons.org/licenses/by-sa/3.0)], via Wikimedia Commons

File:Diatomic phonons.png

2017-03-18T15:24:58Z

2017-03-18T11:36:40Z

Ben14:

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==
Thermal Expansion is the change in shape, area or volume of matter in response to a change in Temperature.<ref>R. A. Schapery, J. Compos. Mater., 1968.</ref> <ref>P. Atkins and J. de Paula, Atkins’ Physical Chemistry, OUP Oxford, 10th edn., 2014.</ref> The coefficient of thermal expansion describes this, at constant pressure<ref>K. V. K. Rao, AIP Conf. Proc., 1974, 17, 219–230.</ref>.

==== Aims of Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of a crystal of MgO. The vibrational energies were used to determine computationally the free energy of the crystal. This was used to predict how the material expands when heated. The Quasi Harmonic (LD) and Molecular Dynamics (MD) approximations were compared for the prediction of thermal expansion.

== Experimental ==

==== Methodology Used ====

==== Tools in Use (computer programmes) ====

== Results and Discussion ==

==== Representing the MgO Structure ====
[[File:BN_MgO_snapshot_of_lattice.PNG|none|thumb|251x251px|Figure to show Snapshot of MgO lattice]]

==== Phonon Dispersion Curves ====
[[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|414x414px|Figure to show the Phonon Dispersion Curve Computed for MgO]]
[[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|413x413px|Figure to show MgO DOS for 1x1x1]]

==== Grid Size Optimisation ====
[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]
<blockquote></blockquote>
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}

==== The Quasi-Harmonic Approximation ====
[[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1455x1455px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==
<references />

===== 1. Calculating the internal energy of an MgO crystal =====

===== 2. Lattice Vibrations - Computing the Phonons =====

===== 3. Calculating the Free Energy in the Harmonic Approximation =====

===== 4. The Thermal Expansion of MgO =====

===== 5. Molecular Dynamics =====

Rep:Mod:MgOben14

2017-03-17T14:21:48Z

Ben14:

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==

==== Aims of Computational Lab ====
A simple model of atomic interactions were used to calculate the energy and vibrations of a crystal of MgO. The vibrational energies were used to determine computationally the free energy of the crystal. This was used to predict how the material expands when heated. The Quasi Harmonic (LD) and Molecular Dynamics (MD) approximations were compared for the prediction of thermal expansion.

== Experimental ==

==== Methodology Used ====

==== Tools in Use (computer programmes) ====

== Results and Discussion ==

==== Representing the MgO Structure ====
[[File:BN_MgO_snapshot_of_lattice.PNG|none|thumb|251x251px|Figure to show Snapshot of MgO lattice]]

==== Phonon Dispersion Curves ====
[[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|414x414px|Figure to show the Phonon Dispersion Curve Computed for MgO]]
[[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|413x413px|Figure to show MgO DOS for 1x1x1]]

==== Grid Size Optimisation ====
[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1142x1142px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]
<blockquote></blockquote>
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}

==== The Quasi-Harmonic Approximation ====
[[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1455x1455px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|744x744px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==

===== 1. Calculating the internal energy of an MgO crystal =====

===== 2. Lattice Vibrations - Computing the Phonons =====

===== 3. Calculating the Free Energy in the Harmonic Approximation =====

===== 4. The Thermal Expansion of MgO =====

===== 5. Molecular Dynamics =====

Rep:Mod:MgOben14

2017-03-16T15:13:45Z

Ben14: /* Phonon Dispersion Curves */

= The Free Energy and Thermal Expansion of MgO =

== Introduction ==

==== Aims of Computational Lab ====

== Experimental ==

==== Methodology Used ====

==== Tools in Use (computer programmes) ====

== Results and Discussion ==

==== Representing the MgO Structure ====
[[File:BN_MgO_snapshot_of_lattice.PNG|none|thumb|249x249px|Figure to show Snapshot of MgO lattice]]

==== Phonon Dispersion Curves ====
[[File:BN calculating the phonon modes of MgO white background.JPG|none|thumb|411x411px|Figure to show the Phonon Dispersion Curve Computed for MgO]]
[[File:BN 1x1x1 calculating the phonon modes of MgO desnity of states 4 peaks white background.JPG|none|thumb|407x407px|Figure to show MgO DOS for 1x1x1]]

==== Grid Size Optimisation ====
[[File:BN Summary of MgO DOS for dimension variations.jpg|none|thumb|1072x1072px|Figure to show DOS Graphs for MgO using various Grid Sizes ]]
<blockquote></blockquote>
{| class="wikitable"
|+Table to show the Optimisation of Grid Size
!Dimensions
!<nowiki> </nowiki>Helmholz Free Energy (eV)
|-
|<nowiki> </nowiki>1x1x1
|<nowiki> </nowiki>-40.330301
|-
|<nowiki> </nowiki>2x2x2
|<nowiki> </nowiki>-40.926609
|-
|<nowiki> </nowiki>3x3x3
|<nowiki> </nowiki>-40.926432
|-
|<nowiki> </nowiki>4x4x4
|<nowiki> </nowiki>-40.926450
|-
|<nowiki> </nowiki>5x5x5
|<nowiki> </nowiki>-40.926463
|-
|<nowiki> </nowiki>6x6x6
|<nowiki> </nowiki>-40.926471
|-
|<nowiki> </nowiki>7x7x7
|<nowiki> </nowiki>-40.926475
|-
|<nowiki> </nowiki>8x8x8
|<nowiki> </nowiki>-40.926478
|-
|<nowiki> </nowiki>16x16x16
|<nowiki> </nowiki>-40.926482
|-
|<nowiki> </nowiki>32x32x32
|<nowiki> </nowiki>-40.926483
|-
|<nowiki> </nowiki>64x64x64
|<nowiki> </nowiki>-40.926483
|}

==== The Quasi-Harmonic Approximation ====
[[File:Summary of graphs to show Free energy, lattice constant and cell volume against temperature.JPG|none|thumb|1448x1448px|(From Left to Right) a) Figure to show Graph of Free Energy of MgO against Temperature, b) Graph of Lattice Constant of MgO against Temperature, c) Graph of Cell Volume against Temperature]]

==== Molecular Dynamics Calculations ====
[[File:BN Graph to compare LD with MD.JPG|none|thumb|734x734px|Figure to show Graph for Comparison of the LD (Quasi-Harmonic) and MD (Molecular Dynamics) Approach]]

== Conclusion ==

== References ==

===== 1. Calculating the internal energy of an MgO crystal =====

===== 2. Lattice Vibrations - Computing the Phonons =====

===== 3. Calculating the Free Energy in the Harmonic Approximation =====

===== 4. The Thermal Expansion of MgO =====

===== 5. Molecular Dynamics =====

File:BN Summary of MgO DOS for dimension variations.jpg

2017-03-16T15:11:56Z

Ben14: