2015-11-20T02:41:27Z

Xl7712:

=Introduction=

Magnesium oxide naturally exists as crystal based on face-centred cubic lattice with the lattice points taken by Mg2+ and the octahedral holes filled with O2-.
[[File:MgO unit cell X.png | right | x300px | 300px | thumb | Fig.1 conventional cell ( in black ) and primitive cell ( in light blue ) of MgO; The volume of primitive cell is 1/4 of that of conventional cell.]]
{| class="wikitable"
! Type of Unit Cell
! Shape
! Parameter
! Internal Angel
! Volume
! Number of MgO
|-
| Conventional
| Cube
| ac = bc = cc = 4.212 ‎Å
| αc = βc = γc = 90•
| Vc = 74.725 Å3
| Nc = 4
|-
| Primitive
| Rhombohedron
| ap = bp = cp = 2.978 ‎Å
| αp = βp = γp = 60•
| Vp = 18.6812 Å3
| Np = 1
|}

Vibrations of a solid system are related to many of its physical properties such as free energy, heat capacity, expansion, phase transition, thermal conductivity and dielectric phenomena at low frequencies. This study compares two methods for simulation of MgO crystal vibrations. '''Quasi-harmonic Approximation''' considers vibrations as phonons representing elementary vibrational modes in which a lattice of particles uniformly oscillates at a single frequency. '''Molecular Dynamics''' allows the particle in the system to interact for a given time period and the coordinates of the particles are numerically solved based on Newton's Laws ''Rfinal = Rinitial + vfinal*dt = Rinitial + vinitial*dt + a*dt2 = Rinitial + vinitial*dt + (F/m)*dt2''. Both methods were conducted on Linux based programme GULP (General Utility Lattice Program) via the user interface for constructing and visualizing provided by DL Visualize.

=Theory=

In statistical mechanics, the physical properties of a system are in Boltzmann Distribution '''''ni / N = exp (-βui) / q''''' where '''''β = 1 / (kbT)''''' and '''''q = Σj=1levels exp (-βuj)'''''. This means that once the partition function q is correctly expressed, the properties of the system can be calculated. 要舉個例子嗎？In this experiment, in accordance with harmonic oscillation model, the vibrational frequency ω must be quantised and summing over the frequencies will lead to the partition function.

MgO crystal is made of repeating unit cells, so it is sensible to start with the simplest model first to see how frequency ω is related to the repeating structure. When 1-dimensional chain of one kind of atom vibrate, they can have several different types of vibrations and each one can be described as a wave with a wavelength equal to the length of the repeating unit (Fig.2) and plotting the vibrational frequencies VS the k-vectors (showing directions and wavelengths of vibrations) gives a graph like Fig.3. If each atom in this chain is superseded by a MgO, there is now a pair of ions in each repeat unit, a' = 2a, hence -π/(2a) < k < π/(2a) and folding branch occurs (Fig.4).

{| class="wikitable"
! [[File:1D vibrations L.png | x400px | 400px | thumb | Fig1.2 1D vibrations; wave vector '''''k = 2π / λ'''''; -π/2 < k < π/2; frequency ω increases as k increases.]]
! [[File:Typical dispersion curve X.png | x250px | 250px | thumb | Fig1.3 typical ω(k) is plotted as a dispersion curve; k = 0 at Γ point ]]
! [[File:1D Diatomic Chain X.png | x280px | 280px | thumb | Fig1.4 1D MgO chain; There exist two frequencies for one k-vector due to different oscillation between Mg2+ and O2-.]]
|}

Both structures mentioned above are limited in 1 dimension. When this is expanded to two dimension, particles can also vibrate up and down with respect to the horizontal axis, hence k-vectors are expressed as (kx, ky) in Cartesian coordinate system, and the ω(k) plot becomes a dispersion surface with frequency ω showed in z-axis. It is now easy to see that for 3-dimensional MgO crystal, k points includes (kx, ky, kz), and there will be four Cartesian axises for a ω(k) plot, which is not able to show in real life. In this case, a certain path in the 3-dimensional solid is set and the coordinates through the path were set as the k points, thus ω(k) can be plotted against the path and it is again back to the dispersion curve.

Once the all the vibrational branches are obtained, sum over them to form the partition function and the vibrational energy levels can be computed.

=Results and Discussion=

[[File:Phonon dispersion X.png | thumb | right | x350px | 350px | Fig.5 Dispersion Curve of MgO lattice vibrations; path W-L-G-W-X-K on horizontal axis with coordinates shown in blue]]

The lattice energy of MgO calculated is -41.075 eV, and this is the potential energy holding the lattice together induced by electrostatic interaction between Mg2+ and O2- ions, which means to move all the ions in the lattice apart to infinity requires an energy of 41.075 eV. Also, this equals to the internal energy of an ideal MgO lattice as perfect crystals have no vibrations.

As mentioned in the Theory part, to understand the variation of frequencies with k, a dispersion curve is essential. To deal with the 3-dimensional MgO infinite lattice, a conventional path in the k-space is used to compute the vibrational modes, and for Fig.5, 50 points along the path was computed and shows all the phonon modes.

The strategy to sum up the phonon modes is to construct the Density of Sate (DOS), indicating the probability of a phonon to be in a certain frequency. It is important to sum up phonons for an adequate number of k points so that the distribution of them can be represent the distribution of phonos of an infinite lattice. The following shows the process of finding the best number of k points for DOS.

{| class="wikitable"
! Shrinking Factor
! Density of State
! Comments
|-
| 1
| [[File:Phonon DOS 1 X.png|thumb|x250px|250px]]
|The DOS for 1x1x1 grid was computed for a single k-point 'L'. There are 4 distinct peaks, the two around 300 and 400 cm-1 is double in intensity compared to the 700 and 800 cm-1, which is corresponding to the two branches across point 'L' shown in the dispersion curve.
summing over more k points to get a distribution
|-
| 2
| [[File:Phonon DOS 2 X.png|thumb|x250px|250px]]
| 7 distinct peaks, highest peak at 400 cm-1
lower density of each peaks compared to the DOS above
A grid size of 2*2*2 is sufficient to get the correct highest value.
summing over more k points to get a distribution
|-
| 4
| [[File:Phonon DOS 4 X.png|thumb|x250px|250px]]
| more peaks, highest peak at 400 cm-1 with lowered density
not a smooth distribution yet
|-
| 8
| [[File:Phonon DOS 8 X.png|thumb|x250px|250px]]
| highest peak at 400 cm-1 with even lowered density
Distribution features appear.
large fluctuations
|-
| 16
| [[File:Phonon DOS 16 X.png|thumb|x250px|250px]]
| a general shape of the distribution with small fluctuations
no much change in the densities
|-
| 32
| [[File:Phonon DOS 32 X.png|thumb|x250px|250px]]
| nice distribution compared to the upper one
acceptable computation time
|-
| 64
| [[File:Phonon DOS 64 X.png|thumb|x250px|250px]]
| There is no much change compare to the DOS computed along grid 32*32*32, however takes a few minutes longer to compute.
The DOS is converging which means it is converging.
|}

As the grid size increases, more possible vibrations are sampled and the distribution is smoothened, nevertheless, the change in DOS decreases each time the grid size is doubled. Computing over more k-points requires more resources and time, which is obvious from grid 32*32*32 to 64*64*64. A compromise can be grid 32*32*32 which can give a good enough distribution which is a close approximation to the infinite lattice economically.

Since there is a way to compute all the phonon modes in MgO infinite lattice, the free energy of it can also be calculated. GULP searches for the minimum free energy with respect to the structure via calculating the internal energy and phonons at a sequence of geometries. Similarly, the computing path is the same as that for computing DOS, so there is also the grid size problem.

{| class="wikitable"
! Shrinking Factor
! Helmholtz Free Energy (eV)
! Accuracy
|-
| 1
| - 40.930301
|
|-
| 2
| - 40.926609
| 1 meV
|-
| 3
| - 40.926432
| 0.5 meV
|-
| 4
| - 40.926450
|
|-
| 5
| - 40.926463
| 0.1 meV
|-
| 6
| - 40.926471
|
|-
| 7
| - 40.926475
| 0.01 meV
|-
| 8
| - 40.926478
|
|-
| 9
| - 40.926479
|
|-
| 10
| - 40.926480
|
|-
| 11
| - 40.926481
| convergence
|-
| 12
| - 40.926481
|
|}

A grid size of 11*11*11 (i.e. 0.01 meV accuracy) was chosen for the following calculations based on Quasi-harmonic Approximation. The Helmholtz free energy and the cell colume were optimised to observe the variations with different temperature. As temperature is raising, the Helmholtz Free Energy becomes more negative, while the cell volume is expanding, Both of the variations can be well expressed by polynomial equations. Calculation failed when temperature is close to the melting point of MgO, 3125 K {http://www.rsc.org/chemistryworld/2014/08/magnesium-oxide-mgo-podcast} The reason can be the vibrations is so large that atoms clashes into each other causing computing errors.

[[File:A VS T X.png | thumb | x500px | 500px | Fig.5 Helmholtz Free Energy VS Temp. Quasi-harmonic Approximation]]
[[File:QHA V VS T X.png | thumb | x500px | 500px | Fig.6 Cell Volume VS Temp. calculated by Quasi-harmonic Approximation]]
{| class="wikitable"
! Temperature (K)
! Helmholtz Free Energy A (eV)
! Lattice Volume (Å3)
|-
| 0
| -40.9019
| 18.8365
|-
| 100
| -40.9024
|18.8383
|-
| 200
| -40.9094
|18.8562
|-
| 300
| -40.9281
|18.8900
|-
| 400
| -40.9586
|18.9325
|-
| 500
| -40.9994
|18.9801
|-
| 600
| -41.0493
|19.0312
|-
| 700
| -41.1071
|19.0851
|-
| 800
| -41.1719
|19.1413
|-
| 900
| -41.2430
|19.1997
|-
| 1000
| -41.3110
|19.2601
|-
| 1200
| -41.4887
|19.3872
|-
| 1400
| -41.6755
|19.5233
|-
| 1600
| -41.8780
|19.6698
|-
| 1800
| -42.0944
|19.8287
|-
| 2000
| -42.3237
|20.0029
|-
| 2300
| -42.6895
|20.3047
|-
| 2600
| -43.0800
|20.6889
|-
| 2900
| -43.4948
|21.3217
|}

In order to compare with literature thermal expansion, several volumes were calculated by substituting some specific temperatures into the trend line equation in Fig.6. to get the predicted cell volumes. The cell volumes were then transferred into molar volume by multiplying Avogadro constant NA, the units were also changed to what is used in the literature.
{| class="wikitable"
! Temperature (K)
! Cell Volume (Å3)
! Molar Volume (cm-3mol-1)
! Literature Molar Volume (cm-3mol-1)
|-
|298
| 18.8851
| 11.3688
|11.2434
|-
|455
| 18.9425
| 11.4034
|11.3004
|-
|710
| 19.0601
| 11.4742
|11.4109
|-
|1096
| 19.2570
| 11.5927
|11.6211
|-
|1527
| 19.4673
| 11.7193
|11.8218
|-
|2106
| 19.7624
| 11.8970
|12.2287
|-
|2703
| 20.3014
| 12.2214
|12.6887
|-
|2986
| 20.7658
| 12.5010
|12.9244
|-
|3015
| 20.8248
| 12.5365
|12.9723
|}

Generally, the trend and the magnitude are consistent between the predicted values and literature values, but the slope of the literature value is steeper than that of the predicted values, so before 多少度 the predicted values are slightly higher, they come together at aldjfhalkdfjhjh, afterwards这里要不要plot个图啊？

[[File:V VS T both X.png | thumb | x400px | 400px | Fig.6 suitable V values from Quasi-harmonic Approx. and Melecular Dynamics plotted in one graph]]

As mentioned in the introduction part, another simulation method called Molecular Dynamics was also used to calculate the equilibrium energy and volume. The MD obtained values before the melting point of MgO are similar compared to those obtained by Quasi-harmonic Approximation, and after the m.p. MD can compensate the failure of Quasi. When around m.p, there is a range where the volume almost constant, indicating phase changing. When T reaches 4000 K, the volume is lifted by more than 10 Angstrom higher, phase change is completed and the volume of liquid phase will continue increase with the raising temperature but with a steeper gradient. Atkins ref?????? If the temperature goes on increasing, the volume will become infinite as the gas phase does not have a volume without any pressure.

The change of cell volume can be describes as thermal expansion coefficient α=(1/V0)(dV/dT). This property can be calculated for both of data sets obtained from both methods, and v0 is the zero-point volume in each case and the (dV/dT) is the gradient of the trend lines.
[[File:Expansion coefficient both X.png|thumb|x400px|400px|Fig.7 Comparison between αV obtained from different methods. It can be concluded that the difference is decreasing from low temperature to 1500 K after which the difference in expansion coefficient tends to keep constant.]]

The properties of materials (solids, liquids, gasses) are a statistical average over the many different energy states of the molecules making up the material. The vibrational free energy of H2 can be computed analytically by summing over the harmonic vibrations of the molecule. This cannot be done by hand for a real material containing many atoms.

In this laboratory you will use a simple model of atomic interactions to calculate the energy and vibrations of a crystal of MgO. These vibrational energy levels will then be used to compute the free energy of the crystal and to predict how the material expands when heated. In the last final stage you will go beyond the harmonic (and quasi-harmonic) approximation and expand the crystal using a technique called molecular dynamics - essentially reproducing the actual vibration motions of the atoms. Fortunately the computer will do most of the work !'''
Vibrations are quantised and can be seen as particles called phonons. The energy of phonon is expressed in XXXX equation.
When the ensemble is large enough, the
Vibrational energies of H2 are quantised with the expression εvib = (n+1/2)ћω

Body of the text
write it like a scientific paper (well-articulated sentences,
NOT a list of two-word answers)
analyse critically obtained data and given answers
round numerical answers to a specific number of decimal
places (i.e. 4)
add literature/web citations whenever a comparison with
experimental data is required
add explicitely every formu

Pictures
max 20
reasonably sized (NOT one-page sized pictures, but still
readable)
white background (follow the instructions given on the
website clicking on the link ’How to save a picture for your
report’)
described in caption or in the text

Graphs
add labels and units
add a critical comment whenever required (NOT a merely
descriptive comment)la used one to obtain results
check spelling
Tables
add labels and units
round numerical answers to a specific number of decimal
places (i.e. 4)
repeat heading if the table cannot fit in a single page

Conclusions
give a general description of your calculations and your
main findings
outline the differences between the methods in use and the
results obtained
analyse critically these differences

2015-11-19T20:58:23Z

Xl7712:

Rep:MgO:XL

2015-11-19T20:49:31Z

Xl7712:

=Introduction=

Magnesium oxide naturally exists as crystal based on face-centred cubic lattice with the lattice points taken by Mg2+ and the octahedral holes filled with O2-.
[[File:MgO unit cell X.png | right | x300px | 300px | thumb | Fig.1 conventional cell ( in black ) and primitive cell ( in light blue ) of MgO; The volume of primitive cell is 1/4 of that of conventional cell.]]
{| class="wikitable"
! Type of Unit Cell
! Shape
! Parameter
! Internal Angel
! Volume
! Number of MgO
|-
| Conventional
| Cube
| ac = bc = cc = 4.212 ‎Å
| αc = βc = γc = 90•
| Vc = 74.725 Å3
| Nc = 4
|-
| Primitive
| Rhombohedron
| ap = bp = cp = 2.978 ‎Å
| αp = βp = γp = 60•
| Vp = 18.6812 Å3
| Np = 1
|}

Vibrations of a solid system are related to many of its physical properties such as free energy, heat capacity, expansion, phase transition, thermal conductivity and dielectric phenomena at low frequencies. This study compares two methods for simulation of MgO crystal vibrations. '''Quasi-harmonic Approximation''' considers vibrations as phonons representing elementary vibrational modes in which a lattice of particles uniformly oscillates at a single frequency. '''Molecular Dynamics''' allows the particle in the system to interact for a given time period and the coordinates of the particles are numerically solved based on Newton's Laws ''Rfinal = Rinitial + vfinal*dt = Rinitial + vinitial*dt + a*dt2 = Rinitial + vinitial*dt + (F/m)*dt2''. Both methods were conducted on Linux based programme GULP (General Utility Lattice Program) via the user interface for constructing and visualizing provided by DL Visualize.

=Theory=

In statistical mechanics, the physical properties of a system are in Boltzmann Distribution '''''ni / N = exp (-βui) / q''''' where '''''β = 1 / (kbT)''''' and '''''q = Σj=1levels exp (-βuj)'''''. This means that once the partition function q is correctly expressed, the properties of the system can be calculated. 要舉個例子嗎？In this experiment, in accordance with harmonic oscillation model, the vibrational frequency ω must be quantised and summing over the frequencies will lead to the partition function.

MgO crystal is made of repeating unit cells, so it is sensible to start with the simplest model first to see how frequency ω is related to the repeating structure. When 1-dimensional chain of one kind of atom vibrate, they can have several different types of vibrations and each one can be described as a wave with a wavelength equal to the length of the repeating unit (Fig.2) and plotting the vibrational frequencies VS the k-vectors (showing directions and wavelengths of vibrations) gives a graph like Fig.3. If each atom in this chain is superseded by a MgO, there is now a pair of ions in each repeat unit, a' = 2a, hence -π/(2a) < k < π/(2a) and folding branch occurs (Fig.4).

{| class="wikitable"
! [[File:1D vibrations L.png | x400px | 400px | thumb | Fig1.2 1D vibrations; wave vector '''''k = 2π / λ'''''; -π/2 < k < π/2; frequency ω increases as k increases.]]
! [[File:Typical dispersion curve X.png | x250px | 250px | thumb | Fig1.3 typical ω(k) is plotted as a dispersion curve; k = 0 at Γ point ]]
! [[File:1D Diatomic Chain X.png | x280px | 280px | thumb | Fig1.4 1D MgO chain; There exist two frequencies for one k-vector due to different oscillation between Mg2+ and O2-.]]
|}

Both structures mentioned above are limited in 1 dimension. When this is expanded to two dimension, particles can also vibrate up and down with respect to the horizontal axis, hence k-vectors are expressed as (kx, ky) in Cartesian coordinate system, and the ω(k) plot becomes a dispersion surface with frequency ω showed in z-axis. It is now easy to see that for 3-dimensional MgO crystal, k points includes (kx, ky, kz), and there will be four Cartesian axises for a ω(k) plot, which is not able to show in real life. In this case, a certain path in the 3-dimensional solid is set and the coordinates through the path were set as the k points, thus ω(k) can be plotted against the path and it is again back to the dispersion curve.

Once the all the vibrational branches are obtained, sum over them to form the partition function and the vibrational energy levels can be computed.

=Results and Discussion=

[[File:Phonon dispersion X.png | thumb | right | x350px | 350px | Fig.5 Dispersion Curve of MgO lattice vibrations; path W-L-G-W-X-K on horizontal axis with coordinates shown in blue]]

The lattice energy of MgO calculated is -41.075 eV, and this is the potential energy holding the lattice together induced by electrostatic interaction between Mg2+ and O2- ions, which means to move all the ions in the lattice apart to infinity requires an energy of 41.075 eV. Also, this equals to the internal energy of an ideal MgO lattice as perfect crystals have no vibrations.

As mentioned in the Theory part, to understand the variation of frequencies with k, a dispersion curve is essential. To deal with the 3-dimensional MgO infinite lattice, a conventional path in the k-space is used to compute the vibrational modes, and for Fig.5, 50 points along the path was computed and shows all the phonon modes.

The strategy to sum up the phonon modes is to construct the Density of Sate (DOS), indicating the probability of a phonon to be in a certain frequency. It is important to sum up phonons for an adequate number of k points so that the distribution of them can be represent the distribution of phonos of an infinite lattice. The following shows the process of finding the best number of k points for DOS.

{| class="wikitable"
! Shrinking Factor
! Density of State
! Comments
|-
| 1
| [[File:Phonon DOS 1 X.png|thumb|x250px|250px]]
|The DOS for 1x1x1 grid was computed for a single k-point 'L'. There are 4 distinct peaks, the two around 300 and 400 cm-1 is double in intensity compared to the 700 and 800 cm-1, which is corresponding to the two branches across point 'L' shown in the dispersion curve.
summing over more k points to get a distribution
|-
| 2
| [[File:Phonon DOS 2 X.png|thumb|x250px|250px]]
| 7 distinct peaks, highest peak at 400 cm-1
lower density of each peaks compared to the DOS above
A grid size of 2*2*2 is sufficient to get the correct highest value.
summing over more k points to get a distribution
|-
| 4
| [[File:Phonon DOS 4 X.png|thumb|x250px|250px]]
| more peaks, highest peak at 400 cm-1 with lowered density
not a smooth distribution yet
|-
| 8
| [[File:Phonon DOS 8 X.png|thumb|x250px|250px]]
| highest peak at 400 cm-1 with even lowered density
Distribution features appear.
large fluctuations
|-
| 16
| [[File:Phonon DOS 16 X.png|thumb|x250px|250px]]
| a general shape of the distribution with small fluctuations
no much change in the densities
|-
| 32
| [[File:Phonon DOS 32 X.png|thumb|x250px|250px]]
| nice distribution compared to the upper one
acceptable computation time
|-
| 64
| [[File:Phonon DOS 64 X.png|thumb|x250px|250px]]
| There is no much change compare to the DOS computed along grid 32*32*32, however takes a few minutes longer to compute.
The DOS is converging which means it is converging.
|}

As the grid size increases, more possible vibrations are sampled and the distribution is smoothened, nevertheless, the change in DOS decreases each time the grid size is doubled. Computing over more k-points requires more resources and time, which is obvious from grid 32*32*32 to 64*64*64. A compromise can be grid 32*32*32 which can give a good enough distribution which is a close approximation to the infinite lattice economically.

Since there is a way to compute all the phonon modes in MgO infinite lattice, the free energy of it can also be calculated. GULP searches for the minimum free energy with respect to the structure via calculating the internal energy and phonons at a sequence of geometries. Similarly, the computing path is the same as that for computing DOS, so there is also the grid size problem.

{| class="wikitable"
! Shrinking Factor
! Helmholtz Free Energy (eV)
! Accuracy
|-
| 1
| - 40.930301
|
|-
| 2
| - 40.926609
| 1 meV
|-
| 3
| - 40.926432
| 0.5 meV
|-
| 4
| - 40.926450
|
|-
| 5
| - 40.926463
| 0.1 meV
|-
| 6
| - 40.926471
|
|-
| 7
| - 40.926475
| 0.01 meV
|-
| 8
| - 40.926478
|
|-
| 9
| - 40.926479
|
|-
| 10
| - 40.926480
|
|-
| 11
| - 40.926481
| convergence
|-
| 12
| - 40.926481
|
|}

A grid size of 11*11*11 (i.e. 0.01 meV accuracy) was chosen for the following calculations based on Quasi-harmonic Approximation.

[[File:A VS T X.png | thumb | x500px | 500px | Fig.5 Helmholtz Free Energy VS Temperature calculated by Quasi-harmonic Approximation. Calculation failed at 3100 K which is close the melting point 3125 K {http://www.rsc.org/chemistryworld/2014/08/magnesium-oxide-mgo-podcast} of MgO]]
{| class="wikitable"
! Temperature (K)
! Helmholtz Free Energy A (eV)
! Lattice Volume (Å3</sip>)
|-
| 0
| -40.9019
| 18.8365
|-
| 100
| -40.9024
|18.8383
|-
| 200
| -40.9094
|18.8562
|-
| 300
| -40.9281
|18.8900
|-
| 400
| -40.9586
|18.9325
|-
| 500
| -40.9994
|18.9801
|-
| 600
| -41.0493
|19.0312
|-
| 700
| -41.1071
|19.0851
|-
| 800
| -41.1719
|19.1413
|-
| 900
| -41.2430
|19.1997
|-
| 1000
| -41.3110
|19.2601
|-
| 1200
| -41.4887
|19.3872
|-
| 1400
| -41.6755
|19.5233
|-
| 1600
| -41.8780
|19.6698
|-
| 1800
| -42.0944
|19.8287
|-
| 2000
| -42.3237
|20.0029
|-
| 2300
| -42.6895
|20.3047
|-
| 2600
| -43.0800
|20.6889
|-
| 2900
| -43.4948
|21.3217
|-
| 3100
| calculation failed
| the reason can be the vibrations is so large that atoms clashes into each other causing computing errors
|}

The properties of materials (solids, liquids, gasses) are a statistical average over the many different energy states of the molecules making up the material. The vibrational free energy of H2 can be computed analytically by summing over the harmonic vibrations of the molecule. This cannot be done by hand for a real material containing many atoms.

In this laboratory you will use a simple model of atomic interactions to calculate the energy and vibrations of a crystal of MgO. These vibrational energy levels will then be used to compute the free energy of the crystal and to predict how the material expands when heated. In the last final stage you will go beyond the harmonic (and quasi-harmonic) approximation and expand the crystal using a technique called molecular dynamics - essentially reproducing the actual vibration motions of the atoms. Fortunately the computer will do most of the work !'''
Vibrations are quantised and can be seen as particles called phonons. The energy of phonon is expressed in XXXX equation.
When the ensemble is large enough, the
Vibrational energies of H2 are quantised with the expression εvib = (n+1/2)ћω

Body of the text
write it like a scientific paper (well-articulated sentences,
NOT a list of two-word answers)
analyse critically obtained data and given answers
round numerical answers to a specific number of decimal
places (i.e. 4)
add literature/web citations whenever a comparison with
experimental data is required
add explicitely every formu

Pictures
max 20
reasonably sized (NOT one-page sized pictures, but still
readable)
white background (follow the instructions given on the
website clicking on the link ’How to save a picture for your
report’)
described in caption or in the text

Graphs
add labels and units
add a critical comment whenever required (NOT a merely
descriptive comment)la used one to obtain results
check spelling
Tables
add labels and units
round numerical answers to a specific number of decimal
places (i.e. 4)
repeat heading if the table cannot fit in a single page

Conclusions
give a general description of your calculations and your
main findings
outline the differences between the methods in use and the
results obtained
analyse critically these differences

Rep:MgO:XL

2015-11-19T20:36:04Z

Xl7712:

=Introduction=

Magnesium oxide naturally exists as crystal based on face-centred cubic lattice with the lattice points taken by Mg2+ and the octahedral holes filled with O2-.
[[File:MgO unit cell X.png | right | x300px | 300px | thumb | Fig.1 conventional cell ( in black ) and primitive cell ( in light blue ) of MgO; The volume of primitive cell is 1/4 of that of conventional cell.]]
{| class="wikitable"
! Type of Unit Cell
! Shape
! Parameter
! Internal Angel
! Volume
! Number of MgO
|-
| Conventional
| Cube
| ac = bc = cc = 4.212 ‎Å
| αc = βc = γc = 90•
| Vc = 74.725 Å3
| Nc = 4
|-
| Primitive
| Rhombohedron
| ap = bp = cp = 2.978 ‎Å
| αp = βp = γp = 60•
| Vp = 18.6812 Å3
| Np = 1
|}

Vibrations of a solid system are related to many of its physical properties such as free energy, heat capacity, expansion, phase transition, thermal conductivity and dielectric phenomena at low frequencies. This study compares two methods for simulation of MgO crystal vibrations. '''Quasi-harmonic Approximation''' considers vibrations as phonons representing elementary vibrational modes in which a lattice of particles uniformly oscillates at a single frequency. '''Molecular Dynamics''' allows the particle in the system to interact for a given time period and the coordinates of the particles are numerically solved based on Newton's Laws ''Rfinal = Rinitial + vfinal*dt = Rinitial + vinitial*dt + a*dt2 = Rinitial + vinitial*dt + (F/m)*dt2''. Both methods were conducted on Linux based programme GULP (General Utility Lattice Program) via the user interface for constructing and visualizing provided by DL Visualize.

=Theory=

In statistical mechanics, the physical properties of a system are in Boltzmann Distribution '''''ni / N = exp (-βui) / q''''' where '''''β = 1 / (kbT)''''' and '''''q = Σj=1levels exp (-βuj)'''''. This means that once the partition function q is correctly expressed, the properties of the system can be calculated. 要舉個例子嗎？In this experiment, in accordance with harmonic oscillation model, the vibrational frequency ω must be quantised and summing over the frequencies will lead to the partition function.

MgO crystal is made of repeating unit cells, so it is sensible to start with the simplest model first to see how frequency ω is related to the repeating structure. When 1-dimensional chain of one kind of atom vibrate, they can have several different types of vibrations and each one can be described as a wave with a wavelength equal to the length of the repeating unit (Fig.2) and plotting the vibrational frequencies VS the k-vectors (showing directions and wavelengths of vibrations) gives a graph like Fig.3. If each atom in this chain is superseded by a MgO, there is now a pair of ions in each repeat unit, a' = 2a, hence -π/(2a) < k < π/(2a) and folding branch occurs (Fig.4).

{| class="wikitable"
! [[File:1D vibrations L.png | x400px | 400px | thumb | Fig1.2 1D vibrations; wave vector '''''k = 2π / λ'''''; -π/2 < k < π/2; frequency ω increases as k increases.]]
! [[File:Typical dispersion curve X.png | x250px | 250px | thumb | Fig1.3 typical ω(k) is plotted as a dispersion curve; k = 0 at Γ point ]]
! [[File:1D Diatomic Chain X.png | x280px | 280px | thumb | Fig1.4 1D MgO chain; There exist two frequencies for one k-vector due to different oscillation between Mg2+ and O2-.]]
|}

Both structures mentioned above are limited in 1 dimension. When this is expanded to two dimension, particles can also vibrate up and down with respect to the horizontal axis, hence k-vectors are expressed as (kx, ky) in Cartesian coordinate system, and the ω(k) plot becomes a dispersion surface with frequency ω showed in z-axis. It is now easy to see that for 3-dimensional MgO crystal, k points includes (kx, ky, kz), and there will be four Cartesian axises for a ω(k) plot, which is not able to show in real life. In this case, a certain path in the 3-dimensional solid is set and the coordinates through the path were set as the k points, thus ω(k) can be plotted against the path and it is again back to the dispersion curve.

Once the all the vibrational branches are obtained, sum over them to form the partition function and the vibrational energy levels can be computed.

=Results and Discussion=

[[File:Phonon dispersion X.png | thumb | right | x350px | 350px | Fig.5 Dispersion Curve of MgO lattice vibrations; path W-L-G-W-X-K on horizontal axis with coordinates shown in blue]]

The lattice energy of MgO calculated is -41.075 eV, and this is the potential energy holding the lattice together induced by electrostatic interaction between Mg2+ and O2- ions, which means to move all the ions in the lattice apart to infinity requires an energy of 41.075 eV. Also, this equals to the internal energy of an ideal MgO lattice as perfect crystals have no vibrations.

As mentioned in the Theory part, to understand the variation of frequencies with k, a dispersion curve is essential. To deal with the 3-dimensional MgO infinite lattice, a conventional path in the k-space is used to compute the vibrational modes, and for Fig.5, 50 points along the path was computed and shows all the phonon modes.

The strategy to sum up the phonon modes is to construct the Density of Sate (DOS), indicating the probability of a phonon to be in a certain frequency. It is important to sum up phonons for an adequate number of k points so that the distribution of them can be represent the distribution of phonos of an infinite lattice. The following shows the process of finding the best number of k points for DOS.

{| class="wikitable"
! Shrinking Factor
! Density of State
! Comments
|-
| 1
| [[File:Phonon DOS 1 X.png|thumb|x250px|250px]]
|The DOS for 1x1x1 grid was computed for a single k-point 'L'. There are 4 distinct peaks, the two around 300 and 400 cm-1 is double in intensity compared to the 700 and 800 cm-1, which is corresponding to the two branches across point 'L' shown in the dispersion curve.
summing over more k points to get a distribution
|-
| 2
| [[File:Phonon DOS 2 X.png|thumb|x250px|250px]]
| 7 distinct peaks, highest peak at 400 cm-1
lower density of each peaks compared to the DOS above
A grid size of 2*2*2 is sufficient to get the correct highest value.
summing over more k points to get a distribution
|-
| 4
| [[File:Phonon DOS 4 X.png|thumb|x250px|250px]]
| more peaks, highest peak at 400 cm-1 with lowered density
not a smooth distribution yet
|-
| 8
| [[File:Phonon DOS 8 X.png|thumb|x250px|250px]]
| highest peak at 400 cm-1 with even lowered density
Distribution features appear.
large fluctuations
|-
| 16
| [[File:Phonon DOS 16 X.png|thumb|x250px|250px]]
| a general shape of the distribution with small fluctuations
no much change in the densities
|-
| 32
| [[File:Phonon DOS 32 X.png|thumb|x250px|250px]]
| nice distribution compared to the upper one
acceptable computation time
|-
| 64
| [[File:Phonon DOS 64 X.png|thumb|x250px|250px]]
| There is no much change compare to the DOS computed along grid 32*32*32, however takes a few minutes longer to compute.
The DOS is converging which means it is converging.
|}

As the grid size increases, more possible vibrations are sampled and the distribution is smoothened, nevertheless, the change in DOS decreases each time the grid size is doubled. Computing over more k-points requires more resources and time, which is obvious from grid 32*32*32 to 64*64*64. A compromise can be grid 32*32*32 which can give a good enough distribution which is a close approximation to the infinite lattice economically.

Since there is a way to compute all the phonon modes in MgO infinite lattice, the free energy of it can also be calculated. GULP searches for the minimum free energy with respect to the structure via calculating the internal energy and phonons at a sequence of geometries. Similarly, the computing path is the same as that for computing DOS, so there is also the grid size problem.

{| class="wikitable"
! Shrinking Factor
! Helmholtz Free Energy (eV)
! Accuracy
|-
| 1
| - 40.930301
|
|-
| 2
| - 40.926609
| 1 meV
|-
| 3
| - 40.926432
| 0.5 meV
|-
| 4
| - 40.926450
|
|-
| 5
| - 40.926463
| 0.1 meV
|-
| 6
| - 40.926471
|
|-
| 7
| - 40.926475
| 0.01 meV
|-
| 8
| - 40.926478
|
|-
| 9
| - 40.926479
|
|-
| 10
| - 40.926480
|
|-
| 11
| - 40.926481
| convergence
|-
| 12
| - 40.926481
|
|}

A grid size of 11*11*11 (i.e. 0.01 meV accuracy) was chosen for the following calculations based on Quasi-harmonic Approximation.

[[File:A VS T X.png | thumb | x400px | 400px | Fig.5 Helmholtz Free Energy VS Temperature calculated by Quasi-harmonic Approximation]]
{| class="wikitable"
! Temperature (K)
! Helmholtz Free Energy A (eV)
|-
| 0
| -40.9019
|-
| 100
| -40.9024
|-
| 200
| -40.9094
|-
| 300
| -40.9281
|-
| 400
| -40.9586
|-
| 500
| -40.9994
|-
| 600
| -41.0493
|-
| 700
| -41.1071
|-
| 800
| -41.1719
|-
| 900
| -41.2430
|-
| 1000
| -41.3110
|-
| 1200
| -41.4887
|-
| 1400
| -41.6755
|-
| 1600
| -41.8780
|-
| 1800
| -42.0944
|-
| 2000
| -42.3237
|-
| 2300
| -42.6895
|-
| 2600
| -43.0800
|-
| 2900
| -43.4948
|-
| 3100
| calculation failed
|}

The properties of materials (solids, liquids, gasses) are a statistical average over the many different energy states of the molecules making up the material. The vibrational free energy of H2 can be computed analytically by summing over the harmonic vibrations of the molecule. This cannot be done by hand for a real material containing many atoms.

In this laboratory you will use a simple model of atomic interactions to calculate the energy and vibrations of a crystal of MgO. These vibrational energy levels will then be used to compute the free energy of the crystal and to predict how the material expands when heated. In the last final stage you will go beyond the harmonic (and quasi-harmonic) approximation and expand the crystal using a technique called molecular dynamics - essentially reproducing the actual vibration motions of the atoms. Fortunately the computer will do most of the work !'''
Vibrations are quantised and can be seen as particles called phonons. The energy of phonon is expressed in XXXX equation.
When the ensemble is large enough, the
Vibrational energies of H2 are quantised with the expression εvib = (n+1/2)ћω

Body of the text
write it like a scientific paper (well-articulated sentences,
NOT a list of two-word answers)
analyse critically obtained data and given answers
round numerical answers to a specific number of decimal
places (i.e. 4)
add literature/web citations whenever a comparison with
experimental data is required
add explicitely every formu

Pictures
max 20
reasonably sized (NOT one-page sized pictures, but still
readable)
white background (follow the instructions given on the
website clicking on the link ’How to save a picture for your
report’)
described in caption or in the text

Graphs
add labels and units
add a critical comment whenever required (NOT a merely
descriptive comment)la used one to obtain results
check spelling
Tables
add labels and units
round numerical answers to a specific number of decimal
places (i.e. 4)
repeat heading if the table cannot fit in a single page

Conclusions
give a general description of your calculations and your
main findings
outline the differences between the methods in use and the
results obtained
analyse critically these differences

2015-11-19T13:25:40Z

Xl7712:

File:Phonon DOS 64 X.png

2015-11-19T13:25:30Z

Xl7712: