Disorder in energy materials and strategies to model it

ABSTRACT

The functionality of the materials used for energy applications is critically determined by the physical properties of small active regions such as dopants, dislocations, interfaces, grain boundaries, etc. The capability to manipulate and utilize the inevitable disorder in materials, whether due to the finite-dimensional defects (such as vacancies, dopants, grain boundaries) or due to the complete atomic randomness (as in amorphous materials), can bring innovation in designing energy materials. With the increase in computational material science capabilities, it is now possible to understand the complexity present in materials due to various degrees of disorder resulting in pathways required for optimizing their efficiencies. This article provides a critical overview of such computational advancements specifically for designing realistic materials with various types of disorders for sustainable energy applications such as catalysts and electrochemical devices. The ultimate goal is to gain a thorough knowledge of the traditional approaches (implemented via tools such as density functional theory, and molecular dynamics) as well as modern approaches such as machine learning that exist for modeling the disorder present in materials, thereby identify the future opportunities for energy materials design and discovery.

Graphical abstract

KEYWORDS:

• Computational material science
• classical and quantum calculations
• disordered materials
• energy materials
• catalysis and electrochemical devices

1. Introduction

While remarkable research efforts have been made to achieve the sustainable-fuel future, built on economically viable and secure energy supply [1–6], the commercialization of associated technologies remains a challenge. The major hurdle in achieving the economic success from the existing sustainable technologies is that the materials used at their core are extremely expensive (such as Pt, Nafion) [7], and the materials that are cheap and abundant (such as glasses, ceramics) are not active enough to be used commercially. Improving the performance of inexpensive and abundant materials is a challenge because of the poor understanding of their fundamental properties arising from their complex structures. Often there exist structural defects in these materials whether punctual, which results in only a small deviation in their properties keeping their overall intrinsic properties mostly intact, or as severe as in amorphous materials, which result in materials having entirely new properties as compared to their crystalline counterparts.

Fully understanding, controlling, and using disorder present in materials is a challenging task. Despite that, evidence exists that disorder in materials is capable of boosting the performances of novel sustainable technologies, significantly [8,9]. Some shreds of evidence include disordered-electrode of rechargeable Li-batteries having both higher capacity and energy density than ordered electrodes [8,9]; high proton conductivity of solid oxide fuel cell (SOFC) electrolytes due to the presence of grain boundaries (GBs) [10]; the high catalytic activity of amorphous silica, mixed-metal oxides (like Cu-Mn mixed oxide), alloys such as Ni-B and Ni-P, and disordered nickel-organic photocatalyst [11–15]. Furthermore, it is also well known that many crystalline materials transform into their amorphous forms while in use. For example, VOHPO₄ · 0.5H₂O which is used to catalyze n-butane to maleic anhydride is introduced into the reactor as a crystalline material [16], but it transforms to a catalytically active disordered material during the reaction process. Another example is crystalline Co₃O₄ which converts to CoO_x(OH)_y on the surface during the oxygen evolution reaction [17].

Understanding and utilizing the disordered materials for energy applications require extensive laboratory testing, which can often become time-consuming and expensive and does not necessarily provide insights into the role that disorder plays in the chemical processes. The standard computational approaches are also not enough to mimic the disordered materials and to study their properties from nano- to electronic-scale [18]. However, with the recent developments in computational materials science, it is now possible to provide the long-awaited progress in our understanding of these complex materials, thereby guiding experiments to synthesize sustainable energy materials. This review, therefore, aims to provide the critical appraisal of the existing advanced computational approaches targeted to design realistic energy materials having disorder due to the co-existence of various types of defects. Specifically, sustainable catalysis and electrochemical devices (such as fuel cells and batteries) are focused in this review as the applications of complex disordered materials.

The structure of this review paper is as follows. In section 2, an overview of the possible disorders present in energy materials will be provided. Additionally, the successes achieved in identifying the role of disorder for sustainable technologies will be discoursed. In section 3, the current computational tools and strategies used for designing realistic energy materials with various degrees of disorder will be reviewed. Lastly, in section 4, future directions and opportunities that exist for the development of computational tools to enable the characterization, control, and utilization of material disorder for energy technologies will be discussed.

2. Types of disorders and their role in energy technologies

The commonly used nanocrystalline energy materials such as metal oxides, ceramics, transition metals, alloys, semiconductors, perovskites, etc. usually contain a partial low-symmetry environment due to one or multiple imperfections. Some energy materials, usually referred as ‘disordered’ materials, also lack long-range order completely and have the same key features as crystalline materials only at medium or short length scales. The disorder present in energy materials can result in minor to major deviation in their properties when compared to their pure crystalline forms. There exist the following five types of defects that can lead to disorder in materials:

1. 0D (dimensional) or point defects (such as vacancies, interstitials, etc.)

2. 1D or linear defects (such as dislocations, disclinations, etc.)

3. 2D or planar defects (such as grain boundaries (GBs), surfaces, etc.)

4. 3D or extended defects (such as pores, cracks, etc.)

5. Random arrangement of atoms (at least) at the long-range scales.

Usually, the materials prepared in laboratories contain a combination of these defects. Therefore, to get a better perspective, one can categorize disorder in materials as 1^st degree where crystalline materials are slightly disordered due to the random distribution of 0D defects; 2^nd degree where crystalline materials are disordered due to the presence of higher-order (1D, 2D, and 3D) defects such as grain boundaries, interfaces, dislocations, disclinations, pores, etc. in addition to 0D defects; and, 3^rd degree where materials are mostly disordered, i.e. have a random network of atoms and retain crystallinity only at short or medium range.

Among these categories, materials with 1^st degree disorder have been investigated the most via atomistic modeling techniques due to the simplicity of designing their structures [19–21]. Consequently, materials with 1^st degree disorder have been used tremendously in enhancing efficiencies of sustainable energy technologies. It has been shown that the existence of substitutional point defects in electrode materials increases the diffusivity [22–24], phase transformation rate [25], and even the cyclability [26] and stability [27] of the cell. Further, for H₂ fuel generation by solar water splitting, the modification of the promising photoanode BiVO₄ by substituting V with W [28], Mo [29], or P [30] or the increment of O vacancies by hydrogen treatment improves the charge transport [31]. Besides, several groups use point-defected catalysts for hydrogen evolution reaction (HER), and oxygen reduction reaction (ORR) in fuel cell technology [32]. Lastly, one of the most significant materials discovered in the 21^st century, ‘graphene,’ does not possess high catalytic properties by itself, either. Its catalytic properties increase considerably when treated with nitrogen to create point defects on the mesoporous carbon surface (Figure 1(a)). Further, the nitrogen-doped graphene can also bind through the pyridinic, and pyrrolic sites with tiny amounts of atomic Cobalt (Co) deposited on the surface (Figure 1(b, c)), resulting in an inexpensive and efficient catalyst for HER [33]. The effect of these substitutional point defects on the catalytic activity of graphene has also been investigated by Density Functional Theory (DFT) study (Figure 1(d)) and explained via charge population analysis (Figure 1(e)) [34].

Figure 1. Atomic cobalt on nitrogen-doped carbon (Co-NG) catalyst. a) STEM image of the Co-NG catalyst showing the defective and disordered nitrogen-doped carbon structure. b) and c) HAADF-STEM images of the Co-NG catalyst showing the substitutional Co defects dispersed in the carbon matrix at 1 nm and 0.5 nm scale, respectively. Adapted from Fei et al. [33] d) 3D representation of the Co-NG catalyst. e) DFT calculation for differential charge densities of Co-NG catalyst after absorbing hydrogen or oxygen. Yellow isosurfaces show the electron gain and cyan the electron loss. Adapted by permission from Springer Nature Customer Service Centre GmbH: Springer-Nature, Nature Materials, Atomic-level tuning of Co–N–C catalyst for high-performance electrochemical H₂O₂ production, Euiyeon Jung et al [34], Copyright 2020

In contrast, very few studies focussed on materials with higher-degree disorders due to the inherent heterogeneity in their composition and structure [35,36]. The realistic materials with higher-order defects usually contain several defects in the close vicinity to each other. For example, dislocations in bi-layer electrolytes used in SOFCs mostly occur to accommodate the lattice misfit of the two materials forming the GB or an interface [37]. Further, there exist several catalysts, such as hollow PtNi nanoparticles supported on carbon [38], in which GBs and vacancies co-exist and are responsible for the enhanced kinetics of the involved chemical reactions. However, these kinds of defected materials are usually modeled with a single defect, due to their complex structures, and large supercell sizes required (to reduce interactions between defects) [39,40], with some studies also providing methods to model nanomaterials having disorder due to their composition (e.g. Nanoalloys) and surface reconstruction (e.g. nanoparticles) [41].

Finally, various experimental studies have shown that even the materials with 3^rd order of disorder can prove exceptional for sustainable technologies. Amorphous metal oxides used as electrode materials for alkali ion batteries are a good example [26]. Lee et al. have recently described that the high energy density and high capacity of the disordered Li_1.211Mo_0.467Cr_0.3O₂ cathode result in higher Li-ion diffusion, and lower mechanical stress and capacity loss in contrast to the crystalline cathodes [8]. Another example is the amorphous molybdenum dioxide electrode (a-MoO₂) reported by Ku et al. [42], where the structural defects stabilize the Li storage allowing an unexpectedly high capacity (four Li per unit of MoO₂) and provide a high rate diffusion path for Li-ions through the opened vacancies and void spaces. Along with amorphous metal oxides [43], amorphous ceramics [44,45], and amorphous silicon [46], are also being studied for developing batteries with better performance. Furthermore, Hasannaeimi et al. studied recently the behavior of Pd- and Pt-based amorphous alloys catalyst in contrast with the pure Pt and Pd for hydrogen oxidation reaction in PEMFCs. The results showed that the disordered alloy structure favored the charge transfer on their surface due to their lower work function compared to pure catalysts, resulting in their higher catalytic activity compared to pure crystalline catalysts [47]. Clearly, the disorder in materials is not always bad. It is just not well understood, resulting in materials with disorder not been used extensively for energy applications. In the next section, we highlight some of the compelling strategies used so far to model disorder in materials from nano- to electronic- scales, thereby encouraging their use for future studies focusing on sustainable technologies.

3. Computational materials design strategies to characterize disorder

The structure of a crystalline material is defined by lattice symmetry, and the lattice constants, which can be efficiently designed at any scale for computational analyses. To study non-periodic disordered materials from nano to electronic level is a more difficult undertaking. Thanks to high computing power and modern techniques, it is now possible to design materials with structural defects. Below we summarized the current computational techniques and strategies capable of capturing disorder in materials reasonably accurately.

3.1. Modeling 1^st degree disorder in materials

A wide variety of semiconductors, metals, and alloys are used for the development of new electrochemical storage and conversion device technologies. In all of them, following three types of point defects can be found in different concentrations depending on the crystal structure of the materials and their electrochemical properties:

1. Vacancies, where an atom is missing from its regular crystal lattice.

2. Substitutional, where an atom occupies a site that does not correspond to it usually, also named as antisite.

3. Interstitial, where an atom from the lattice (self-interstitial) or impurities (interstitial impurity) occupies a lattice position that usually is unoccupied.

When a point defect is created, the potential energy and the configuration of the lattice changes, and knowledge of this change becomes necessary for determining the material behavior. To model point defects in materials, one of the most relevant and extended approaches is DFT, which is based on quantum mechanics. DFT has emerged as the standard worldwide first-principles tool for modeling point defects in materials at the ground state. DFT allows the description of many-body electronic structure in terms of a single particle equation, giving accurate results for defect formation energy and electrochemical properties. Nevertheless, the exact exchange and correlation functionals are not known. Therefore, several approximations are used to solve the Kohn-Sham equations. The detailed information of the theory behind first principle calculations used to study point defects, its challenges, and the physical properties calculated for different kinds of materials (metals, semiconductors, or alloys) can be found in the book of Drabold&Estreicher [19] and the reviews of Freysoldt et al. [20] and Spitaler et al. [21]

In the last decades, the development of an efficient approximation to solve Kohn-Sham equations [48] has been one of the main issues for the theoretical community. Although the local functional of density (LDA) [49] and the generalized gradient approximation (GGA) [50] fulfill their function for a wide range of studies, they tend to fail in the description of electronic properties related to the bandgap. Hybrid functionals, DFT with Hubbards correction (DFT+U), and quasiparticle Green’s function (GF) method are some of the latest advances to overcome this deficiency. Also, Quantum Monte Carlo (QMC) methods represent an alternative to DFT for the study of defect formation energies [51]. That being said, the aim of QMC is not to replace but to complement DFT calculations offering a higher level of accuracy at a higher cost [52].

Furthermore, two different models are usually used to study point defects in solids via DFT: cluster model [53,54] and periodic supercells [55–57]. Cluster models represent the perturbed local structure due to the point defects. A representative defect and its close surrounding are selected and treated explicitly through molecular techniques [58]. These techniques are well known for the description of the properties, such as vibrational modes of local defects. When the cluster model is constructed for plane waves basis set calculations, the size of the vacuum in the supercell should be large enough to avoid the cluster–cluster interaction that causes artifacts associated with the chosen boundary conditions [19,20,58]. These artifacts are mostly created by the interactions of overlapping wave functions, which are elastic, magnetic, and electrostatic [20]. On the other hand, when the cluster model is constructed for localized atomic orbitals basis set calculations, these artifacts disappear because the periodic background is lost. The chemistry of the distortion remains; therefore, it is possible to study the electronic properties of the selected area [53,59]. However, the band structure is lost, and related information is not possible to be obtained through this method. For representing surfaces, in either case, specific rules must be followed when clusters are used. It is necessary to use saturators in order to deal with the dangling bonds of the surface atoms that have a negative effect on the entire model. Commonly hydrogens are used as saturators as they work as an electron supplier and logical termination that reduce the adverse effects mentioned [60,61]. However, the interaction of the hydrogen edges can provide misleading information about the electron distribution or conduction, creating a polarized surface because of the different electron affinity between crystal areas and the termination edges. This phenomenon impedes the proper description of charged defects when electrons are added or removed in semiconductors [21]. Figure 2 represents an example of such a cluster model used to represent the bulk and the surface of the experimentally prepared In₂O₃ nanoparticles [62].

Figure 2. The absorption spectra of (a) bulk In₂O₃ obtained from cluster model, (b) bulk In₂O₃ obtained from periodic model (c) (111) In₂O₃ surface obtained from cluster model and (c) experimental samples obtained from steady-state experiments. Insets of (a), (b), (c) and (d) show the cluster model of bulk In₂O₃, periodic model of bulk In₂O₃, cluster model of (111) In₂O₃ surface and HR-TEM image of experimental In₂O_3-x(OH)_y sample, respectively. In, O and H atoms are represented by blue, red, and white spheres, respectively. Adapted with permission from reference [61], Copyright 2016 American Chemical Society

Periodic supercells are the best option for describing defects in either 3D or 2D (slabs) models. Periodic supercell models consist of a host material containing a defect surrounded by dozen to hundreds of atoms repeated through space periodically [55]. Due to the development of much more powerful computers, the high cost necessary to simulate periodic models is now affordable, and most studies opt for this choice. One has to consider that the periodic boundary conditions convert the study of point defects into the study of an array of periodically repeated point defects. The repetition of the defect increases the concentration of them in the material overcoming the usual amount. Thus, artifacts created by the defect–defect interactions affect the calculations’ results. As in the cluster model approach, the effect of these artifacts can be reduced by increasing the supercell size. However, when point defects are studied in materials such as nanoparticles, nanotubes, surfaces, etc. using periodic boundary conditions, it is necessary to construct the supercell model with broad enough vacuum regions and thick enough material. The reason behind this is to avoid the defect–surface interaction from adjacent repeated cells or between both surfaces (inside the same slab model) [63]. The convergence of the adequate vacuum region size must be carried out carefully, always considering the critical effect that an increment of vacuum space will have on the calculation time and cost. The supercell size required to cancel these artifacts would be too large [64,65]; hence, when the computational power was not as high as currently is, cluster models with localized atomic basis set were popular approaches. A different solution to solve the long-term interaction between defects in adjacent cells could be the post-treatment [64,66] of the calculated results, as provided by Varvenne et al. [67]This reported correction can be used for any ab initio study of point defects, including charged defects. The proposed approach only requires to know the elastic constant of the perfect crystal.

Furthermore, to study the electronic properties of semiconductors, charged defects must be considered by introducing a compensating background to avoid electrostatic divergence [68]. However, as previously mentioned, in addition to the artifacts created by the finite-size, the GGA frequently underestimates the bandgap. Hence, post-DFT methods, such as DFT+U, and hybrid DFT are good choices [69,70]. Recently, several groups have made the effort of creating tools to automate the setup and post-analysis of charged point defects in semiconductors and insulators [71–73]. The latest one reported by Broberg et al. consists of a Python Charged Defects Toolkit (PyCDT) [74] in collaboration with the Material Projects database [75].

To simplify the study of point defects for a broad range of materials, various techniques used in conjunction with the DFT method have been reported over the past decade. The most significant ones are GF techniques and Machine Learning (ML). The GF approach serves as a valuable companion to quantum mechanics to solve a complicated system. The system is represented by features such as the geometry, type of material, loading condition, etc. Then, for solving a problem related to the system, the GF is derived by numerically handling the changes in these features [76]. The GF techniques for the study of defects in materials were used initially as an approximation to provide a reasonable solution with lower demand for computational power [77–81]. However, this technique evolved with the development of computer technology, overcoming the challenges faced by the original GF technique from a programming point of view [20]. Korzhavyi et al. [82] used a combination of ab initio calculations and locally self-consistent Green’s function (LSGF) techniques for the first time to study the vacancy formation energy and the lattice volume change. They explicitly treated the region containing the defect considering the lattice crystal as a finite space and employed the first-principles to calculate the physical properties of the supercell. Then, with the LSGF method, formation energies were calculated. In this study, 30 different lattices close to 30 different metal compositions were studied, where the Green’s function eliminated the vacancy–vacancy interaction reducing the computational effort considerably. Further, Yang et al. [76] reformulated the multiscale Green’s function (MSGF) method for modeling point defects in nanostructures based on Tewary’s previous work using LSFG [83]. The methodology utilizes a hybrid LSFG, which is obtained by solving the boundary problem under a Continuum Green’s Function (CGF). It can be applied to describe a point defect and a free surface in a semi-infinite solid, covering the range between the atomistic nanoscale to the continuum macroscale. Thus, making it possible to relate the physical properties at the atomistic level to the observable macroscopic measurements [76,83].

DFT works well for small unit cells; however, when large cells are treated, the computational power required increases substantially. Furthermore, different approaches and modeling have to be done for every new material that needs to be studied. The combination of high throughput DFT techniques with ML provides the solution. By combining DFT with ML, one can cover a wide range of compounds, enabling the acceleration of materials design simultaneously, reducing the cost, time, and effort involved. ML is the branch of artificial intelligence that attends to develop algorithms that use available data to carry out tasks or predictions accurately. In other words, the ML algorithms use the provided data as valuable knowledge to predict the results; then, those results are used to further develop the method. Hence, ML offers the possibility to predict materials properties using existing data taken from the DFT calculation; then, the predictions provide more data that can be used to improve the ML models [84]. Concerning the study of point defects in energy materials, Medasani et al. reported, for the first time, a combinative methodology for the study of cubic B2 intermetallic structures [85]. First, they performed DFT thermodynamic calculation using periodic supercells for a range of 100 intermetallic compounds containing four types of defects. The setup and post-analysis were done by using PyCDT method mentioned previously [74]. Comparing the energies with the optimized perfect crystal structure obtained from the Materials Project database, the dominant defect for each compound was determined. Later on, an ML model based on decision trees were used over the calculated range of structures that worked as a database. The concluding methodology resulted in a fit accuracy of 98% and a prediction accuracy of 75%. After realizing the power of this predictive technique, the ML model was applied to classify the defects in 846 B2 compounds. This approach, if extended to semiconductors, insulators, and multicomponent materials, could result in a robust methodology for the study of point defects due to the almost inexistent computational cost and the self-improving capabilities. However, as the automating procedures for material design are not based on physical principles but on experimental and theoretical databases, they should be carried out with great care and caution, particularly for studying difficult problems such as associated with defected materials.

Overall, point defects are the most common defects present in materials, as well as the easiest to model. There exist several methods to analyze the point defects in energy materials. Depending on the available computational resources and problems at hand, one can choose the appropriate technique. That being said, we still need to develop better approximations that serve as the standard formalism for the study of point defects in any material. Hence, we will not have to rely on computationally expensive post-DFT methods to adjust the mistreatment of long-range interactions, and the inadequate description of bandgaps and charged point defects. Together, the improvement in computational power and automating procedures can lead to significant breakthroughs for the study of 1^st degree of disordered materials.

3.2. Modeling 2^nd degree disorder in materials

2D defects such as interfaces (i.e two grains with zero orientation) and grain boundaries (i.e two grains with non-zero orientations) are the exigent part of solids used in energy technologies. Additionally, 1D defects such as dislocations (defects originating in crystal structure from the displacement jump) and disclinations (defects originating in crystal structure from the rotation incompatibility) arise naturally in the idealized description of interfaces and grain boundaries from a microscopic view. There have been several attempts to analyze these complex defects in the past decades in order to enhance material performances. In particular, for designing novel materials for batteries, fuel cells, and catalysis, the understanding of the transport of ions, electrons, and phonons across 1D-defects (dislocations, disclinations) and 2D-defects (interfaces, grain boundaries, etc.) is crucial.

When two grains meet at a boundary at a certain angle, the atoms at the boundary rearrange to obtain the lowest energy structure between the two crystal lattice directions. This results in a disordered region of a few atomic layers between the two grains, the magnitude of which depends on the angle between the two lattices. Furthermore, the meeting of two grains not only result in disorder due to GB but also due to other defects such as vacancies, dislocations, etc. which arise due to the stress and strains created by different lattices of the two grains (Figure 3). It is therefore expected that any phenomenon (such as ion transport), which depends upon the degree of disorder between the two grains would change based on the angle between them. As there can be an infinite number of GB configurations, depending upon the orientation between the adjacent planes, the relative position of two crystals, and resulting local rearrangement of the atoms near the GB, efficient search strategies in atomistic modeling are needed to generate a realistic material with GB structure along with other defects.

Figure 3. CeO₂ and Yttria Stabilized Zirconia (YSZ) interface system is commonly used as a SOFC electrolyte. a) Morphology of the film is shown by the STEM image. b) Misfit dislocation was found at the interface (indicated with a white arrow). (a) and (b) are adapted from reference [37]: Song, K.; Schmid, H.; Srot, V.; Gilardi, E.; Gregori, G.; Du, K.; Maier, J. and van Aken, P. A., Cerium reduction at the interface between ceria and yttria-stabilized zirconia and implications for interfacial oxygen non-stoichiometry, APL Mater, 2, 10.1063/1.4867556, 2014; licensed under a Creative Commons Attribution (CC BY) license. c) A computational model of CeO₂/YSZ interface (prepared by Amorphization and Recrystallization method), containing 0D, 1D and 2D defects [100]. In (c) O^2-, Ce⁴⁺, Zr⁴⁺, Ce³⁺, and Y³⁺ are represented by red, white, light yellow, bright yellow, and purple spheres, respectively. Adapted with permission from K. K. Ghuman, E. Gilardi, D. Pergolesi, J. A. Kilner, T. K. Lippert, J. Phys. Chem. C 2020, DOI 10.1021/acs.jpcc.0c01797. Copyright (2020) American Chemical Society

There exist several techniques to generate GBs. The easiest one is to use just the geometrical tools which can generate simple GBs. In this method, two crystals are overlayed with a specified rotation angle, and the atoms that don’t fit in the new coinciding lattice formed can be considered as excess, and are deleted manually [86–88]. Another set of methods includes mimicking the natural process of GB formation. These methods start with two crystal slabs with a certain angle between them and separated by randomly placed atoms imitating the liquid phase. This is followed by cooling the system using the Molecular Dynamics (MD) technique until the solid GB is formed [89]. One of the major disadvantages of these methods is that they can generate spurious results if the system gets trapped in local minima while relaxation. The third set of methods generate a large number of possible initial GB configurations and select the one which is most stable [90,91]. The last set of methods include those that perform ‘smart’ search, followed by a global structure optimization via Genetic Algorithm targeting the prediction of GB structures [92,93].

Due to the complexity involved in the materials having GBs, all the above-mentioned computational methods focus solely on modeling GBs in the materials. Other defects such as dislocations/disclinations, point defects, nanopores, etc. that exist along with the grain boundaries and result due to the strain/stress released during GB formation are ignored. As in real materials, GBs co-exist with 0D, 1D, 2D, and 3D defects, from hereon we focus on elaborating only those methods that are capable of designing 2^nd order disorder in materials with multiple defects.

The ideal way of capturing all the microstructural details of experimental samples is to simulate the crystallization process itself. Indeed, there are many elegant studies that explore the crystallization process via expensive computation [94,95]. However, an outstanding approach that stands out is provided by Dean Syale [96,97]. In this approach, known as ‘Amorphization and Recrystallisation (A&R),’ the two crystals forming the thin film are first forced to undergo a transformation to an amorphous state, followed by the recrystallization process. One of the advantages of the A&R method is that the thin film structure prepared by this method is determined by the recrystallization process rather than the initial structure. Furthermore, as recrystallization is started from the high-energy amorphous structure, this methodology allows a thorough exploration of the configurational space resulting in thin films with realistic distribution of point defects and line defects.

Furthermore, A&R technique is sometimes utilized with Near-coincidence-site lattice (NCSL) theory to model realistic bi-layer materials, such as Gd-doped Ceria (GDC)/Yttria stabilized Zirconia (YSZ) used as SOFC electrolytes, by mimicking experimental samples closely with all the microstructural defects. Sayle et al. described the basic principles of NCSL theory by applying it to cubic and hexagonal systems (BaO-MgO and CeO₂-Al₂O₃) [98,99]. Fisher and Matsubara [100] later extended the theory to rectangular systems (NiO-ZrO₂). Ghuman et al. [101] has recently applied this technique, with slight modifications, to design realistic interfaces present in CeO₂/YSZ electrolyte, commonly used in SOFCs. In brief, in this method first, a large supercell of an appropriate surface interface is obtained by aligning near-coincidence lattice points by satisfying the following relationships: A1a1≈ A2a2 and B1b1≈ B2b2, where a1, a2, b1, and b2 are the surface lattice parameters and A1, A2, B1, and B2 are adjustable parameters. Because the lattice parameters do not have the same dimensions, they do not match perfectly. Therefore, an adjustment from one of the lattice or both of them is necessary to obtain a common lattice fit, resulting in a system with strain. The magnitude of strain in the initial structure can be reduced by using large supercell sizes. To define the degree of fit between the two surfaces and indicate the supercell size variables such as ‘misfit percentage’ and ‘planar coincidence density,’ respectively, are used. This initial configuration provided by NCSL is then used to search a lower energy configuration, by increasing the temperature of the system close to the melting point. The MD simulation carried out at a certain temperature allows the slabs to overcome energy barriers that might prevent structure optimization. Subsequently, the system is recrystallized by applying high pressure while maintaining a high temperature. Once the system is recrystallized, the pressure is released and the temperature is reduced gradually until it reaches close to 0 K. Clearly, the initial configurations based on NCSL theory cannot automatically capture the 0D and 1D defects such as dislocations, vacancies, interstitials, etc. present in real materials. Further, it is also not practical to introduce these defects ‘by hand’ as the experimental data regarding the detailed structure of these defects is limited. However, NCSL with the A&R technique can explore the synergy that exists in ‘real’ systems and therefore provides a solution to designing materials with 2^nd order disorder.

Another noteworthy technique used to design GBs, specifically their segregation structures, is the one that uses hybrid molecular dynamics/Monte Carlo (MD/MC) algorithm. As discussed by Sadigh et al. this technique is ideal for simulating multiphase materials having several components, thereby forming large supercells with millions of atoms, due to the use of variance-constrained semi-grand canonical ensemble (VC-SGC-MC) implemented via MC. Further, whether VC-SGC-MC method is combined with MD, one can also analyze the complex phenomenon such as chemical mixing and precipitation, phase segregation, structural relaxations, and thermal vibrations. The starting model for the simulation using this method is prepared by distributing crystalline grains evenly in the supercell with random orientations, and subsequently filling them with the Voronoi volumes around each grain. The simulation itself is carried out by an MD run disrupted every n^th step to execute a certain number of MC trials. The simulation cell achieves the structural as well as chemical equilibrium much faster via this combined MD/MC technique as compared to when only MD is used. It should be noted that this technique is currently applicable to binary alloys only and needs to be tested and adapted for a vast range of other energy materials [102].

Furthermore, another important technique that can model the amorphous nanoporous metals and semiconductors, that are currently sought after for applications such as hydrogen storage, is called ‘expanding lattice method’ [103]. The simulations using this technique is started with a large crystalline supercell having a density close to the real value, which is further lowered during the simulation by increasing its volume. Subsequently, these supercells are subjected to the ab initio or parameterized-Tersoff-based-MD simulations at temperatures below the corresponding bulk melting points. Finally, the supercell is optimized, resulting in the amorphous materials with pores along with the crystallographic directions.

Eventually, while there exist several studies that explore properties of the complex materials having higher-order defects using ML and artificial neural networks [104–107], to our knowledge these techniques are not yet applied for designing realistic computational models of the materials having various higher-order defects, therefore, are excluded from this review.

3.3. Modeling 3^rd degree disorder in materials

Materials with 3^rd degree disorder, i.e. having no long-range order, usually referred to as amorphous materials, have proven useful for various sustainable energy applications. Conceptually an ideal amorphous material, defined as a material having perfect random networks, can be constructed by using certain rules as stated in a review article by Zbigniew et al. [108] Assuming a negligible effect of gravitation, friction, or vibrations, one can create amorphous material by a simple technique where first a sphere, representing an atom, is created with a sufficiently large number of equally spaced points on its surface. From them, a certain number of points are selected at random such that when new spheres are placed at them they do not overlap and only could touch the starting sphere. This step is followed by the successive addition of the new layers of the spheres onto sites formed by the adjacent spheres below. This process can be continued, layer by layer, until the aimed size of the amorphous structure is reached.

However, to construct a realistic amorphous material is not simple. This is because the structure of the amorphous material is not just non-repetitive, but also contain lots of irregularities due to structural defects. Hence, it is not feasible to construct it ‘by hand’ or by using the Bloch theorem as in case of crystalline materials. To represent the amorphous materials appropriately, their supercell should be large enough to represent the extended structure, which is so far impossible. Therefore, currently, the supercells representing amorphous materials are modeled with a much smaller number of atoms representing realistic information of the short and middle-range structural ordering only. These models are more-appropriately referred as ‘disordered’ materials. In crystallography, a disordered material is described as a material with defects relative to its crystalline counterpart [109]. In principle, we should be able to restore the disordered materials to their ordered counterparts by removing the defects present in them. As for disordered materials that have a short and middle-range order, a smaller supercell can give realistic information, and therefore, one can use first-principles computational simulation. It should be noted that amorphous materials are the materials with ‘random’ rather than ‘disordered’ atomic arrangement and not necessarily can be restored to the perfect crystalline state.

The most common and frequently used technique to prepare the disordered materials is the one demonstrated by Car-Parrinello and is commonly called as ‘melting and quenching’ technique [110]. In this technique, the initial configuration on which heating and quenching are performed is usually chosen to be the crystalline counterpart of the amorphous materials, with a rescaled lattice parameter to have a density equal to that of amorphous material. For example, in order to prepare amorphous titania, the initial configuration is usually that of rutile or anatase with scaled parameters to have a density equivalent to that of amorphous titania, i.e., 3.8 gcm⁻³ [111,112]. Then by using a well-tested force field for crystalline counterpart, melting and quenching are performed. For simulations, periodic boundary conditions are used. The crystalline sample is first subjected to geometric optimization in microcanonical ensemble followed by melting it at the temperature greater than its melting point in order to eliminate its initial structure memory [113]. Thereafter, sequential heating and quenching are performed followed by equilibration. Finally, the samples quenched at the room temperature are relaxed at 0 K to attain the ground state. Throughout the sample preparation process, careful checks are performed to make sure that no dangling bonds are present [114]. Finally, the amorphicity of the prepared models is confirmed by comparing their structural properties with the available experimental data. The major drawback of this method is that one cannot obtain the electronic-level information of the prepared amorphous material. Therefore, if the problem at hand needs electronic-level information the final samples obtained from MD simulations at 0 K temperature can be utilized for DFT calculations if they are small enough, i.e., a couple of hundred atoms [115–118]. In Figure 4 we represent a model mimicking amorphous TiO₂/Au nanorods photocatalyst reported by Li et al. [119] Amorphous TiO₂ model is prepared by melting and quenching technique, and is further decorated with gold, either in the dispersed-ion form (Figure 4(e)) or in cluster form (Figure 4(f)) to represent the experimental TiO₂/Au samples.

Figure 4. a) 3D representation of crystalline-TiO₂/Au/amorphous-TiO₂ photoelectrodes. b), c), and d) represent the low-magnification TEM, SEM and HRTEM images of TiO₂/Au/a-TiO₂, respectively. a), b), c), and d) are adapted with permission from the references [118], Copyright 2018 Wiley-VCH. e) and f) represent the models mimicking the part of the crystalline-TiO₂/Au/amorphous-TiO₂ photoelectrode. Amorphous TiO₂ is prepared by melting and quenching technique [115] and Au metal is represented e) as ions atomically dispersed on amorphous supports, or f) as metal clusters interacting with amorphous TiO₂ supports. The dotted box represents the core C/A interface area

Another approach called the ‘undermelt-quench’ approach is reported by Valladares et al. [120] recently to design disordered materials [121,122]. In this approach, first, a supercell of the crystalline counterpart of the amorphous material is constructed with the same density as the amorphous phase. In the first simulation step, the temperature of the supercell is increased just below the melting temperature of the material in a small number of steps with each time step three to four times the default value. This heating step is followed by an instantaneous cooling of the supercell to 0 K carried in the number of steps same as used for heating the system. This results in the formation of the amorphous material due to the random atomic diffusion during the melting stage. After the quenching step, the stress generated in the supercell is released by subjecting it to several annealing cycles at the temperature determined by the experiment. This is followed by the geometry optimization which provides the final amorphous structure having structure only partially resembling the local structure of the disordered material. Further, it should be noted that the difference between the ‘melting and quenching’ approach and the ‘undermelt-quench’ approach resides in the heating procedure. In the undermelt approach, the crystalline supercells are not melted but are heated just below the melting temperature, which helps in avoiding the over-coordination of the atoms of the supercell.

To overcome the timescales limitation of MD, and the problem of freezing in too much liquid-like character via melting and quenching technique, a technique called ‘decorate and relax’ was introduced by Tafen et al. [123] The peculiarity of this approach is that the starting structure used in this method includes basic a priori information about the chemical order and the coordination number of the amorphous material under preparation, which advances the simulation in the right configuration space. Further, it has been reported that the peaks of the static structure factor of the starting structure resemble the first sharp diffraction peak of glasses indicating that the starting models considered in this approach already exhibit the intermediate-range order, resulting in a much faster method to prepare disordered material.

In contrast to the ‘direct’ methods that use an interatomic potential to perform the simulations to find the disordered structure that is a local minimum of the energy functional, there also exist ‘inverse’ or ‘reverse’ MC [124–126] methods which do not need any interaction potential for simulations. Instead, these methods use structural information from the pair correlation function, coordination number, and the local bonding environment, obtained from the experiments (such as x-ray and neutron diffraction) to construct the suitable starting model. Thus, to start the simulations reverse MC methods use a configuration that satisfies a certain set of constraints appropriate for the material under study. During the simulations, atoms of the starting model are randomly moved under the periodic boundary condition until one obtains a structure that satisfies the input experimental data and the pre-defined constraints. Since this method is used to obtain a material configuration satisfying the desired experimental properties, it is called ‘inverse’ or ‘reverse’ MC. The biggest challenge in using this method is the limited availability of the (structural, electronic, and spectroscopic) experimental data for the material under study.

Despite the significant advances, in the ‘direct’ and the ‘inverse’ approaches, it is difficult to realize realistic models for amorphous materials. With inverse approaches, only indirect statistical knowledge can be gained about the local atomic structure of the material. While direct approaches either describes a system with electronic-scale accuracy (via DFT) which allows only limited system sizes and time scales or with the atomic-scale accuracy (via MD/classical force field), which although require much less computational effort and can provide structural information at the nanoscale, are rarely enough to fully describe the amorphous material structure.

To overcome the limitations of these decades-old approaches a solution is recently provided by a data-based ML approach [127]. One of the areas where the ML approach is proving very useful is in the generation of the interatomic potentials for atomistic simulations of amorphous materials [128–132]. In order to generate the interatomic potentials, the ML approach uses the available datasets of atomic forces and energies generated by the quantum-mechanical approaches. ML-based interatomic potentials are reported to enable simulations with an accuracy comparable to DFT, but with a much smaller computational cost. Amorphous carbon is one of the many reported materials for which ML potential has been reported [133]. This potential is based on the Gaussian approximation potential (GAP) framework [130] and the Smooth Overlap of Atomic Positions (SOAP) atomic similarity kernel [134]. This potential has been reported to captures not only the complicated structural, mechanical, and surface properties of amorphous carbon but also its growth mechanism [127,133]. Amorphous LixSi phases which are used in batteries are also recently simulated via neural network potentials [135,136]. Although ML-based interatomic potentials have resulted in unprecedented accuracy in the models of the structure of the amorphous materials, they are also critically dependent on the data available from the quantum-mechanical calculations. As a result, we are still far away from having a readily available database of ML-based potentials for designing commonly used energy materials with the structural disorder.

4. Future perspectives

In this review, we highlighted the vital role that disorder play in energy technologies and the current strategies to model it at the nano- and electronic-scales. We specifically elaborated those techniques that model realistic materials by taking the entire complexity of the system into account, rather than just focusing on the partial material characteristics that originated due to a single type of defect.

Despite all the progress in implementing disordered materials in energy applications, it is clear that a thorough understanding of the fundamental mechanisms associated with the disorder in materials is still lacking. Figure 5 represents the overview of the challenges (i.e. bridging the gaps between theoretical predictions and experimental observations or between techniques used at different scales) that exists for utilizing disorder for sustainable energy technologies, thereby, emphasizes the need of the development of a robust multi-scale computational framework through rigorous validation at each scale with carefully designed experiments.

Figure 5. Representation of the current challenges (i.e. bridging the gaps between models at different scales, and between theoretical predictions and experimental observations- highlighted with yellow arrows) for designing realistic materials for energy applications. The yellow box illustrate the schematic design of an aspirational high-performance energy material having non-symmetrical grains, with each grain highlighting the potential applications of materials with disorder in the sustainable energy sector. Abbreviations: MD: Molecular dynamics, DFT: Density Functional Theory, MC: Monte Carlo, TEM: Transmission Electron Microscopy, PLD: Pulse Laser Deposition, HR-STEM: High-Resolution Scanning Transmission Electron Microscopy, OER: Oxygen Evolution Reaction, ORR: Oxygen Reduction Reaction, HER: Hydrogen Evolution Reaction. The top left (computational image) is adapted with permission from K. K. Ghuman, E. Gilardi, D. Pergolesi, J. A. Kilner, T. K. Lippert, J. Phys. Chem. C 2020, DOI 10.1021/acs.jpcc.0c01797. Copyright (2020) American Chemical Society. The bottom left (experimental image) is adapted from ref [37]. Song, K.; Schmid, H.; Srot, V.; Gilardi, E.; Gregori, G.; Du, K.; Maier, J. and van Aken, P. A., Cerium reduction at the interface between ceria and yttria-stabilised zirconia and implications for interfacial oxygen non-stoichiometry, APL Mater, 2, 10.1063/1.4867556, 2014; licensed under a Creative Commons Attribution (CC BY) license

Clearly, new techniques are needed to understand the synergetic effect of multiple defects usually present in materials prepared in laboratories. Specifically, new strategies are required for understanding and manipulating disorder at the electronic level. Most of the strategies for designing realistic disordered material use atomic-level techniques such as MD and MC. However, to understand the underlying mechanisms we must find a way to use realistic material models prepared by MD/MC for electronic-level analysis via techniques such as DFT. Currently, the large sizes of the models prepared via MD/MC restrict their use in DFT simulations. The only way to understand electronic-level details of the nanometer-sized model is by analyzing a part of it at a time, which is obviously not the most robust way of analyzing materials. The development of ML techniques that enables the simulation of macroscale models with an accuracy close to DFT, in our opinion, could be one of the most promising technologies to overcome the scale bridging necessary for understanding the disordered materials at electronic scale as well as designing macroscopic materials. Finally, to make significant advances, the role of disorder in materials predicted by computational methods should be confirmed with sophisticated dynamic experimental measurements.

We hope that by highlighting the current computational methods used to design disorder in materials will inspire new opportunities for cross-disciplinary research and development of advanced tools for understanding and thereby utilizing disorder in materials for sustainable energy technologies.

Acknowledgments

The authors acknowledge the financial support from the Natural Sciences and Engineering Research Council of Canada (NSERC) Discovery grant program, [Funding Reference number RGPIN-2020-05924].

Disclosure statement

The authors declare no conflict of interest.

Additional information

Funding

This work was supported by the Natural Sciences and Engineering Research Council of Canada (NSERC).