The Phonon Thermodynamics of Iron and Ce

The Phonon Thermodynamics of Iron and Cementite - CaltechTHESIS
CaltechTHESIS
A Caltech Library Service
About
Browse
Deposit an Item
Instructions for Students
The Phonon Thermodynamics of Iron and Cementite
Citation
Mauger, Lisa Mary
(2015)
The Phonon Thermodynamics of Iron and Cementite.
Dissertation (Ph.D.), California Institute of Technology.
doi:10.7907/Z9TQ5ZH3.
Abstract
The vibrational properties of materials are essential to understanding material stability and thermodynamics. In this thesis I outline vibrational thermodynamic models and the experimental tools that provide evidence on phonon behavior. The introductory section discusses the history of metallurgy and thermodynamic theory, with an emphasis on the role of iron and cementite, two important components of steels. The thermodynamic framework for understanding vibrational material behavior is provided alongside the growing body of experimental and computational tools that provide physical insight on vibrational properties. The high temperature vibrational behavior of iron and cementite are explored within this context in the final chapters.
Body-centered-cubic iron exhibits decreasing phonon energies at elevated temperatures. The observed energy change in not uniform across phonon modes in iron, and specific phonon modes show significant decreases in energy that are not explained by simple vibrational models. This anomalously
energy decrease is linked to the second-nearest-neighbor interactions in the bcc structure, through examination of fitted interatomic force constants. The large changes in phonon energy result in a significant increase in the vibrational entropy, called the nonharmonic vibrational entropy, which emulates the temperature behavior of the magnetic entropy across the Curie temperature. The nonharmonic vibrational entropy is attributed to interactions between the vibrations and state of magnetic disorder in the material, which persists above the magnetic transitions and extends the stability region of the bcc phase.
Orthorombic cementite, Fe
C, exhibits anisotropic magneto-volume behavior in the ferromagnetic phase including regions very low thermal expansion. The phonon modes of cementite show anomalous temperature dependence, with low energy phonon modes increasing their energy at elevated temperatures in the ferromagnetic phase. This behavior is reversed after the magnetic transition and these same phonon modes lower their energies with temperature, consistent with observed thermal expansion. This atypical phonon behavior lowers the vibrational entropy of cementite up to the Curie temperature. The experimentally observed increase in low energy acoustic phonons affects the elastic behavior of Fe
C, increasing the isotropy of elastic response. First principles calculations link the observed phonon energy increases to specific vibrational modes that are polarized along the b-axis, which aligns with the closest Fe-Fe bonding direction. The nonharmonic behavior of the vibrational modes are discussed in the context of other observations of anomalous anisotropic magneto-volume behavior in Fe
C.
Item Type:
Thesis (Dissertation (Ph.D.))
Subject Keywords:
Iron, Cementite, Fe3C, Phonons, Thermodynamics, Magnetism
Degree Grantor:
California Institute of Technology
Division:
Engineering and Applied Science
Major Option:
Applied Physics
Awards:
The Lucy Guernsey Service Award, 2012; Graduate Deans’ Award For Outstanding Community Service, 2015
Thesis Availability:
Public (worldwide access)
Research Advisor(s):
Fultz, Brent T.
Thesis Committee:
Fultz, Brent T. (chair)
Jackson, Jennifer M.
Johnson, William Lewis
Kochmann, Dennis M.
Minnich, Austin J.
Schwab, Keith C.
Defense Date:
28 May 2015
Non-Caltech Author Email:
Lisa.Mauger (AT) gmail.com
Funders:
Funding Agency
Grant Number
Carnegie DOE Alliance Center (CDAC)
UNSPECIFIED
Record Number:
CaltechTHESIS:06092015-103051110
Persistent URL:
DOI:
10.7907/Z9TQ5ZH3
Related URLs:
URL
URL Type
Description
DOI
Publication related to chapter 6
Default Usage Policy:
No commercial reproduction, distribution, display or performance rights in this work are provided.
ID Code:
9013
Collection:
CaltechTHESIS
Deposited By:
Lisa Mauger
Deposited On:
09 Nov 2016 19:33
Last Modified:
04 Oct 2019 00:09
Thesis Files
Preview
PDF
- Final Version
See Usage Policy.
49MB
Repository Staff Only:
item control page
CaltechTHESIS is powered by
EPrints 3.3
which is developed by the
School of Electronics and Computer Science
at the University of Southampton.
More information and software credits
Phonon thermodynamics of iron and cementite

Thesis by

Lisa Mary Mauger

In Partial Fulfillment of the Requirements
for the Degree of
Doctor of Philosophy

California Institute of Technology
Pasadena, California

2015
(Defended May 28, 2015)

c 2015
Lisa Mary Mauger

iii

To my mother, whose love and dedication I will always cherish.

Acknowledgments
This thesis was possible in no small part because of my adviser, Prof. Brent Fultz. It is with
deep gratitude that I acknowledge all the help he provided me along the way. Brent was patient
enough to let me explore various corners of our field, and watchful enough to reel me in when I
drifted too far. I have learned to appreciate the deep wisdom in his comments and suggestions,
even when they are too profound for my initial digestion. His direction, support and encouragement
were essential to the work presented here, and I am extremely grateful for his guidance in shaping
my scientific development. Perhaps most importantly, Brent has put together a fantastic team of
budding scientists that I have been lucky enough to work with over the years.
Science is necessarily collaborative, but in our group beamtime is very collaborative and I can’t
say enough words of thanks to everyone who worked through the night to make so many parts
of my PhD possible. These people have been my colleagues and collaborators over the years, but
also my closest friends. Sally June Tracy, Jorge Muñoz and Matthew Lucas spent more evenings,
weekends, and holidays alongside me in synchrotron experimental halls than any of us would care
to enumerate. Sally June Tracy was an invaluable collaborator, a steady source of scientific insight,
experimental dedication and unique humor. Jorge Muñoz provided enlightened conversation about
science and many other facets of life. Matthew Steven Lucas is an excellent experimentalist, and
I am so thankful that I had the opportunity to learn so much from him over the years. I want
to thank so many other current and former Fultz group members, for their help and thoughtful
insights along the way: Olivier Delaire, Jiao Lin, Chen Li and Xiaoli Tang, for mentoring me as I
was getting started, and for their thoughtful advice on avoiding common graduate student pitfalls
– more of which I wish I had heeded. To current members Tian Lan, Hillary Smith, Dennis Kim,
and Max Murialdo, who have helped me perform the numerous essential lab activities, I appreciate
your patience, consideration, and insight. Jane Herriman and Olle Hellman have bravely supported
my recent adventures in computational material science. Their guidance has been enlightening, and
I only wish I had more time to spend learning from them both. And to the future of the Fultz
group, Nick Weadock, Fred Yang, Xiaoli Xu, and Heng Yang, it has been a pleasure getting to know
you, I look forward to seeing great things in the years to come. I want to acknowledge thoughtful
contributions from the many staff that have been part of the Fultz group over the years. Channing

Ahn, Michael McKerns, and John McCorquodale all dedicated their time to teaching me the ropes
when I arrived. Mike Vondrus played a vital role, not only in machining the components we needed
on my less than forgiving timelines, but also by providing design insight that greatly improved each
project. I would also like to acknowledge the many administrative staff who are so vital to the
interworkings of Caltech. Pamela Alberston gave me a thorough education on department history
and protocals. Joan Sullivan helped our Keck laboratories weather the debris, noise, flooding and
repairs that come with a sixty year old building. Christy Jenstad, our fantastic APh/MS option
representative, who joyfully helped me navigate degree requirements while becoming a dear friend.
I would like thank my thesis committee, for overseeing the final stages of my development here
at Caltech – Prof. Jennifer Jackson, Prof. Bill Johnson, Prof. Dennis Kochmann, Prof. Austin
Minnich and Prof. Keith Schwab, I appreciate the time and consideration you put into reading
this document, and I look forward to thoughtful discussions yet to come. Prof. Jackson and Prof.
Johnson have played especially important roles in the execution of this work, through valuable
ongoing collaborations between our groups.
My graduate research provided me with opportunities to interact with wonderful scientists around
the country and I greatly appreciate the time each of them took to train and mentor me. The staff
of the Carnegie Institute of Washington were extremely important to my development as a scientist.
They supported my research through the Carnegie - Department of Energy Alliance Center (CDAC),
which provided a diverse scientific community with a wealth of experience and insight. I would like to
specifically thank Dr. Steve Gramsch, Dr. Somayazulu and Dr. Russ Hemley for their support and
collaboration. I would also like to thank the staff at the High Pressure Collaborative Access Team
at Argonne National Lab. Their innovation and dedication has produced a unique and excellent
facility. Dr. Yuming Xiao and Dr. Paul Chow provided countless hours of user support. Curtis
Kenney-Benson and Eric Rod also played an essential role in making sure all the pieces fit together.
My experimental achievements were built on the foundation provided by these dedicated scientists
and engineers, and I remember our time together fondly. I would also like to acknowledge the staff
of Sector 3 at the Advanced Photon Source for their support and contributions. Drs. Ercan Alp,
Tom Tollner, Michael Hu, Jiyong Zhao and Wolfgang Sturhahn. The insight, support, dedication
and humor lead to many productive beamtimes and other fun memories.
My friends have played an enormous role in keeping me happy and productive, greatly enriching
time spent in graduate school. I was lucky enough to find and befriend the most wonderful collection
of people here at Caltech, whose friendship I hope to maintain for years to come. My wonderful
roommates, Jessica Pfeilsticker, Katie Casgrain, and Gabriela Venturini – fate threw us together
but sharing a life with the three of you has brought an uncountable number of cherished memories.
Gabriela was a vital source of love and support through so many challenges and adventures over
the years. My wonderful friends Samantha Johnson and Sunita Darbe have been there to push

vi
me forward - at school, in the gym, on the dance floor - greatly enriching my graduate experience.
Finally I would like to acknowledge the unwaivering support of my boyfriend Chris, his appearance
in my life has been a gift that I treasure.
My family has been a wonderful source of encouragement during my education, and I appreciate
them more and more every day. I want to thank my mother Mary, for her strength and encouragement; I can only hope to emulate her selfless dedication that provided my siblings and I with so
many wonderful opportunities and experiences. I also want to thank her partner Dan, for taking
care of her when I can’t be there, and joining her on so many road trip adventures that ended here
in southern California. I want to acknowledge the support of my sister Jenny and my brother David.
They have grown into exciting adults with rich lives since I left for graduate school, and they provide
me counsel on a number of topics. Over the course of my PhD I been delighted to get much closer
with my grandmother Shirley. She has been a source of support and encouragement, bringing joy
and wisdom to my life journey.

vii

Abstract
The vibrational properties of materials are essential to understanding material stability and thermodynamics. In this thesis I outline vibrational thermodynamic models and the experimental tools that
provide evidence on phonon behavior. The introductory section discusses the history of metallurgy
and thermodynamic theory, with an emphasis on the role of iron and cementite, two important components of steels. The thermodynamic framework for understanding vibrational material behavior is
provided alongside the growing body of experimental and computational tools that provide physical
insight on vibrational properties. The high temperature vibrational behavior of iron and cementite
are explored within this context in the final chapters.
Body-centered-cubic iron exhibits decreasing phonon energies at elevated temperatures. The
observed energy change is not uniform across phonon modes in iron, and specific phonon modes show
significant decreases in energy that are not explained by simple vibrational models. This anomalous
energy decrease is linked to the second-nearest-neighbor interactions in the bcc structure, through
examination of fitted interatomic force constants. The large changes in phonon energy result in a
significant increase in the vibrational entropy, called the nonharmonic vibrational entropy, which
emulates the temperature behavior of the magnetic entropy across the Curie temperature. The
nonharmonic vibrational entropy is attributed to interactions between the vibrations and state of
magnetic disorder in the material, which persists above the magnetic transition and extends the
stability region of the bcc phase.
Orthorombic cementite, Fe3 C, exhibits anisotropic magneto-volume behavior in the ferromagnetic phase including regions of very low thermal expansion. The phonon modes of cementite
show anomalous temperature dependence, with low energy phonon modes increasing their energy at
elevated temperatures in the ferromagnetic phase. This behavior is reversed after the magnetic transition temperature and these same phonon modes lower their energies with temperature, consistent
with observed thermal expansion. This atypical phonon behavior lowers the vibrational entropy of
cementite up to the Curie temperature. The experimentally observed increase in low energy acoustic
phonons affects the elastic behavior of Fe3 C, increasing the isotropy of elastic response. First principles calculations link the observed phonon energy increases to specific vibrational modes that are
polarized along the b-axis, which aligns with the closest Fe-Fe bonding direction. The nonharmonic

viii
behavior of the vibrational modes are discussed in the context of other observations of anomalous
anisotropic magneto-volume behavior in Fe3 C.

Contents
Acknowledgments

Abstract

vii

1 Introduction

1.1

Iron and Steel - Ancient History . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1.2

The Dawn of Thermodynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2 Crystals and Phonons

2.1

Crystal Lattices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2.2

Phonons . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2.3

Observations of Phonons . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3 Thermodynamics

3.1

Thermodynamic Relations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.2

Debye Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.3

Quasi-Harmonic Models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.4

Anharmonic Effects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4 Experimental Methods
4.1

4.2

The Mössbauer Effect . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.1.1

Nuclear Forward Scattering . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Nuclear Resonant Inelastic X-ray Scattering . . . . . . . . . . . . . . . . . . . . . . .

5 Computational Methods

5.1

Computational Quantum Mechanics . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.2

Born von-Kármán Fitting of Phonon Spectra . . . . . . . . . . . . . . . . . . . . . .

5.2.1

Fitting Phonon DOS - Practical Matters . . . . . . . . . . . . . . . . . . . . .

6 Anharmoncity in BCC Fe At Elevated Temperatures
6.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

39
39

6.2

Experimental . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.3

Force Constant Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.4

Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.4.1

Phonons . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.4.2

Quasiharmonic Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.4.3

Vibrational Entropy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.4.4

Born-von Kármán Fits . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.5.1

Phonons and Born-von Kármán Model Dispersions . . . . . . . . . . . . . . .

6.5.2

Vibrational Entropy and Free Energy . . . . . . . . . . . . . . . . . . . . . .

6.6

Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.7

Supplemental Material . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.5

7 Cementite

7.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.2

Experimental . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.3

Computational . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.4

Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.4.1

Structure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.4.2

Phonons . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.4.3

Vibrational Entropy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.4.4

Low Energy Phonon Modes . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.4.5

Elastic Constants . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.4.6

Electronic DOS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.5

8 Summary

List of Figures
1.1

Fe-C Phase Diagram . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1.2

Free energy of iron polymorphs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2.1

Phonon Dispersions and Density of States for α-Fe . . . . . . . . . . . . . . . . . . . .

3.1

Debye model for Fe thermodynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.2

Thermal expansion of Fe and Fe3 C . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.3

Quasi-Harmonic Debye model for Fe thermodynamics . . . . . . . . . . . . . . . . . .

3.4

Temperature-dependent γT (T ) for α-Fe . . . . . . . . . . . . . . . . . . . . . . . . . .

3.5

Early thermodynamic assessment of Fe . . . . . . . . . . . . . . . . . . . . . . . . . .

4.1

Mossbauer isomer shift and electric field gradient transition energies . . . . . . . . . .

4.2

Mossbauer isomer shift and hyperfine magnetic field transition energies . . . . . . . .

4.3

NFS Timing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.4

NRIXS Excitations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.5

NRIXS Spectra S(E) of α Fe . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4.6

Density of States (DOS) of α-Fe . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.1

Phonon Effects from 1NN Longitudinal . . . . . . . . . . . . . . . . . . . . . . . . . .

5.2

Phonon Effects from 2NN Longitudinal . . . . . . . . . . . . . . . . . . . . . . . . . .

5.3

Qualities of Fit for BvK DOS Test . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.4

5NN Fits . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.1

Fe phonon DOS Spectra . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.2

Fe phonon modes vs temperature . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.3

Fe Entropy Contributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.4

DOS Fits and Dispersion Comparisons . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.5

BvK Fitted Dispersions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.6

Thermal dispersion behavior divided by branch . . . . . . . . . . . . . . . . . . . . . .

6.7

BvK Fitted Force Constants . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

xii
6.8

Thermal Behavior of High Symmetry Phonon Modes . . . . . . . . . . . . . . . . . . .

6.9

Nonharmonic Vibrational Entropy Comparison . . . . . . . . . . . . . . . . . . . . . .

6.10

Nonharmonic Heat Capacity Contribution . . . . . . . . . . . . . . . . . . . . . . . . .

6.11

Fe FCC Phonon DOS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.12

Comparison of anharmonic phonon behavior with magnetic entropy . . . . . . . . . .

6.13

Thermal behavior of Lamb-Mössbauer Factor . . . . . . . . . . . . . . . . . . . . . . .

7.1

The structure of Fe3 C . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.2

NRIXS phonon partial DOS of Fe3 C . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.3

Comparison of experimental and calculated phonon pDOS . . . . . . . . . . . . . . . .

7.4

Mean phonon energy of the Fe pDOS . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.5

Vibrational entropy of the Fe pDOS . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.6

Sound velocities of Fe3 C with temperature . . . . . . . . . . . . . . . . . . . . . . . .

7.7

Comparison of NRIXS pDOS in the ferromagnetic and paramagentic phases

. . . . .

7.8

Calculated phonon Fe pDOS at high temperatures . . . . . . . . . . . . . . . . . . . .

7.9

Illustrations of phonons that stiffen with temperature at X and Z . . . . . . . . . . . .

7.10

Illustrations of Y-point phonon modes . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.11

Illustrations of Γ point low energy phonon modes . . . . . . . . . . . . . . . . . . . . .

7.12

Calculated Fe pDOS resolved onto distinct lattice sites . . . . . . . . . . . . . . . . . .

7.13

Thermal trends of elastic modulie extracted from DFT

. . . . . . . . . . . . . . . . .

7.14

Electronic density of states of Fe3 C . . . . . . . . . . . . . . . . . . . . . . . . . . . .

xiii

List of Tables
5.1

Interatomic force constants of 300K α-Fe . . . . . . . . . . . . . . . . . . . . . . . . .

6.1

The vibrational entropy of α-Fe

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6.2

BvK Force Constants for α-Fe including 2NN . . . . . . . . . . . . . . . . . . . . . . .

6.3

BvK Force Constants for α-Fe including 4NN . . . . . . . . . . . . . . . . . . . . . . .

7.1

Fe3 C experimental and computed lattice parameters . . . . . . . . . . . . . . . . . . .

7.2

Fe3 C experimental and computer lattice sites . . . . . . . . . . . . . . . . . . . . . . .

7.3

Total Vibrational Entropy of Fe3 C . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

7.4

Calculated elastic moduli of Fe3 C . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Chapter 1

Introduction
1.1

Iron and Steel - Ancient History

Great advances in human civilizations have come from the development of metalworking technology. However, early man lacked a sophisticated understanding of thermodynamics, so metalworking
advances were based on carefully refined recipes rather than physical principles. Despite their lack
of physical understanding, early man learned to isolate metals with heat and alloy substances to
improve mechanical strength. These methodologies were largely dependent on the quality of raw
materials and this knowledge was occasionally lost to time, only to be rediscovered generations later
by different routes.
The discovery of iron and steels transformed civilizations, but this technology was hard won by
early blacksmiths. While a steel composition can be as simple as iron with a few atomic percent
carbon, the development of steel took centuries. We now understand that the mechanical benefits of
steels that early man sought are enabled by the temperature-driven polymorphism in iron. The iron
carbon phase diagram shown in Fig. 1.1 contains three thermodynamically stable regions of pure iron
that form stable solid solutions with carbon. With increasing temperature pure iron transforms from
a ferromagnetic bcc structure to a paramagnetic bcc structure (1043K) to a fcc structure (1185K)
and then back to a bcc structure (1667K) before melting. This re-entrant temperature dependent
polymorphism is not common among the elemental metals, and gives iron alloys a number of unique
properties. The variability in mechanical properties of steels can be controlled to some extent by
composition, but is more dramatically altered by metastable microstructure induced by controlled
temperature cycling through the polymorphic transitions of iron.
Humans discovered how to smelt elemental iron around 2000 BCE, near present-day Turkey [2].
However early iron smelting furnaces could only reach temperatures around 1400K, notably below
the 1811K melting point of iron. Instead the carbon monoxide created by burning charcoal would
reduce iron oxide based ores to a mixture of metallic iron, charcoal, and silicate inclusions. The
impurity phases were removed by high temperature mechanical working, which occupied much time

Liquid
BCC*

1811K
1667K

FCC

1185K

paramagnetic
TC

1043K

ferromagnetic

BCC

⇐ increasing C content (at. %)

Pure Fe

Figure 1.1: The Fe-C phase diagram, maps all stable phases with blue regions. This figure was
constructed using the 1992 phase diagram developed by Okamoto, et. al. [1]. The phase regions
labeled rt refer to the body-centered-cubic (bcc) solid solution phase, while the phase regions labeled
ht refer to the face-centered-cubic (fcc) solid solution phase.
for early blacksmiths. Pure iron is also a relatively soft and ductile metal, which limited the early
use of iron to small domestic applications. The discovery of steel occurred nearly 1000 years later,
providing much greater strength and revolutionizing weaponry [2]. But producing quality steels
was quite difficult for early man, who could not appreciate the thermodynamic forces driving his
processes.
The methodologies that lead early blacksmiths to high strength steels were tedious, and provided
many opportunities for failure with small changes in temperature or chemical composition. When
smiths reheated their pure worked iron in a charcoal fire, they were unknowingly diffusing small
quantities of carbon into their material. Without realizing it, early smiths were heating iron through
the polymorphic transition into the face-centered-cubic structure. This phase has many interstitial
locations for carbon atoms to occupy, and will support carbon concentrations up to 8 atomic %.
Iron is strengthened by adding only 1-2 atomic % carbon, but if larger quantities are diffused at
high temperatures, the mixture will become a hard, brittle, unworkable material which was referred
to as pig iron. Longer heating times increased carbon concentration in steels at a rate which varied

with temperature. But when the temperature of the iron carbon solid solution is lowered back to
the body-centered-cubic phase stability range, the carbons have fewer favorable interstitial sites and
attempt to diffuse out to form Fe3 C shown on the left side of the phase diagram in Fig 1.1. However
if this diffusion limited process is interrupted, by rapidly lowering the temperature of the alloy,
the carbon atoms get stuck in unpreferable interstitial sites and they distort they crystal structure
locally, resulting in a significantly harder material called martensite.
This was the true birth of steels. The rapid quenching of hot Fe-C alloys in water can provide
nearly a five fold increase in strength, making a product much harder than bronze, a Cu-Sn alloy,
which was the technological standard of the time [2]. Once it was realized that the strength of steels
could be adjusted by thermocycling alone, the technological potential was quickly realized. Today
martensitic steels still provide the best strength per cost per unit volume of modern engineered
materials [3].

1.2

The Dawn of Thermodynamics

The physical phenomena early man explored in manipulating metal alloys is now encompassed in
the fields of modern thermodynamics and metallurgy. They refined ore into metals using heat in
early furnaces. These metals were alloyed, heat treated, and worked by hand, adding a variety of
properties to the components that are now understood to come from the atomic and microstructural
arrangements these processes induced.
The efficient and effective processing of metals was a great driver of thermodynamic understanding. Thus modern thermodynamics grew up beside the industrial revolution, when the old
methodologies of metal working were traded for modern processes of great scale. The pioneering
work of Josiah Willard Gibbs laid the groundwork for our modern understanding of phase diagrams
and material equilibria [4]. Gibbs was a mathematician by training and focused on the geometries
of early phase maps to draw connections between phase stability, energy and entropy. The Gibbs
free energy, G,
G = H − TS

(1.1)

is related to the enthalpy, H, temperature, T , and entropy, S, of a substance in a given state.
The early work of Gibbs emphasized the generality of thermodynamic principles to include material
systems of all kinds [4]. His pioneering insights slowly brought scientific unity to the practices of
chemical metallurgy, by uncovering the underlying principles in the centuries of collected physical
observations and material preparation methodologies.
Around the same time, great scientific minds were beginning to quantify the nature of heat in
solids. The Einstein model of solids found a quantum mechanical description of lattice vibrations as
quantum harmonic oscillators [5]. Peter Debye improved on this model to encompass the observed

low temperature heat capacity behavior [6]. This laid the groundwork for physical interpretations
of observed thermodynamic properties in solids.
It was quickly understood that the unique polymorphism of iron provided for the great diversity
of technological properties of iron-based alloys. Improving iron-based steels was of intense industrial
interest, and thus improved thermodynamic understandings were readily applied to iron alloys.
Improvements in calorimetry around the same time produced experimental measures of free energy
derivatives. Attempts to quantify the free energies driving the diverse microstructures of steels
quickly followed. An early assessment of the free energies of iron was provided by Austin as shown
in Fig. [7].

Figure 1.2: The free energy of bcc (α-Fe) and fcc (γ-Fe) extracted from calorimetry measurements
[7].
This was the beginning of physical metallurgy, where detailed methodologies could finally be
understood in terms of physical consequence. The chemical consequences of iron’s polymorphic
phase transitions were mapped on phase diagrams where thermodynamic insight was used to explore new compositions. The equilibrium positions of carbon in iron changes substantially with
crystal structure, and cementite or Fe3 C will precipitate for many compositions at modest temperatures. The mechanical effects of stable carbides cementite, and other less stable configurations
were heavily studied to understand the mechanical responses of steels with varied thermal histories. Phase transformation kinetics were understood in terms of diffusion limited formation of stable
carbides and metastable structures. And understanding the microstructure-related mechanical properties that evolved under various thermocycling conditions became the focus of decades of modern

steel-metallurgy.
Structural steels weren’t the only technological materials to highlight the interesting physical
behavior of Fe during the industrial revolution. The unique physical properties of Fe also participate
in materials that are famous for something other than their strength. The 1920 Nobel prize in
physics was awarded to Charles Guillaume “in recognition of the service he has rendered to precision
measurements in physics by his discovery of anomalies in nickel steel alloys” [8]. The anomaly that
so captivated the Nobel committee was a thermo-volume behavior of FeNi alloys resulting in a near
vanishing thermal expansion for a wide temperature range at a specific composition. The Fe64 Ni36
Invar composition takes a disordered face centered cubic structure, and exhibits a thermal expansion
coefficient which is more than 10 times lower than the thermal expansion of Fe, Ni, and other
distant FeNi compositions. Invar materials were quickly deployed to improve the accuracy of clocks
and many other sensitive measurement techniques, even though the physical origin of Invar eluded
material physicists nearly another century. Since Guillaume’s initial discovery a number of other
iron transition metal alloys have been found to exhibit Invar behavior, including FePt and FePd.
The anomalously low thermal expansion exhibited by these alloys is caused by temperature-induced
magnetic transitions which exhibit magneto-volume behavior.

Chapter 2

Crystals and Phonons
2.1

Crystal Lattices

Both the models of Einstein and Debye relied on the notion of crystalline solids. Which was verified
by the early work of W. H. and W. L. Bragg [9]. The x-ray pattern observed with Bragg diffraction
is the result of regularly repeating arrays of atoms in crystalline materials. The arrangement of the
atoms in a crystal can be reduced to a set of translational symmetry operations that relate every
atomic position in a perfect crystal onto all equivalent positions, thus defining the lattice symmetry
of crystalline solids, and simplifying their structure to a primitive unit cell that may be tessellated
to map out every lattice position. One may then define lattice translation vectors, r, as
r = x1 a1 + x2 a2 + x3 a3 ,

(2.1)

where xj are integers and a1 , a2 , and a3 are the primitive lattice translation vectors, where a1 ·a2 ×a3
defines the primitive unit cell volume. A reciprocal lattice can then be defined for each type of lattice,
based on the constraints of Bragg’s law, which observed the relationship between the wavevector
of the incoming radiation and the structure of the crystal lattice. The reciprocal lattice has a
complementary set of vectors, q, defined as
q = y1 b1 + y2 b2 + y3 b3 ,

(2.2)

where the prefactors yi are again integers and b1 , b2 , b3 are defined as the primitive vectors of the
reciprocal lattice. These special reciprocal space vectors (also referred to here as q-space vectors)
can be constructed from the real space lattice vectors as
b1 = 2π

a2 × a3
a3 × a1
a1 × a2
, b2 = 2π
, b3 = 2π
a1 · a2 × a3
a1 · a2 × a3
a1 · a2 × a3

The primitive cell of reciprocal space is commonly referred to as the first Brillouin zone.

(2.3)

The principles of lattice symmetry discussed here can easily be extended to include compounds
with multiple atomic species. In this case the primitive unit cell is still the smallest volume that
can be tessilated to map all space, defined by three primitive lattice vectors. The only distinction
is that the lattice has a basis, which is typically described by vectors that connect the positions of
unique atoms in the primitive cell. It is within this symmetry that we can now begin to define the
lattice modes, or quanta of vibrations called phonons [10].

2.2

Phonons

The equilibrium interatomic distances offer an energy efficient packing that optimizes the interactions
of the atomic electrons. The positions geometrically optimize the interatomic forces to set the
equilibrium distances. We will describe this potential energy of the interatomic interaction as φ(R),
where R is the distance between a pair of atoms. The potential energy of a crystal, U , can then be
described by summing over all the pairs of atoms in a crystal,
U=

1X
φ(ri − rj ).
2 i,j

(2.4)

If an atom indexed by i is perturbed a small distance from its equilibrium position, ri , to a new
position, ri + ui , the neighboring unperturbed atoms would exert a force on the displaced atom
to return it to its equilibrium position. A Taylor expansion of the potential energy U for small
displacements of the atoms, ui , from their equilibrium positions gives
U=

1X
1X
1X
φ(ri − rj ) +
(ui − uj ) · ∇φ(ri − rj ) +
[(ui − uj ) · ∇]2 φ(ri − rj ) + O(u3i ). (2.5)
2 i,j
2 i,j
4 i,j

Here the first term is the constant equilibrium crystal potential, the linear second term is the restoring
forces that sum to zero over the crystal, and the remaining term is the harmonic term of the original
potential φ(R), which can be written as
Uh =

(ri − rj )µ

i,j&µ,ν=x,y,z

∂ 2 φ(ri − rj )
(ri − rj )ν .
∂(ri − rj )µ ∂(ri − rj )ν

(2.6)

This simplification is called the harmonic approximation because it neglects higher order terms in
the potential. The gradient term at the heart of this expression gives the force constant between
two atoms in a specific Cartesian direction. This expression can be generalized by defining a force
constant matrix, K, that encompasses all the atomic interactions represented in the derivatives of
φ:
Uh =

1X
ui K(ri − rj )uj .
2 i,j

(2.7)

With the proper symmetry constraints and Born-von Kármán periodic boundary conditions,
we cleverly select a general description of the displacements that is appropriate for translational
symmetry:
ui (t) = ei(k·ri −ωt) ,

(2.8)

where is the polarization of the atomic displacement in Cartesian coordinates. We can solve for
the equations of motion given our harmonic potential, where
M üi = −

∂U h
=−
K(ri − rj )uj .
∂ui

(2.9)

Our system of N atoms has 3N discrete vibrational modes that can be supported by the crystal,
and our equation of motion becomes
M ω2 =

K(rj )e−ik·rj = D(k).

(2.10)

This expression is called the dynamical matrix expression, where D(k) is the dynamical matrix.
The dynamical matrix contains the force constants (or potential derivatives) for every pair of atomic
interactions in the crystal. For an atom in a real crystal we know that the largest contribution
to the restoring forces will come from the atoms immediately around it, so we can truncate the
dynamical matrix to include only the most pertinent nearest-neighbor restoring forces, often without
a significant loss of accuracy. The dynamical matrix expression can be solved for the normal mode
frequencies ω and mode wavevectors at every k in reciprocal space.
This methodology connects the symmetry of the lattice with the allowed normal modes and
the interatomic forces driving them; however, quantum mechanical considerations are required to
extend this description of lattice modes to quantized vibrational excitations called “phonons”. These
considerations ensure that vibrations are properly counted a low temperatures, where their discrete
nature becomes apparent. The energy of a crystal is described by 3N quantum harmonic oscillators
with frequencies from the dynamical matrix expression, but governed by Bose-Einstein occupation
statistics. Thus the Hamiltonian for a crystal transitions from its classical definition,
H=

1 X 2 1X
P +
ui D(ri − rj )uj
2M i i
2 i,j

(2.11)

to its quantum representation,
H=

k,s

h̄ωs (k)(αks
αks + ) =
(nks + )h̄ωs (k),
k,s

(2.12)

where the phonon creation operator αks
and phonon annihilation operator αks are defined as

αks
=√

e−ik·ri s (k)

and
1 X −ik·ri
αks = √
s (k)
N i

M ωs (k)
ui − i
pi
2h̄
2h̄M ωs (k)

(2.13)

M ωs (k)
ui + i
pi
2h̄
2h̄M ωs (k)

(2.14)

and
h̄ωs (k)

nks = (e kB T − 1)−1 .

(2.15)

The important distinction between these two models is readily observed in experimental heat capacities. The heat capacity, C, is the temperature derivative of the internal energy of a material, U . In
the classical harmonic oscillator formalism, U (T ) is a linear function of T , so the heat capacity
C=

∂U
(U eq + 3N kB T ) = 3N kB
∂T
∂T

(2.16)

is constant for all temperatures. However, the heat capacity of a set of quantum harmonic oscillators
contains a temperature-specific term
C=

X ∂
X ∂ h̄ωs (k)
∂U
(U eq +
(nks + )h̄ωs (k)) =
(nks )h̄ωs (k) =
∂T
∂T
∂T
∂T e h̄ωkBs (k)
−1
ks

(2.17)

which recovers the experimentally-observed temperature dependence of the heat capacity at very
low temperatures. Using quantum harmonic oscillator formalism in the calculation of lattice thermodynamic variables will include the zero-point vibrational energy of the solid [11].

2.3

Observations of Phonons

Direct measurements of phonons can provide valuable insight into the physical basis for thermal
behavior. The first maps of phonons in solids were completed using reciprocal-space-resolved methods like triple-axis inelastic neutron scattering from single crystals. The geometrical nature of these
measurements permit fine control of the incoming energy and momentum of the neutron, and similar
control of the energy and momentum of the outgoing neutron may be detected. Shull and Brockhouse pioneered early work on neutron scattering techniques that were capable of resolving phonons
in solids. Measurements of the phonons of bcc α-Fe followed immediately after the development
of these instruments, with two separate papers reporting the phonon dispersions of iron published
in 1967 [12, 13]. The early phonon dispersion measurements of iron are shown in Fig.2.1, where
triple-axis neutron measurements generated the set of points that are resolved in q-space and en-

10
ergy. These points are then analyzed with the harmonic model developed in the previous section.

Figure 2.1: Phonon dispersions from triple-axis inelastic neutron scatter (points), overlaid with
Born-von Kármán model fits and the resulting phonon density of states shown on the right [12].
The phonon dispersion points are fit to force constants in a dynamical matrix, typically using a least
squares optimization that seeks to minimize the system of equations against all the observed phonon
measurements. This requires truncating the dynamical matrix to a subset of all the interactions in
the material, typically limiting restoring forces to the closest neighboring atoms. Minkewicz utilized
the atomic interactions for the first through fifth nearest neighbors, obtaining the fits shown in
Fig. 2.1. The force constants can be used to describe phonon behavior in other portions of q-space
where measurements were not collected. They can also be integrated over all of q-space to provide
a phonon density of states (DOS), which is also shown in Fig. 2.1. The phonon density of states can
be readily used to describe the phonon thermodynamics of materials.
Measurements of phonon dispersions using triple-axis neutron spectrometers are routinely conducted today at research reactors. Additionally, methods have been developed to measure the
complete phonon DOS of polycrystalline materials using neutron time-of-flight techniques, and specialized inelastic x-ray scattering methods that are discussed in the next chapter. Phonon DOS
measurements are typically collected much more quickly than dispersion curves, and can be directly
used in thermodynamic expressions without fitting. These methods offer the advantage of measuring
every phonon state available to a material, but they do so without the q-space resolution provided by
phonon dispersion measurements. Therefore, mapping the interatomic forces to the phonon density
of states is less direct because the reciprocal space information has been lost. Sharp features in the
phonon density of states called van Hove singularities may be matched to some points on phonon
dispersion, but direct mapping is complicated. Because the phonon DOS is integrated over all of

11
q-space, there is no way to distinguish between degenerate features. It was therefore believed that
in practice, the interatomic force constants cannot be directly extracted from phonon DOS, since
the relationship is not invertible. However fits may still be accomplished by iteratively exploring
various force constant configurations and comparing the density of phonon states they generate with
experimental observations. This has been demonstrated on relatively simple lattices in elemental
solids and binary alloys when sufficient information is available, [14–18] and will be discussed in
Section 5.2.

Chapter 3

Thermodynamics
3.1

Thermodynamic Relations

Thermodynamics was useful to 20th century metallurgists insofar as it can be used to predict and
extrapolate properties beyond those directly measured. A functional thermodynamic understanding
of iron was readily applied to understanding the technologically-relevant properties of iron alloys.
The Gibbs free energy of a solid, G, can be divided into enthalpy and entropy terms. Under
constant volume conditions, the enthalpy of a solid, H, is largely determined by the internal energy,
U , which can be characterized as the energy involved in assembling a set of atoms into their solid
configuration. The entropy of a solid, S, enumerates the way heat is stored in a material. Both enthalpy and entropy can be extracted by integrating the measured heat capacity at constant pressure,
CP using the expressions

Z T
H(T ) =

CP (T 0 )dT 0

(3.1)

CP (T 0 ) 0
dT .
T0

(3.2)

and
Z T
S(T ) =

In solids at finite temperatures (above ambient conditions) the free energy contribution from entropy, T S, changes more rapidly than the enthalpy, H, dominating the thermal effects on the free
energy. For ordered crystalline solids, lattice vibrations make the largest entropic contribution, Svib .
Electronic excitations also create entropy, Selec , though noticeably smaller than vibrational entropy.
Solids that exhibit magnetic ordering will also have magnetic excitations that perturb spins from
their ground state orientation, providing magnetic entropy, Smag . These three contributions, vibrational, electronic, and magnetic excitations, enumerate the ways that heat can be stored in iron and
comprise its entropy. These contributions are hoped to be adiabatically separable, providing the
expression
S = Svib + Selec + Smag .

(3.3)

13
Understanding the entropy of a solid will greatly inform the temperature-dependent behavior
of its free energy, which is essential to developing the theoretical basis for high temperature phase
diagrams. Early thermodynamic research sought to reconcile these physical excitations with the
aggregate thermodynamics that drives phase transitions. Experimental heat capacities provided a
basis for comparing observed thermodynamic properties with theoretical models. But direct measurement of the phonon density of states of a material can also provide complementary information
that may be used to assess the phonon-specific contributions to thermodynamics.
The phonon contribution to the entropy, Svib , may be calculated directly from the phonon DOS,
g(E), at the temperature that the phonon DOS was acquired
Svib (T ) = 3kB

gT (E){(n + 1) ln(n + 1) − n ln(n)}dE.

(3.4)

Where the integral goes over all phonon energies, and the Planck function n is a function of energy

and temperature only, simplified from Eqn 2.15 to n = (e kB T − 1)−1 . Additionally the phonon
contribution to the heat capacity, Cpvib , may be calculated from the phonon density of states,
Cpvib =

∂n
gT (E)EdE.
∂T

(3.5)

These expressions provide another route to the thermodynamic behavior of materials that focuses
on the phonon contributions alone, by using the phonon density of states. Since the phonon contribution to thermodynamics is almost always the largest thermodynamic contribution at finite
temperatures, early thermodynamic models focused on quantifying the phonon behavior through
various formalisms.

3.2

Debye Model

The Debye model for the vibrational response of the solid makes use of the quantum mechanical
nature of phonons, but largely ignores the details of how phonons relate to the symmetry of the
structure. Debye simplified the normal mode relationships of a crystal considerably by assuming
that the phonon frequencies ω obeyed a linear relationship with respect to the reciprocal lattice
vector, k, ω = c|k|, where c is the sound velocity of the phonon. The assumption of linear phonon
branches only applies rigorously in the long wavelength limit (at very low |k|). The complicated
mathematical formulations of the previous section are simplified to three isotropic acoustic phonon
branches. By selecting an isotropic cutoff for the outer limit of reciprocal space, kD = 3 6π 2 ρ,
where ρ is the atomic density of the material, the number of normal modes available to the material
is properly set at 3N. The Debye model for vibrations has only one free parameter, the Debye
temperature ΘD = h̄ckD , which can be obtained with the slope of the linear acoustic branches, c.

14
The identical linear isotropic acoustic branches can be integrated over the spherical Brillioun zone
to provide a phonon density of states, which enumerates the available phonon modes in the crystal
by their energy level. The Debye heat capacity and phonon density of states are plotted in Fig. 3.1
for ΘD = 420, which is the Debye temperature commonly used for α-Fe.

Energy HmeVL

Cp HtotL

CP HhL

10
Svib HkB atomL

Cp HkB atomL

0.06
0.04
0.02
0.00

DOS H1meVL

0.08

QD = 420

0 200 400 600 800 10001200
Temperature HKL

STotal

Svib HhL

200 400 600 800 10001200
Temperature HKL

research papers

Figure 3.1: Debye model phonon density of state for ΘD = 420, a typical value used for α-Fe [19].
The corresponding vibrational heat capacity compared with α-Fe experimental
heat# capacity [20].
also strongly temperature dependent, at least for the
RT the SGTE database
The Debye vibrational entropy compared with the total entropy
of α-Fe
from
T
xp
T
dT ;
1
Tr
temperature interval covered in this experiment. Similar
Tr
[21]
behaviour is observed, for example, in -iron (Kohlhaas et al.,

Volume HÅ3 L

1967), associated with the ferromagnetic phase transition at
where VTr is the volume at a chosen reference temperature, Tr,
about 1045 K. This competition between the thermal expanand
(T) is the thermal expansion coef®cient, having the form
capacity
from
the Debye
sion arisingThe
from heat
vibrations
of thederived
atoms and
the increase
in model is capable of reproducing both the empirical
volumeDulong-Petit
with decreasing high
temperature
arising from
spontaneous
2
ao aobserved
1 T:
temperature
limit,
and also the low temperature T 3 T
behavior
in measured
magnetostriction provides the basis for the family of `Invar'
Ê , atemperature
= 154.8to
(1) low
This ®t gave
values of by
VTr fitting
alloys (such
Fe with 36 The
wt% Debye
Ni), which
have very low
heat as
capacities.
temperature
is commonly
determined
o = ÿ4 (2)
10ÿ5 Kÿ1 and a1 = 1.6 (3) 10ÿ7 Kÿ2, for a chosen Tr of 300 K.
thermal expansion around room temperature.
heatstudies
capacity
Once a Debye temperatureThe
is obtained,
theuncertainties
full vibrational
thermodynamics
a1 arise from the
large standard
in ao and
Current
of Fedata.
3C are mainly concerned with its
limited
temperature
range
and
the
small
number of data
possibleofoccurrence
in
the
Earth's
core;
it
is,
therefore,
its
a crystalline solid are mathematically accessible. This model is, however, a strictly harmonic
points available. The form of these parameters, with ao
physical properties in a high-pressure `non-magnetic' state
negativeof and
positive,
re¯ects
the strong
(with Pauli
paramagnetism
butoften
in which
are for
no local
1 strongly
approach,
which is
toothere
simple
the behavior
the aphonon
modes
in real
crystals.
The
temperature dependence of seen in the data, but probably
magnetic moments on the atoms) that are of most interest to
harmonicSuch
formalism
explain
phenomena
likethat
thermal
solids
thermal
also implies
it wouldexpansion
be unwise tofuse
themand
to extrapolate
Earth scientists.
a phase fails
only to
exists
abovenatural
60 GPa
to much higher temperatures. Expressed in this way, takes
(VocÏadlo, Brodholt et al., 2002) and is therefore not readily
resistivity of materials. These effect arise from other
interactions
that
ourexpression
ÿ5 are
Kÿ1 attruncated
460 K and in
6 (3)
10ÿ5 Kÿ1 at
values
of 3 (2) 10
accessible experimentally. The accessible high-temperature
600Debye
K, in good
agreement
results
of Reed
& Root
paramagnetic
(in which
there areThe
localinability
magneticof the
for thephase
interatomic
potential.
model
to dealwith
withthethe
physical
effects
of
(1998) from which a temperature-independent value of 5.2 (2)
moments but they are randomly disordered) will, however, be
ÿ5
ÿ1
thermal expansion
to severalphase
modifications
“quasi-harmonic”
can be
derived [de®ned asmodels.
= (1/V)(dV/dT)].
10areKcalled
more representative
of the led
core-forming
than the that
ferromagnetically ordered material (for further discussion see
below). To determine the thermal expansion coef®cient of this
phase in the 25.5
form tabulated by Fei (1995), the eight data points
shown in Fig. 2 for which T 460 K were ®tted to
25.0
Α-Fe
Γ-Fe
24.5
24.0
23.5

500

1000

Temperature HKL

1500

Figure 3.2: The experimental volume expansion of pure iron from low temperatures until melting
[22]. The observed thermal volume expansion of Fe3 C through the 460K magnetic transition [23].

3.3

Quasi-Harmonic Models

Quasi-harmonic models attempt to rectify the shortcomings of the harmonic models by taking into
account the effects of thermal expansion. In a strictly harmonic model, the vibrational degrees of
freedom have no dependence on the volume of a solid or its temperature. In 1926, Eduard Grüneisen
proposed a thermodynamic equation of state for matter that incorporates the vibrational effects from
changes in volume at finite temperatures. He defined the phonon mode Grüneisen parameter, γj ,
γj = −

∂lnωj
V ∆ωj
|T ≃ −
∂lnV
ωj ∆V

(3.6)

a unitless scaling parameter defined in terms of the phonon frequency, ωj , and the volume, V ,
of the solid [24]. Grüneisen used this description to develop a thermodynamic equation of state
which incorporates quantum mechanical lattice contributions. However experimental data on γj for
individual phonons is extremely rare.
More often an average bulk Grüneisen parameter γT is constructed to model the bulk material
behavior
γT = V

∂P
αKT
|V =
∂U
CV ρ

(3.7)

where V is the volume, P is pressure, α is the thermal volume expansion, KT is the isothermal
bulk modulus, CV is the heat capacity at constant volume, and ρ is the atomic density [24]. While
the microscopic Grüneisen parameter, γj , is an exact thermodynamic definition, models that use
the thermal Grüneisen parameter, γT , often include approximations such as an isotropic crystalline
response. However, γT can be readily obtained from ambient measured bulk properties of a material,
yielding values typically lie between 1 and 2 for most well-behaved materials [24]. The thermal
Grüneisen parameter for α-Fe is 1.81 [24], and the thermal Grüneisen parameter for Fe3 C is between
2.0 and 2.4 depending on which values for the bulk modulus you trust. This γT can then be used in
the microscopic definition to scale observed phonon frequencies with temperature
ω(T ) = ω0 (1 − γT

V − V0
).
V0

(3.8)

If the Debye model provides the phonon DOS, then the Debye DOS can be scaled with thermal
expansion to provide the QH vibrational entropy of a material. However, the Debye DOS can
readily be replaced with an experimentally determined phonon density of states without altering the
nature of the model. A quasi-harmonic model that utilizes an average thermal γT has no frequency
dependence; all phonons shift in energy the same way. The effect of the quasi-harmonic model is
completely independent of the vibrational spectra being scaled.
The thermodynamics effects of thermal expansion based phonon softening (phonons shift to lower

16
energy) can be ascertained by calculating the heat capacity and the vibrational entropy from the energy shifted phonon DOS. The quasi-harmonic vibrational entropy can be calculated using Eqn. 3.4,
and the heat capacity can be calculated using Eqn. 3.5. Under the quasiharmonic approximation
the heat capacity can also be re-written to directly include thermal expansion using a Grüneisen
parameter without using the phonon DOS. This expression is called the Nerst-Grüneisen expression,
CP = CV (1 + 3γT α(T )T ),

(3.9)

where α(T ) is the linear thermal expansion. The thermodynamic contributions from the quasiharmonic model are compared with the harmonic Debye model curves of the previous section in Fig.
3.3.
10

0.96
0.94
0.92

ΓT HTL
0 200 400 600 800 1000

Cp HtotL

CP HqhL

CP HhL

Temperature HKL

10
Svib HkB atomL

0.98

Γ = ΓT H300L

Cp HkB atomL

ΩQHA Ω0

1.00

0 200 400 600 800 10001200
Temperature HKL

STotal

Svib HqhL

Svib HhL

200 400 600 800 10001200
Temperature HKL

Figure 3.3: Debye model phonon density of state for ΘD = 420, a typical value used for α-Fe [19].
The corresponding vibrational heat capacity compared with α-Fe experimental heat capacity [20].
The Debye vibrational entropy compared with the total entropy of α-Fe from the SGTE database
[21].
High temperature calculations with the QHA often employ a constant Grüneisen parameter,
though the bulk thermal properties encompassed in the thermal Grüneisen parameter have been
observed to vary with temperature in a number of materials. In an attempt to rectify this problem,
a temperature-dependent Grüneisen parameter can be used. For the case of pure iron, a wealth of
information on temperature-dependent properties is available. Multiple assessments of temperature
dependent thermal expansion [22, 25–27] and bulk modulus [28–30] can be used. The heat capacity
at constant volume, CV , can also be calculated in a number of ways with different approximations.
We examined the temperature dependence of the Grüneisen parameter by constructing two
separate functions for γT (T ) that sampled the full range of observed values in the literature. The
results of these assessments are shown in Fig. 3.4, and show variations in γT (T ) between 1.5 and
2.1 over the temperature range of interest. Phonon frequencies can be scaled according to the
temperature-dependent Grüneisen parameter using the following expression:

ω(T ) = ω0

TY
i =T
Ti =1

[1 − γT (Ti )(

V (Ti )
− 1)],
V (Ti−1 )

(3.10)

ΓT HTL: KT -RayneDever
Α-TPRC, Cv -Debye, Ρ-TPRC

ΓT HTL

2.5

2.0
ΓT HTL: KT -AdamsDever
Α-Liu, Cv -NRIXS, Ρ-Seki

1.5

200

400

600

800

Temperature HKL

1000

1200

Figure 3.4: The temperature-dependent Grüneisen parameter γT (T ) assembled from various literature values. Black γT (T ) from KT [28, 30], α [22], CV from a Debye model (ΘD =420K), and ρ [22].
Red γT (T ) from KT [28, 29], α [25], CV from 14K NRIXS measurements, and ρ [26].
which reduces to Eqn 3.8 when γT is a constant. The result of this more careful calculation was nearly
the same as calculated quasi-harmonic change in frequency for a constant γT = γT (300K), as shown
in the left panel of Fig. 3.3. The resulting changes in the vibrational entropy from including γT (T )
were always below 0.5% at 1180K. Therefore, in the case of bcc iron, the addition of a temperaturedependent Grüneisen parameter has only a small effect on the calculated phonon energies. The
temperature dependence is influenced much more by the selected thermal expansion values.
Experimental heat capacities provided some of the earliest data for comparing observed thermodynamic properties with theoretical models. Experimental heat capacities and the thermodynamic
connections to the free energy provided a theoretical framework for the exploration of many metal
alloys, including iron. Efforts to exploit the predictive powers of thermodynamics of elemental metals to understand alloys drove many early (and modern) thermodynamic studies. Early studies of
thermodynamic contributions from various physical excitations like those conducted by Weiss employed simple Debye-like model [19]. Weiss sought to categorize the thermodynamics contributions
for pure iron by applying theoretical frameworks to the experimentally determined heat capacity. He
used a quasi-harmonic Debye model to account for the vibrational contribution and linear electronic
Grüneisen parameter model to account for the electronic contribution. Assuming the remainder
to be magnetic, he quantified the thermodynamic contributions to the free energy of pure iron in
Fig 3.5. This assessment suggested that in the absence of magnetism, FCC iron should be the
entropically-stable configuration. However, even this early assessment by Weiss pointed out that
a large quantity of the heat capacity, and an even greater share of the entropy, evolved from the
lattice vibrations as shown in Fig. 3.5. The model he used, however, is a quite simplified version,

18
Phonons
Electrons
Magnons

600
400

TS HmeVatomL

H HmeVatomL

800

200
400

600

800

1000

Temperature HKL

1200

-200
-400
-600
-800
400

Phonons
Electrons
Magnons
600

800

1000

Temperature HKL

1200

Figure 3.5: The enthalpic and entropic contributions to free energy [19].
and captures only aggregate vibrational behavior. His assumption that magnetic entropy must be
the remaining unassigned entropy points out a shortcoming that has long thwarted thermodynamic
assessments of magnetic materials. Few models exist that can accurately describe finite temperature
magnetic disorder, and there are even fewer experimental studies for comparison. The state of spin
disorder in magnetic materials at finite temperatures remains an active field of study today [31].

3.4

Anharmonic Effects

The harmonic approximation of the interatomic potential simplifies many aspects of the physics of
vibrations, and this approximation is normally valid for low temperatures. However there are several
important physical phenomena that cannot be resolved by harmonic descriptions. The harmonic
model assumes that phonons are quantum harmonic oscillators, which are non-interacting. However,
we know that phonons do interact in real systems; phonons may interact with each other and
also with other excitations that occur in real materials. Thermal expansion is the most obvious
anharmonic thermal effect, and while the quasi-harmonic approximation may improve the accuracy
of thermodynamic models to deal with observed thermal expansion – it is by no means a physically
rigorous approach. Nonharmonic interactions also result in routinely-observed phenomena like finite
thermal conductivities, which result from phonon scattering events that indicate real phonons have
finite lifetimes. Phonons are also known to play a role in electrical resistivity, where electronic
carriers scatter, creating phonons that can increase the temperature of a material.
Models that go beyond the harmonic approximation typically do so using perturbation theory.
Additional terms from the Taylor series in Eqn. 2.5 are retained for a more physical description of
the interionic potential. Additional terms in the periodic potential are often thought of in terms of
the quantum interactions they might produce. The cubic term encompasses three phonon processes,
such as two phonons combining to produce a new one, or a phonon decaying into two others. The
quartic term, by extension, enumerates the four phonon processes, including a phonon decaying into

19
three new phonons, or two phonons interacting to create two new phonon states. These processes are
governed by the laws of energy and momentum conservation. Thus not all combinations are possible;
new states must have the proper energy and momenta, which are governed by the allowed crystal
modes and described by the phonon dispersions. While the probability of multiphonon processes are
low at low temperatures, they do have appreciable effects at finite temperatures when large phonon
populations are present in the material.
In perturbation theory, it is often possible to keep only the next highest order term to improve
on the physical description. While this could be done with the interionic potential, there are many
physical arguments for keeping both the cubic and quartic term. The cubic term is asymmetric in
nature, creating physical situations where the potential may become unstable if only this term is
applied. There are also many crystalline symmetry constraints on the phonon processes produced
by the cubic interaction term that limits the number of anharmonic interactions that are described
by this formalism. Incorporation of the quartic term improves the limiting behavior of the net
potential (since a Hamiltonian retaining only the cubic term may be unstable) [11]. Further, observations of high temperature phonon behavior suggest that in many instances quartic interactions
may contribute comparable thermodynamic effects to those from cubic term interactions.
Experimental characterizations of phonons in materials at high temperatures can begin to quantify the importance of these effects. Anharmonic phonon-phonon interactions affect the observed
phonon spectra by both shifting their absolute energies, and also broadening their energy and qspace signatures. This is apparent in measurements of phonon dispersions and DOS measurements
when the thermal broadening of specific features overcomes the instrument resolution, resulting in
a broadening that scales with temperature. In high temperature phonon measurements of bcc Ti,
Zr, and Hf, very broad phonon signatures have been resolved in specific q-space directions. This
anomalously large broadening extends over a significant energy range and has been implicated as
a dynamic precursor of the first-order martensitic transformation between the bcc and hexagonal
crystal structures [32–34].
The anharmonic phonon broadening caused by more typical cubic and quartic interactions are
expected to be Lorentzian in nature [35]. In phonon DOS measurements the anharmonic broadening
of DOS features can be modeled using a modified Lorentzian function [14, 36]. This construction
provides a route for estimating mean phonon lifetimes from spectra broadenings. Results from
investigations of relatively isotropic elemental metals have shown that anharmonic effects are often
linear with temperature.
A detailed mathematical derivation of the effects of the cubic perturbation theory terms in the
calculation of thermodynamic variables is performed by Wallace [35]. He states that to leading order in perturbation theory, anharmonic effects on the entropy are likely to be linear in temperature.
Anharmonic free energy contributions have quadratic temperature dependence accordingly, both in

20
contributions from the entropy and the internal energy. When considering only lattice contributions,
anharmonic effects should cause the heat capacity to increase linearly in temperature, rather than
approaching the Dulong-Petit limit. Wallace also cautions that since phonon measurements inherently capture the anharmonic phonon frequencies, the expression for calculating thermodynamic
contributions given in Section 3.1 should be used with some caution. Using measured frequencies
in calculating the vibrational entropy as given in Eqn 3.4 is accurate to first order in perturbation
theory. But the anharmonic calculation of the heat capacity using Eqn 3.5 or internal energy from
measured anharmonic vibrational frequencies can lead to errors in accounting for large anharmonic
contributions.

Chapter 4

Experimental Methods
4.1

The Mössbauer Effect

Rudolf Mössbauer won the 1961 Nobel Prize in Physics for “his researches concerning the resonance
absorption of gamma radiation and his discovery in this connection of the effect which bears his
name” [8]. The effect which bears his name is the recoilless nuclear resonance absorption of gamma
rays by nuclei. The concept of resonant nuclear absorption and flourescence preceded Mössbauer’s
initial work; however, the phenomena had not been observed efficiently because it was argued that
the nuclear recoil should alter the energy of the photon emitted from the decay of the nuclear excited
state. When Mössbauer observed this phenomena in 1958 he realized that recoilless nuclear resonant
absorption was possible because the recoil was not confined to a single nucleus, but rather the entire
crystal in which that atom was embedded. The recoil momentum is taken up by the entire crystal,
whose mass is much much greater than a single nucleus, and accordingly the energy shift associated
with the recoil upon gamma ray emission is negligibly small [37]. The efficient nuclear resonance
comes from a finite probability that the energy transferred to the crystal lattice during a nuclear
decay occurs without the excitation of a vibration in the lattice. These processes produce gamma
rays capable of re-exciting other resonant nuclei in the lattice. There are many nuclei that exhibit
recoil-free resonance, but the ease of which a nuclear excitation can be induced and observed have
limited Mossbauer spectroscopy to a few more practical isotopes.
Mössbauer observed his effect by selecting the appropriate source and absorber pair. This pair
must contain the same nuclear isotope, but the source contains the Mössbauer isotope in an excited
state, while the absorber contains the same isotope in its ground state [38]. The recoilless nuclear
resonant gamma transition can then be observed by varying the energy of the gamma rays from the
source by the small energy of a Doppler shift. If a sufficiently thin absorber is employed, the transmitted gamma ray intensity exhibits a Lorentzian lineshape with a width of twice the Lorentzian
energy width of the nuclear excited state [37]. This measurement can provide inherent information
about the nuclear excited state, but the crystal environment often alters the energy of the resonant

22
absorption of resonant nuclei. The small energy alterations from hyperfine interactions produce
unique signals in Mössbauer spectroscopy about the local electronic environment of the resonant
nucleus. To date the Mössbauer effect has been demonstrated on more than 100 different isotopes,
but only a few are practical with conventional Mössbauer spectroscopy [38]. Rudolf Mössbauer’s
original work used the 191 isotope of iridium, but today 57 Fe is the most commonly used isotope for
Mössbauer spectroscopy. The ubiquity of iron in many geological, technological, and biological materials makes 57 Fe Mössbauer spectroscopy the most commonly employed nuclear resonant isotope.
Natural iron contains about 2% of the 57 Fe isotope, which is sufficient for many studies of natural
materials. This thesis work is on iron alloys, so the following discussion of Mössbauer spectra will
focus on the 57 Fe nuclear excitations.
Mössbauer spectroscopy utilizes a Doppler shift to sweep through energies around the nuclear
resonant energy. Depending on the local environment of the nuclei, a number of hyperfine interactions can be identified in the spectra. These features appear as the splitting of the nuclear resonant
energy levels as a result of the local electronic environment at the nucleus. In addition to the recoil
free fraction, the number of nuclear excitations that produce a resonant recoil-free gamma ray, there
are three other quantities that are commonly measured with Mössbauer spectroscopy. These three
quantities are from the hyperfine interactions, and their schematic effects on the resonant nucleus
are illustrated in Fig. 4.1 and Fig. 4.2.
The isomer shift is a hyperfine interaction that is present in every Mössbauer spectrum. The
isomer shift is a measure of the electron density in the nucleus. It arises from size differences between
of the nuclear excited state and the ground state, which interact with the electron wave function at
the nucleus. The isomer shift results in a small change of the nuclear resonant excitation energy as
shown in Fig. 4.1. In the absence of other hyperfine effects, the isomer shift can be used to track
changes in occupation of the s-electronic states that have finite density at the nucleus.
The electric quadrupole splitting is another hyperfine interaction observed in the Mössbauer
spectra of many materials. This effect splits the excited state of the nucleus as part of an interaction
with the electric field gradient at the nuclear position. The electric quadrupole splitting arises from
the electric quadrupole moment of the nucleus, which has a different orientation in the nuclear
ground and excited states. In the presence of an electric field gradient, which occur at crystal sites
with less than cubic symmetry, the nuclear resonant feature in the spectra is split symmetrically
about the centroid dictated by the isomer shift.
The final hyperfine parameter is the magnetic hyperfine splitting, which further divides the
nuclear excited state in response to magnetic fields present in the material. A hyperfine magnetic
field at the nucleus lifts all degeneracies of the nuclear excited states, as the spin of the nucleus
interacts with the local magnetic field as shown in Fig 4.2. This creates a Mössbauer pattern with
6 features in ferromagnetic bcc iron, whose spacing can be directly related to the strength of the

23
In = ± 3/2

In = ± 1/2
In = ± 3/2

∆E ~ neV

∆E 57Fe =
14.41 keV

In= ± 1/2

In= ± 1/2
Isomer
Shift

Isolated
Nucleus

Electric Field
Gradient

Figure 4.1: Changes in the nuclear excitation that correspond with an isomer shift and an electric
field gradient.

In = + 3/2

In = ± 3/2

In = + 1/2

∆E ~ neV

In = - 1/2
In = - 3/2

∆E 57Fe =
14.41 keV

In= ± 1/2

In = - 1/2
In = + 1/2

Isolated
Nucleus

Isomer
Shift

∆E ~ neV

Hyperfine
Magnetic Field

Figure 4.2: Changes in the nuclear excitation that correspond with an isomer shift and a hyperfine
magnetic field.
magnetic field at the nucleus.
The hyperfine interactions of many resonant isotopes in different materials have been thoroughly
studied throughout the years, providing a near fingerprint recognition quality to this method. The
method is routinely used to characterize the different atomic environments present in both natural
and engineered materials. There is a 57 Fe Mössbauer spectrometer on the Mars rover [39]. Additionally, the technique has been used to study dynamic phenomena such as valence fluctuations and
melting [40, 41].

4.1.1

Nuclear Forward Scattering

New high-brilliance synchrotron sources enable a different approach to Mössbauer spectrometry.
The synchrotron offers a new way to excite resonant nuclei without radioactive sources. Ruby first
proposed this experimental arrangement in 1974 [42], though it was several years before synchrotron
Mössbauer spectroscopy (SMS), also known as nuclear forward scattering (NFS), was realized experimentally. In addition to the high flux of third-generation synchrotron sources in the energy range
of several resonant nuclei, single crystal silicon monochromators played an important role in the
technique by refining the incoming x-ray spectrum [43].
The incoming synchrotron pulse interacts with the resonant nuclei, producing a coherent beam
in the forward direction that is modulated by the hyperfine splittings of the nuclear resonant states
[44](Smirnov for more detail). Unlike conventional Mössbauer spectrometry, which relies on a continuous photon signal modulated in energy, SMS uses a temporally-resolved photon pulse to excite
the nuclei. The synchrotron pulse arriving at the sample excites both the nuclear resonance and
many other atomic excitations. However the relatively long lifetime of the nuclear excited state
makes it possible to make high fidelity measurements of the temporal dependence of nuclear decays
4.3. This

Log[Intensity]

long after electronic processes have abated. This effect is shown schematically in Fig.

153
Time Delay (ns)

306

Figure 4.3: Synchrotron pulses arrive separated by 153ns. The signal from electronic excitations
(blue) occurs within a new nanosecond, while the signal from nuclear resonant excitations (orange)
can be monitored between pulses, when the background signal from other processes is extremely
low.
measurement technique is enabled by the synchrotron timing structure, which produces short pulses
of radiation followed by relatively long periods where no radiation arrives. When NFS spectra are
collected, the detector is gated electronically to reject photons emitted during the initial pulse arrival
and shortly thereafter. The NFS signal from nuclear resonant decay is then collected before the next
pulse arrives. Fast, low background avalanche photodiode detectors are essential for this work.
The results presented in this thesis were collected at the Advanced Photon Source at Argonne
National Laboratory. The Advanced Photon Source is the only facility in the country that routinely
runs NFS experiments. The phonon pulses at the Advanced Photon Source are approximately 34ps
in duration, and arrive at 153ns intervals under standard operating mode. This aligns well with

25
the nuclear resonant lifetime of 57 Fe, which is 141ns, permitting characterization of the interfering
temporal signal from the nuclear resonant excited states over approximately one nuclear lifetime.
NFS can be described as the Fourier transform of the more traditional Mössbauer spectroscopy
energy domain signal. The incoming pulse simultaneously excites all available resonant nuclei, and
the temporal decay of these excitations includes their interference, also known as quantum beats
[43]. Quantum beats develop in a transmitted NFS spectrum, providing information on the energy
splitting of the nuclear resonant levels. If only two nuclear resonant levels are present in a material,
the energy splitting can be inferred directly from the temporal width of the beat pattern T ≈ ∆E

[44]. However if more levels are present, or the spectra are significantly modulated by thickness
effects, the spectra are more complicated. Magnetic transitions that occur as a result of temperature
adjustments are particularly vivid in NFS spectra. NFS has been used as an in-situ monitor of
the state of magnetic order in a material. The abrupt loss of features as a magnetic material is
heated a valuable indication of the magnetic transition. These spectra are typically interpreted with
information from previous conventional Mössbauer studies and the use of specialized fitting tools
like CONUSS [45], which refines the hyperfine material parameters through iterative comparison
with the measured hyperfine spectra.

4.2

Nuclear Resonant Inelastic X-ray Scattering

The nuclear resonant inelastic x-ray scattering (NRIXS) technique makes use of nuclear excitations
to enable the study of vibrations in samples containing Mössbauer nuclei. The premise of acquiring
vibrational information of a material sample from Mössbauer spectroscopy was suggested by Rudolf
Mössbauer in his early work [37]. However an effective methodology for obtaining the vibrational
spectrum of a material was not demonstrated until decades later [46]. The interactions between
vibrations and nuclear resonant excitations was well understood as the method was developed, but
it took the high photon flux and timing structure of the synchrotron pulse to make these interactions
measurable [42]. The technique differs from many other measures of vibrational properties in that
NRIXS is incoherent scattering and inherently samples the full q-dependence of vibrations. The
NRIXS technique makes use of the relatively long lifetime of the nuclear excited state and the
timing of the synchrotron source to create a very low background technique that samples the full
phonon density of states of the nuclear resonant species. However, the technique is only sensitive
to the motions of resonant nuclei, which can be viewed as both a limitation and an advantage. The
limitation of this selectivity is that the resulting phonon DOS is only a partial density of states, i.e.,
the phonon spectrum of atomic motions of the resonant nuclei. The phonon DOS acquired from
NRIXS is therefore defined as a partial density of states (pDOS) for every alloy or compound, but is
a true phonon DOS for elemental materials. An obvious advantage of this selectivity is an absence

26
of background signal. Since nuclear resonant isotopes are the only source of the time gated signal,
NRIXS does not suffer from background signals acquired from auxiliary equipment or materials in the
beam path. This is a clear advantage over neutron techniques where scattering signals from auxiliary
equipment, such as a furnace or cyrostat, or impurities (especially hydrogen-containing impurities
like water) can significantly obscure the signal from the sample. The selectivity of NRIXS vibrational
methods also permits focusing on specific portions of a material by restricting Mössbauer isotopes
to localized areas. Recently, several layered systems have been probed selectively by isotopically
enriching specific regions (layers or nanoclusters), which can be probed with high accuracy [47–50].
The atomic selectivity of NRIXS also makes it a complementary tool to total DOS spectra measured
by other means. The vibrational dynamics of binary alloys can be probed systematically, by collect
NRIXS on two separate resonant species, in the case where both elements have a resonant isotope
[51]. Binary alloys can also more closely analyzed by comparing a NRIXS pDOS with the total DOS
Phonon
Annihilation
Resonance
obtained
by eutron
techniques, as has been demonstrated
in several binary alloysPhonon
[52–61].Creation
In =energy
3/2 of the
In = 3/2probes the phonon spectrum
In = 3/2
The NRIXS technique
of a material by scanning the

incoming phonon. However, unlike traditional Mössbauer spectroscopy which uses the Doppler shift
∆E 57Fe by tens of neV, the high resolution
Fe
∆E 57Fe monochromators alter ∆E
to adjust the incident energies
the57incoming
14.41 keV
14.41 keV
14.41 keV
energy by tens of meV. Detuning the incoming x-ray beam away from the nuclear resonance results
In= 1/2 nuclear
In= 1/2x-ray can only excite the narrow
in a strong suppressionIn=of1/2
scattered photons. The incoming

resonance if it obtains additional energy through interactions with the sample. This is done by
creating or annihilating a phonon in the sample to obtain the nuclear resonant energy, as shown
schematically in Fig.

4.4. As the energy of the incoming photons is varied, only phonons (or

Phonon Annihilation

Resonance

Phonon Creation

In = ± 3/2

∆E 57Fe
14.41 keV

In= ± 1/2

Figure 4.4: An illustration of the NRIXS excitations processes as the incoming x-ray energy is tuned
through the resonant energy. A photon with an energy below the resonant energy (red) must absorb
a phonon of the appropriate energy to excite the nuclear resonance. A photon with an energy above
the resonant energy must give up energy by exciting a phonon before it can excite the nuclear
resonance.
combinations of phonons) that precisely match the difference between the incident energy and the
nuclear resonant excitation energy will create scattered photons. The signal is integrated across the
full time window, without regard to the hyperfine structure modulations. Therefore the observed

27
scattering as a function of energy can be directly mapped as shown in Fig 4.5. At low temperatures

1180K
1081K
981K
870K
740K
647K
510K
298K
130K
30K
-40

-20

20
Energy HmeVL

Figure 4.5: Raw NRIXS scattering spectra that shown vibrational excitations as a function of scattering energy. Spectra are collected on a polycrystalline sample of α-57 Fe over a range of temperatures.
the number of phonons available for annihillation is quite low, so scattering is suppressed when the
incoming photon energy is detuned below the nuclear resonance. However, phonons can still be
readily excited at low temperatures by the incoming beam so scattering above the nuclear resonant
peak still occurs. At higher temperatures the creation and annihilation components of NRIXS
scattering spectra start to be comparable in size.
The partial phonon DOS can be extracted directly from NRIXS scans using the software package
Phoenix [45]. This analysis of the NRIXS scattering spectra, S(E), requires removing the resonant
peak, accounting for the Debye-Waller scattering factor, and generating a phonon DOS consistent
with the remaining scattering signal. The phonon DOS is generated self consistently by utilizing
both the thermally-weighted phonon creation and annihilation signals to generate a DOS function.
This involves more than single-phonon processes, and includes a contribution from two-phonon

28
processes, and n-phonon processes which scale by the approximate probability of those processes at
the temperature of observation. Several phonon DOS spectra from NRIXS measurements of α-57 Fe
are shown in Fig 4.6.

1180K
1081K
981K
870K
740K
647K
510K
298K
130K
30K

20
30
Energy HmeVL

Figure 4.6: Phonon DOS of α-57 Fe at a range of temperatures.
Analysis of the phonon spectra show that, as expected, the vibrational modes of α-Fe do change
with temperature. All the phonon frequencies shift to lower energies at elevated temperatures, which
can also be called“phonon softening”. The observed thermal expansion of Fe from Fig 3.2 accounts
for part of the phonon softening, as the volume of the unit cell increases the energy required to
excite lattice vibrations is slightly reduced. This is evidence that a strictly harmonic model will
not account for the vibrational contributions to thermodynamics for this system, which obviously
exhibits nonharmonic effects.

Chapter 5

Computational Methods
5.1

Computational Quantum Mechanics

The development of quantum mechanics was arguably the most significant advance in understanding
the physical behavior of materials in the 20th century. Quantum mechanics provides a methodology
for predicting how atoms and electrons behave in different spatial configurations. But while quantum
mechanics is an extremely accurate formalism for describing the behavior of matter in simple systems
like the quantum harmonic oscillator, applying quantum mechanics robustly to less idealized manybody systems requires much mathematical rigor and is often intractable [62]. For systems with
several atoms, one must consider the full three dimensional wavefunction character of every nucleus,
every electron, and the interactions between them. The time-independent Schrödinger equation
provides a physical description of matter in the ground state through the expression,
ĤΨ = EΨ,

(5.1)

where E is the scalar ground state energy, Ĥ, the Hamiltonian operator, and, Ψ, the system wavefunction that is an eigenvector of the Hamiltonian. When one considers the interactions of Ni
ionic cores with Ne electrons, the Hamiltonian operator is recast to enumerate the relevant physical
interactions as
ĤΨ = {−

Ne
Ni
Ne
Ni
Ne X
Ni
h̄2 X
h̄2 X
∇2
e2
Zu Zv
eZu
}Ψ = EΨ
∇2 −
2me j=1
2 u=1 Mu
|rj − rk | u,v=1 |ru − rv | j=1 u=1 |rj − ru |
j,k=1
k

(5.2)
where the terms account for the kinetic energy of the electrons, kinetic energy of the ions, potential
energy of electron interactions, potential energy of ion interactions, and the potential energy of
electron-ion interactions, in that order [63]. This many body problem is a formidable mathematical
challenge for small systems and so the predictive power of quantum mechanics for real material

30
systems went unexploited for many decades.
Walter Kohn was awarded the 1998 Nobel prize in chemistry “for his development of the densityfunctional theory” [8]. Density-functional theory solves the Schrödinger equation in 5.2 by making a
few key assumptions and several clever observations to provide computationally-tractable predictions
of real material systems [62, 63]. These points are briefly summarized here. The first inherent
assumption is the Born-Oppenheimer approximation, which states that the ions can be assumed
stationary on the time scale of electronic interactions, which is justified by the large mass mismatch
between the two [62]. This eliminates the second term of Eq. 5.2, and provides constant values for the
fourth term. The Kohn-Hohenberg theorems state that the ground state energy of the Schrödinger
equation is a unique functional of the electron density, and the energy density that minimizes the
functional is the true ground state electron density of the Schrödinger equation [63, 64]. The energy
of the configuration can be then computed by finding the proper electron density (a 3-dimensional
function), rather than the electron wavefunction (a 3Ne dimensional function). Kohn and Sham
approached the problem by assuming a separable set of electron wave functions that satisfy individual
electron Hamiltonians, called the Kohn-Sham equations (which are dependent on knowledge of the
total electron density), and wrapping the unknown bits of physics into an exchange-correlation
potential Vxc which is necessarily approximated [65]. The single electron Hamiltonian can then
be solved self-consistently until the electron density calculated from the solution wave functions is
sufficiently similar to the electron density used to compute those wavefunctions.
This is a very brief summary of the merits of DFT, but the creation of this methodology has
spawned a number of first principles methods that provided useful predictions of real materials
behavior. DFT simply requires information on the ions in a material and their locations, then the
ground state energy of this configuration can be calculated. From these electron energy calculations,
electronic band structures provide information on electronic transitions. Additionally, the structural
arrangements of atoms input into DFT can be refined to find the lowest energy configuration for
the symmetry of the system. Multiple structural symmetry arrangements may be calculated to find
the relative energy cost of different crystal structures, and thus the energy competition of different
material phases can be analyzed at 0K (the ground state).
Modern density functional codes like the Vienna Ab Initio Package (VASP), have built on this
core DFT capability to provide a range of computed material properties [66]. For instance, the
calculation of phonons is now becoming somewhat routine, even though DFT requires fixed ion
positions to calculate electron densities [67, 68]. This is accomplished by cleverly choosing atomic
displacements in the relaxed atomic structure. Displaced atoms provide electron configurations that
are considerably higher than the ground state configuration, which can be related to interatomic
forces which can be used in the dynamical matrix formalism to calculate phonon spectra [69]. Similarly elastic material response can be extracted by examining the electron response to spatial strains

31
of the equilibrium unit cell atomic structure. These adaptations are more computationally costly
because they require large supercells and many additional steps, but they are routinely performed
as part of computational DFT studies today.
Density functional theory is an extremely powerful tool and has fueled a revolution in materials
physics, however it is an inherently ground state methodology, providing properties at 0K. Real materials, however, are very rarely prepared or used under such conditions. More recent developments
have focused on adapting this methodology to explain properties of materials at the ground state,
and non-equilibrium behavior. For instance if we are interested in phase diagrams of materials (which
I am) we need to consider the energy and entropy of materials at finite temperatures. This returns
us to dealing with non-harmonic effects like thermal expansion, which plays an important role in
material phase stability and entropy at finite temperatures. There is a DFT approach to thermal
expansion with a quasi-harmonic model, although typically different from the formalism developed
in Chapter 3. Typical DFT quasi-harmonic thermal expansion is calculated by minimizing the free
energy of the material (including electronic and vibrational entropic contributions) with respect to
the volume of the structure. This is computationally expensive because it requires iteratively calculating vibrational frequencies with the methods described above. Nevertheless, it does provide a
completely first principles approach to thermal expansion that does not require experimental input
on thermal behavior.
While the strength of DFT is ground state phenomena, molecular dynamics (MD) is an inherently high-temperature computational approach to materials behavior. Molecular dynamics (MD)
provides a classical computational approach for calculating atomic motions in materials [70]. The
atoms in the system (with their electrons included) are described by an interaction potential which is
dependent on the distances between the atoms involved. Newton’s equations of motion can then be
solved for the forces on each atom, determined from the interaction potential,to bring the atoms into
their next spatial configuration. The methodology can be very accurate when very small time steps
are used between calculating the forces and propagating their motion, and, more importantly, that
the interatomic potential correctly reproduces the relevant atomic interactions at every length scale.
The first requirement can be handled readily by rapidly expanding computational power, which permits many time steps on relatively large systems compared with DFT. The second requirement is
somewhat more challenging, since some inherently quantum phenomena in atomic interactions may
not be readily reproduced by simple interatomic potentials. These methods have been advanced by
working in conjunction with first principles simulations, which can provide sophisticated parameterizations of interatomic potentials. However, MD may succeed when DFT fails in capturing high
temperature and non-equilibrium behavior of materials. This is achieved by constraining the energy
in the system using thermostats that accurately reproduce finite temperature effects (that are within
the scope of the interatomic potential).

32
Ab Initio Molecular Dynamics (AIMD) is a combination approach that uses DFT within each MD
time step to realize quantum mechanical behavior at finite temperatures and time scales [66]. This
hybrid approach is more computationally expensive than either method in isolation, but it provides
benefits that could not be realized by either individually. AIMD has been increasingly utilized
to examine the finite temperature behavior of materials, including thermal expansion and high
temperature phonon behavior. Methods for extracting vibrational information from AIMD are still
an active field of study, with very promising early results [71, 72]. They present a unique opportunity
to examine computed vibrational material behavior beyond the limitations of the harmonic and
quasi-harmonic models without many-body perturbation theory.

5.2

Born von-Kármán Fitting of Phonon Spectra

The harmonic model for lattice vibrations constructed in Section 2.2 can be readily applied to experimental measurements given the proper set of assumptions, as noted in Section 2.3. Often called
the Born von-Kármán (BvK) model, this formalism relates interatomic force constants of crystalline
solids to phonon frequencies [11]. To solve for phonon frequencies, the number of interatomic interactions in the dynamical matrix is typically limited to the first few nearest-neighbor interactions.
This truncation of interaction forces limits the size of the dynamical matrix (which might otherwise
be as large as the square of the number of atoms involved). Short-range forces are quite physically
reasonable for well behaved materials. The dynamical matrix can then be solved for specific phonon
frequencies along the direction of interest in k-space, assuming one has the interatomic force constants in the dynamical matrix D(k). The interatomic force constant elements of the dynamical
matrix are defined as derivatives of the crystal potential, but this potential, even if harmonic, cannot
be derived analytically for the response of real solid material.
When experimental phonon dispersion are available, the elements of the dynamical matrix can
be fit to observed phonon energies and momenta by a least squares approach, by constructing a
linear set of equations for each observed frequency,

ω12



1 (k1 )

D(k1 )



1 (k1 )

 


 2 
 

 ω2   2 (k2 )   D(k2 )   2 (k2 ) 
.
M  . 
..
..
..
=

 ..  

 

ωN
N (kN )
D(kN )
N (kN )
The force constant tensors for the first five nearest neighbors of the bcc structure from phonon

33
dispersion fits [12] are

1XX

K1NN = 1XY
1XY

1XY
1XX
1XY

2XX
1XY K2NN =  0
1XX
1XY

2YY

3XX

3XY
3NN
2YY

3XY
3XX

0 
3ZZ
(5.3)

4XX
K4NN = 4XY
4XY

4XY
4YY
4XZ

5XX
4XZ  K5NN = 5XY
5XY
4YY

4XY

5XY
5XX
5XY

5XY

5XY .
5XX

(5.4)

The bcc Fe phonon dispersion data collected by neutron triple-axis measurements shown in Fig 2.1
were fit to a Born von-Kármán model. The interatomic force constants that resulted from the
least squares fit are shown in Table 5.1. The force constant matrices can be projected onto their
Table 5.1: Interatomic force constants for bcc Fe at 300K from neutron triple-axis measurements
[12].
1NN

2NN

3NN

4NN

5NN

Kij
(N/m)

2 (111)

2 (200)

2 (220)

2 (311)

2 (222)

1XX = 16.88
1XY = 15.01

2XX = 14.63
2YY = 0.55

3XX = 0.92
3XY = 0.69
3ZZ = -0.57

5XX = -0.29
5XY = 0.32

Longitudinal
Ave. Transverse

46.9 N/m
1.87 N/m

14.63 N/m
0.55 N/m

1.61 N/m
-0.17 N/m

4XX = -0.12
4YY = 0.03
4XZ = 0.0007
4YZ = 0.52
0.0026 N/m
-0.313 N/m

0.035 N/m
-0.61 N/m

bonding directions to provide longitudinal and transverse force constants, which describe the phonon
vibrations parallel and perpendicular to the bonding axis, respectively. The two transverse force
constants perpendicular to the bonding direction have been averaged here. These longitudinal and
average transverse force constants provide a more intuitive comparison of effects across nearestneighbor interactions. The values of the longitudinal force constants are much higher than the
average transverse phonon modes. Additionally we see that the magnitude of the force constants
drop off rapidly for greater interatomic distances.
In the bcc structure, the different bond force constants have notably different effects on the
phonon spectra as demonstrated by Fig. 5.1 and Fig. 5.2. When a single longitudinal force constant
is varied, the vibrational effects vary throughout the phonon dispersions and the DOS. Scaling of
the first nearest-neighbor longitudinal force constant adjusts the energies of every phonon in the
dispersions except those in the lowest energy transverse phonon branch between the Γ and N points.
Scaling the second nearest-neighbor longitudinal force constant leaves some phonon branches fixed,
like the lowest energy modes between Γ and H, but strongly influences others, like the low energy Γ

DOS H1meVL

Energy HmeVL

40
30
20

140% L1
100% L1
60% L1

Energy HmeVL

Figure 5.1: The vibrational effects of varying the first nearest-neighbor longitudinal force constant
from the 300K observed values [12].

DOS H1meVL

Energy HmeVL

20
140% L2
100% L2

60% L2

Energy HmeVL

Figure 5.2: The vibrational effects of varying the second nearest-neighbor longitudinal force constant
from the 300K observed values [12].
to N branch, which was unaffected by the first nearest-neighbor longitudinal force constant.
While the dynamical equation used here is the result of a harmonic approximation of the interatomic potential, the force constants fits can be used to identify nonharmonic behavior from how
the force constants vary as a result of temperature or applied pressure. If the dynamical matrix
equation is fit to elevated temperature or pressure spectra, the deviation of the resulting fits from
their low-temperature values describes the non-harmonic nature of the interatomic potential. Information on the interatomic force constants can be quite valuable for interpreting how each set of
nearest-neighbor interactions contributes to the dynamical stability of the structure and the energies
of various phonon modes, and interpreting the nature of the bonding between each nearest-neighbor
shell of atoms.
In direct phonon DOS measurements like NRIXS, extracting interatomic force constants is less
traditional and less straightforward than for phonon dispersions. Phonon DOS measurements lack
k-space resolution, and many distinct k-space features are degenerate in energy. Furthermore,
an expression for the interatomic force constants cannot be directly compared with phonon DOS
spectra, because the DOS requires evaluation of the dynamical matrix over the whole Brillouin zone.

35
This step cannot be inverted, so solving the dynamical matrix for the phonon DOS is inherently
a forward model. However, modern computational optimization methodologies do offer routes for
tackling “univertable problems”. These methods simply explore large regions of parameters space
looking for the “best” answer by running many many iterations of the forward model.
To extract interatomic force constants from observed phonon DOS the method developed in our
group uses a genetic algorithm (named for its similarity to the principle of natural selection [73]).
The genetic algorithm optimization of interatomic force constants to fit a phonon density of states
is implemented through the open-source package mystic [74]. To fit a phonon DOS, a candidate
set of force constants is randomly generated within a set of reasonable force constant bounds. The
candidate dynamical matrix is used to generate a phonon DOS that can be compared with the
experimentally obtained DOS. The genetic algorithm optimization simultaneously generates many
candidate solutions called a “population”. The candidate solutions that best reproduce the experimentally observed phonon DOS (minimizing the mean squared error with respect to the observed
DOS) are selected as the parents of the next generation. The parent solutions seed the subsequent
population of candidate force constants by selecting random combinations of parent parameters and
also introducing random “mutations”. The next generation is then evaluated and the optimization
continues until the population converges on a set of force constants that provide the “best fit” to
the experimental density of states. These optimized force constants can then be used to generate
phonon dispersions with q-space resolved information.
Force constant optimizations will reproduce the phonon DOS by design, but the existence of
a single physically meaningful solution is not guaranteed. While the phonons generated by force
constant optimization are necessarily confined to the crystal symmetry embedded in the construction
of the dynamical matrix, nonphysical results are still possible. The most obvious instance of these are
force constants that generate negative phonon frequencies, which result from imaginary eigenvalues
of the dynamical matrix problem. Negative phonon frequencies suggest dynamical instability of
the model in specific directions in reciprocal space, suggesting that phonons initiated along these
wavevectors would proceed without restoring forces and destabilize the system. Therefore candidate
solutions that generate negative frequencies are discarded by our algorithm.
Force constant optimizations also cannot guarantee the uniqueness of a solution for a given
phonon DOS. The optimization will eventually conclude with a set of interatomic force constants
that reproduce the experimental DOS, but there is no guarantee that another optimized fit could
not be found with a notably different set of force constants. The issue of uniqueness is inherent in
phonon DOS optimizations, since distinct features in the dispersions are often degenerate in energy.
This problem is significantly exacerbated by the experimental resolution of the measurements, which
blur features together, adding ambiguity.

5.2.1

Fitting Phonon DOS - Practical Matters

There are several practical matters to consider when attempting to fit a phonon DOS using the
BvK formalism described above. Probably the most important is well defined phonon DOS features.
Density of states measurements with well-defined van Hove singularities have features that increase
the likelihood of finding unique and reliable solutions for the interatomic force constants. This
requires phonon energy dispersion features at distinct energies. Phonon dispersions with many
features in a small energy range are unlikely to yield high-fidelity BvK fits, especially with imperfect
experimental resolution. These feature-heavy energy regions are more common in low-symmetry
crystals with large unit cells, so simple elementary and binary systems are obvious candidates for
BvK fitting. The instrumental resolution of the experimentally-observed phonon DOS is another
critical input for the BvK fitting. While the highest resolution measurements are always preferred,
an accurate resolution function is necessary for comparing resolution-convolved calculated phonon
DOS with experimental spectrum.
One must also consider the number of variables that can accurately be fit to a phonon DOS
spectra. Born-von Kármán fits to phonon dispersions commonly employ five nearest neighbors, or
13 variables for bcc Fe. However, this many variables may not be well constrained in phonon DOS
fitting. To test the limits of the method with a known case, the phonon DOS derived from the
Minkiewicz force constants in Table 5.1 was calculated and convolved with a typical NRIXS instrument resolution function (a near-Gaussian function with a 2.3meV FWHM). This “experimental”
phonon DOS (with known force constants) was then fit using our BvK optimization. We varied the
number of parameters used in the optimization by including increasingly distant nearest-neighbor
interactions. This resulted in five different optimizations, from fitting with only the first nearest
neighbor (1NN) interactions and two force constants, to fitting with five nearest neighbors (5NN)
and 13 force constants. The χ2 values in Fig 5.3 show the quality of fit for each different fitting
model. The quality of fit improves with additional nearest-neighbor interactions, although notably
less beyond the 3NN configuration. The 5NN configuration provides the lowest χ2 value, indicating the most exact reproduction of the DOS. However, the quality of DOS fit does not necessarily
correspond to accurately reproducing the interatomic force constants.
For each fitting model, the errors in the force constants were assessed to see how well this
information is reproduced by the fitting in the presence of experimental resolution. The sum of
the absolute difference between the fitted force constants (fcs) and the known values are divided
by the number of fcs variables used to provide the quality-of-fit metric displayed in Fig 5.3. The
force constant fit metric suggests that including the 5NN force constants as variables decreases the
accuracy compared to the 4NN model. It is also important to note that the quality-of-fit from the
2NN and 3NN models are not significantly worse than the 4NN fits, suggesting that using only two
nearest neighbors may be sufficient for reproducing the phonon DOS.

SÈDfcsÈðfcs HNmL

Χ2 fit value

0.1
0.01
0.001
-4

1NN

2NN

3NN

4NN

5NN

2.0
1.0
0.5
0.2

1NN

2NN

3NN

4NN

5NN

Figure 5.3: Qualities of fit from various configurations of the BvK 300K Fe DOS optimization. The
left figure shows the χ2 values for the fit DOS. The right panel shows the errors between fitted force
constants and the known force constants used to generate the DOS.
The phonon DOS and dispersions from the fitted force constants are compared with the known
phonon dispersions, and the experimentally convolved DOS in Fig 5.4. It is important to note that
the phonon dispersion data were not used in the force constant optimization; only the experimental
resolution convolved phonon DOS (shown in black) was used. Therefore, agreement between the
calculated dispersions and experimental dispersions speaks to the merits of our fitting approach.
All models which incorporate more than one nearest-neighbor interaction reproduce the phonon
DOS and the density of states well. The 2NN-4NN neighbor fits are very similar in the quality of
dispersions and force constants, which suggests any of these might be a reasonable choice. While the
5NN model obviously contains the best phonon DOS fit, it does not reproduce the force constants
notably better than lower neighbor models, suggesting that 13 variables may be over fitting our
phonon DOS.
The inclusion of more neighbors for fitting larger data sets in Chapter 6 have also led to more
nonphysical behavior and the appearance of multiple basins of solutions. The most persistent
non-physical behavior observed is the appearance of non-linearity in the low-k phonon dispersion
branches. The fitted force constants dispersions showed a slight positive curvature of the low-k
branches that correspond to physically unreasonable group velocities for these modes.
The inclusion of more variables also tended to produce multiple basins of solutions within a
temperature series. These basins of solutions were distinct groups of fitted force constants which
produced comparable phonon DOS, but different dispersions. The most common solution basins
observed for bcc Fe were characterized by a difference in the high energy phonon dispersion modes at
the N and H points. One basin of solutions would over-estimate the H point modes and underestimate
the N point modes (like the 1NN fits from Fig. 5.4), while the other would show the opposite
behavior. Including fewer force constant variables decreased the likelihood of finding multiple basins
of solutions.

DOS H1meVL

Energy HmeVL
30

Energy HmeVL

DOS H1meVL

Energy HmeVL
30

Energy HmeVL

Energy HmeVL
10

Energy HmeVL

DFCS HNmL

Longitudinal Force Constant

1NN 2NN 3NN 4NN 5NN

2NN Fit

35
30
25
20
15
10

-1
-2

Ave. Transverse Force Constant

-1
-2

Longitudinal Force Constant

1NN 2NN 3NN 4NN 5NN

3NN Fit

35
30
25
20
15
10

-1
-2

Ave. Transverse Force Constant

-1
-2

Longitudinal Force Constant

1NN 2NN 3NN 4NN 5NN

4NN Fit

30
20

-1
-2

Ave. Transverse Force Constant

-1
-2

Longitudinal Force Constant

1NN 2NN 3NN 4NN 5NN

5NN Fit

-1
-2

1NN 2NN 3NN 4NN 5NN

5NN Fit
Χ2 =0.000066

DOS H1meVL

35
30
25
20
15
10

1NN 2NN 3NN 4NN 5NN

4NN Fit
Χ2 =0.000363

Ave. Transverse Force Constant

1NN 2NN 3NN 4NN 5NN

3NN Fit
Χ2 =0.000161

DFCS HNmL

DOS H1meVL

Energy HmeVL
30

DFCS HNmL

Energy HmeVL

DFCS HNmL

-1
-2

1NN 2NN 3NN 4NN 5NN

2NN Fit
Χ2 =0.002232

DFCS HNmL

Energy HmeVL

DFCS HNmL

Longitudinal Force Constant

1NN 2NN 3NN 4NN 5NN

1NN Fit

-1
-2

20
DFCS HNmL

Energy HmeVL

DOS H1meVL

DFCS HNmL

1NN Fit
Χ2 =0.039209

-1
-2

Ave. Transverse Force Constant

1NN 2NN 3NN 4NN 5NN

Figure 5.4: Fitting Details for 5NN Fits. Fitted phonon DOS and dispersions are shown in red.
Errors in the longitudinal and transverse fitted force constants are shown in the charts on the right.

Chapter 6

Anharmoncity in BCC Fe At
Elevated Temperatures
Phonon densities of states (DOS) of bcc α-57 Fe were measured from room temperature through the
1044K Curie transition and the 1185K fcc γ-Fe phase transition using nuclear resonant inelastic
x-ray scattering. At higher temperatures all phonons shift to lower energies (soften) with thermal
expansion, but the low transverse modes soften especially rapidly above 700K, showing strongly
nonharmonic behavior that persists through the magnetic transition. Interatomic force constants
for the bcc phase were obtained by iteratively fitting a Born-von Kármán model to the experimental
phonon spectra using a genetic algorithm optimization. The second-nearest-neighbor fitted axial
force constants weakened significantly at elevated temperatures. An unusually large nonharmonic
behavior is reported, which increases the vibrational entropy and accounts for a contribution of 35
meV/atom in the free energy at high temperatures. The nonharmonic contribution to the vibrational entropy follows the thermal trend of the magnetic entropy, and may be coupled to magnetic
excitations. A small change in vibrational entropy across the α-γ structural phase transformation
is also reported.

6.1

Introduction

In its metallic form, iron exhibits fascinating physics, plays a central role in geophysics, and is
of paramount importance to metallurgy. Iron is polymorphic under temperature, pressure, and
alloying, and both its magnetic properties and its mechanical properties undergo major changes
with crystal structure. The thermodynamics of the temperature-induced polymorphism of iron have
been of interest for many years. A proper thermodynamic treatment of metallic iron must consider
the energetics as well as the degrees of freedom of electrons, phonons, and spins, and the couplings
between them. Although this is a complex problem, it has received longstanding interest both for
its own sake, and for predicting the phases of iron alloys with an eye to controlling them [75, 76].

40
There have been a large number of heat capacity [20, 22] and elastic constant [28, 30] measurements of iron at various temperatures, and the thermodynamic entropy of iron is sufficiently reliable
to be used in Calphad-type calculations of free energy [21, 77]. There have been a number of efforts to create predictive thermodynamics models by resolving the entropy into contributions from
phonons, spins, and electrons [19, 31, 77]. Phonons make the largest contribution to the entropy
at elevated temperatures, and therefore the accuracy of the phonon entropy is critical. A harmonic
model can account for most of the vibrational entropy of elemental solids. The vibrational entropy
of iron is quite large, however, exceeding 6 kB /atom at 1000K, so even errors of a few percent are
thermodynamically important. The quasiharmonic model of vibrational entropy incorporates the
phonon frequency shifts that result from finite temperature thermal expansion, but it neglects many
other nonharmonic physical interactions. Phonons interact through the cubic and quartic parts of
the interatomic potential [35]. These anharmonic phonon-phonon effects further change the phonon
frequencies and shorten their lifetimes, resulting in thermal broadening of phonon spectra [14–16].
Thermal excitations of electrons and magnons also affect the phonon frequencies through adiabatic
electron-phonon and magnon-phonon interactions. The impact of these physical interactions on the
vibrational entropy and free energy has been shown to be important in many materials [78, 79], and
their role in the vibrational thermodynamics of iron warrants further investigation.
Inelastic neutron scattering studies of phonon dispersions in iron provide essential information
on the phonon contribution to entropy, and how it changes with temperature [12, 80–83]. High temperature phonon dispersions show significantly decreased phonon frequencies with thermodynamic
implications [83]. These measurements also provide insight into the mechanism of the polymorphic
transitions, and correlate with the inherent weaknesses of the bcc structure [33, 34, 84, 85]. However,
the existing experimental results are somewhat sparse in temperature.
Ab initio investigations have attempted to identify individual contributions to the free energy of
Fe and its alloys, but earlier studies relying on quasiharmonic approximations at high temperatures
had limited success. Electronic structure calculations on iron have advanced considerably in the
past few years, and recent work has carefully considered the different contributions of magnetism
and vibrations to the thermodynamics of the bcc phase. Only recently have computational developments permitted DFT calculations to reproduce the observed high temperature phonon behavior
by including the finite temperature magnetic configurations and electron-phonon coupling [86–89].
These developments suggest we may soon sort out the complex interactions in polymorphic iron,
but experimental validation is still needed. Measurements of phonon dynamics through the Curie
point at 1043K up to the γ-Fe phase transformation can provide further insight into the physical
interactions and thermodynamics governing the complex behavior of iron.
Here we report results of an experimental study of the vibrational properties of bcc α-Fe at
elevated temperatures, and an analysis of its interatomic interactions and thermodynamic functions.

41
Nuclear resonant inelastic x-ray scattering (NRIXS) was used to measure vibrational spectra of
bcc Fe and obtain reliable phonon densities of states (DOS) which, unlike most phonon dispersion
measurements, can be used directly in thermodynamic functions. A methodology was developed for
reliably extracting the temperature-dependent interatomic force constants, and consequently phonon
dispersions, from the phonon DOS spectra. Much of the high temperature nonharmonic phonon
dynamics depends on the rapid softening of the 2NN interatomic forces, and the resulting softening
of transverse phonons in Γ-N direction. The vibrational entropy is assessed with different models for
predicting high temperature thermodynamics. We report a large nonharmonic contribution to the
phonon entropy, and suggest that it originates with effects of magnetic excitations on the phonon
spectra.

6.2

Experimental

Nuclear resonant inelastic x-ray scattering (NRIXS) measurements were performed on bcc α-Fe at
high temperatures. NRIXS is a low background technique that provides direct access to the full
phonon density of states (DOS) [42, 46]. NRIXS spectra were collected from a 25 µm thick Fe foil of
99.9% purity and 95% 57 Fe isotopic enrichment. Measurements were performed at beamline 16ID-D
of the Advanced Photon Source at Argonne National Laboratory using a radiative heating furnace.
This NRIXS vacuum furnace used a narrow kapton window to permit the x-rays to access the sample.
The Fe foil was either held by two Ta heat shields adjacent to a thermocouple, or mounted directly
on the thermocouple. The NRIXS measurements performed below room temperature employed a
He flow Be-dome cryostat. The temperatures were accurate to within ±20K, where ambiguity comes
from comparing the furnace thermocouple measurements to in-situ nuclear forward scattering and
the NRIXS-derived detailed balance temperature calculations following the procedures described in
the literature [45, 90].
An avalanche photodiode was positioned at approximately 90◦ from the incident beam to collect
re-radiated photons beginning approximately 20 ns after the synchrotron pulse. The incident photon
energy was tuned to 14.413 keV using a high-resolution silicon crystal monochromator to provide a
narrow distribution of energies with a FWHM of 2.3 meV. The incident photon energy was scanned
through a range of ±120 meV, centered on the nuclear resonant energy. The Phoenix reduction
package was used to extract phonon DOS spectra from the collected spectra [45]. Lamb-Mössbauer
factors from this reduction are compared with literature values in Fig. 6.13 of the Supplemental
Material.

6.3

Force Constant Analysis

Many thermal properties of crystalline solids can be explained by a simple model of the crystal
as a set of massive nuclei whose interactions act like springs, providing a restoring force against
displacements. This model was developed by Born and von Kármán (BvK), and transforms the
real space interatomic forces into a dynamical matrix [35]. While this model is commonly employed
for fitting phonon dispersions, its utility for fitting phonon DOS spectra is less straightforward.
Phonon DOS spectra are an aggregate of all phonon modes in reciprocal space; therefore, fitting
force constants to a phonon DOS spectrum is more challenging than modeling phonon dispersions.
To model our phonon DOS spectra, trial force constants were used to construct a dynamical matrix,
D(q), which was diagonalized at a randomly-distributed set of 3.375 million q-points in the first
Brillouin zone to collect the spectrum of phonon frequencies, ω 2 ,
M ω 2 = D(q),

(6.1)

where M is the mass of the atom and is the polarization of the phonon mode corresponding to
reciprocal space vector q. This BvK model was embedded in a genetic algorithm global optimization
framework, where trial sets of force constants were generated randomly according to the differential
evolution algorithm and the resulting DOS are compared with experimental data [73, 74]. Each
NRIXS DOS was fit independently to obtain a force constant tensor that minimized the sum of
squared differences between the model and the experimental phonon DOS. The optimizations used
populations of 50 members which “evolve” until they converge, typically after a few hundred generations, on a set of force constants that gave the best fit to the experimental NRIXS DOS. These
optimizations were repeated several times to ensure convergence. For the optimization process, the
highest energy feature of each phonon DOS spectrum was fit to a Gaussian distribution, the distribution was then used to replace the high energy tails of each DOS used for fitting. This was done to
standardize the phonon cutoff energy across the data set and suppress fitting to higher energy noise
(which results from the data reduction and would not contribute meaningfully to the optimizations).
The BvK optimization was tried in several different configurations, each permitting a different
number of nearest-neighbor (NN) force constants to vary. The largest optimizations included atomic
interactions through 5NN (13 independent force constants), which is consistent with the number of
variables commonly used to fit neutron triple-axis dispersion data in this system [12, 80, 91]. The
fitting process was also performed with fewer variable nearest-neighbors forces (leaving more distant
force constants fixed to 300K tensorial force constants from literature [12]). The most restrictive
case limited the dynamics to interactions through 2NN (4 tensorial force constants).
To test our methodology, we calculated the phonon DOS using force constants from Minkiewicz,
et al.,[12] and convolved it with our NRIXS experimental resolution function. This DOS was opti-

43
mized using the genetic algorithm DOS fitting method, with the number of variable nearest-neighbor
force constants ranging from 2NN to 5NN. The optimizations that included interactions through the
4NN shell reproduced the known force constants accurately across several test cases, although only
slightly better than optimizations that included interactions through only the 2NN shells. It was
found that allowing variations through the 5NN shell noticeably increased the error in the tensorial
force constants found by the algorithm. Accordingly, the results presented here are from the optimizations where only the first 2NN force constants were allowed to vary, except in Fig. 6.7, where
values for fits through 4NN are included for comparison.

6.4

Results

6.4.1

Phonons

The 57 Fe phonon DOS spectra for α-Fe from 30K to just below the γ-Fe transition at 1185K are
shown in Fig. 6.1. All phonon modes shift to lower energies (soften) with increasing temperature,
although some soften noticeably more than others. The phonon DOS of bcc Fe has three features
corresponding to the Van Hove singularities of the longitudinal and two transverse acoustic phonon
branches. The mean energies of the three features were obtained by simultaneously fitting three
Lorentzian curves to the measured DOS spectra. The temperature dependence of these DOS features,
displayed in Fig. 6.2, show that the low transverse phonons soften more than other modes. While
phonon softening with increasing temperature is ubiquitous in most materials, the large preferential
softening of certain modes suggests strongly nonharmonic behavior of the lattice vibrations in α-Fe
at temperatures above 700K.

6.4.2

Quasiharmonic Model

The quasiharmonic model for predicting phonon frequencies employs the measured thermal expansion and a Grüneisen parameter to account for how phonon frequencies deviate from the harmonic
model at elevated temperatures. The quasiharmonic phonon frequencies, ωiqh (T ), are
ωiqh (T ) = ωi300K (1 − γ th

VT − V300K
),
V300K

(6.2)

where ωi300K is the measured value of the ith phonon frequency at 300K, γ th is the thermal Grüneisen
parameter and VT is the observed volume of the system at temperature T . This expression comes
ωi
V ∆ωi
from the definition of the microscopic mode Grüneisen parameter, γi = (− ∂∂ ln
ln V )T ≃ − ωi ∆V , where

a thermal Grüneisen parameter [24], γ th , is commonly used in the absence of detailed experimental
observations of the mode Grüneisen parameters, γi [36]. The thermal Grüneisen parameter can
be calculated from observed bulk material properties; in the following analysis Anderson’s value

1180K
1158K
1150K
1135K
1123K
1100K
1081K
1040K
1039K
1020K
1016K
1010K
1003K
1000K
990K
980K
DOS 1meV

970K
962K
952K
933K
900K
860K
830K
805K
771K
740K
700K
647K
592K
517K
510K
300K
299K
298K
292K
130K
50K
30K

Energy meV

Figure 6.1: The 57 Fe phonon DOS extracted from NRIXS measurements at various temperatures.
The spectra are normalized and offset for comparison.

Mean Energy meV

Longitudinal
30
High Transverse

Low Transverse

200

400

600

800

Temperature K

1000

1200

Figure 6.2: The measured NRIXS DOS were fit with three Lorentzian curves to find a characteristic
mode energy for each phonon branch. The softening of these mode features is compared with a
quasiharmonic model prediction from low temperature measurements (dashed lines). The Curie
temperature at 1044K is marked by a vertical line.
of 1.81 for α-Fe is used [24]. The quasiharmonic prediction from ambient temperature is shown by
dashed lines in Fig. 6.2 for each acoustic mode feature in the phonon DOS spectra. At temperatures
beyond 800K the mean phonon energies for each acoustic branch soften more rapidly than predicted
by the quasiharmonic model. At 1180K, just before the γ-Fe structural transition, the average
phonon energy has decreased by 19% of its low temperature value, more than twice the 8% softening
predicted by the quasiharmonic model.
A temperature-dependent thermal Grüneisen parameter was also calculated from the observed
bulk properties of α-Fe using the expression
γ th (T ) =

α(T )BT (T )ν(T )
CV (T )

(6.3)

where BT (T ) is the bulk modulus [28–30], α(T ) is the linear thermal expansion [22, 25], ν(T ) is

46
the crystalline volume per atom [22], and CV (T ) is calculated by integrating the low temperature
phonon spectra [12, 80] or a Debye model [19]. The temperature-dependent Grüneisen parameters
that can be created by various combinations of physical constants from the literature range from
1.7 to 2.2 over the temperatures of interest. However, including temperature-dependent parameters
in our analysis did not significantly alter the quantitative results provided by the quasiharmonic
model. Our further analysis therefore used the simpler approach with a constant thermal Grüneisen
parameter of 1.81 [24].

6.4.3

Vibrational Entropy

S SGTE

Entropy kB atom

Svib

Sqh
Sh

400

600

800

Temperature K

1000

1200

Figure 6.3: Our measured vibrational entropy (points) compared with the Scientific Group Thermodata Europe total entropy (thick line), the quasiharmonic vibrational entropy estimate (dashed
line), and the harmonic vibrational entropy (thin line).
The total entropy of iron is often split into component entropies attributed to different physical
phenomena,
S(T ) = Svib (T ) + Sel (T ) + Smag (T ),

(6.4)

where Svib (T ) is the vibrational entropy, Sel (T ) is the electronic entropy, and Smag (T ) is the magnetic entropy. While this division neglects the complex interplay of excitations in real materials at

47
elevated temperatures, it can still be useful for deconstructing thermodynamic models, especially
when experimental observations focus on a subset of the physical interactions. Accurate values of
the vibrational entropy, Svib (T ), can be obtained directly from the experimentally measured phonon
DOS as
Svib (T ) = 3kB

gT (E){(n + 1) ln(n + 1) − n ln(n)}dE,

(6.5)

where kB is the Boltzmann constant, gT (E) is the measured DOS at temperature, T , and n is a
Planck distribution evaluated at T , for a given energy, E. When experimental DOS spectra are
available for a given temperature, this expression provides accurate entropy values that include
both quasiharmonic effects and also nonharmonic effects (to first order) [35]. The total vibrational
entropies from NRIXS DOS spectra, Svib , are compared with the total entropy from the SGTE
database, S, in Fig. 6.3, together with the entropies of the harmonic and quasiharmonic models.
The total calculated entropies are also shown in Table 6.1. The total vibrational entropy, Svib (T ),
can be divided into component entropies as
Svib (T ) = Sh (T ) + ∆Sqh (T ) + ∆Snh (T ),

(6.6)

where Sh (T ) is harmonic vibrational entropy, ∆Sqh (T ) ≡ Sqh (T )−Sh (T ) is the purely quasiharmonic
contribution, and ∆Snh (T ) ≡ Svib (T ) − Sqh (T ), is the nonharmonic contribution. Figure 6.3 shows
that both the harmonic model and the quasiharmonic model significantly underestimate the vibrational entropy obtained from NRIXS measurements. Above 1000K, the nonharmonic vibrational
entropy, ∆Snh , is larger than the quasiharmonic contribution, ∆Sqh . At the highest temperatures
the nonharmonic vibrational entropy, ∆Snh , results in a 0.35kB /atom (5%) increase over the vibrational entropy provided by the quasiharmonic model.

6.4.4

Born-von Kármán Fits

Tensorial force constants were optimized to fit each NRIXS DOS spectrum using the genetic evolution
fitting algorithm, permitting the calculation of phonon dispersions at each temperature. Typical fits
and dispersions are shown in Fig. 6.4. The fitting procedure reproduces the DOS spectra quite well,
and also generates phonon dispersions consistent with previous triple-axis neutron measurements
[12, 83]. The calculated phonon dispersions corresponding to each NRIXS DOS measurement are
displayed along the high symmetry directions and the [ξξ2ξ] direction in Fig. 6.5. The phonon
dispersions change monotonically with temperature, exhibiting significant softening at the highest
temperatures, consistent with the phonon DOS spectra. These calculated dispersions show the trends
identified by prior neutron scattering measurements [83], and elastic moduli extracted by fitting the
low-q portions of the dispersion branches are in good agreement with measured elastic moduli [28, 30].

Table 6.1: The total vibrational entropy, Svib , calculated from the phonon density of states at each
temperature.
T (K)
30
50
130
298
299
300
300
510
523
595
647
687
740
773
805
842
870
900
923
934
962
963
981
1000
1003
1010
1023
1024
1033
1034
1039
1081
1090
1123
1135
1150
1158
1180

Svib (kB /atom)
0.020
0.114
1.122
3.115
3.145
3.160
3.135
4.747
4.781
5.206
5.489
5.702
5.916
6.080
6.240
6.400
6.492
6.627
6.713
6.777
6.888
6.877
6.962
7.007
7.034
7.066
7.104
7.108
7.146
7.232
7.208
7.383
7.432
7.585
7.557
7.612
7.666
7.751

Figure 6.4: NRIXS DOS fits and corresponding BvK dispersions are compared with literature.
Panels A and C compare the BvK fits (dashed lines) to the phonon DOS spectra at 300K and
1158K, respectively. Panel B compares our calculated 300K dispersions with the 295K neutron
triple-axis measurements of Minkiewicz, et al., (dots) [12]. Panel D compares our calculated 1158K
dispersions with the 1173K neutron triple-axis measurements of Neuhaus, et al. (dots) [83].
At elevated temperatures the optimized BvK fits began to segregate into two distinct solution basins.

Energy meV

30K

30
20
1180K

ΞΞ2Ξ 0.5ΞΞ0 0

Ξ00

1 ΞΞΞ 0.5 ΞΞΞ 0

Figure 6.5: Phonon dispersions resulting from force constant optimizations which permitted only
2NN force constants to vary (4 variables). Colors correspond to temperatures as labeled in Fig 6.1.
The basins had similar energy spectra, but quite different phonon dispersions and force constants.
The second basin of fits was characterized by H-point phonon energies that were significantly higher
than the N-point longitudinal phonons. They usually had lower qualities of fit than the primary
basin, especially when the optimizations included higher nearest-neighbor interactions. Fits from
this second basin were easily identified as erroneous from their discontinuous changes in tensorial
force constants with temperature, and their departure from measured dispersions [83]. Thus, BvK
fits from the second solution basin were excluded from further analysis.
The calculated dispersions demonstrate that some phonons soften significantly more than others.
Figure 6.6 shows the changes in the different dispersion branches with temperature relative to the 30K
frequencies. The measured phonon DOS spectra exhibited a 19% decrease in phonon frequencies at
the highest temperatures. Figure 6.6 makes clear that this average decrease is not evenly distributed
across all the phonon modes. Most of the nonharmonic softening of the phonon DOS appears
to originate in a few regions of the Brillouin zone. Most notably, anomalously large softening is
observed in the low transverse modes along the Γ-N direction, where all phonons soften by more than
twice the average decrease observed in the phonon DOS. The low transverse phonon branch T2 [ξξ0]
corresponds to the [110] phonon polarization direction and softens significantly with temperature,
consistent with the limited number of phonon dispersion studies on Fe at elevated temperature [83].
Large softening also occurs for the [ξξ2ξ] branch and between the H and P high symmetry points
at the 2/3 L [ξ,ξ,ξ] mode. Thermal softening seems to increase near the Γ point on several high
symmetry branches, but this increased softening at low-q may be an artifact of the fitting method.
Extracting a phonon DOS from the NRIXS spectra includes a removal of the elastic peak centered
at zero energy transfer, requiring that an extrapolation be used at energies below 4 meV (marked

Frequency Decrease %

100
80
60
40
20

Longitudinal

ΞΞ2Ξ

ΞΞ0

Ξ00

ΞΞΞ

High Transverse

ΞΞ2Ξ

ΞΞ0

Ξ00

ΞΞΞ

Low Transverse

ΞΞ2Ξ

ΞΞ0

Ξ00

ΞΞΞ

Figure 6.6: Percentage change in frequency with temperature of dispersions of Fig. 6.5. The
black lines denote the low energy dispersion regions below the dashed line in Fig. 6.5. The colors
correspond to temperatures as labeled in Fig. 6.1.
by the horizontal dotted line in Fig 6.5). Accordingly, Fig. 6.6 has black lines that delimit the low-q
region corresponding to the elastic peak extrapolation, below which our fits are less reliable.
By projecting each nearest-neighbor tensorial component along the NN bond direction, axial
and transverse force constants for bcc Fe were obtained for each nearest-neighbor pair as a function
of temperature, as shown in Fig. 6.7. With increasing temperature there is a large reduction in
the first-nearest-neighbor (1NN) and second-nearest-neighbor (2NN) longitudinal force constants.
The 2NN longitudinal force constant undergoes an especially strong softening. Above the Curie
temperature the magnitude of the 2NN longitudinal force constant is reduced to 60% (2NN fits) or
40% (4NN fits) of its low temperature value. The 1NN longitudinal force constant decreases by only
20% in the same temperature range.

Figure 6.7: Extracted force constants from BvK fits. The force constant results using up through
4NN (11 variables) are shown in gray squares and those resulting from using up through 2NN (4
variables) are shown in black points.
The transverse force constants are calculated as an average of the tensorial force constants projected onto two vectors orthogonal to the bond direction. The values of these force constants are
small and their trends are less reliable than for the larger longitudinal forces. The magnitude of
the 1NN average transverse force constant decreases rapidly at elevated temperatures and becomes
negative beyond 800K, indicating a weakness of the bcc structure to shear stress [33, 34, 85]. The
2NN transverse force constant appears to increase modestly with temperature, but this could be a
compensation for changes in longer-range interatomic forces that were not varied in the fitting procedure. The inclusion of additional variables produced the same general trends as those displayed
in Fig. 6.7, however, but with considerably more scatter. The large decrease in the 1NN and 2NN
longitudinal force constants with temperature occurred for every optimization configuration that
was tried.

ìì ì

L ΞΞ0

ìì ì
ì ìììììì
ììì
ìì
ìì
ìì
ìì

Energy meV

35
30
25
20

T1 ΞΞ0

òò ò

òò ò òò òò
òòòòòòò
òòòòòò
òò
òòòò
òò

àà à
23 L ΞΞΞ

ôô ô

T2 ΞΞ0

à àà
à àà
àààààà
ààà
àà
àà
àààà
ô ô
àà

ôô ôô
ôôô ô
ôôô
ôôô
ôô
ôô
ôôôô
ôô

15
10

200

400

600

800

Temperature K

1000

1200

Figure 6.8: Energies of specific phonon dispersion modes at the temperatures measured. The modes
are compared with their quasiharmonic estimates (grey dashed lines). The colors correspond to
temperatures as labeled in Fig. 6.1.

6.5

Discussion

6.5.1

Phonons and Born-von Kármán Model Dispersions

There is significant phonon softening in bcc Fe at elevated temperatures. A quasiharmonic model
accounts for some of the measured phonon softening, but underestimates the thermal trends. Direct
analysis of the phonon spectra shows that all the phonon DOS features exhibit softening beyond
the prediction of the quasiharmonic model, and this excess softening is most notable in the low
transverse modes. Both the departure from the quasiharmonic model at moderate temperatures
and the differential mode softening are indicative of strongly nonharmonic behavior. However, the

54
thermal broadening of the DOS features associated with phonon anharmonicity is small.
Phonon dispersions calculated from fitted force constants show that most of the nonharmonic
softening occurs in low energy phonon branches, while most of the higher energy longitudinal phonons
soften by an amount closer to that predicted by the quasiharmonic model. The temperature dependence of several phonon dispersion modes are shown in Fig. 6.8. The largest thermal softening
is found for the low transverse T2 [ξξ0] branch. Anomalous softening of these phonons has been
associated with dynamical precursors toward the fcc transition [33, 34, 84, 85]. A combination of a
displacement along the T2 [ξξ0] phonon mode, coupled with low-q shearing consistent with T2 [ξξ0]
and low-q shearing along the T2 [ξξ2ξ] branch is a possible path for the structural transformation
[33]. All these modes soften anomalously with temperature. Softening of the modes on the Brillioun
zone face between the H and P high symmetry points, most noticeably at 2/3 L[ξξξ], have been
associated with the structural instability of the bcc lattice under pressure towards the hexagonal ωphase. The dynamical precursors to the α-γ transition in Fe seem to originate with the softening of
the [ξ,ξ,0] branch, which is much larger than the softening of the [ξ,ξ,ξ] branch that is characteristic
of the structural ω-phase transition in the Group 4 bcc metals (Ti, Zr, Hf) [33, 34] and Cr [92] at
elevated pressures. A large decrease in 2NN longitudinal forces was reported in bcc chromium at
high temperatures, but Cr melts before the 2NN longitudinal force constant reaches the low values
seen here for Fe [92]. The soft phonons shown in Fig. 6.8 begin to deviate from quasiharmonic
behavior several hundred degrees below the magnetic transition, and continue to soften above the
Curie temperature. This anomalous phonon softening occurs in the same temperature range as
the rapid decrease in the magnetization of α-Fe [93]. Magnetic short range order has long been
suspected of being important for the phonon thermodynamics of the paramagnetic phase [94], and
recent DFT calculations that account for paramagnetic interactions have successfully predicted the
phonon dynamics at these temperatures [86].
Simulations were performed to vary the longitudinal and traverse force constants individually for
each nearest-neighbor pair. Adjusting the 1NN longitudinal force constant with the others fixed had
no effect on the T2 [ξξ0] branch. Decreasing the 2NN longitudinal force constant relative to the others
resulted in a rapidly softening T2 [ξξ0] branch that made the system dynamically unstable (imaginary
phonon frequencies) when this force constant dropped below zero. When the 1NN transverse force
constant decreased below -5 N/m, the phonons of the Γ-N branch also became dynamically unstable.
A strong decrease of the 2NN longitudinal force constant with temperature can produce the
significant non-harmonic softening observed in the Γ-N phonon branch [84]. Similarly, anomalous
behavior was seen in the 2NN magnetic exchange interaction parameters in a detailed study of
magnon-phonon coupling in bcc Fe [95]. It was found that including vibrational effects (including local volume and orientation) had a strong effect on the magnetic exchange interaction of the
nearest-neighbor pairs, most notably for the second nearest-neighbors. The abnormal 2NN magnetic

55
exchange behavior may be linked to the anomalous 2NN mode softening seen here.

6.5.2

Vibrational Entropy and Free Energy

Our measurements of the phonon DOS spectra over a range of temperatures permits the direct
assessment of vibrational entropy and vibrational free energy of bcc α-Fe. The vibrational entropy
from NRIXS measurements increases faster than predicted by the quasiharmonic model, ∆Sqh + Sh .
Any linear trends extracted from our force constants do not coincide with the volume normalized
values provided by Klotz, et al., for bcc Fe under pressure [91]. The purely volume-dependent
(quasiharmonic) effect from measurements at elevated pressure are quite different from our measured
nonharmonic effects at high temperature. Furthermore, linear fits to our tensorial force constants
are not capable of accurately reproducing the temperature dependence of the measured vibrational
entropy of Fe. There is a noticeable disagreement on either side of the magnetic transition; the force
constants have a nonlinear thermal trend through the Curie temperature.
The discrepancy between the quasiharmonic vibrational entropy, ∆Sqh + Sh , and the measured
vibrational entropy, Svib , is the nonharmonic entropy contribution, ∆Snh . Figure 6.9 shows that the

Entropy kB atom

0.6
ììììììì
ììììììììììììììì
ììììììì
ììììì
ìììì
0.5
ììì
ìì
0.4
ì Magnetization
ææ
ìì
Smag
0.3
ìæ æ
ìæ
ææ
0.2
æ DSnh
ææ
ææ
ææ æ
0.1
æææ
æ ææ
0.0 æ
400

600

800

Temperature K

1000

MM0

vibrational entropy of α-Fe has a significant nonharmonic contribution, ∆Snh . The nonharmonic

1200

Figure 6.9: The nonharmonic vibrational entropy, ∆Snh from measured phonon DOS spectra,
compared to the magnetization of bcc Fe [93], and the magnetic vibrational entropy, Smag , obtained
by subtracting Svib and Sel [77] from the SGTE total entropy, S [21].
vibrational entropy, ∆Snh , is compared to the magnetic entropy, Smag , and the magnetization also
shown in Fig. 6.9. Smag was calculated by subtracting our vibrational entropy, Svib , and also the
electronic contribution, Sel , as described by Jacobs, et al., [77] from the total entropy of the SGTE

56
database, S [21]. At temperatures just below the α-γ phase transition, ∆Snh changes the free energy
by about 35meV/atom.

Cp tot

Cp kB atom

T ¶Svib ¶T

CP qh

CP h

200

400

600

800

Temperature K

1000

1200

Figure 6.10: The heat capacity calculated from fits to our vibrational entropy measurements (solid
line), compared with the measured heat capacity from White, et al.[20] in points, heat capacity from
the harmonic model (thin line), and heat capacity from the quasiharmonic estimate (dashed line).
Measurements of heat capacity at constant pressure, CP , have provided some of the most important experimental information on the thermodynamics of α-Fe [19, 20, 31, 77, 95], and the
contributions from vibrational models are compared in Fig. 6.10. Because heat capacity is obtained
as a derivative quantity of the phonon entropy, we present heat capacity curves obtained from polynomial fits to our experimental ∆Snh results. As such, the curves in Fig. 6.10 should be reliable for
gradual trends, but possible features near the Curie temperature may be missing. Nevertheless, it
is clear that the ∆Snh of Fig. 6.9 makes a significant contribution to the heat capacity, larger than
the usual quasiharmonic contribution also shown in Fig. 6.10.
The nonharmonic vibrational entropy can be written as
∆Snh (T ) = Sppi (T ) + Sepi (T ) + Smpi (T ),

(6.7)

which includes the vibrational entropy from phonon-phonon interactions, Sppi (T ), vibrational entropy from electron-phonon interactions, Sepi (T ), and vibrational entropy from magnon-phonon
interactions, Smpi (T ). Experimental measurements of phonon DOS spectra cannot alone be used
to identify the individual terms Sppi , Sepi , or Smpi of ∆Snh . Nevertheless, the thermal trends are
suggestive. The Sppi contribution from phonon-phonon interactions (often called the “anharmonic”

57
contribution), arises from both the cubic and quartic components of the phonon potential [35].
With the cubic contribution comes a lifetime broadening from the imaginary part of the phonon
self energy. Even at the highest temperature, the lifetime broadening of the phonon DOS in α-Fe
is small compared to other systems [14–16, 96, 97]. For a damped harmonic oscillator, the shift
in phonon frequencies, ∆, associated with lifetime broadening is ∆ = −Γ2 /(2ε), where Γ is the
linewidth broadening and ε is the oscillator mode energy [36]. The broadening of the DOS features
was assessed by examining the widths of the Lorentzian fits used to create Fig. 6.2. From this
measured broadening, the classical anharmonic shift, ∆, is at least an order of magnitude smaller
than the observed high temperature shifts. Finally, phonon-phonon interactions from both cubic
and quartic perturbations increase linearly with temperature, and the ∆Snh in Fig. 6.9 does not
follow a linear trend. It appears that phonon-phonon interactions are not the main contribution to
∆Snh at high temperatures.
Electron-phonon coupling has been investigated by spin-polarized DFT calculations [89], and
effects were found to be modest. These calculations did find large differences in the electron-phonon
interactions for the majority and minority spin electrons, but did not consider disordered spin configurations. Second-nearest-neighbor magnetic exchange interactions were reported to be anomalously
sensitive to local atomic configurations [95], and we found that 2NN force constants decrease significantly at temperatures where the spin order was decreasing rapidly. The ∆Snh curve has a strikingly
similar shape to the magnetic entropy curve in Fig. 6.9. The magnon dispersions in iron have a
maximum energy approximately an order of magnitude higher than the phonon dispersions [98], but
perhaps 5% of the magnons are in the energy range of phonons in Fe. More processes involving two
magnons may affect the phonon self energies. A detailed analysis of phonon-magnon interactions is
required for further progress, but it seems plausible that Smpi is large.
Two phonon DOS spectra were acquired when the sample was in the fcc phase above 1185K, and
these are shown in Fig. 6.11. From these measurements, the change in vibrational entropy across the
α-γ phase transition at 1185K was found to be 0.05 kB /atom. This is notably smaller than previous
literature values of 0.091 and 0.14 kB /atom [83], which are similar to the SGTE recommended value
of 0.103 kB /atom [21], but the latter also includes magnetic and electronic contributions. The fcc γFe DOS spectra of Fig. 6.11 are compared with a 1428K γ-Fe DOS, calculated using force constants
from the literature [82] and convolved with our experimental NRIXS resolution for comparison.
The mean energies of the three DOS spectra are quite similar, but the two NRIXS DOS spectra
measured here are significantly broader than the DOS calculated from dispersion measurements.
High temperature phonon measurements of several fcc metals exhibited significant phonon lifetime
broading effects owing to phonon-phonon interactions, so these effects may make a more significant
contribution to the thermodynamics of fcc γ-Fe [14–16, 96, 97]. A more systematic study of the shape
of the phonon DOS in the fcc phase should help determine if there is a large lifetime broadening,

DOS 1meV

1428K
1265K

1210K

Energy meV

Figure 6.11: The phonon DOS at two temperatures in the fcc phase, compared with a spectrum
from neutron triple-axis measurements at 1428K [82].
but if so we may have underestimated the vibrational entropy for the fcc phase [99].

6.6

Conclusions

Nuclear resonant inelastic x-ray scattering was used to measure the phonon DOS of bcc α-Fe from
low temperature up through the α-γ transition. The vibrational entropy deviated significantly
from predictions of quasiharmonic theory by as much as 0.35kB /atom (or a free energy contribution
35 meV/atom) at 1150K. The nonharmonic contribution ∆Snh was distinctly nonlinear with temperature, and occurred without significant broadening of the phonon lineshape, unlike typical behavior
with phonon-phonon interactions. The temperature-dependence of ∆Snh followed the magnetic entropy, however, suggesting that the change of magnon-phonon interactions with temperature makes
a significant contribution to the nonharmonic phonon softening of α-Fe. The vibrational entropy of
the bcc-fcc Fe transition at 1185 K was found to be smaller than the assessed thermodynamic value.
A Born–von Kármán model was fit to the experimental phonon DOS spectra, and used to extract interatomic force constants. Full phonon dispersions were then calculated from the Born–von
Kármán force constants. These dispersions showed that the anomalous softening originates primarly
from low transverse modes along the Γ-N high symmetry direction, in agreement with single crystal triple-axis neutron studies. The anomalous softening originates with the large softening of the
2NN longitudinal force constant, which may be consistent with the atypical sensitivity of the 2NN

59
exchange interaction to local atomic configurations.

6.7

Supplemental Material

A Born-von Kármán (BvK) model was used to extract tensorial force constants from phonon DOS
measurements at 39 temperatures in the bcc phase. The BvK model was embedded in a genetic evolution algorithm to optimize tensorial force constants for each NRIXS measurement independently.
The optimizations were run in several configurations, where different numbers of nearest-neighbor
force constants were permitted to vary, with more distant neighbor force constants fixed to the
ambient temperature values [12]. The optimizations that varied the first four nearest neighbors (11
variables) most accurately reproduced known force constant values in our test case, but only slightly
better than the optimizations that included only two nearest neighbors (4 variables). Klotz, et al.,
also found that varying the force constants for only two nearest neighbors could reproduce phonon
trends in bcc Fe at elevated pressure [91]. Table 6.2 contains the force constants from the fits which
varied only the first and second nearest neighbors. Table 6.3 contains the force constants from fits
which permitted the first four nearest neighbors to vary.
The force constants from the BvK fits permit calculations of phonon frequencies at any k-point
in the Brillouin zone. Some phonon dispersions along crystallographic directions of high symmetry
were presented in figures in the manuscript. As shown in Fig. 6, the T2 [ξξ0] branch from Γ to N
has a particularly strong temperature dependence. Figure 8 shows the temperature dependence of
the energy for a particular phonon on this branch, the 1/2 T2 [110] phonon at the N-point. In Fig. 1
below we replot this curve by subtracting the quasiharmonic contribution to the softening, showing
only the nonharmonic contribution to the temperature dependence. The data in Fig. 1 are distinctly
nonlinear with temperature, unlike expectations from anharmonic phonon perturbation theory [36].
Figure 1 also shows the magnetic entropy, Smag (T ), obtained as described in the manuscript. The
similarities of the two curves in Fig. 1 suggest the importance of magnon-phonon interactions for
the thermal softening of the T2 [ξξ0] phonon branch.
The Lamb-Mössbauer factor (LMF) describes the ratio of recoil free to total nuclear resonant
absorption. The LMF extracted from the large NRIXS data set collected in this study is provided
in Fig. 6.13. The temperature dependent behavior of the 57 Fe LMF is in agreement with previous
experimental assessments.

Table 6.2: The fitted tensorial force constants for each temperature when only two nearest neighbors
(4 variables) are permitted to vary. The remaining force constants are fixed to room temperature
values from Minkiewicz et. at., [12] (3XX = 0.92, 3XY = 0.69, 3ZZ = -0.57, 4XX = -0.12, 4XZ =
0.0007, 4YY = 0.03,4YZ = 0.520, 5XX = -0.29, 5XY = 0.32). Force constants are reported in units
of N/m.
T (K)
30
50
130
298
299
300
300
510
523
595
647
687
740
773
805
842
870
900
923
934
962
963
981
1000
1003
1010
1023
1024
1033
1034
1039
1081
1090
1123
1135
1150
1158
1180

1XX
17.79
17.53
17.36
16.34
16.37
15.78
16.62
15.34
15.92
14.86
14.84
14.05
14.82
13.52
13.45
13.27
13.19
12.67
12.88
12.53
12.35
12.33
12.27
12.34
12.21
12.27
12.10
12.15
11.82
11.52
11.91
12.07
11.52
10.83
11.47
11.23
10.96
10.83

1XY
15.57
15.66
15.54
15.60
15.37
15.63
14.89
14.77
15.17
15.09
14.68
15.05
14.73
15.17
14.81
14.78
14.91
14.63
14.94
14.82
14.80
14.82
14.69
14.69
14.53
14.75
14.90
14.71
14.99
14.60
14.37
14.44
14.45
13.89
14.47
14.53
14.27
14.51

2XX
16.87
16.79
17.07
17.49
16.58
17.20
15.67
14.91
15.25
15.74
13.51
14.97
14.34
15.06
14.04
14.25
14.12
12.99
13.93
13.83
13.52
13.56
13.22
13.32
11.20
12.55
12.37
12.70
12.62
11.84
11.28
10.88
11.03
11.06
11.07
11.42
10.59
11.56

2YY
0.42
0.30
0.44
1.50
1.15
2.31
0.63
1.29
1.30
1.76
1.63
2.05
1.09
2.39
1.73
2.06
2.02
2.55
2.08
2.25
1.85
2.50
2.19
2.31
2.25
1.98
2.43
2.24
2.35
2.07
2.10
1.10
2.00
2.18
1.93
1.75
2.03
1.98

Table 6.3: The fitted tensorial force constants for each temperature when four nearest neighbors
(11 variables) are permitted to vary. The remaining force constants are fixed to room temperature
values from Minkiewicz et. at., [12] (5XX = -0.29, 5XY = 0.32). Force constants are reported in
units of N/m.
T (K)
30
50
130
298
299
300
300
510
523
595
647
687
740
773
805
842
870
900
923
934
962
963
981
1000
1003
1010
1023
1024
1033
1034
1039
1081
1090
1123
1135
1150
1158
1180

1XX
17.71
17.69
17.49
17.50
17.00
16.67
16.89
16.43
16.55
16.15
14.13
14.34
15.73
14.72
14.04
14.41
14.56
14.40
13.90
13.77
14.77
14.27
11.88
12.33
12.94
13.57
12.35
11.64
12.18
13.82
13.13
13.55
13.32
10.34
13.15
14.18
11.26
10.02

1XY
15.36
15.47
15.40
15.04
15.15
14.96
15.37
14.45
14.92
14.72
15.37
15.40
14.47
14.97
14.36
14.43
14.61
14.34
14.80
14.62
14.08
14.29
15.23
15.00
15.17
14.86
15.62
15.36
15.42
14.17
14.31
14.32
13.99
14.41
14.33
13.87
14.63
14.93

2XX
16.92
16.08
16.84
16.18
14.69
14.46
14.79
12.12
14.38
13.23
10.24
11.97
13.06
13.12
15.15
13.58
12.36
8.15
11.97
12.05
9.50
11.19
9.96
10.41
7.88
8.38
7.62
7.80
9.59
6.52
6.17
8.58
9.75
8.54
8.55
6.43
7.53
8.50

2YY
-1.14
-0.12
0.32
0.26
1.01
1.33
1.02
1.78
1.00
0.89
2.31
2.85
0.76
2.22
-0.15
1.15
0.80
4.38
1.60
2.07
-0.25
2.34
3.55
3.93
3.52
3.60
4.16
4.44
3.46
4.05
4.10
1.10
1.76
3.47
1.74
2.59
3.15
3.90

3XX
2.070
1.130
1.100
1.160
0.707
0.505
1.100
-0.003
0.965
0.752
0.390
-0.207
0.777
0.041
1.560
1.010
0.983
-1.270
0.804
0.399
1.050
0.088
-0.519
-0.857
-0.649
-0.976
-1.020
-1.060
-0.845
-1.660
-1.490
0.271
0.246
-0.479
0.370
-0.921
-0.464
-1.110

3XY
1.090
0.369
0.723
0.355
0.348
0.496
0.547
0.271
0.720
0.029
2.050
1.150
0.110
1.690
2.830
0.007
1.250
-0.239
0.260
0.207
0.069
-0.744
1.370
1.060
0.186
-0.529
0.773
1.320
-0.675
-1.910
-1.250
0.136
1.380
1.270
-0.554
-1.850
-0.543
1.300

3ZZ
-0.50
-0.37
-0.55
-0.44
0.20
-0.03
0.71
1.02
0.22
0.46
1.98
1.73
0.29
0.97
0.05
-0.05
1.18
3.13
0.79
1.02
1.57
1.10
2.27
1.56
3.23
3.02
3.85
3.36
2.43
3.39
3.39
1.74
1.52
2.30
2.12
2.78
2.49
1.70

4XX
-0.270
-0.451
-0.653
-0.799
-0.645
-0.474
-0.565
-0.752
-0.393
-0.373
0.572
0.459
-0.144
0.081
0.574
-0.086
-0.279
-0.761
-0.370
-0.206
-1.090
-0.660
1.350
0.757
0.469
0.116
0.741
1.060
0.863
-0.875
-0.268
-0.149
-0.401
1.320
-0.326
-1.230
0.352
1.440

4XZ
-0.137
-0.054
-0.029
-0.187
-0.029
0.074
-0.101
-0.342
-0.202
-0.022
0.113
0.003
-0.107
-0.564
-0.683
-0.133
-0.521
-0.609
-0.026
-0.263
-0.380
-0.445
0.180
0.024
-0.363
-0.271
-0.219
0.123
0.250
-0.495
-0.057
-0.423
-0.779
0.006
-0.349
-0.695
0.114
0.130

4YY
-0.045
0.146
0.151
0.087
0.164
0.216
-0.030
0.130
-0.155
-0.078
0.025
-0.139
-0.261
-0.201
-0.607
-0.362
-0.347
-0.054
-0.223
-0.327
-0.142
-0.216
-0.311
0.019
-0.337
-0.355
-0.208
-0.042
-0.094
-0.039
0.035
-0.446
-0.536
-0.397
-0.541
-0.148
-0.139
0.128

4YZ
0.650
0.505
0.320
0.534
0.418
0.367
0.569
0.170
0.287
0.653
0.392
0.377
0.560
-0.229
-0.206
0.556
0.130
0.013
0.537
0.427
0.683
0.537
0.624
0.484
1.130
0.913
0.778
0.758
2.000
1.080
1.180
0.495
-0.179
0.754
0.772
0.757
1.350
0.733

0.5

ôô 0.4
ôôô 0.3
ôô
0.2
ôô
0.1
ôô ôô

0 ôô ô

ôô ô

200

400

Smag HkB atomL

ÈDEnh È 12 T2 @110D HmeVL

0.0

600

800

Temperature HKL

1000 1200

Lamb Mossbauer Factor

Figure 6.12: The deviation of the 1/2 T2 [110] from the low temperature quasiharmonic approximation shown in Fig. 8. The nonharmonic deviation of this phonon mode tracks the magnetic entropy
of bcc Fe from Fig. 9.

0.8
0.6
0.4
0.2

200

400

600

800

Temperature K

1000

1200

Figure 6.13: Lamb-Mössbauer factors calculated from measured NRIXS spectra. The experimental
data from this study (presented in black) are compared with literature values in open circles [100],
open squares [101], open up triangles [102], and open down triangles [103].

Chapter 7

Cementite
Phonon densities of states (pDOS) of 57 Fe3 C were measured from low temperatures through the
Curie transition using nuclear resonant inelastic x-ray scattering. The cementite pDOS reveal that
the low energy acoustic phonons move to higher energies (stiffen) with temperature before the
magnetic transition. Such behavior is unusual because phonon frequencies typically move to lower
energies (soften) in conjunction with the quasi-harmonic approximation and finite-temperature thermal expansion. The unexpected stiffening observed in cementite suggested that the quasi-harmonic
model is not sufficient to describe the experimentally observed properties of cementite and that
previous claims of extreme elastic anisotropy in cementite may not hold at moderate temperatures.
Computational results were obtained using density functional theory (DFT) to complement our experimental results. These show that the unexpected stiffening observed experimentally in cementite
is reproduced by accounting for finite temperature phonon-phonon interactions. The anomalous
temperature response of phonons in Fe3 C is linked to the low energy acoustic phonon branches with
polarizations along the [010] direction. The effect was further localized to the motions of the FeII site
within the orthorhombic structure, which participates disproportionately in the anomalous phonon
stiffening.

7.1

Introduction

Cementite, Fe3 C, is the most common carbide observed in steels. Accordingly it has gathered much
interest for its significant role in thermomechanical processing and how its presence in different
microstructures affects observed material properties. Cementite has also recently been considered as
a candidate light-element-containing phase for the Earth’s inner core, motivating a variety of studies
at elevated pressures.
Cementite, Fe3 C, has an orthorhombic crystal structure that is ferromagnetically ordered at
ambient pressures below 460K and at ambient temperatures below 8GPa [104, 105]. The material properties of cementite change notably across the magnetic transition, which has complicated

64
extrapolations of ambient behavior to the high temperatures or high pressures of interest. Cementite is metastable at ambient conditions, and has been synthesized in its pure form only recently
[106]. Despite the longstanding metallurgical interest in Fe3 C, many of its physical properties are
still poorly constrained. This is due to the large variety of materials synthesis routes (alloyed steel
sub-components, mechanical milling, and high pressure synthesis) which produce varied cementite
products. The temperature-dependent heat capacity of cementite remains poorly constrained by a
wide range of experimental values reported on samples of varying origin and purity, and the high
temperature stability Fe3 C is still being debated [23, 107–110]. Only recently have detailed studies
on large, high-purity, cementite samples begun to illustrate the atypical physical properties of this
material.
The temperature-dependent thermal expansion of Fe3 C, evaluated by neutron diffraction illustrated a very anomalous magneto-volume behavior, including ranges of anomalously small and
anisotropic thermal expansion in the ferromagnetic phase [23]. Similarly, high pressure diffraction
studies have shown that the ferromagnetic phase is much less compressible than paramangetic or nonmagnetic orthorhombic high pressure phases. Many first principles computational studies have made
predictions of the material behavior of cementite under a range of physical conditions, but in the
absence of concrete experimental details these calculations are difficult to validate. Quasi-harmonic
density functional theory calculations of Fe3 C fail to capture the anomalous thermal expansion, but
provide a computational thermodynamic examination that falls well within the range of observed
thermodynamic behavior [110, 111].
Understanding the mechanical behavior of cementite is critical to a fundamental understanding
of its role in the hardening of steels, but comprehensive experimental studies are unavailable. First
principles studies have suggested that ferromagnetic cementite has an unusual level of elastic and
shear anisotropy [109, 112–114]. This has been extended to computational predictions of strainstiffening in cementite under extreme deformations in specific crystal orientations [113]. The limited
number of experimental measures of elastic behavior report that cementite does appear anisotropic,
but not to the extent predicted by density functional theory. Additionally, a low temperature study
on the ultrasonic sound velocities reported anomalous behavior at low temperatures [115]. The highly
anisotropic elastic behavior of cementite is likely linked to the complex anisotropic magneto-volume
behavior. A more concrete understanding of these relationships would have important consequences
for understanding the role of cementite in the mechanical behavior of steels.
The small thermal expansion of cementite has encouraged several groups to draw comparision
with invar materials that exhibit high-spin to low-spin magnetic transitions [116–120]. These are
largely based on the existence of magneto-volume anomalies, and the suggestion that the valence
electron per transition metal atom (e/a) ratio in Fe3 C is quite similar to the invar FeNi composition
if one assumes the carbon atoms donate electrons to the valence band. Under applied pressure

65
Fe3 C exhibits a number of phenomena associated with magnetic transitions at different pressures.
Synchrotron Mössbauer spectra have shown a loss of magnetic beats around 7GPa [105, 121] while
X-ray Magnetic Circular Dichronism (XMCD) indicates a change at 10GPa [118]; x-ray emission
spectroscopy (XES) shows a loss of magnetic character at 25GPa [122], and high pressure x-ray
diffraction (XRD) studies show changes in the volume compression at both 10GPa and 70GPa [123].
The combination of these pressure effects has been interpreted as a ferromagnetic to paramagnetic
phase transition around 8GPa, and also a paramagnetic to non-magnetic transition around 25GPa.
However, few observations have been made of electronic and magnetic behavior across the ambient
pressure transition at 460K.
Vibrational spectra of Fe3 C from low temperatures through the magnetic transition will help
elucidate the underlying physics in this material. They provide direct access to the vibrational
entropy of the material, and how it changes through the magnetic transition. Phonons are quite
sensitive to changes in bonding and magnetic configurations, so the thermal vibrational behavior may
improve the understanding of underlying magnetic phenomena. The elastic behavior of a material is
closely linked to its low-q vibrational modes, so the high temperature phonon behavior can provide
information on the thermal trends in elastic moduli.

7.2

Experimental

High pressure, high temperature synthesis was used to prepare 57 Fe3 C at stoichiometric composition
with 95% 57 Fe isotopic enrichment. The material was prepared by placing 95% isotopically enriched
57

Fe powder and graphite powder inside a MgO crucible. The crucible was placed inside a large

volume press where it was held at 2GPa and 1273K for 24 hours as in [124]. Pieces of the MgO
crucible were visually and magnetically separated from the 57 Fe3 C crystals which were ground into
powder under acetone.
Nuclear resonant inelastic x-ray scattering (NRIXS) measurements were performed on 57 Fe3 C
at high temperatures. NRIXS is a low background technique that provides direct access to the full
phonon density of states (DOS) [42, 46]. Measurements were performed at beamline 16ID-D of the
Advanced Photon Source at Argonne National Laboratory using a radiative heating furnace. The
powder sample was mounted on a thermocouple using Cotronics Cermabond 7020 aluminia-based
ceramic compound. The NRIXS measurements performed below room temperature employed a He
flow Be-dome cryostat, with the powder mounted in cryogenic vacuum grease. The temperatures
were accurate to within ±20K, where ambiguity comes from comparing the furnace thermocouple
measurements to in-situ nuclear forward scattering and the NRIXS-derived detailed balance temperature calculations following the procedures described in literature [45, 90].
An avalanche photodiode was positioned at approximately 90◦ from the incident beam to collect

66
re-radiated photons in a reflection geometry, beginning approximately 20 ns after the synchrotron
pulse. The incident photon energy was tuned to 14.413 keV using a high-resolution silicon crystal
monochromator to provide a narrow distribution of energies with a FWHM of 2.3 meV. The incident
photon energy was scanned through a range of ±170 meV, centered on the nuclear resonant energy.
The Phoenix reduction package was used to extract phonon DOS spectra from the collected spectra
[45].

7.3

Computational

Vibrational spectra were calculated using Density Functional Theory (DFT) across the range of
temperatures at which experimental data were taken. The vibrational spectra of Fe3 C were assessed
by two separate methods; quasi-harmonic scaled-volume calculations at 0K, and constant-volume
calculations at finite temperatures. This allowed us to separate and distinguish the effects of quasiharmonic thermal expansion and the effects of phonon-phonon interactions through finite temperature constant volume calculations. These vibrational spectra were, in turn, used to calculate the
elastic constants of cementite across a range of temperatures.
We used the Vienna Ab initio Simulation Package (VASP) [125–128] with a generalized gradient
approximation (GGA) exchange correlation functional as parameterized by Perdew, Burke, and
Ernzerhof [129] to calculate the vibrational spectra and elastic constants of cementite. For phonon
calculations, we modeled cementite as a 2 × 2 × 3 supercell of 192 atoms; elastic constant calculations
were performed on a 16-atom unit cell. With a 2 × 2 × 3 supercell, we achieved convergence with
respect to the system’s total energy and vibrational spectra using a Monkhorst-Pack [130] generated
k-point mesh of 3 × 3 × 3 and a plane wave energy cutoff of 800 eV; for unit cell calculations, we
employed a 11 × 9 × 13 k-points. All computational data reported herein relies on spin-polarized
calculations to model the system as ferromagnetic.
An Fe3 C supercell was fully relaxed to find the theoretical equilibrium lattice parameters. To
calculate vibrational spectra, we introduced a random set of displacements {u~i |1 ≤ i ≤ N } characteristic of temperature T to the relaxed supercell of N atoms. The displacements {u~i } were generated
as a linear combination of plane waves
~ui =

X i ck p
−2 ln ξ1 ei2πξ2 .
√k
mi

(7.1)

Here are the normal mode eigenvectors for the modes commensurate with the supercell, and mi the
mass of atom i. The amplitudes ck are derived from the same normal modes, with the displacement

67
length per mode given by Errea, et al. [131]
ck =

h̄
coth
2ωk

h̄ωk
2kB T

(7.2)

where ωk are the normal mode frequencies. The numbers ξ1 and ξ2 are uniformly distributed
random number between 0 and 1 to produce the standard Box-Muller transform for generating
normally distributed amplitudes. This distribution approximates the inclusion of zero-point motion,
and as such has nonzero displacements at 0K, and connects seamlessly to the classical limit at high
temperature.
We thereby generated sets of displacements for a random and representative selection of points
from the system’s phase space, as described by a Bose-Einstein distribution. Static DFT calculations
on these structures and post-processing by the TDEP method [71, 132, 133] yielded the interatomic
force constants and phonon DOS of Fe3 C, both at zero and finite temperatures.
First, we calculated vibrational spectra with the quasi-harmonic approximation (QHA) by accounting for temperature-induced volume changes to the system while using the 0K potential energy surface. To do this, the ground state lattice parameters were scaled by the experimentally
observed thermal expansion [23] for temperatures between 0K and 600K in 50K steps. Second,
the constant-volume finite-temperature (CVFT) calculations were performed using the equilibrium
lattice parameters, but the temperature in Eqn. 7.2 was increased to 200, 400, and 800K, adjusting
the amplitudes of displacements in the vibrational calculations to provide finite temperature effects.
We calculated changes to the elastic constants using low energy phonon group velocities from
the phonon dispersion relations and the Christoffel equations. This was performed for both the
scaled-volume QHA series, and the finite temperature constant volume calculations. These phononderived elastic constants were also compared to elastic constants calculated by the density functional
perturbation theory (DFPT) [134] as implemented in VASP, using unit cells scaled by the observed
experimental volume changes.

7.4

Results

7.4.1

Structure

Cementite takes the orthorhombic structure with the Pnma space group, containing 12 iron atoms
and 4 carbon atoms per unit cell as shown in Fig. 7.1. The experimental and computed equilibrium
lattice parameters from this study are compared with other literature values in Table 7.1. The DFT
calculated lattice parameters are in good agreement with the exception of the b lattice parameter,
which appears slightly lower than previously reported results. The experimental lattice parameters
were confirmed by powder x-ray diffraction (XRD) using a laboratory Cu Kα source, and Reitveld

Figure 7.1: (Color online) The structure of cementite. Carbon atoms are black, FeI atoms are gold,
while FeII atoms are red.
Table 7.1: Unit cell parameters of cementite from various experimental and computational studies.
Source
this study
Jiang [112]
Dick [111]
Nikolussi [114]
Haeglund [135]
this study
Wood [23]
Gao [105]
Wood [23]

5.0429
5.04
5.035
5.036
5.089
5.086
5.081
5.0814
5.082

6.7028
6.72
6.716
6.724
6.743
6.754
6.753
6.751
6.733

4.4816
4.48
4.480
4.480
4.523
4.520
4.515
4.516
4.512

Functional
PAW/PBE
PAW/PBE
PAW/PBE
PAW/PBE
LMTO
Experimental (300K)
Experimental (300K)
Experimental (300K)
Experimental (4K)

refinement, giving interplanar spacings very similar to synchrotron XRD measurements on other
samples of Fe3 C generated by the same high pressure synthesis methods [136].
In cementite, the Fe atoms occupy two distinct sites on the lattice, with 8 FeII atoms occupying
the general site (8d) shown in red, and four FeI occupying the special site (4c) shown in gold. The
carbon atoms also occupy the 4c site, with prismatically coordinated positions between the Fe atoms.
The lattice positions of each of the unique crystal sites are shown in Table 7.2. The calculated and
experimentally-determined lattice sites are also in close agreement with previous results.
Both Fe sites in the orthogonal cementite structure have carbon first nearest neighbors, although
the Fe-C bond is somewhat shorter for FeI atoms [137]. The FeII sites have 11 Fe neighbors in the
second nearest neighbor shell, while the FeI site have 12. On average the Fe neighbors of the FeII
site are closer than the FeI site, with the closest Fe-Fe distance existing between FeII sites aligned

69
Table 7.2: Unit cell parameters of the unique crystallographic sites in cementite from various experimental and computational studies.
Source
this study
Jiang [112]
this study
Wood [23]

FeI (x,y,z)
(0.036,0.025,0.837)
(0.036,0.025,0.837)
(0.035,0.250,0.838)
(0.034,0.250,0.841)

FeII (x,y,z)
(0.176,0.068,0.332)
(0.176,0.068,0.332)
(0.185,0.059,0.334)
(0.184,0.057,0.333)

C (x,y,z)
(0.877,0.250,0.438)
(0.876,0.250,0.438)
(0.898,0.250,0.447)
(0.894,0.250,0.450)

Method
PAW/PBE
PAW/PBE
300K Exp.
300K Exp.

along the b-axis.

7.4.2

Phonons

Nuclear Resonant Inelastic X-ray Scattering (NRIXS) spectra were measured at 17 temperatures
from 14K through the magnetic transition at 460K, up to the reported limit of material stability
at 600K [23]. The NRIXS spectra were reduced to provide the 57 Fe partial phonon DOS curves in
Fig. 7.2. The phonon partial DOS show only small changes through the temperature range. The
most apparent change is softening and broadening of the feature near 35 meV. The lowest temperature phonon partial DOS is compared with the 0K DFT calculated partial DOS in Fig. 7.3 at the
calculated equilibrium unit cell volume reported in Table 7.1. The lowest temperature measurement
is in excellent agreement with our DFT calculations when the experimental resolution is taken into
account and the energy axis is scaled by 5.5% to align the mean energies of the experimental and
calculated phonon spectra. Our phonon dispersions calculated using DFT are consistent with dispersions from previous computational studies of cementite [112]. The phonon partial DOS curves
from DFT show that displacements of the Fe atoms dominate the low-energy phonons, while carbon
atom motions dominate the higher-energy phonon branches. This segregation of the phonon DOS
can be related to the atomic mass mismatch of the two species. While most phonon modes include
motions of both Fe and C atoms, the carbon motions dominate the high frequency phonons, and
the iron motions dominate the lower energy phonons. All Fe pDOS curves are normalized to 67
meV, because experimental spectra were too noisy to make out the subtle DOS features beyond this
energy. While the DFT calculations find some intensity in the Fe pDOS above 70meV, this intensity
makes only a small contribution to integrated quantities.
The average phonon energies from the Fe partial DOS (pDOS) are plotted in Fig. 7.4. The experimental mean energies are compared with their QH DFT counterparts, and also a quasi-harmonic
model with a thermal Grüneisen parameter calculated from the 300K bulk properties of Fe3 C using
the expression
γ th =

αKT ν
= 2.24,
CV

(7.3)

where KT is the bulk modulus [122], α(T ) is the volume thermal expansion [23], ν(T ) is the volume

600K
550K
500K
480K
463K
DOS H1meVL

444K
423K
410K
375K
365K
330K
297K
295K
294K
232K
105K
14K

Energy HmeVL

Figure 7.2: (Color online) The 57 Fe partial phonon DOS extracted from NRIXS measurements
at various temperatures. The spectra are normalized and offset for comparison. Experimentally
determined errors are shown as partially shaded regions along each DOS.

14K 57 Fe pDOS
DFT Fe pDOS
DOS H1meVL

DFT C pDOS

Energy HmeVL

Figure 7.3: (Color online) The 57 Fe phonon DOS from NRIXS at 14K, compared with partial
density of states from DFT. The experimental DOS is shown in purple, where the shaded region is
representative of the experimental uncertainties in the phonon pDOS measurement. The DFT Fe
pDOS is shown in black, along with the C pDOS (gray dashed).
per atom [23], and CV (T ) is the heat capacity at constant volume calculated by integrating the
DFT total phonon DOS. The experimental average phonon energies are nearly constant (within
experimental scatter) to the Curie temperature at 460K. At higher temperatures, where cementite
becomes paramagnetic, the phonon energies begin to decrease (or soften) with temperature. While
the observed thermal expansion in the ferromagnetic temperature region is quite small, the QH estimate from DFT predicts a 0.2 meV decrease in mean phonon energy up to the magnetic transition,
while the experimental data are constant, or perhaps undergo the opposite trend. The near constant
behavior of the mean phonon energy with temperature suggests there must be nonharmonic phonon
behavior approximately equal and opposite to the thermal expansion driven phonon softening.

7.4.3

Vibrational Entropy

The total entropy of a material is often separated into component entropies attributed to different
physical phenomena,
S(T ) = Svib (T ) + Sel (T ) + Smag (T ),

(7.4)

where Svib (T ) is the vibrational entropy, Sel (T ) is the electronic entropy, and Smag (T ) is the magnetic entropy. Vibrational entropy makes the largest contribution at finite temperatures, so small

26.4
NRIXS

HmeVL

26.2
26.0
25.8

QH DFT

25.6

ΓT QH Model

25.4
25.2

100

200

300

400

Temperature HKL

500

600

DSNRIXS-Harmonic HkB atomL

Figure 7.4: (Color online) The mean energy of the 57 Fe partial phonon DOS extracted from NRIXS
measurements at various temperatures. The thick dashed line is the mean energy calculated from
a Grüneisen parameter quasi-harmonic model. The thin dashed line is the Fe partial DOS mean
energy from DFT calculations.

0.08
¶S

0.06
0.04
0.02

ΓT QH Model

¶V

*VHTL

QH DFT

0.00 ì

-0.02 NRIXS
CVFT DFT
-0.04
100 200 300

TC
400

Temperature HKL

500

600

Figure 7.5: (Color online) The change in vibrational entropy with temperature as it departs from
the harmonic model.
changes in vibrational entropy can have a large impact on material stability. Accurate values of the
vibrational entropy, Svib (T ), can be obtained directly from the phonon DOS as
Svib (T ) = 3kB

gT (E){(n + 1) ln(n + 1) − n ln(n)}dE,

(7.5)

73
where kB is the Boltzmann constant, gT (E) is the DOS at temperature T , and n is a Planck
distribution evaluated at T , for a given energy E. When experimental DOS curves are available for a
given temperature, this expression provides accurate entropy values that include both quasiharmonic
effects and also nonharmonic effects (to first order) [35].
The vibrational entropy, Svib (T ), can be further divided into component entropies as
Svib (T ) = Sh (T ) + ∆Sqh (T ) + ∆Sah (T ),

(7.6)

where Sh (T ) is harmonic vibrational entropy, ∆Sqh (T ) ≡ Sqh (T )−Sh (T ) is the purely quasiharmonic
contribution, and ∆Sah (T ) ≡ Svib (T ) − Sqh (T ), is the anharmonic contribution. For the 57 Fe partial
DOS of cementite from NRIXS, we obtain a partial vibrational entropy for the Fe atoms. The Fe
partial entropy calculated from NRIXS pDOS is compared to quasi-harmonic models for ∆Sqh (T ) in
Fig 7.5. For clarity, the harmonic contribution to the vibrational entropy, Sh (T ), has been subtracted
so the non-harmonic behavior can be more closely assessed. The quasi-harmonic contribution,
∆Sqh (T ) adjusts the vibrational entropy to account for the thermal expansion of a material at finite
temperatures. While the thermal expansion of cementite is small in the magnetic phase, it still has
an appreciable effect on the Fe partial vibrational entropy. Several quasi-harmonic contributions
are compared in Fig. 7.5. The vibrational entropy derived from QH DFT and the Grüeneisen
parameter model are compared with a quasi-harmonic entropy calculated from the entropic volume
dependence from high pressure NRIXS measurements [105]. All the methods examined diverge
notably from the measured Fe partial vibrational entropy at temperatures just below TC , by more
the 0.03kB /atom. The measured partial vibrational entropy is noticeably lower than predicted by
quasi-harmonic models, suggesting an anharmonic contribution that is approximately equal and
opposite to the quasi-harmonic contribution. The constant volume finite temperature (CVFT)
calculations do agree with the measured Fe partial vibrational entropy. This suggests that the
elevated temperature deviations from the quasi-harmonic estimate below the magnetic transition
originate from phonon-phonon interactions that are well reproduced by the CVFT calculations.
Total entropy estimates are provided in Table 7.3, which are calculated by adding a harmonic
C partial vibrational contribution calculated from the 0K DFT C pDOS. The errors include both
the calculated quasiharmonic and anharmonic C phonon shifts, which contribute approximately
0.01 kB per atom at TC and oppose each other in sign. The deviation from the ab initio calculated
vibrational entropies are quite small, provided that the calculated phonon spectra are scaled to match
the observed NRIXS mean energies. However, they do differ from the DFT calculated vibrational
entropy of previous studies [111], by about 0.17 kB /atom.

74
Table 7.3: The total vibrational entropy calculated from NRIXS Fe pDOS and DFT C pDOS.
Temperature (K)
14
105
232
294
295
297
330
365
375
410
423
444
463
480
500
550
600

7.4.4

Total Svib kB /atom
0.003 ± 0.000
0.744 ± 0.001
2.205 ± 0.002
2.786 ± 0.004
2.783 ± 0.004
2.794 ± 0.004
3.088 ± 0.005
3.351 ± 0.007
3.429 ± 0.008
3.658 ± 0.009
3.745 ± 0.010
3.880 ± 0.011
4.007 ± 0.012
4.113 ± 0.013
4.238 ± 0.014
4.520 ± 0.016
4.792 ± 0.019

Low Energy Phonon Modes

The Debye sound velocities were extracted from the low energy region of the phonon partial DOS
curves following the methodology of Hu, et al., [138] using energies below 11 meV. The room temperature sound velocities are in reasonable agreement with the powder sample measurements of Gao,
et al., at 300K [139] and are presented in Fig. 7.6. We see a noticeable increase in the Debye sound

VD HkmsL

3.05

3.00

2.95

ææ

2.90

2.85

100

200

300

400

Temperature HKL

500

600

Figure 7.6: The Debye sound velocities calculated from the NRIXS pDOS. A reference calculated
from the NRIXS pDOS of the measured powder sample in Gao et. al.,[136] is shown as a gray
diamond.
velocity in the ferromagnetic phase with increasing temperature. The increase in the sound velocity
extracted from our phonon pDOS measurements reverses itself once the magnetism is lost, and the
Debye velocity decreases with temperature in the paramagnetic phase.

Ferromagnetic Fe3 C

463K
DOS H1meVL

14K

8 10 12

Energy HmeVL

Paramagnetic Fe3 C

DOS H1meVL

600K
463K

8 10 12

Energy HmeVL

Full Temperature Range

600K
DOS H1meVL

14K

Energy HmeVL

8 10 12

Figure 7.7: (Color online) The 57 Fe phonon partial DOS are compared at various temperatures,
shaded areas show uncertainty in measured spectra. The calculated phonon difference spectra are
shown below the phonon DOS and the inset shows the low energy region. Panel A highlights changes
over the ferromagnetic temperature region. Panel B shows changes that occur after the magnetic
transition. Panel C compares the highest and lowest phonon spectra.
This is consistent with early ultrasonic measurements that showed a small decrease in the Debye
sound velocity with temperature below 300K [115], but our decrease in Debye velocity below 300K

76
is nearly an order of magnitude larger. This may be a result of the material quality of the ultrasonic
study, or perhaps the anomalous stiffening we observe at finite energy has a smaller effect on phonons
in the long wavelength limit sampled by ultrasonic measurements. Dodd, et. al., found that the
decrease in the sound velocity on decreasing temperature was due to changes in the shear wave
velocities; the longitudinal sound velocity showed a normal increase with decreasing temperature.
Prior NRIXS data on 57 Fe3 C at elevated pressures show that the low-energy phonons stiffened
normally with pressure [105, 115, 139]. Nevertheless, the shear velocity increased much more rapidly
with pressure in the ferromagnetic phase, compared to the paramagnetic phase.
Closer examination of the low energy phonon spectra in Fig. 7.7 shows that phonons near 9meV
actually stiffen (or increase in energy) as the material is heated. This stiffening of the low energy
phonon modes accounts for the increase in Debye sound velocities observed with increasing temperature. This stiffening of phonon energies with temperature is opposite to what might be expected
in a material with a net volume expansion over the same temperature region. This anomalous stiffening of the low energy Fe phonon modes persists until the magnetic transition is reached, then
the low-energy phonons soften with temperature. This anomalous low-energy phonon stiffening was
not reproduced by QH DFT calculations performed at volumes scaled to match the observed thermal expansion, explaining the discrepancy between the QH models and the NRIXS observed values
for mean energies and vibrational entropies at finite temperatures in Fig. 7.4. The volume-scaled
DFT phonon dispersions soften monotonically at all energies, indicating that the anomalous lowenergy phonon stiffening in the NRIXS pDOS cannot be solely attributed to the anisotropic thermal
expansion.
The CVFT calculations captured the anomalous low-energy stiffening behavior at finite temperatures, as shown in Fig 7.8. The CVFT Fe pDOS shows qualitative agreement with the change in
low energy phonon modes observed in the NRIXS pDOS at 400K. These constant volume calculations show the low energy phonons stiffening continuously with increasing temperature. The C
pDOS makes a very small contribution to the total DOS at low energies, but it also shows behavior
consistent with the low energy stiffening seen in the Fe pDOS. The low temperature phonon dispersions show which phonon branches comprise the 6-12meV shoulder shown in the Fe pDOS. This
energy range encompasses low-energy transverse acoustic phonon modes, and also the lowest energy
optical modes near the Γ point. At elevated temperatures, these modes move to higher energies as
shown in Fig 7.8, resulting in the loss of the low energy shoulder in the Fe partial phonon DOS.
The longitudinal acoustic phonon modes have higher energies, above the low energy shoulder, and
remain fairly constant between 0 and 400K in the CVFT phonon dispersions. In addition to recreating our thermal trends in the Fe pDOS, these calculated phonon dispersions are consistent with
the ultrasonic measurements, which show different temperature behaviors of longitudinal and shear
sound waves, relating to the longitudinal and transverse acoustic branches, respectively [115].

DOS H1meVL

CVFT DFT

400K
0K

@000D

Energy HmeVL

@100D @110D

@010D @000D

8 10 12

@001D

@101D

Energy HmeVL

60
50
40
30
20
10

Figure 7.8: (Color online) The Fe partial phonon DOS from CVFT calculations are convolved with
an experimental resolution function. The calculated phonon difference spectra are shown below the
phonon DOS, and the inset highlights the changes in the low energy region over the ferromagnetic
temperature region. Panel B shows the corresponding phonon dispersions at the high symmetry
points at 0K (black) and 400K (red).
The low-energy transverse acoustic phonon branch stiffens along many of the high-symmetry
directions of the Brillouin zone except the [010] direction between Γ and Y. The transverse acoustic
phonons have their lowest energies along the Γ - Z direction, but this branch stiffens noticeably

78
with temperature, even at low-q. This stiffening explains the observed change in the Debye sound
velocity with increasing temperature. This stiffening is particularly large at the X and Z points,
where the lowest energy transverse acoustic phonons stiffen by 26% and 40% respectively. The Y
point shows almost no stiffening with temperature in the CVFT calculations. Characteristic phonon
displacements for the low-energy Z and X point modes are shown in Fig. 7.9. Both sets of modes are
characterized by large displacements of Fe atoms along the b-axis. In contrast characteristic atom
displacements for the low-energy Y point modes are shown in Fig 7.10. These Y-point modes which
do no stiffen significantly with temperature, characterized by large displacements of the Fe atoms
perpendicular to the b-axis.

Several other k-points in high symmetry directions were examined,

Figure 7.9: (Color online) Phonon mode illustrations from the X-point (top) and the Z-point
(bottom), which undergo anomalous thermal stiffening at moderate temperatures.
and the thermal stiffening of these modes were correlated with phonon polarization along the [010]
direction.

Figure 7.10: (Color online) Phonon mode illustrations from the Y-point. These phonon modes do
not undergo anomalous thermal stiffening.
The two lowest-energy optical modes also show large stiffening effects near the Γ point, with a
16% phonon energy increase for the lowest point, and a 18% energy increase for the highest point.
These phonon modes are shown in Fig. 7.11. The lowest energy optical phonon mode at Γ shows
somewhat different behavior than the other stiffening modes, with significant phonon displacements
of the Fe atoms in the c-axis direction, but with components along the a and b directions. The
highest energy optical phonon is characterized by large displacements of all atoms along the b-axis,
and this mode exhibited large thermal stiffening.

Figure 7.11: (Color online) Phonon mode illustrations for the lowest energy Γ-point phonons. These
phono modes undergo significant thermal stiffening.

80
The Fe partial DOS can be further divided onto the unique Fe sites in the orthorhombic cementite
lattice as shown in Fig. 7.12. The FeII pDOS contains much greater intensity in the low-energy

DOS H1meVL

FeI pDOS

400K

8 10 12

Energy HmeVL

DOS H1meVL

FeII pDOS

400K
0K

Energy HmeVL

8 10 12

Figure 7.12: (Color online) The Fe phonon partial DOS from CVFT calculations, convolved with an
experimental resolution function, projected onto the two distinct iron lattice sites. The calculated
phonon difference spectra are shown below the phonon DOS. The inset highlights the changes in
the low energy region over the ferromagnetic temperature region. Panel A shows the phonon DOS
of the FeI sites (4 atoms/unit cell). Panel B shows the phonon DOS of the FeII sites (8 atoms/unit
cell).

81
phonon shoulder that exhibits anomalous stiffening with temperature. The FeI sites participate in
these phonon modes to a lesser degree, and exhibit less stiffening with temperature. This suggests
that the anomalous stiffening seen in the Fe partial DOS at finite temperatures is directly attributable
to the FeII lattice site. Additionally the FeII pDOS contains the DOS feature at 35meV, which softens
significantly with temperature in agreement with the observed NRIXS pDOS behavior. Figures 7.9
and 7.11 show that the FeII atoms have the largest displacements, often along the b-axis, or in other
patterns that strongly distort the alignment of the layers of trigonal prisms stacked along the b-axis.

7.4.5

Elastic Constants

Cementite has been identified as having a particularly high elastic anisotropy, which can be directly
linked to the low-q behavior of the acoustic phonons [112]. Computational DFT studies of the
elastic moduli of cementite have consistently reported an anomalously low C44 value near 20GPa
[112, 114]. The nine single crystal elastic constants of Fe3 C were calculated with DFPT as described
in Section 7.3. Results are compared with values from the literature in Table 7.4, which includes
elastic constants extracted from the DFT low-q phonon group velocities.
Table 7.4: Single crystal elastic moduli extracted from DFPT and from calculated phonon dispersions.
Source
this study (phonons)
this study (dfpt)
Jiang [112] (phonon)
Jiang [112] (energy)
Nikolussi [114] (energy)

C11
315
380
384
388
385

C22
315
345
325
345
341

C33
292
302
283
322
316

C12
...
160
...
156
157

C23
...
152
...
162
162

C13
...
160
...
164
167

C44
26
20
26
15
13

C55
137
132
134
134
131

C66
137
135
125
134
131

Our calculated elastic moduli at 0K compare well with similar approaches in the literature, with a
C44 elastic constant that is quite low. Our temperature-dependent observations of the Fe pDOS suggest the low temperature acoustic modes stiffen significantly with temperature, up to the magnetic
transition. This implies that the elastic constants should change significantly with temperature.
Elastic constants that are strongly temperature dependent would explain the discrepancies between
the highly anisotropic elastic constants calculated at 0K and less anisotropic room temperature
experimental observations of elastic behavior [112, 114, 140–142]. Since our CVFT phonon calculations show such excellent agreement with our experimental results, we will interpret the changes
in those spectra, and their implications for finite temperature elastic constants in Fe3 C. The values
presented here are the first calculations of finite temperature elastic constants for cementite, which
depart noticeably from the low temperature calculated values.
The single crystal elastic constants can be related to the low energy phonon group velocities
through the Christoffel equations with the appropriate symmetry considerations. For orthorhombic

82
cementite the relationships between the acoustic phonon branches and the elastic constants are

C11
[100]
VT [010] =

C22
[010]
VT [100] =

C33
[001]
VT [100] =

[100]
VL[100] =

[010]
VL[010] =

[001]
VL[001] =

C66
[100]
VT [001] =

C66
[010]
VT [001] =

C55
[001]
VT [010] =

C55

(7.7)

C44

(7.8)

C44

(7.9)

[hkl]

where V[uvw] is the phonon group velocity of the longitudinal (L) or transverse (T) phonon in the
[hkl] direction with the [uvw] polarization, and ρ is the temperature-dependent density of cementite
[23]. The C44 elastic constant is directly related to the group velocities of low-q phonons in the
lowest transverse phonon branch along the Γ-Z and Γ-Y directions. The high-temperature phonon
dispersions in Fig. 7.8 show a large thermal stiffening of the low-q branch in the Γ-Z direction.
This suggests that the anomalously low elastic constants calculated for C44 should increase with
temperature. This increase from 0K to near ambient temperatures may be the reason that elastic
moduli measurements at 300K have not reproduced the extreme shear anisotropy reported from 0K
DFT calculations [113].
The temperature dependence of the elastic constants calculated from the CVFT dispersions are
shown in Fig 7.13. The anomalously low C44 undergoes a 55% increase between 0 and 400K when
the thermal expansion of Fe3 C is included as a temperature dependent density. The C55 and C66
decrease slightly as a result of modest softening with temperature.
If we calculate the Zener shear anisotropy ratios using these temperature-dependent terms we find
that all three shear ratios become more isotropic with increasing temperature. We were unable to
reliably extract temperature dependent information on the remaining off-symmetry elastic constants
(C12 ,C13 ,C23 ) from the phonon spectra, so these were assumed constant in the calculation of shear
anisotropy. However, these three off-axis elastic moduli would need to change by between 1825% with temperature (nearly double the observed changes of all elastic moduli except C44 ), to
maintain the low-temperature shear anisotropy ratios. The low-energy phonon results suggest that
cementite may exhibit significantly less anisotropy at 400K than 0K calculations suggest. This is in
agreement with recent experimental studies including nano indentation of single crystals of Fe3 C that
suggest the elastic response is less dependent on crystallographic orientation [140]. Experiments that
calculated the Young’s modulus from bending small oriented single crystals of cementite extracted
from an Fe matrix also showed a slightly less anisotropic response than predicted by first principles
at 0K [141].

350

300 õ

Cii HGPaL

250
á C11
ó C22
õ C33

200
150

200

à C44
ç C55
í C66

100
50

100

300

DCii H%L

400

40
20
0 ç

100

200

300

Temperature HKL

400

Figure 7.13: The thermal trends of elastic moduli extracted from the phono dispersions at low q,
calculated by CVFT.

7.4.6

Electronic DOS

The CVFT calculated electron density of states (eDOS) for cementite is shown in Fig 7.14, with the
electronic contributions resolved into their orbital characters. The eDOS projected onto the carbon
site lacks d-electron character, as expected. The carbon site electron DOS shows a large deeply
bound concentration of p-electrons below the Fermi surface, similar to findings from Khmelevskyi,
et. al. [112, 119]. This invalidates the simple suggestion that the thermal expansion anomaly in Fe3 C
might be explained by arguments regarding how carbon atoms donate their conduction electrons,
increasing the Fe valence to levels comparable with the FeNi invar composition. The FeII majority
site eDOS has a larger concentration of electrons at the Fermi level than the minority FeI site. The
FeII majority site has a large d-electron feature at energies just below the Fermi level. The high
temperature calculated density of states at 400K shows no major change in features aside from a
general thermal smoothing. At 400K, the FeII site has a larger increase in electrons at the Fermi
level than the FeI site. The calculated electronic DOS at the Fe sites both undergo 7 % increases
in their d-electron occupations at the Fermi level, but the FeII site also shows a 4% increase in
p-electron levels. While the calculated magnetic moments on both Fe sites remain nearly constant

84
in high temperature CVFT calculations, their differing electronic character is suggestive of selective
bond interaction behavior.

Electron DOS

100
50

FeI H4cL

s electrons
p electrons
d electrons

-50
-100
-10

Electron DOS

100
50

-5

Energy HRyL

FeII H8dL

s electrons
p electrons
d electrons

-50
-100
-10

Electron DOS

100

-5

Energy HRyL

C H4cL

s electrons
p electrons

50
-50
-100
-10

-5

Energy HRyL

Figure 7.14: (Color online) The calculated electronic DOS, resolved into orbital contributions at
the three distinct lattice sites.

7.5

Conclusions

The NRIXS Fe pDOS have been measured from low temperature through the magnetic Curie transition at 460K and into the paramagnetic phase. The phonon spectra show an unusual thermal
stiffening of low-energy phonon frequencies that is well reproduced by CVFT calculations. The
low-energy phonon stiffening counteracts the softening expected from thermal expansion, keeping
the mean energy nearly constant with temperature and lowering the vibrational entropy by 0.03
kB /atom compared to estimates based on quasi-harmonic models. The low energy stiffening is particularly notable in the Γ - Z direction in the phonon dispersions, where the whole branch shifts to
higher energies with increasing temperature in the ferromagnetic phase. The stiffening of several
low-energy acoustic branches has important implications for the single crystal elastic moduli of Fe3 C,
which become more isotropic with temperature in the ferromagnetic phase. This low-energy phonon
stiffening is evident in the dispersions throughout the Brillouin zone, where the low-energy transverse
acoustic modes have polarizations along the [010] direction. Motions of Fe atoms at the FeII sites
are strongly correlated with this anomalous stiffening behavior, as shown by the site projected FeII
pDOS. The FeII sites have the closest Fe-Fe distance in the orthorhombic unit cell along the [010]
b-direction. The FeII site also has a greater concentration of electrons at the Fermi surface, specifically more d-electrons which may be related to its role in the anomalous high temperature phonon
stiffening. The anomalous low-energy phonon stiffening stops abruptly at the magnetic transition,
and experimental spectra show these modes soften with temperature in the paramagnetic phase.
The magnetic ordering of cementite has a strong effect of the vibrational behavior, and also the
elastic constants. These thermal effects influence the anisotopy of the material structure and the
interatomic interactions.

Chapter 8

Summary
The vibrations of materials play a vital role in defining high-temperature material thermodynamics.
The formalism for understanding material vibrations is well developed; however, finite temperature
effects are often underestimated or ignored in thermodynamic assessments. Nonharmonic effects
are critical to understanding and quantifying the Gibbs free energy at temperatures of material
phase transformations, and they can also provide information on atomic bonding. Interatomic force
constants provide an especially advantageous avenue for exploring the specifics of material dynamic
stability, specifically the precursors to, and mechanisms of, diffusionless structural transformations.
Assumptions of harmonic interatomic potentials are embedded in many theoretical expressions of
material properties and they often pass without notice. Remembering or removing these assumptions
is essential for in-depth understanding of materials, especially materials under extreme conditions
of temperature and pressure. The role of nonharmonic effects in elastic response was considered
here, but many other physical effects including thermal conductivity and diffusion also rely on the
accurate descriptions interatomic interactions and phonon dynamics.
Recent computational innovations have lead to exciting new tools for understanding atom interactions in hot materials, but experimental verification of computed properties are still extremely
important. Experimental measurements of phonon spectra are a robust benchmark for assessing
the validity of computational material models. Development of advanced predictive computational
methodology in material thermodynamics will necessarily include nonharmonic effects.
Understanding phonon dynamics is essential for a complete physical picture of materials at high
temperatures. However, phonons are also strongly influenced by other physical phenomena such
as magnetic ordering and magnetic excitations. The high temperature vibrational studies of bcc
α-Fe and Fe3 C described here show that magnetic order strongly affects the interatomic interactions
underlying vibrational frequencies. Additionally, these studies highlight the vast range of effects
resulting from phonon interactions with magnetic order. Temperature-induced magnetic disorder in
bcc Fe reduced the energy of specific phonon modes, increasing the entropy, both before and after
the magnetic Curie transition. In orthorhombic cementite, the ferromagnetic phase showed thermal

87
increases in phonon energy that lower the vibrational entropy, an effect that is promptly reversed at
the magnetic transition. Understanding magnetic phenomena in materials remains an active frontier
in materials research, and the interactions between magnetic ordering and vibrational response are
likely important in a large number of materials.

Bibliography
[1] H. Okamoto. The C-Fe (carbon-iron) system. Journal of Phase Equilibria, 13(5):543–565,
October 1992.
[2] Stephen L. Sass. The substance of civilization: materials and human history from the stone
age to the age of silicon. Arcade Pub. : Distributed by Little, Brown and Co, New York, 1st
ed edition, 1998.
[3] Joseph R. Davis. ASM Specialty Handbook: Heat-Resistant Materials. ASM International,
January 1997.
[4] Josiah Willard Gibbs. The Scientific Papers of J. Willard Gibbs. Longmans, Green and
Company, 1906.
[5] A. Einstein. Die Plancksche Theorie der Strahlung und die Theorie der spezifischen Wrme.
Annalen der Physik, 327(1):180–190, January 1906.
[6] P. Debye. Zur Theorie der spezifischen Wrmen. Annalen der Physik, 344(14):789–839, January
1912.
[7] J. B. Austin. Entropy, heat Content, and Free Energy of Iron. Industrial & Engineering
Chemistry, 24(12):1388–1391, December 1932.
[8] Nobelstiftelsen, editor. Physics. Nobel lectures, including presentation speeches and laureates’
biographies. Published for the Nobel Foundation by Elsevier Pub. Co. etc, Amsterdam ; New
York, etc, 1967.
[9] W. L. Bragg. The Specular Reflection of X-rays. Nature, 90(2250):410–410, December 1912.
[10] Charles Kittel. Introduction to solid state physics. Wiley, Hoboken, NJ, 8th ed edition, 2005.
[11] Neil W. Ashcroft and N. David Mermin. Solid state physics. Saunders College, Philadelphia,
1976.
[12] V. J. Minkiewicz, G. Shirane, and R. Nathans. Phonon dispersion relation for iron. Physical
Review, 162(3):528, October 1967.

89
[13] B.N. Brockhouse, H.E. Abou-Helal, and E.D. Hallman. Lattice vibrations in iron at 296K.
Solid State Communications, 5(4):211–216, April 1967.
[14] Max G Kresch. Temperature Dependence of Phonons in Elemental Cubic Metals Studied by
Inelastic Scattering of Neutrons and X-Rays. PhD. thesis. California Institute of Technology,
2009.
[15] M. Kresch, O. Delaire, R. Stevens, J. Y. Y. Lin, and B. Fultz. Neutron scattering measurements
of phonons in nickel at elevated temperatures. Physical Review B, 75(10):104301, March 2007.
[16] M. Kresch, M. Lucas, O. Delaire, J. Y. Y. Lin, and B. Fultz. Phonons in aluminum at
high temperatures studied by inelastic neutron scattering. Physical Review B, 77(2):024301–9,
January 2008.
[17] M. L. Winterrose, L. Mauger, I. Halevy, A. F. Yue, M. S. Lucas, J. A. Muñoz, H. Tan, Y. Xiao,
P. Chow, W. Sturhahn, T. S. Toellner, E. E. Alp, and B. Fultz. Dynamics of iron atoms across
the pressure-induced Invar transition in Pd3 Fe. Physical Review B, 83(13):134304, April 2011.
[18] L. Mauger, M. S. Lucas, J. A. Muñoz, S. J. Tracy, M. Kresch, Yuming Xiao, Paul Chow,
and B. Fultz. Nonharmonic phonons in α-iron at high temperatures. Physical Review B,
90(6):064303, August 2014.
[19] R. J. Weiss and K. J. Tauer. Components of the thermodynamic functions of iron. Physical
Review, 102(6):1490, June 1956.
[20] G. K. White and M. L. Minges. Thermophysical properties of some key solids: An update.
International Journal of Thermophysics, 18(5):1269–1327, September 1997.
[21] A. T. Dinsdale. SGTE data for pure elements. Calphad, 15(4):317–425, October 1991.
[22] Y. S Touloukian. Thermophysical Properties of Matter. IFI/Plenum, New York, 1970.
[23] I. G. Wood, Lidunka Voadlo, K. S. Knight, David P. Dobson, W. G. Marshall, G. David
Price, and John Brodholt. Thermal expansion and crystal structure of cementite, Fe 3 C,
between 4 and 600K determined by time-of-flight neutron powder diffraction. Journal of
Applied Crystallography, 37(1):82–90, January 2004.
[24] Orson L. Anderson. The grueneisen ratio for the last 30 years. Geophysical Journal International, 143(2):279–294, 2000.
[25] Y.C Liu, F Sommer, and E.J Mittemeijer. Calibration of the differential dilatometric measurement signal upon heating and cooling; thermal expansion of pure iron. Thermochimica
Acta, 413(12):215–225, April 2004.

90
[26] Ichiro Seki and Kazuhiro Nagata. Lattice Constant of Iron and Austenite Including Its Supersaturation Phase of Carbon. ISIJ International, 45(12):1789–1794, 2005.
[27] N Ridley and H Stuart. Lattice parameter anomalies at the Curie point of pure iron. Journal
of Physics D: Applied Physics, 1(10):1291–1295, October 1968.
[28] D. J. Dever. Temperature dependence of the elastic constants in iron single crystals: relationship to spin order and diffusion anomalies. Journal of Applied Physics, 43(8):3293–3301,
August 1972.
[29] J. J. Adams, D. S. Agosta, R. G. Leisure, and H. Ledbetter. Elastic constants of monocrystal
iron from 3 to 500k. Journal of Applied Physics, 100(11):113530–113530–7, December 2006.
[30] J. A. Rayne and B. S. Chandrasekhar. Elastic constants of iron from 4.2 to 300K. Physical
Review, 122(6):1714, June 1961.
[31] F. Körmann, A. Dick, B. Grabowski, B. Hallstedt, T. Hickel, and J. Neugebauer. Free energy of
bcc iron: Integrated ab initio derivation of vibrational, electronic, and magnetic contributions.
Physical Review B, 78(3):033102, July 2008.
[32] J. Trampenau, A. Heiming, W. Petry, M. Alba, C. Herzig, W. Miekeley, and H. R. Schober.
Phonon dispersion of the bcc phase of group-IV metals. III. bcc hafnium. Physical Review B,
43(13):10963–10969, May 1991.
[33] W. Petry, A. Heiming, and J. Trampenaux. Phonons at martensitic phase transitions of bcc-ti,
bcc-zr and bcc-hf. MRS Online Proceedings Library, 166:C215–C228, 1989.
[34] W. Petry.

Dynamical precursors of martensitic transitions.

Journal de Physique IV,

05(C2):C2–15–C2–28, February 1995.
[35] Duane C. Wallace. Thermodynamics of Crystals. Dover Publications, January 1998.
[36] Brent Fultz.

Vibrational thermodynamics of materials.

Progress in Materials Science,

55(4):247–352, May 2010.
[37] R. L. Mossbauer. Recoilless Nuclear Resonance Absorption. Annual Review of Nuclear Science,
12(1):123–152, 1962.
[38] Dominic P. E. Dickson and Frank J. Berry, editors. Mössbauer spectroscopy. Cambridge
University Press, Cambridge ; New York, 1986.
[39] G. Klingelhöfer, B. Bernhardt, J. Foh, U. Bonnes, D. Rodionov, P. A. De Souza, Ch Schröder,
R. Gellert, S. Kane, P. Gtlich, and E. Kankeleit. The Miniaturized Mössbauer Spectrometer

91
MIMOS II for Extraterrestrial and Outdoor Terrestrial Applications: A Status Report. In
P. Gtlich, B. W. Fitzsimmons, R. Rffer, and H. Spiering, editors, Mössbauer Spectroscopy,
pages 371–379. Springer Netherlands, 2003.
[40] Jennifer M. Jackson, Wolfgang Sturhahn, Michael Lerche, Jiyong Zhao, Thomas S. Toellner,
E. Ercan Alp, Stanislav V. Sinogeikin, Jay D. Bass, Caitlin A. Murphy, and June K. Wicks.
Melting of compressed iron by monitoring atomic dynamics. Earth and Planetary Science
Letters, 362:143–150, January 2013.
[41] S. J. Tracy, L. Mauger, H. J. Tan, J. A. Muñoz, Yuming Xiao, and B. Fultz. Polaron-ion
correlations in Lix FePO4 studied by x-ray nuclear resonant forward scattering at elevated
pressure and temperature. Physical Review B, 90(9):094303, September 2014.
[42] W. Sturhahn and V.G. Kohn. Theoretical aspects of incoherent nuclear resonant scattering.
Hyperfine Interactions, 123-124(1):367–399, March 1999.
[43] U. van Bürck, D. P. Siddons, J. B. Hastings, U. Bergmann, and R. Hollatz. Nuclear forward
scattering of synchrotron radiation. Physical Review B, 46(10):6207, 1992.
[44] G. V. Smirnov. General properties of nuclear resonant scattering. Hyperfine Interactions,
123-124(1-4):31–77, March 1999.
[45] W. Sturhahn. CONUSS and PHOENIX: evaluation of nuclear resonant scattering data. Hyperfine Interactions, 125(1):149–172, March 2000.
[46] A.I. Chumakov and W. Sturhahn. Experimental aspects of inelastic nuclear resonance scattering. Hyperfine Interactions, 123-124(1):781–808, March 1999.
[47] B. Roldan Cuenya, W. Keune, R. Peters, E. Schuster, B. Sahoo, U. von Hörsten, W. Sturhahn,
J. Zhao, T. S. Toellner, E. E. Alp, and S. D. Bader. High-energy phonon confinement in
nanoscale metallic multilayers. Physical Review B, 77(16):165410, April 2008.
[48] B. Roldan Cuenya, A. Naitabdi, J. Croy, W. Sturhahn, J. Y. Zhao, E. E. Alp, R. Meyer,
D. Sudfeld, E. Schuster, and W. Keune. Atomic vibrations in iron nanoclusters: Nuclear
resonant inelastic x-ray scattering and molecular dynamics simulations. Physical Review B,
76(19):195422, November 2007. Copyright (C) 2010 The American Physical Society; Please
report any problems to prola@aps.org.
[49] T. Ruckert, W. Keune, B. Sahoo, W. Sturhahn, T. S. Toellner, E. E. Alp, and R. Röhlsberger.
Atomic Vibrational Dynamics of Thin Films Studied by Nuclear Resonant Inelastic X-Ray
Scattering: Amorphous Tb1−x Fex Alloys. Hyperfine interactions, 144(1):65–76, 2002.

92
[50] B. Sahoo, W. Keune, W. Sturhahn, T.S. Toellner, and E.E. Alp. Atomic vibrational dynamics
of amorphous Fe-Mg alloy thin films. Journal of Physics and Chemistry of Solids, 66(12):2263–
2270, December 2005.
[51] E. A Tanis. Phonon Density of States of 57-Iron and 161-Dysprosium in DyFe3 by Nuclear
Resonant Inelastic X-Ray Scattering Under High Pressure. Master’s thesis, University of
Nevada, 2010.
[52] Hillary L. Smith, B. C. Hornbuckle, L. Mauger, B. Fu, S. J. Tracy, G. B. Thompson, M. S.
Lucas, Y. Xiao, M. Y. Hu, J. Zhao, E. Ercan Alp, and B. Fultz. Changes in vibrational entropy
during the early stages of chemical unmixing in fcc Cu6% Fe. Acta Materialia, 61(19):7466–
7472, November 2013.
[53] J. A. Muñoz, M. S. Lucas, L. Mauger, I. Halevy, J. Horwath, S. L. Semiatin, Yuming Xiao,
Paul Chow, M. B. Stone, D. L. Abernathy, and B. Fultz. Electronic structure and vibrational
entropies of fcc Au-Fe alloys. Physical Review B, 87(1):014301, January 2013.
[54] M. S. Lucas, L. Mauger, J. A. Muñoz, I. Halevy, J. Horwath, S. L. Semiatin, S. O. Leontsev,
M. B. Stone, D. L. Abernathy, Yuming Xiao, Paul Chow, and B. Fultz. Phonon densities of
states of face-centered-cubic Ni-Fe alloys. Journal of Applied Physics, 113(17):17A308, 2013.
[55] J. A. Muñoz, M. S. Lucas, O. Delaire, M. L. Winterrose, L. Mauger, Chen W. Li, A. O. Sheets,
M. B. Stone, D. L. Abernathy, Yuming Xiao, Paul Chow, and B. Fultz. Positive Vibrational
Entropy of Chemical Ordering in FeV. Physical Review Letters, 107(11):115501, 2011.
[56] M. S. Lucas, J. A. Muoz, O. Delaire, N. D. Markovskiy, M. B. Stone, D. L. Abernathy,
I. Halevy, L. Mauger, J. B. Keith, M. L. Winterrose, Yuming Xiao, M. Lerche, and B. Fultz.
Effects of composition, temperature, and magnetism on phonons in bcc Fe-V alloys. Physical
Review B, 82(14):144306, October 2010.
[57] M. S. Lucas, J. A. Muñoz, L. Mauger, Chen W. Li, A. O. Sheets, Z. Turgut, J. Horwath, D. L.
Abernathy, M. B. Stone, O. Delaire, Yuming Xiao, and B. Fultz. Effects of chemical composition and B2 order on phonons in bcc FeCo alloys. Journal of Applied Physics, 108(2):023519,
July 2010.
[58] M. S. Lucas, O. Delaire, M. L. Winterrose, T. Swan-Wood, M. Kresch, I. Halevy, B. Fultz,
Jingzhu Hu, M. Lerche, M. Y. Hu, and M. Somayazulu. Effects of vacancies on phonon entropy
of B2 FeAl. Physical Review B, 80(21):214303, December 2009.
[59] A. D. Christianson, M. D. Lumsden, O. Delaire, M. B. Stone, D. L. Abernathy, M. A. McGuire,
A. S. Sefat, R. Jin, B. C. Sales, D. Mandrus, E. D. Mun, P. C. Canfield, J. Y. Y. Lin, M. Lucas,

93
M. Kresch, J. B. Keith, B. Fultz, E. A. Goremychkin, and R. J. McQueeney. Phonon Density
of States of LaFeAsO1−x Fx . Physical Review Letters, 101(15):157004–4, October 2008.
[60] O. Delaire, M. Kresch, and B. Fultz. Vibrational entropy of the gamma-alpha martensitic
transformation in Fe71Ni29. Philos. Mag., 85(30):3567–3583, 2005.
[61] M. S. Lucas, A. Papandrew, B. Fultz, and M. Y. Hu. Partial phonon densities of states of
Fe-57 in Fe-Cr: Analysis by a local-order cluster expansion. Phys. Rev. B, 75(5):8, 2007.
[62] David S. Sholl and Janice A. Steckel. Density functional theory: a practical introduction.
Wiley, Hoboken, N.J, 2009.
[63] Richard M. Martin. Electronic structure: basic theory and practical methods. Cambridge
University Press, Cambridge, UK ; New York, 2004.
[64] P. Hohenberg and W. Kohn. Inhomogeneous Electron Gas. Physical Review, 136(3B):B864–
B871, November 1964.
[65] W. Kohn and L. J. Sham. Self-Consistent Equations Including Exchange and Correlation
Effects. Physical Review, 140(4A):A1133–A1138, November 1965.
[66] R. Car and M. Parrinello. Unified Approach for Molecular Dynamics and Density-Functional
Theory. Physical Review Letters, 55(22):2471–2474, November 1985.
[67] Atsushi Togo, Fumiyasu Oba, and Isao Tanaka. First-principles calculations of the ferroelastic transition between rutile-type and cacl2 -type sio2 at high pressures. Physical Review B,
78(13):134106, October 2008.
[68] A. van de Walle, M. Asta, and G. Ceder. The Alloy Theoretic Automated Toolkit: A User
Guide. Calphad, 26(4):539–553, December 2002. arXiv: cond-mat/0212159.
[69] R. P. Feynman. Forces in Molecules. Physical Review, 56(4):340–343, August 1939.
[70] B. J. Alder and T. E. Wainwright. Studies in Molecular Dynamics. I. General Method. The
Journal of Chemical Physics, 31(2):459–466, August 1959.
[71] Olle Hellman, Peter Steneteg, I. A. Abrikosov, and S. I. Simak. Temperature dependent
effective potential method for accurate free energy calculations of solids. Physical Review B,
87(10):104111, March 2013.
[72] C. Z. Wang, E. Tosatti, and A. Fasolino. Structure, Transition, and Dynamics of a Displacive
Incommensurate Surface Reconstruction. Physical Review Letters, 60(25):2661–2664, June
1988.

94
[73] Rainer Storn and Kenneth Price. Differential evolution a simple and efficient heuristic for
global optimization over continuous spaces. Journal of Global Optimization, 11(4):341–359,
December 1997.
[74] Michael M. McKerns, Leif Strand, Tim Sullivan, Alta Fang, and Michael A. G. Aivazis. Building a framework for predictive science. arXiv:1202.1056 [cs], February 2012.
[75] Werner Pepperhoff and Mehmet Acet. Constitution and Magnetism of Iron and its Alloys.
Springer, 1st edition. edition, November 2010.
[76] Fritz Körmann, Abed Al Hasan Breidi, Sergei L. Dudarev, Nathalie Dupin, Gautam Ghosh,
Tilmann Hickel, Pavel Korzhavyi, Jorge A. Muñoz, and Ikuo Ohnuma. Lambda transitions
in materials science: Recent advances in CALPHAD and first-principles modelling. Physica
Status Solidi B, 251(1):53–80, January 2014.
[77] Michel H. G. Jacobs and Rainer Schmid-Fetzer. Thermodynamic properties and equation of
state of fcc aluminum and bcc iron, derived from a lattice vibrational method. Physics and
Chemistry of Minerals, 37:721, May 2010.
[78] O. Delaire, M. Kresch, J. A. Muñoz, M. S. Lucas, J. Y. Y. Lin, and B. Fultz. Electron-phonon
interactions and high-temperature thermodynamics of vanadium and its alloys. Physical Review B, 77(21):214112–12, June 2008.
[79] O. Delaire, M. S. Lucas, J. A. Muñoz, M. Kresch, and B. Fultz. Adiabatic electron-phonon
interaction and high-temperature thermodynamics of a15 compounds. Physical Review Letters,
101(10):105504–4, 2008.
[80] C Van Dijk and J Bergsma. In Neutron Inelastic Scattering, volume 1, page 157, Vienna, 1968.
IAEA.
[81] S. K. Satija, R. P. Comes, and G. Shirane. Neutron scattering measurements of phonons in
iron above and below tc. Physical Review B, 32(5):3309–3311, September 1985.
[82] J. Zarestky and C. Stassis. Lattice dynamics of γ-Fe. Physical Review B, 35(9):4500, March
1987.
[83] J. Neuhaus, W. Petry, and A. Krimmel. Phonon softening and martensitic transformation in
α-fe. Physica B: Condensed Matter, 234-236:897–899, June 1997.
[84] H Hasegawa, M W Finnis, and D G Pettifor. Phonon softening in ferromagnetic BCC iron.
Journal of Physics F: Metal Physics, 17(10):2049–2055, October 1987.

95
[85] A. Heiming, W. Petry, and J. Trampenau. Are Martensitic Phase Transitions In Pure Group
3 and 4 Metals Driven By Lattice Vibrations ? Journal de Physique IV, 01(C4):C4–83–C4–88,
November 1991.
[86] F. Körmann, A. Dick, B. Grabowski, T. Hickel, and J. Neugebauer. Atomic forces at finite
magnetic temperatures: Phonons in paramagnetic iron. Physical Review B, 85(12):125104,
March 2012.
[87] F. Körmann, A. Dick, T. Hickel, and J. Neugebauer. Rescaled monte carlo approach for
magnetic systems: Ab initio thermodynamics of bcc iron. Physical Review B, 81(13):134425,
April 2010.
[88] I. Leonov, A. I. Poteryaev, V. I. Anisimov, and D. Vollhardt. Calculated phonon spectra of
paramagnetic iron at the - phase transition. Physical Review B, 85(2):020401, January 2012.
[89] Matthieu J. Verstraete. Ab initio calculation of spin-dependent electronphonon coupling in
iron and cobalt. Journal of Physics: Condensed Matter, 25(13):136001, April 2013.
[90] Wolfgang Sturhahn and Jennifer M. Jackson. Geophysical applications of nuclear resonant
spectroscopy. Geophysical Society Of America Special Papers, 421:157, 2007.
[91] S. Klotz and M. Braden. Phonon dispersion of bcc iron to 10 GPa. Physical Review Letters,
85(15):3209–3212, October 2000.
[92] J. Trampenau, W. Petry, and C. Herzig. Temperature dependence of the lattice dynamics of
chromium. Physical Review B, 47(6):3132, February 1993.
[93] J. Crangle and G. M. Goodman. The magnetization of pure iron and nickel. Proceedings of
the Royal Society of London A, 321(1547):477–491, March 1971.
[94] A J Holden, V Heine, and J H Samson. Magnetic contributions to thermal expansion of
transition metals: implications for local moments above t c. Journal of Physics F: Metal
Physics, 14(4):1005–1020, April 1984.
[95] Junqi Yin, Markus Eisenbach, Don M. Nicholson, and Aurelian Rusanu. Effect of lattice vibrations on magnetic phase transition in bcc iron. Physical Review B, 86(21):214423, December
2012.
[96] Xiaoli Tang, Chen W. Li, and B. Fultz. Anharmonicity-induced phonon broadening in aluminum at high temperatures. Physical Review B, 82(18):184301, November 2010.
[97] Xiaoli Tang and B. Fultz. First-principles study of phonon linewidths in noble metals. Physical
Review B, 84(5):054303, August 2011.

96
[98] A. T. Boothroyd, T. G. Perring, A. D. Taylor, D. McK. Paul, and H. A. Mook. High energy spin
waves in iron measured by neutron scattering. Journal of Magnetism and Magnetic Materials,
104–107:713–714, February 1992.
[99] M. Palumbo, B. Burton, A. Costa e Silva, B. Fultz, B. Grabowski, G. Grimvall, B. Hallstedt,
O. Hellman, B. Lindahl, A. Schneider, P. E. A. Turchi, and W. Xiong. Thermodynamic
modelling of crystalline unary phases. Physica Status Solidi B, 251(1):14–32, January 2014.
[100] T. A. Kovats and J. C. Walker. Mssbauer absorption in fe57 in metallic iron from the curie
point to the γ-δ transition. Physical Review, 181(2):610–618, May 1969.
[101] U. Bergmann, S. D. Shastri, D. P. Siddons, B. W. Batterman, and J. B. Hastings. Temperature
dependence of nuclear forward scattering of synchrotron radiation in α-57 Fe. Physical Review
B, 50(9):5957–5961, September 1994.
[102] A. I. Chumakov, R. Rüffer, A. Q. R. Baron, H. Grünsteudel, and H. F. Grünsteudel. Temperature dependence of nuclear inelastic absorption of synchrotron radiation in α-57 Fe. Physical
Review B, 54(14):R9596–R9599, October 1996.
[103] B. Kolk, A. L. Bleloch, and D. B. Hall. Recoilless fraction studies of iron near the curie
temperature. Hyperfine Interactions, 29(1-4):1377–1380, February 1986.
[104] Jie Li Walker, D. and S. M. Clark B. Kalkan. Thermal, compositional, and compressional
demagnetization of cementite. American Mineralogist, Submitted.
[105] Lili Gao, Bin Chen, Jingyun Wang, Esen E. Alp, Jiyong Zhao, Michael Lerche, Wolfgang
Sturhahn, Henry P. Scott, Fang Huang, Yang Ding, Stanislav V. Sinogeikin, Craig C. Lundstrom, Jay D. Bass, and Jie Li. Pressure-induced magnetic transition and sound velocities
of Fe3c: Implications for carbon in the Earth’s inner core. Geophysical Research Letters,
35(17):L17306, September 2008.
[106] M. Ron and Z. Mathalone. Hyperfine Interactions of 57 Fe in Fe3 C. Physical Review B, 4(3):774–
777, August 1971.
[107] X. L. Dong, Z. D. Zhang, Q. F. Xiao, X. G. Zhao, Y. C. Chuang, S. R. Jin, W. M. Sun, Z. J. Li,
Z. X. Zheng, and H. Yang. Characterization of ultrafine γ-Fe(C), α-Fe(C) and Fe3 C particles
synthesized by arc-discharge in methane. Journal of Materials Science, 33(7):1915–1919, April
1998.
[108] Takahiro Miki and Koutarou Ishii. Decomposition Behavior of Fe3 C under Ar Atmosphere.
ISIJ International, 54(1):29–31, 2014.

97
[109] K. O. E. Henriksson and K. Nordlund. Simulations of cementite: An analytical potential for
the Fe-C system. Physical Review B, 79(14):144107, April 2009.
[110] Bengt Hallstedt, Dejan Djurovic, Jörg von Appen, Richard Dronskowski, Alexey Dick, Fritz
Körmann, Tilmann Hickel, and Jörg Neugebauer. Thermodynamic properties of cementite ().
Calphad, 34(1):129–133, March 2010.
[111] A. Dick, F. Körmann, T. Hickel, and J. Neugebauer. Ab initio based determination of thermodynamic properties of cementite including vibronic, magnetic, and electronic excitations.
Physical Review B, 84(12):125101, 2011.
[112] C. Jiang, S. G Srinivasan, A. Caro, and S. A Maloy. Structural, elastic, and electronic properties of Fe3 C from first principles. Journal of Applied Physics, 103(4):043502–043502–8,
February 2008.
[113] Chao Jiang and Srivilliputhur G. Srinivasan. Unexpected strain-stiffening in crystalline solids.
Nature, 496(7445):339–342, April 2013.
[114] M. Nikolussi, S. L. Shang, T. Gressmann, A. Leineweber, E. J. Mittemeijer, Y. Wang, and
Z. K. Liu. Extreme elastic anisotropy of cementite, Fe3c: First-principles calculations and
experimental evidence. Scripta Materialia, 59(8):814–817, October 2008.
[115] S. P. Dodd, G. A. Saunders, M. Cankurtaran, B. James, and M. Acet. Ultrasonic study of the
temperature and hydrostatic-pressure dependences of the elastic properties of polycrystalline
cementite (Fe3C). physica status solidi (a), 198(2):272–281, August 2003.
[116] Toshihiko Shigematsu. Invar Properties of Cementite (Fe1−x Mex C, Me=Cr, Mn, Ni. Journal
of the Physical Society of Japan, 39(4):915–920, 1975.
[117] E. Duman, M. Acet, T. Hülser, E. F. Wassermann, B. Rellinghaus, J. P. Itié, and P. Munsch.
Large spontaneous magnetostrictive softening below the Curie temperature of Fe3C Invar
particles. Journal of Applied Physics, 96(10):5668–5672, November 2004.
[118] E. Duman, M. Acet, E. F. Wassermann, J. P. Itié, F. Baudelet, O. Mathon, and S. Pascarelli.
Magnetic Instabilities in Fe3C Cementite Particles Observed with Fe K-Edge X-Ray Circular
Dichroism under Pressure. Physical Review Letters, 94(7):075502, February 2005.
[119] Sergii Khmelevskyi, Andrei V Ruban, and Peter Mohn. Electronic structure analysis of the
pressure induced metamagnetic transition and magnetovolume anomaly in Fe3 Ccementite.
Journal of Physics: Condensed Matter, 17(46):7345–7352, November 2005.
[120] Guillaume Fiquet, James Badro, Eugene Gregoryanz, Yingwei Fei, and Florent Occelli. Sound
velocity in iron carbide (Fe3c) at high pressure: Implications for the carbon content of the

98
Earth’s inner core. Physics of the Earth and Planetary Interiors, 172(12):125–129, January
2009.
[121] C. Prescher, L. Dubrovinsky, C. McCammon, K. Glazyrin, Y. Nakajima, A. Kantor, M. Merlini, and M. Hanfland. Structurally hidden magnetic transitions in Fe3c at high pressures.
Physical Review B, 85(14):140402, April 2012.
[122] Jung-Fu Lin, Viktor V. Struzhkin, Ho-kwang Mao, Russell J. Hemley, Paul Chow, Michael Y.
Hu, and Jie Li. Magnetic transition in compressed Fe3c from x-ray emission spectroscopy.
Physical Review B, 70(21):212405, December 2004.
[123] Lili Gao. Density, magnetic properties and sound velocities of iron-rich materials at high temperature and high pressure. PhD thesis, University of Illinois at Urbana-Champaign, December
2010.
[124] J. Li, H. K. Mao, Y. Fei, E. Gregoryanz, M. Eremets, and C. S. Zha. Compression of Fe3C to
30 GPa at room temperature. Physics and Chemistry of Minerals, 29(3):166–169, April 2002.
[125] G. Kresse and J. Hafner. Ab initio molecular dynamics for open-shell transition metals. Physical Review B, 48(17):13115–13118, November 1993.
[126] G. Kresse and J. Furthmller. Efficiency of ab-initio total energy calculations for metals and
semiconductors using a plane-wave basis set. Computational Materials Science, 6(1):15–50,
July 1996.
[127] G. Kresse and J. Furthmüller. Efficient iterative schemes for ab initio total-energy calculations
using a plane-wave basis set. Physical Review B, 54(16):11169–11186, October 1996.
[128] G. Kresse and D. Joubert. From ultrasoft pseudopotentials to the projector augmented-wave
method. Physical Review B, 59(3):1758–1775, January 1999.
[129] John P. Perdew, Kieron Burke, and Matthias Ernzerhof. Generalized Gradient Approximation
Made Simple. Physical Review Letters, 77(18):3865–3868, October 1996.
[130] Hendrik J. Monkhorst and James D. Pack. Special points for Brillouin-zone integrations.
Physical Review B, 13(12):5188–5192, June 1976.
[131] Ion Errea, Matteo Calandra, and Francesco Mauri. Anharmonic free energies and phonon dispersions from the stochastic self-consistent harmonic approximation: Application to platinum
and palladium hydrides. Physical Review B, 89(6):064302, February 2014.
[132] O. Hellman, I. A. Abrikosov, and S. I. Simak. Lattice dynamics of anharmonic solids from
first principles. Physical Review B, 84(18):180301, November 2011.

99
[133] Olle Hellman and I. A. Abrikosov. Temperature-dependent effective third-order interatomic
force constants from first principles. Physical Review B, 88(14):144301, October 2013.
[134] Stefano Baroni, Stefano de Gironcoli, Andrea Dal Corso, and Paolo Giannozzi. Phonons and
related crystal properties from density-functional perturbation theory. Reviews of Modern
Physics, 73(2):515–562, July 2001.
[135] J. Häglund, G. Grimvall, and T. Jarlborg. Electronic structure, x-ray photoemission spectra,
and transport properties of Fe3C (cementite). Physical Review B, 44(7):2914–2919, 1991.
[136] Lili Gao, Bin Chen, Jiyong Zhao, Esen E. Alp, Wolfgang Sturhahn, and Jie Li. Effect of
temperature on sound velocities of compressed Fe3C, a candidate component of the Earth’s
inner core. Earth and Planetary Science Letters, 309(34):213–220, September 2011.
[137] E. J. Fasiska and G. A. Jeffrey. On the cementite structure. Acta Crystallographica, 19(3):463–
471, September 1965.
[138] Michael Y. Hu, Wolfgang Sturhahn, Thomas S. Toellner, Philip D. Mannheim, Dennis
E. Brown, Jiyong Zhao, and E. Ercan Alp. Measuring velocity of sound with nuclear resonant inelastic x-ray scattering. Physical Review B, 67(9):094304, March 2003.
[139] Lili Gao, Bin Chen, Michael Lerche, Esen E. Alp, Wolfgang Sturhahn, Jiyong Zhao, Hasan
Yava, and Jie Li. Sound velocities of compressed Fe3 C from simultaneous synchrotron X-ray
diffraction and nuclear resonant scattering measurements. Journal of Synchrotron Radiation,
16(6):714–722, November 2009.
[140] Jon Alkorta and Javier Gil Sevillano. Assessment of elastic anisotropy and incipient plasticity
in Fe3C by nanoindentation. Journal of Materials Research, 27(01):45–52, 2012.
[141] Bon-Woong Koo, Young Jin Chang, Seung Pyo Hong, Chan Soon Kang, Shin Woong Jeong,
Won-Jong Nam, Il-Jeong Park, Young-Kook Lee, Kyu Hwan Oh, and Young-Woon Kim.
Experimental measurement of Youngs modulus from a single crystalline cementite. Scripta
Materialia, 82:25–28, July 2014.
[142] A. Leineweber. Anisotropic microstrain broadening in cementite, Fe3 C, caused by thermal
microstress: comparison between prediction and results from diffraction-line profile analysis.
Journal of Applied Crystallography, 45(5):944–949, October 2012.