Synergy of Binary Substitutions for Improving the Cycle Performance in LiNiO2 Revealed by Ab Initio Materials Informatics.

We explore LiNiO2-based cathode materials with two-element substitutions by an ab initio simulation-based materials informatics (AIMI) approach. According to our previous study, a higher cycle performance strongly correlates with less structural change during the charge-discharge cycles; the latter can be used for evaluating the former. However, if we target the full substitution space, full simulations are infeasible even for all binary combinations. To circumvent such an exhaustive search, we rely on Bayesian optimization. Actually, by searching only 4% of all of the combinations, our AIMI approach discovered two promising combinations, Cr-Mg and Cr-Re, whereas each atom itself never improved the performance. We conclude that the synergy never emerges from a common strategy restricted to combinations of "good" elements that individually improve the performance. In addition, we propose a guideline for the binary substitutions by elucidating the mechanism of the crystal structure change.

Introduction

Lithium-ion batteries (LIBs) have a high-voltage and high-capacity and hence are widely used as secondary batteries for mobile devices and hybrid/electric cars.[1,2] The cathode material in LIBs is one of the most important factors that determine the battery performance. Thus, its development has recently attracted much attention from the viewpoint of industrial applications. A practical solution for improving the performance is called “atomic substitution”.[3−18] For example, consider LiNiO2 (LNO),[19−54] which is known to have a higher capacity and lower cost than LiCoO2 (LCO). However, LNO has fewer cycle characteristics than LCO. Indeed, the substitution of Ni sites in LNO with Co and Al, Li(NiCoAl)O2, a commercially used cathode material, has prolonged the cycle life compared to LNO itself. Note that Co and Al are, respectively, known to improve the rate performance and thermal stability.[22−25] This indicates that “individually good” elements come together to tune the corresponding characteristics. This is a commonly used strategy for improving the battery performances. A situation may occur, however, where a combination of “no-good” elements has the potential of improving the performance. However, it is really unclear if such a synergistic effect emerges. Unfortunately, however, we cannot straightforwardly establish an experimental verification even for binary substitutions having a huge number of possible combinations. This is because such an exhaustive search based on experiments requires enormous amounts of time and high costs for the synthesis of candidate materials.

One of the most promising solutions to the above exhaustive search problem is materials informatics (MI),[26−28] a recently emerging paradigm in materials science combined with information and data science. There are a number of successful MI studies that have explored new materials with desirable properties.[29−35] Since a sufficiently large amount of experimental data on material properties is generally unavailable unlike in other research areas (e.g., bioinformatics), computational approaches are useful for generating the data. In most cases regarding electronic properties, ab initio simulations can satisfactorily generate data to construct machine learning-based prediction models, thereby tackling the exhaustive search problem. This may be called “ab initio materials informatics” (AIMI), although computationally proposed candidates should be verified experimentally. In battery materials, AIMI approaches have been used successfully to explore cathode coating materials,[33] and new cathodes with better capacity and thermal stability.[34]

Focusing on the cycle performance of LIBs, an AIMI approach based on a high-throughput screening combined with density functional theory (DFT) calculations has been applied to co-substituted LiFePO4 cathode materials, and it successfully found element combinations that could prolong the cycle characteristics.[35] The volume change during the charge–discharge cycle and planer mismatch adopted in their study worked as the screening criterion. This success can be attributed to the fact that capacity fading is caused by microcracks during the charge–discharge cycles,[36] although ab initio simulations cannot directly evaluate the cycle performance. Very recently, we have investigated the unary substitutions that improve the LNO cathode by computational screening.[37] Our findings are as follows: (i) DFT simulations must incorporate van der Waals (vdW) corrections in exchange–correlation functionals adopted to describe the structural changes during electrochemical cycling and (ii) our descriptor analysis based on sparse modeling elucidates important descriptors correlating with the contraction. Finally, our AIMI approach indicated that Nb is the most promising candidate for the substitute element.

Our previous study provides the basis for the further exploration of new cathode materials beyond unary substitutions. As we mentioned before, however, binary substitutions have a much larger search space, so we employ a Bayesian optimization technique[38,39] to efficiently discover the best binary combinations from the huge search space. In this study, we efficiently found synergetic binary substitutions by computing only about 4% of the search space, and we elucidate the mechanism of structural changes in a MI manner, which will be helpful for further explorations by AIMI approaches.

Results and Discussion

According to our previous study, the cycle performance is dominated by the c-axis contraction, Δdave (for detail, see the Methodology section). Δdave values for all of the unary substitutions (65 elements) are summarized in the Supporting Information. Their comparisons with LNO value (Δdave = 0.156) classify the binary substitutions into “positive element” and “negative element” substitutions depending on whether or not they suppress a structural change, i.e., the positive/negative has a Δdave value less/greater than 0.156. Furthermore, the negative elements can be divided into the following two groups: those that are dissolved into the Ni site and the others that are moved from the Ni site to the Li site by structural relaxation (see Table ).

Table 1

Classification of the Presence/Absence of the Δdave Suppression Effect in Unary Substitutionsa

positive element	negative element
V, Nb, Ge, Ir, Au	Pb, Pa, Hf, Rh, Fe, Cu
Al, Ti, Ta, Mg, Ga	Po, W, Zn, Cr, Zr, Pd
Ru, Bi, Dy, Os, Sb	Mo, Ca, Yb, Cd, Ag
Mn, Tc, Y, Tb, Re	Gd*, Ce*, Eu*, Tl*, Pm*, Sm*
Tm, Lu, Pt, Ho, In	Hg*, Na*, Th*, Nd*, Pr*, Sr*
Er, Sc, Sn	La*, Ac*, Ba*, Ra*, K*, Rb*
	Fr*, Cs*

In the negative element, the element marked with * moves from the Ni site to the Li site.

V is the best dopant among the unary substitutions, having Δdave = 0.124. We note that Zr and Na are classified into negative elements by the ab initio simulations, but these elements are known to improve the cycle performance as follows: Zr stabilizes the cation ordering of Li and Ni;[16] Na stabilizes solid solutions at the Li site (pillar effect).[18] These effects are beyond the scope of the present simulations, so we no longer consider these elements, which prevents us from obtaining inconsistent results.

As mentioned above, V is the best element for unary substitutions. Here, we investigate if binary substitutions further suppress Δdave using Bayesian optimization.[38,39] First, we randomly select 31 candidates and evaluate their Δdave values by ab initio DFT simulations. Next, using the descriptors shown in Table , we construct a prediction model of Δdave based on the Gaussian process regression learned from the data set. We finally conduct Bayesian optimization for exploring the best binary combinations.

Table 2

List of Descriptors Considered in the Present Study

elemental information	crystal structure information
atomic number Z	lattice volume V
atomic mass m	lattice constant a, b, c
electronegativity EN	angle formed by basic vector α, β, γ
covalent radius CR	substituted positions of X1 and X2x, y, z
first-ionization energy IE
maximum oxidation state O
coefficient of van der Waals C6

Figure shows how quickly the Bayesian optimization and random search found the smallest Δdave. For Bayesian optimization, we considered two sets of descriptors including all and the “best” selected ones (explained later). The average number of observations (Nave) required for finding the optimized solution using Bayesian optimization with all of the descriptors and the random search were 13 and 15, respectively. As can be seen in Figure , the Bayesian optimization rapidly reaches the smallest Δdave as compared with the random search, but its efficiency is not so high. In particular, when the number of observations is less than 15, the success probability of Bayesian optimization with a set of all of the descriptors is comparable to that of random sampling. This is probably because the set of all of the descriptors shown in Table is inefficient for predicting Δdave accurately. Therefore, as explained below, we conducted the descriptor analysis by LASSO regression[40] and then extracted the most important ones in terms of the Δdve prediction; we anticipate that this reduced set of descriptors can improve the efficiency of Bayesian optimization.

Figure 1

Success probability for the Bayesian optimization with all of the descriptors (blue) and the best descriptors (red). Random sampling (black) is also shown for comparison.

In our descriptor analysis, we first divide the data into training/test sets with the ratio of 8/2. The prediction model of Δdave is based on the LASSO regression learned from the training set. With the use of this trained model, we predict Δdave values for the test set and compare them to the corresponding ab initio values. Its prediction accuracy is evaluated in terms of the mean square error. The LASSO is known to automatically select important descriptors for prediction. As a result, we obtained six descriptors (V, α, m × m, EN × EN, a, z × z), and using them, Bayesian optimization was implemented again, which is denoted by the “best descriptors”. As expected, Nave was 10 for the best descriptors, and the success probability improved (see Figure ). We thus expect that these descriptors are efficient even with further searches, and then, we proceed with the Bayesian optimization for the remaining 1626 candidates.

Thirty-nine observations for the Bayesian optimization with the best descriptors were carried out after 31 observations for random sampling, where their search history is shown in Figure . We show all of the numerical data obtained in the present study in the Supporting Information. Meanwhile, the random sampling (observations 1–31) chooses Δdave values distributed from large to small values, and the Bayesian optimization (32–70) selects smaller Δdave values (N.B., our acquisition function is the maximum probability of improvement; see the Methodology section).

Figure 2

Search history for the binary substitutions. Black-dashed line shows the smallest Δdave for the unary substitutions.

From the above observations by the Bayesian optimization, we discovered V–Ga, Ti–Fe, Al–Mn, Al–Cr, Mg–Cr, and Cr–Re as the combinations that suppress Δdave better than the best unary substitution. Although V–Ga and Al–Mn are composed of positive elements, the other combinations include negative elements, Fe or Cr. Nevertheless, for the latter, their Δdave values showed more suppression as compared to the best unary substitution, which can be interpreted as being synergistic effects.

To elucidate why the synergy emerges, the suppression mechanism in the binary substitutions is investigated here. We construct a regression model of Δdave predictions learned from 27 combinations where dopants do not move to the Li sites. Changes in interlayer distances d during cycling are anticipated to be determined by the “Coulomb repulsive interaction between oxygen in the NiO6 layer”, “vdW attractive interactions between the NiO6 layers”, “ionic radii of doping cations”, etc. Accordingly, the following quantities are considered as our descriptors entering the regression model: “a product of the averaged Bader charges[41,42] over the facing oxygen layers”, “a sum of the averaged vdW coefficients over the facing Ni–O slabs”, “a sum of the averaged ionic radii over the facing Ni layers”, and “the difference in the averaged Bader charges within the oxygen layer between charged and discharged states”. In total, we obtained 10 descriptors for the Δdave prediction. At the charge rate used in this study, Ni–O slabs can be classified into the following two cases (see Figure a): (A) Ni–O slabs sandwiching a Li-discharged layer and (B) Ni–O slabs sandwiching a Li-charged layer. The above descriptors are calculated for each case. Note that two independent O layers exist in the (B) Ni–O slab, and the corresponding descriptors distinguished by their magnitudes are treated independently.

Figure 3

(a) Descriptor classification. (A) Ni–O slabs sandwiching a Li-discharged layer and (B) Ni–O slabs sandwiching a Li-charged layer. (b) Heat map of regression coefficients arranged according to the models’ coefficient of determination R2. The largest R2 was 0.77. From left, “a product of Badar charge on oxygen (A, B small, B large)”, “a sum of vdW coefficients in the Ni–O slab (A, B small, B large)”, “a sum of ionic radii of cations (A, B small, B large)”, and “a difference in Bader charge on oxygen between charged and discharge states”. See text for more details.

We consider all possible combinations of 10 descriptors to construct the regression models. Figure b shows their regression coefficients in a heat map form. Looking at models with higher R2 (coefficient of determination) values, it can be seen that the charged product of the oxygen layer (A)/the sum of the vdW coefficient has a negative/positive contribution. This is because the larger Coulomb repulsion and smaller vdW attraction suppress the shrinkage in the c-axis direction. These two descriptors are the most important factors describing the change along the c-axis. In other words, to improve the cycle characteristics, it is preferable to substitute Ni atoms with those that have a low vdW coefficient, holding large Bader charges on oxygen layers. We also found that a large change in the Bader charge on oxygen accompanying the charge–discharge cycle significantly contributes to the suppression of structural change. This is because a smaller change in the electronic structure during the charge–discharge cycle suppresses the structural change. Although the ionic radius of the cation is less important than the product of the oxygen charge and the sum of the vdW coefficient, it is required for constructing an accurate model.

It should be noted that cations substituted with Ni require a small ionic radius in (A) and inversely a large one in (B). Consequently, any unary substitution cannot satisfy both the requirements simultaneously; hence, its Δdave suppression is smaller than that of the binary substitution. In other words, these inverse requirements are the origins of the synergistic effects. Indeed, Cr–Re (Δdave = 0.08) shows the most suppression of Δdave change, where a large Cr cation and a small Re cation are doped into (B) and (A), respectively. It is concluded that understanding why synergistic effects occur is important to design the best performance material.

Conclusions

We investigated unary and binary substitutions in LiNiO2 to improve its cycle performance. Owing to their strong correlation, the cycle characteristics were evaluated from the perspective of structural change along the c-axis during cycling, and this was investigated for various doping elements using ab initio simulations. For the unary substitutions, an exhaustive search based on the ab initio simulations revealed that V is the best for unary substitution. As for the binary substitutions, our Bayesian optimization explored the candidate combinations more efficiently than random sampling, and the results revealed synergistic effects for combinations of doping elements that never improved the cycle performance individually (Cr–Mg and Cr–Re). We also proposed guidelines for battery material design. Our ab initio materials informatics approach presented here can be applied to other cathode materials such as LiPFeO4 in which the cycle performance is determined by the volume changes,[35] and can be regarded as a promising fundamental technique for further exploration of high-performance battery materials.

Methodology

All our ab initio DFT simulations were carried out using the Vienna ab initio simulation package (VASP)[43,44] with the projector-augmented wave (PAW) potentials. The plane-wave basis set cutoff energy and k-point mesh were set to be 650 eV and 3 × 3 × 1, respectively, throughout this study. We employ the van der Waals exchange–correlation (vdW-XC) functional[45,46] to accurately describe the crystal structural change during the charge–discharge cycle. In general, the DFT + U method is known to reproduce lattice constants compared to experimental values, but we assume U = 0 because a previous study demonstrates that the c-axis contraction can be well-described even by DFT without U, which avoids high computational costs.[37]Figure shows the LNO crystal structure with a rhombohedral R3̅m unit cell.[21] The supercell used in the study is 2 × 2 × 1 of the conventional cell. We assume that doping elements are dissolved in the Ni sites. For unary substitutions, we replace Ni with X (Ni/X = 92:8); for binary substitutions, each of the two Ni is replaced by X1 and X2 (Ni/X1/X2 = 84:8:8), respectively. We may expect that our findings shown in the present study still remain valid even for lower levels of doping, since lattice constants a/c linearly decrease/increase as the doping ratio increases.[23]

Figure 4

Crystal structure of LiNiO2. The supercell consists of four conventional cells. The structural change is measured as the average over interlayer distances, d (see text).

To design a high-capacity and high-cycle performance battery, we must limit the candidates. The valence of Ni changes during the charge–discharge cycle. When nine Li ions are removed from the LNO unit cell (83% charged state), the valence of 10 Ni changes from +3 to +4 and that of two Ni remains at +3. Assuming that only Ni redox could compensate for charge neutrality in doping system, the sum of the valence of X1 and X2 must exceed 3 + 3 = 6. Owing to the above limitation, the number of possible candidates amounts to 1657 when considering the third and lower rows in the periodic table up to Pa (atomic number 91) except for the elements for which pseudopotentials are unavailable. The structural changes are evaluated by the change in the crystal structures between 0 and 83% charged states. Although there are several possibilities for doping positions and a delithiated structure for each substitution, only the most stable structure is chosen by comparing their energies. According to our previous study,[37] the cycle performance strongly correlates with the c-axis contraction, which is measured by the following quantitywhere dd (dc) denotes the interlayer distances for the discharged (charged) state (see Figure ).

We performed ab initio DFT simulations for all of the unary substitutions (65 elements) to evaluate Δdave. In the binary substitutions, all of the possible candidates amounted to ∼105; however, in reality, such a huge number of simulations is infeasible. To efficiently discover the best binary combination(s), we thus adopted Bayesian optimization, which is an experimental design algorithm associated with Gaussian process regression.[38,39] It has been reported that the design and choice of descriptors appropriate for describing material features uniquely are crucial for conducting Bayesian optimization.[29] The interlayer distances, d (i = 1, 2, 3), are determined as a balance between the “Coulomb repulsive interaction between oxygen in the NiO6 layer” and “vdW attractive interactions between the NiO6” layers.[47−52] Although those interactions would be good descriptors, it is quite hard to directly evaluate them. Instead, as listed in Table , we collected “elemental information of dopants” and “structural information of 0% charged state” as descriptor candidates; all of the quantities were symmetrized with respect to the exchange between two dopants, e.g., m + m and m × m for atomic mass m; we eventually selected those quantities as descriptors by random sampling. Based on Gaussian regression with the descriptors selected above, we constructed a prediction model of Δdave and then performed Bayesian optimization under the maximum probability of improvement (as an acquisition function). We employed the Common Bayesian Optimization Library (COMBO)[53] to implement our Bayesian optimization.