Computational Modelling of Ceria-Based Solid Oxide Fuel Cell Electrolyte Materials
Authors: Ahmed Mohamed , Rodley David , Jones Thomas , Abdolvand Amin , Lightfoot Alison , Fruchtl Herbert , Baker Richard 
Published: July 2021

The significant advances in materials development for solid oxide fuel cells (SOFCs) and modelling techniques have challenged researchers to develop streamlined methodologies aiming at faster and more reliable fabrication of highly efficient cells. Instead of adopting a trial-and-error technique to discover the optimal composition of fuel cell component materials, alternative tools and techniques should be developed for predicting the likely best compositions in a cost-effective and timely manner. Many formulations have been studied for application as SOFC electrolytes over the years and doubly-doped ceria materials have been of considerable recent interest.
Here, a simulation methodology to calculate the lattice parameter and the oxygen ion migration energy of ceria-based electrolyte formulations was developed, and the results were analysed and benchmarked against corresponding results obtained experimentally in order to verify the efficacy of the simulation methodology. To this end, a supercell consisting of 2 x 2 x 2 ceria unit cells doubly-doped with samarium and gadolinium was modelled and simulated using molecular mechanics force field methods using the CP2K computational chemistry suite. Two compositional series were modelled with the general formula Ce32-(x+y)SmxGdyO64-(x+y)/2 where x and y are the number of samarium and gadolinium dopant atoms per supercell respectively. The first compositional series, Series One, had equal amounts of dopants but different total dopant concentration with values of x = y = 1, 2, 8, or 16. Series Two had a fixed total dopant concentration (x + y = 8) but different amounts of each dopant in each structure such that x had integer values from 0 to 8 (x = 0 ... 8) and y = 8 - x. Moreover, two doping strategies were simulated for each compositional series depending on how the dopant atoms and resulting vacancies were distributed throughout the crystal. The first doping strategy places one dopant atom in two or more of the eight minicubes comprising the 2 x 2 x 2 supercell for Series One structures and, for Series Two structures, one dopant atom is placed in each minicube starting with one gadolinium atom in each minicube and then progressively replacing one gadolinium atom with one samarium atom in each of the subsequent structures of this series until there is one samarium atom in each minicube. The second doping strategy places the dopants exclusively in the four minicubes that let the dopants be as far as possible from each other. In all cases, the oxygen vacancies are chosen to be the nearest neighbours, i.e., in the nearest lattice positions, to the dopant atoms.
For each structure, the lattice parameter was calculated using a cell optimisation technique utilizing the conjugate gradient method where the soft particle mesh Ewald (SPME) sums method and the DIPole Polarizable Ion Model (DIPPIM) for the long- and short-range potentials respectively. To estimate the activation energies associated with hopping of an oxygen ion from a specific site to a neighbouring site in the structures studied, the following steps were applied:

The positions of the dopant atoms in the pure crystal were assigned.

The at a distance within 3.0 Å was generated.

The climbing image nudged elastic band (CI-NEB) method was applied to determine all the activation barriers of the doped crystal.

The average energy barrier was calculated.

The results obtained for the lattice parameter were comparable to those observed experimentally and exhibited the same linear trends satisfying Vegard's law, albeit the trend for the fixed dopant concentration series had a slight negative slope. Simulation-derived values of the activation energy showed the same trends as those observed experimentally. The simulations also produced transition state curves showing the saddle point and energy minima for each structure; information that could not be obtained experimentally and that can shed some light on the mechanism of ionic conduction. It was found that the position of oxygen vacancies within the crystal structure play a key role in determining the activation energy of ion hopping, and that the different patterns of doping do not affect the overall trends in lattice parameter and activation energy. They do however produce different patterns in the energy profiles contributing to the energy barriers to ion hopping within the crystal structure.
This work demonstrates that computational simulations of SOFC electrolyte materials provide results comparable with those obtained experimentally and could also provide insights into the mechanisms of ionic conduction. The incongruence of the computational and experimental results was attributed to the limitations of the molecular mechanics force field methodology utilised, with the expectation that an ab initio density functional theory (DFT) calculation would yield closer conformance.

Figure 1