One of the reasons why we want to have the fully optimized unit cell of a crystal, is that the corresponding energy is needed to construct a computed phase diagram. Let’s see how that works.
Now comes the icing on your Fe-Al exercise cake. Up to now, you have found the volume, shape and internal positions that a FeAl-crystal would adopt if you would start from the particular hcp-based arrangement that was described at the start of this chapter. Let me say that again: you started from a particular arrangement of atoms, that was an ideal packing as if Fe-atoms and Al-atoms were hard spheres of the same size, ordered in a particular way. In reality, they are not hard spheres and they are not exactly equal in size. Therefore, you obtained a lowest-energy geometry that deviates somewhat from the starting structure. However, this does not guarantee that this lowest-energy solution will be the crystal structure adopted by the FeAl-alloy. Indeed, the experimental crystal structure is this one: a CsCl-type crystal. If you do the same geometry optimization process for this type of crystal, using the PBE functional, the same basis set size and the same pseudopotentials, you will find an energy that is lower than the one you found for the hcp-based crystal (for the same number of atoms). You don’t need to do this calculation explicitly if you don’t have the time right now (the resulting energy is in the table further down anyway), but it’s instructive to do it once, sooner or later.
We can conclude from this that our geometry optimization procedure will give us the quantitative geometrical details of the crystal that is qualitatively similar to our starting point, yet it will not lead to a crystal that is drastically different from the starting point. In some sense, DFT brings you to the local minimum in ‘crystal structure space’ that is nearest to your initial guess. In order to predict what would be the crystal structure really adopted by FeAl — the global minimum in crystal structure space — we would have to run geometry optimizations for thousands or millions of such initial guesses (ab initio crystal structure prediction is a field in its own right, that is under rapid development — a not very recent yet instructive starting point to read on this is available here).
Will Fe-atoms and Al-atoms combine into a Fe-Al alloy at all, or would they rather prefer to stay in a pure Fe-crystal separated from a pure Al-crystal? DFT can predict this as well. Here is a table with total energies for 4 fully geometry-optimized crystals: bcc-Fe, fcc-Al, FeAl in the orthorhombic cell (FeAl-ortho), and FeAl in the experimentally observed CsCl-type structure (FeAl-CsCl):
From this information you can calculate the formation energy (per atom) with respect to the elemental ground state crystals. It is the energy gain per atom if two solids in their elemental ground state (here: fcc-Al and bcc-Fe) react to form an alloy (here: FeAl in two different crystal structures). Use the data in this table to find the formation energies for all four crystals, and put them on the picture underneath (you will know you did it right if you find the same value for FeAl-CsCl as the one that is already on the picture).
These formation energies per atom are given in this picture for several different alloys:
(reproduced from http://dx.doi.org/10.1016/j.surfcoat.2017.05.091)
From this picture, you see there is indeed an energy gain when fcc-Al and bcc-Fe combine to form FeAl. That means the alloy is in principle possible. And the crystal structure that will occur, is the one for which the energy gain is maximal.
It is interesting to determine the formation energies of the crystals we found as intermediate steps during the geometry optimization of orthorhombic FeAl. Put on the picture the formation energy for the initial structure we started from (E=-737.51933681 Ry/cell), the one after volume optimization only (E=-737.53661703 Ry/cell), after b/a optimization (E=-737.55172056 Ry/cell) and after optimization of the atomic positions (E=-737.55692408 Ry/cell). It gives you a feeling which types of structural changes affect the energy gain during alloy formation most.
A picture with formation energy as a function of alloy composition is called a phase diagram. It shows us that if you take one dose of bcc iron for two doses of fcc aluminum, you can form a FeAl2-alloy. If you take 5 doses of iron for 8 doses of aluminum, you can form Fe5Al8. The latter has a favorable formation energy (energy gain with respect to the elemental phases), yet it is not a stable crystal. Why not: you can lower the energy of Fe5Al8 even further if you decompose it in 3 doses of FeAl2 and 2 doses of FeAl. You can’t apply that argument to FeAl2 — any decomposition into two or more other alloys will lead to a higher energy instead.
Graphically, it is straightforward to identify on the above picture all crystals that cannot lower their energy further by decomposition: these are the crystals that lie on the convex hull, which is the collection of straight lines that is on or below all data points (the green lines in the picture). Every crystal that is not on the convex hull can be decomposed into a weighted sum of the two crystals at its left and right on the convex hull (we just did that for Fe5Al8). Find out yourself how the decomposition of FeAl3 and Fe3Al would look like.
DFT-computed phase diagrams as in this picture are convenient to find out which alloys are likely to exist in experiment and which are not. And you often don’t need to do DFT calculations yourself to access such information. Indeed, the MaterialsProject database has them for you (take apps/phase diagrams, type Al-Fe, generate). Compare the picture you get there with the one given above. Any similarities or differences?
The Open Quantum Materials Database (OQMD) is very convenient for phase diagrams and decompositions too: select Analysis, Phase Diagrams, and type Al-Fe for the region of phase space. You’ll get a similar picture as before. Simultaneously, this illstrates one of the shortcomings of this kind of databases: they are only a smart as the amount of data they contain. According to OQMD (data in 2017), Al9Fe2 is a stable crystal — it lies on the convex hull. But that is because OQMD has no data on Fe4Al13 or FeAl2. Fortunately, you can feed the phase diagram tool of OQMD with your own data (‘use custom phase data file’).
The GLCP tool in OQMD is particularly useful for decomposition analysis: if you give it ‘Fe5Al8’, it will reply in a blink that this decomposes in given amounts of AlFe and Fe2Al9. No need for you to figure out the algebra. Well, OQMD does not know about FeAl2, hence its answer was not complete, always be aware of that. You can modify the phase data under the decomposition to include other phases or modify their formation energies, if you have more or better information yourself. Good to know: GLCP/OQMD works not only for binary phase diagrams, but for ternary, quaternary, … (no limit) diagrams as well. The algebra becomes rapidly cumbersome when working with many elements, hence this is a convenient tool.
Optional reading: a review paper on thermodynamic stability, phase diagrams, convex hull,… also generalized to different thermodynamic potentials.
So much on phase diagrams…
If you ran into a difficulty that you cannot solve, please report it in the forum hereunder. If you ran into a difficulty that you could solve, please share the solution with us too — it may help others. Feel free to post any other comment or thought about this exercise. And if you can answer questions posted by fellow students, please help them out.