The article information
 Dan Wu, Gongdong Chen, Chaoyi Ge, Zhenpeng Hu, Xuehao He, Xingang Li
 吴丹, 陈功东, 葛超逸, 胡振芃, 何学浩, 李新刚
 DFT+U Analysis on Stability of LowIndex Facets in Hexagonal LaCoO_{3} Perovskite:Effect of Co^{3+} Spin States
 Co^{3+}离子自旋态对六方相LaCoO_{3}的低指数晶面稳定性影响的密度泛函理论研究
 Chinese Journal of Chemical Physics, 2017, 30(3): 295302
 化学物理学报, 2017, 30(3): 295302
 http://dx.doi.org/10.1063/16740068/30/cjcp1703035

Article history
 Received on: March 16, 2017
 Accepted on: March 21, 2017
b. Collaborative Innovation Center for Chemical Science & Engineering(Tianjin), Tianjin 300072, China;
c. School of Physics, Nankai University, Tianjin 300071, China;
d. Department of Chemistry, School of Science, Tianjin University, Tianjin 300072, China
In recent decades, perovskite oxides (ABO$_3$) have aroused considerable interests as the cathode materials in solid oxide fuel cells for use in many energy storage and conversion technologies, such as the oxygen evolution reaction (OER) and oxygen reduction reaction (ORR) [15]. These oxides are also catalyticallyactive and costeffective catalysts for many important chemical reactions, such as NO$_x$ (including NO and NO$_2$) storage reduction [6, 7], CO oxidation [8], and methane combustion [9]. Since those catalytic reactions are surface processes, studies on the surface properties of those materials play an important role in deeply understanding the catalytic reactions. Usually, the experimental observations believed that (001) facets, particularly (001) Bterminated surface, in perovskite oxides (ABO$_3$) are stable and play a pivotal role in these catalytic chemical reactions [1014]. However, to the best of our knowledge, the surface stabilities of perovskite oxides, especially for LaCoO$_3$, are still under debate from theoretical calculations. For example, using the cubicphase structure of LaCoO$_3$ in the calculations, Kahn et al. [15] predicted the Coterminated (111) surface and LaCoOterminated (110) surface to be the most stable; Chen et al. [16] reported that the ground termination was the (001) LaO surface, transitioning to the (111) LaO$_3$ facet in oxygenrich condition; while Zhang et al. [17] found that the LaOand CoO$_2$terminated (001) surfaces were stable.
It is noticed that LaCoO$_3$ has a rhombohedral primitive cell with $R\overline{3c}$ symmetry [1823], which belongs to the hexagonal phase, but the cubic phase was used in the previous calculations [1517]. As well known, results from the firstprinciples calculations are very sensitive to the geometric structure of materials. Therefore, the complexity of previous studies may come from the difference on crystal structures of LaCoO$_3$. Recently, Liu et al. explored the stability of lowindex facets on the hexagonalphase LaCoO$_3$ by using the firstprinciples analysis [24], and they found that the (1102) LaO, (1104) O$_2$ and (0001) LaO$_3$ terminations were thermodynamically stable with a nonspin (NS) state. Though a more accurate crystal structure was employed in their calculations, the results were still in conflict with the common sense from the experimentalists that the (001) CoO$_2$ was stable and active in the reactions [1114]. There might be some key factors ignoredin previous calculations, which can reconcile the disagreement between computational result and experimental observation.
Besides the crystal structure, the spin state of Co$^{3+}$ ion is another important factor which can affect the results of calculations, but was not well considered before. As reported in Ref.[2529], LaCoO$_3$ exhibits nonmagnetic behavior at low temperature (below 90 K) where the Co$^{3+}$ ions are in the low spin (LS) states, and undergoes a transition from semiconductor to metal over 500 K, coincident to the spin state transition from intermediate spin (IS) to high spin (HS). It has also been reported that the spin states of Co$^{3+}$ ions in LaCoO$_3$ can greatly affect oxygen vacancy formation on the surface [30] and their migration in the bulk [31]. However, there is a lack of report on surface stability of LaCoO$_3$ with consideration of the effects of spin states of Co$^{3+}$ ions. To this point, the reason for the above disagreement may arise from the spin states of the Co$^{3+}$ ions in LaCoO$_3$. It is necessary to perform calculations with appropriate Co$^{3+}$ spin states to clarify the properties of LaCoO$_3$ surface systems for the indepth researches.
In this work, we compared the surface grand potentials ($\Omega$) of lowindex facets of hexagonalphase LaCoO$_3$, where the Co$^{3+}$ ions' spin states are in LS, IS (ferromagnetic (FM), and antiferromagnetic (AFM)) and HS (FM, and AFM) states, to get a comprehensive understanding of LaCoO$_3$ surface systems. We determined the stable surfaces of LaCoO$_3$ in special chemical conditions by comparing the surface grand potentials ($\Omega$). As the spin states of Co$^{3+}$ ions change, the region of stability diagram undergoes a significant variation. In addition, the $E_{\textrm{O}_\textrm{v}}$ of stable facets of LaCoO$_3$ were computed in different spin states. It indicates that the energy of surface oxygen vacancy links strongly to the spin moment of Co$^{3+}$ ions. Our simulations indicate that the (1102) CoO$_2$terminated facet in hexagonal LaCoO$_3$ can be stabilized in most spin states when the Co$^{3+}$ ions' spin states are taken into consideration. And this agrees well with the experimental observation.
Ⅱ. Computational details A. Models of bulk and slabThe space group of crystalline LaCoO$_3$ is $R\overline{3c}$ (No.167), which can be constructed with either a rhombohedral or a hexagonal model (FIG. 1). The JahnTeller effect causes a slight distortion of CoO$_6$ octahedron in LaCoO$_3$, resulting in a primitive cell of rhombohedral structure [18] with lattice parameters of $a=b=c$=5.344 Å, $\alpha=\beta=\gamma$=61.01° (FIG. 1(a)). As shown in FIG. 1(b), the experimental lattice parameters of hexagonal model [18] are $a=b$=5.426 Å, $c$=12.991 Å, $\alpha=\beta$=90°, and $\gamma$=120°.
In order to discuss the thermodynamic stability of the lowindex surfaces, we constructed six slab models. Adapted from the models in Ref.[24], a sevenlayer nonstoichiometric slab of LaO (FIG. 2(a)) or CoO$_2$ (FIG. 2(b)) termination was used to simulate the facet along the (1102) direction, a thirteenlayer nonstoichiometric slab of LaO$_3$ (FIG. 2(c)) or Co (FIG. 2(d)) termination was used to simulate the facet along the (0001) direction; and an elevenlayer nonstoichiometric slab of LaCoO (FIG. 2(e)) or O$_2$ (FIG. 2(f)) termination was used to simulate the facet along the (1104) direction. The (1102), (0001) and (1104) directions in the hexagonal phase correspond to the (001), (111) and (110) directions in a pseudo cubic phase. The black frames in FIG. 2 represent the periodic boundaries.
B. General setup for computationDensity functional theory (DFT) calculations were performed with the Vienna ab initio Simulation Package (VASP) [32, 33]. The nuclei and core electrons were treated with the projector augmented wave (PAW) [34] method. Generalized gradient approximation (GGA) with the PerdewBurkeErnzerhof (PBE) [35] form was employed to describe the electron exchange and correlation. For relaxation of bulk LaCoO$_3$ (hexagonal phase), a plane wave basis set with cutoff of 400 eV and a 7×7×3 MonkhorstPack [36] $k$point mesh were used to get the optimal lattice parameters. The optimal parameters of hexagonal phase are $a=b$=5.485 Å and $c$=13.031 Å. For the slab calculations, the plane wave energy cutoff of 400 eV and a 3×3×1 of $k$point mesh were used to get the properties of different surfaces. All the atoms in the bulk and slab were allowed to relax until the maximum force on each atom was smaller than 0.05 eV/Å. In each slab model, a separation over 15 Å in vacuum was introduced to minimize interactions between periodic images.
The oxygen vacancy formation energy ($E_{\textrm{O}_\textrm{v}}$) is calculated using the following equation:
$\begin{eqnarray}E_{\textrm{O}_\textrm{v}}=E(\textrm{defect})  E(\textrm{perfect}) + \frac{1}{2} E(\textrm{O}_2) \end{eqnarray}$  (1) 
where $E\textrm{(defect)}$ is the total energy of the slab with one surface oxygen vacancy, $E\textrm{(perfect)}$ is the energy of the ideal slab and $E(\textrm{O}_2) $ is the energy of the O$_2$ molecule in the gas phase, respectively.
The different Co3d occupations result in three spin states of Co$^{3+}$ ions in LaCoO$_3$, including LS of t$_{2\textrm{g}}$$^6$e$_\textrm{g}$$^0$ (S=0), IS of t$_{2\textrm{g}}$$^5$e$_\textrm{g}$$^1$ (S=1) and HS of t$_{2\textrm{g}}$$^4$e$_\textrm{g}$$^2$ (S=2). Each spin state was constructed by changing the initial magnetic moment of Co$^{3+}$ ion, i.e. LS=0, IS=$\pm$2 and HS=$\pm$4, where the quantity is the difference between alpha and beta electrons in a certain Co$^{3+}$ ion. The initial configurations of all Co$^{3+}$ ions were set by MAGMOM, and then NUPDOWN was chosen to control the total spin of the slab. The selection rule of spin was also applied to get a reasonable energy when oxygen vacancy was generated. For example, two total spin numbers of the defect slab (+2 or 2 from the total spin number of initial slab) were calculated and the lower energy was accepted. The spin state of calculated O$_2$ molecule is always triplet. The partially filled d states in the Co$^{3+}$ ions are not well described by the standard DFT calculations, where the normal GGA methods give a zero band gap of LaCoO$_3$ on the contrast to the experimental value of about 0.6 eV [37]. We thus performed the DFT+$U$ approach, where the onsite $U$ and $J$ were treated as a single effective parameter among the d electrons on the Co$^{3+}$ ions, $U_{\textrm{eff}}=UJ$ [3840]. We fixed $J$=0.49 eV [41] in all our GGA+$U$ calculations, and varied the value of $U$ in determining an optimal parameter $U$ (3.4 eV), which will be discussed later.
C. The surface grand potentialTo compare the stability of different terminations, the surface grand potential ($\Omega$) is implemented in the calculations, as defined in Ref.[42]
$\begin{eqnarray} \Omega = \frac{1}{{2S}}\left( {E_{\textrm{slab}}  N_{\textrm{La}} \mu _{\textrm{La}}  N_{\textrm{Co}} \mu _{\textrm{Co}}  N_\textrm{O} \mu _\textrm{O} } \right)\end{eqnarray}$  (2) 
where $N_{\textrm{La}}$, $N_{\textrm{Co}}$, and $N_\textrm{O}$ are the numbers of La, Co, and O atoms in the slab, respectively, and $\mu$ represents the chemical potential of La, Co, and O atoms species. Since the surface is in thermodynamic equilibrium with the bulk LaCoO$_3$ ($\mu_{\textrm{LaCoO}_3}=E$$_{\textrm{LaCoO}_3}^{\textrm{bulk}}$) and the chemical potentials of each species are related to the chemical potential of the bulk crystal, we thus have the following constraint:
$\begin{eqnarray} E_{\textrm{LaCoO}_3 }^{\textrm{bulk}} =\mu _{\textrm{LaCoO}_3 } = \mu _{\textrm{La}} + \mu _{\textrm{Co}} + 3\mu _\textrm{O} \end{eqnarray}$  (3) 
When $\mu_{\textrm{La}}$ is eliminated by Eq.(2) and (3), we introduce:
$\begin{eqnarray} \Delta \mu _{\textrm{Co}}\hspace{0.15cm}& =& \hspace{0.15cm}\mu _{\textrm{Co}}  E_{\textrm{Co}}^{\textrm{bulk}}\end{eqnarray}$  (4) 
$\begin{eqnarray} \Delta \mu _\textrm{O } \hspace{0.15cm}&=&\hspace{0.15cm} \mu _\textrm{O}  \frac{1}{2}E_{\textrm{O}_2 }^{\textrm{mol}}\end{eqnarray}$  (5) 
Subsequently, the surface grand potential can be determined as:
$\begin{eqnarray}\Omega \hspace{0.15cm}& =&\hspace{0.15cm} \phi \hspace{0.1cm} \hspace{0.1cm} \frac{1}{{2S}}\left[{\Delta \mu _\textrm{O} (N_\textrm{O} \hspace{0.05cm}  \hspace{0.05cm}3N_{\textrm{La}} ) \hspace{0.05cm}+\hspace{0.05cm} \Delta \mu _{\textrm{Co}} (N_{\textrm{Co}} \hspace{0.05cm} \hspace{0.05cm} N_{\textrm{La}} )} \right] \end{eqnarray}$  (6) 
$\begin{array}{l} \phi = \frac{1}{{2S}}[{E_{{\rm{slab}}}}  {N_{{\rm{La}}}}E_{{\rm{LaCo}}{{\rm{O}}_3}}^{{\rm{bulk}}}  \frac{{E_{{{\rm{O}}_2}}^{{\rm{mol}}}}}{2}({N_{\rm{O}}}  3{N_{{\rm{La}}}})  \\ \quad \quad E_{{\rm{Co}}}^{{\rm{bulk}}}({N_{{\rm{Co}}}}  {N_{{\rm{La}}}})] \end{array}$  (7) 
The surface La, Co, and O atoms are assumed to form no condensate on the surface, and the chemical potential of each species must be lower than the energy of an atom in the stable phase. We thus obtain the following upper limits of the chemical potentials:
$\begin{eqnarray} \Delta \mu _\textrm{O} = \mu _\textrm{O}  \frac{1}{2}E_{\textrm{O}_2 }^{\textrm{mol}} < 0\end{eqnarray}$  (8) 
$\begin{eqnarray}\Delta \mu _{\textrm{La}} = \mu _{\textrm{La}}  E_{\textrm{La}}^{\textrm{bulk}} < 0\end{eqnarray}$  (9) 
$\begin{eqnarray}\Delta \mu _{\textrm{Co}} = \mu _{\textrm{Co}}  E_{\textrm{Co}}^{\textrm{bulk}} < 0\end{eqnarray}$  (10) 
By combining Eq.(3) and (9), the lower limits can be determined as follows,
$\Delta \mu _{\textrm{Co}} + 3\Delta \mu _\textrm{O} >  E_{\textrm{LaCoO}_3 }^\textrm{f}$  (11) 
$  E_{{\rm{LaCo}}{{\rm{O}}_3}}^{\rm{f}} = E_{{\rm{LaCo}}{{\rm{O}}_3}}^{{\rm{bulk}}}  E_{{\rm{La}}}^{{\rm{bulk}}}  E_{{\rm{Co}}}^{{\rm{bulk}}}  \frac{3}{2}E_{{{\rm{O}}_2}}^{{\rm{mol}}}$  (12) 
where $E_{\textrm{LaCoO}_3 }^\textrm{f}$ is the formation energy of LaCoO$_3$ bulk crystal with respect to the La and Co atoms in their bulk phases, and the O atom in the gas phase. Once $\Delta \mu_{\textrm{Co}}$ and $\Delta\mu_\textrm{O}$ are determined in the effective ranges, the accessible values of $\Omega$ are thus obtained.
Ⅲ. Results and discussion A. Determining the optimized $U$ parameterIt is well known that $U$ parameter plays an important role in the determination of the valence structure and the crystal structure. Since the band gap is generally related to the electronic structure of the bulk LaCoO$_3$, it is essential to determine an appropriate $U$ value for Co element to correct the mistake. Previous studies indicate that the ground state LaCoO$_3$ is a nonmagnetic semiconductor with Co$^{3+}$ ions in the LS state [28, 29]. Therefore, band gap scan calculations were performed with the parameter $U$ varying from 2.9 eV to 3.9 eV with an interval of 0.1 eV on the LS ground state. The band gap as a function of the exchange parameter $U$ is shown in FIG. 3. With the increase of parameter $U$, the value of the band gap linearly increases. When $U$=3.4 eV is employed, a band gap of 0.61 eV is obtained, which well agrees with the experimental value [37]. The $U$ parameter of 3.4 eV was thus used for the GGA+$U$ calculations to investigate surface stability of LaCoO$_3$ perovskite.
B. Thermodynamic stability in different spin states with GGA+$U$To examine the influence of spin states of Co$^{3+}$ ions on surface stability, the surface grand potentials ($\Omega$) of lowindex facets have been calculated (Table S1S3 in supplementary materials), where the Co$^{3+}$ ions' spin states are in LS, IS (FM, AFM) and HS (FM, AFM) states, respectively. The FM and AFM configurations can be considered as two limits in energy with respect to the spin state of LaCoO$_3$. To determine the accessible values of $\Omega$, we should calculate the effective intervals of $\Delta\mu_{\textrm{Co}}$ and $\Delta \mu_\textrm{O}$ according to the formation energy of LaCoO$_3$ bulk crystal, as defined in Ref.[42]. The obtained formation energy of hexagonal LaCoO$_3$ is 9.31 eV. On the basis of Eq.(11), $\Delta \mu_{\textrm{Co}}$ and $\Delta \mu_\textrm{O}$ are thus restricted within the ranges of (9.31 eV < $\Delta \mu_{\textrm{Co}}$ < 0.00 eV) and (3.10 eV < $\Delta \mu_{\textrm{O}}$ < 0.00 eV), respectively.
The stability diagram of the surface grand potential $\Omega$ of different terminations within the allowed area is displayed in FIG. 4(a) when the Co$^{3+}$ ions are in LS state. Only three terminations out of six lowindex facets are found to be stable: (1102) LaO, (0001) LaO$_3$ and (1102) CoO$_2$. The calculated results indicate that the (1102) LaO termination is stable in low O chemical potential and low Co chemical potential, while its complementary, the (1102) CoO$_2$ termination is in rich O chemical potential and rich Co chemical potential. Additionally, the (0001) LaO$_3$ facet shows a stability domain in moderate Co and rich O environment. FIG. 4(b) reveals the $\Omega$ of different facets of the low spin Co$^{3+}$ ions as a function of $\Delta \mu_{\textrm{Co}}$ with $\Delta \mu_\textrm{O}$=0 eV. It can be noted that the (1102) LaO termination is favored in a large interval (9.31 eV < $\Delta \mu_{\textrm{Co}}$ < 3.29 eV) and the (1102) CoO$_2$ termination becomes stable in a small interval (0.57 eV < $\Delta \mu_{\textrm{Co}}$ < 0 eV). When $\Delta \mu_{\textrm{Co}}$ is restricted to (3.29 eV < $\Delta \mu_{\textrm{Co}}$ < 0.57 eV), only the (0001) LaO$_3$ termination can be observed.
FIG. 5(a) shows the $\Omega$ of different terminations as functions of $\Delta \mu_{\textrm{Co}}$ and $\Delta \mu_\textrm{O}$ in IS state with a FM configuration. Our calculations suggest that (1102) CoO$_2$, (0001) Co and (1104) LaCoO are unstable because their surface grand potentials are always larger than that of at least one stable facet within the allowed region. According to FIG. 5(a), the (1102) LaO termination is the most stable one in Oand Copoor conditions. The (0001) LaO$_3$ and (1104) O$_2$ terminations are favored in a small region corresponding to Oand Corich environments, respectively. FIG. 5(b) shows the $\Omega$ of different facets as a function of $\Delta \mu_{\textrm{Co}}$ in Orich condition ($\Delta \mu_\textrm{O}$=0 eV) in IS state with a FM configuration. It appears that the (0001) LaO$_3$ and (1104) O$_2$ terminations can be stabilized in a very limited area (2.15 eV < $\Delta \mu_{\textrm{Co}}$ < 0.93 eV and 0.93 eV$ < $$\Delta \mu_{\textrm{Co}}$ < 0 eV), respectively. When $\Delta \mu_{\textrm{Co}}$ is lower than 2.15 eV, only the (1102) LaO termination can be observed.
As displayed in FIG. 5(c), the (1102) LaO, (0001) LaO$_3$, (1104) O$_2$ and (1102) CoO$_2$ are thermodynamically most favorable when the IS state is in an AFM configuration. Compared to the stability diagram of the IS state with a FM configuration, the area of the (1102) LaO and (0001) LaO$_3$ facets hardly change, while the region of (1104) O$_2$ termination only exists in a very small domain. The (1102) CoO$_2$ termination emerges and becomes most stable in rich O and Co environment. FIG. 5(d) plots the surface grand potential $\Omega$ of different facets as a function of $\Delta \mu_{\textrm{Co}}$ with $\Delta \mu_{\textrm{O}}$=0 eV in IS state with an AFM configuration. The (1102) LaO termination is favorable in a large interval (9.31 eV < $\Delta \mu_{\textrm{Co}}$ < 2.36 eV). For the complementary termination, the (1102) CoO$_2$ can be stabilized in a small interval (0.46 eV < $\Delta \mu_{\textrm{Co}}$ < 0 eV). When $\Delta \mu_{\textrm{Co}}$ is restricted to (2.36 eV < $\Delta \mu_{\textrm{Co}}$ < 0.46 eV), two stable terminations, the (0001) LaO$_3$ and (1104) O$_2$ terminations, can be obtained.
FIG. 6(a) plots the surface grand potential $\Omega$ of different terminations as functions of $\Delta \mu_{\textrm{Co}}$ and $\Delta \mu_{\textrm{O}}$ in HS state with a FM configuration. The calculated results indicate that only four terminations are thermodynamically stable: (1102) LaO, (0001) LaO$_3$, (1104) O$_2$ and (1102) CoO$_2$ as shown in FIG. 6(a). The (1102) LaO termination is stable at low O chemical potential (Opoor limit) and low Co chemical potential (Copoor limit), as the (0001) LaO$_3$ facet shows a stability domain in moderate Co and rich O environment. In general, the (1104) O$_2$ and (1102) CoO$_2$ terminations become the stable facets in a small domain corresponding to Coand Orich conditions. FIG. 6(b) plots the $\Omega$ of different facets of the high spin Co$^{3+}$ ions with a FM configuration as a function of $\Delta \mu_{\textrm{Co}}$ with $\Delta \mu_{\textrm{O}}$=0 eV (Orich condition). The results indicate that when $\Delta \mu_{\textrm{Co}}$ ranges from 9.31 eV to 4.60 eV, the (1102) LaO termination is favored. When $\Delta \mu_{\textrm{Co}}$ is restricted to (4.60 eV < $\Delta \mu_{\textrm{Co}}$ < 0.98 eV), only the (0001) LaO$_3$ surface can be observed. When $\Delta \mu_{\textrm{Co}}$ increases to 0.45 eV, the (1102) O$_2$ termination disappears and the (1102) CoO$_2$ surface emerges. Notably, these (0001) Co and (1104) LaCoO terminations cannot be stabilized because they are outside of the accessible region.
As shown in FIG. 6(c), when the HS state is in the AFM configuration, the (1102) LaO, (0001) LaO$_3$ and (1102) CoO$_2$ terminations are the most favorable facets. The stable region of (1102) LaO termination of the AFM configuration is larger than that of the FM configuration. However, the area of (0001) LaO$_3$ becomes smaller. And the (1102) CoO$_2$ facet almost does not change. FIG. 6(d) reveals the $\Omega$ of the lowestenergy facets of the high spin Co$^{3+}$ ions with an AFM configuration as a function of $\Delta \mu_{\textrm{Co}}$ in Orich condition ($\Delta$$\mu_\textrm{O}$=0 eV). According to FIG. 6(d), the (1104) LaCoO, (0001) Co and (1104) O$_2$ facets cannot be obtained; however, the (1102) LaO, (0001) LaO$_3$ and (1102) CoO$_2$ terminations can be stabilized in some special Co environment, as 9.31 eV < $\Delta \mu_{\textrm{Co}}$ < 2.36 eV, 2.36 eV < $\Delta \mu_{\textrm{Co}}$ < 0.58 eV, and 0.58 eV < $\Delta \mu_{\textrm{Co}}$ < 0 eV, respectively.
According to the above calculations, as the Co$^{3+}$ ions' spin states change, the stability diagram of lowindex facets in hexagonalphase LaCoO$_3$ undergoes a significant variation. Considering different spin states of Co$^{3+}$ ions and different magnetic configurations, one can find that (1102) LaO termination is stable in a large region with oxygen poor condition. Moreover, the (0001) LaO$_3$ and (1102) CoO$_2$ terminations can be stabilized in an oxygen rich condition, which is a typical experimental environment for LaCoO$_3$. And tuning the chemical potential of Co can tune the final surface exposition of LaCoO$_3$ nanoparticle. The (1102) CoO$_2$ termination is favorable under suitable chemical potential regions, which can reach an agreement on theoretical results and experimental observation. This corresponds to the spin states of Co$^{3+}$ ions in perovskite LaCoO$_3$. The Co$^{3+}$ ions occupying octahedral sites surrounded by oxygen ions can introduce complex magnetic property in LaCoO$_3$, which has influence on the CoO and OO bond strength on the surface.
C. Oxygen vacancy formation energy of the stable facetsSurface oxygen vacancy plays a key role in the catalytic oxidation reactions in perovskite LaCoO$_3$. We therefore analyze how the spin states of Co$^{3+}$ ions change the fundamental properties of oxygen vacancy formation. In order to explore the properties of these possible exposed facets, the $E_{\textrm{O}_\textrm{v}}$ of the surfaces have been calculated. As shown in Table Ⅰ, it can be found that higher Co$^{3+}$ magnetic moments lead to lower oxygen vacancy formation energies, in good agreement with the previous work [30]. They attribute the tendencies in the surface oxygen vacancy formation energies to variations in the O pband center, which can be described as the CoO bond strength. This trend suggests that oxygen vacancy energetics link strongly to the spin moment of Co$^{3+}$ ions. Furthermore, whether the Co$^{3+}$ ions are in a LS, an ISFM, an ISAFM, a HSFM or a HSAFM configuration, the (1102) LaO termination has the highest $E_{\textrm{O}_\textrm{v}}$ among these surfaces. It indicates that the (1102) LaO facet may be not active, since O atoms are strongly bound. In the HSFM configuration, the $E_{\textrm{O}_\textrm{v}}$ of (1102) CoO$_2$ termination is close to zero (0.08 eV), which is lower than that of (0001) LaO$_3$ termination (0.20 eV). When the spin state of Co$^{3+}$ ions is in a HSAFM configuration, the $E_{\textrm{O}_\textrm{v}}$ of (1102) CoO$_2$ facet is 1.05 eV, which is more than the (0001) LaO$_3$ termination (0.69 eV). As a result, it can be predicted that the (1102) CoO$_2$and (0001) LaO$_3$terminated surfaces might have good activity of O atoms, which can play a critical role in the surface reaction processes.
For comparison, we also calculated surface oxygen vacancy formation energies with NS polarized method as used in Ref.[24]. The NS oxygen vacancy formation energy of the (1102) LaO termination is still the maximum. Moreover, the $E_{\textrm{O}_\textrm{v}}$ of (1102) CoO$_2$ termination is greater than that of the (0001) LaO$_3$ surface. Comparing the $E_{\textrm{O}_\textrm{v}}$ in NS state of (1102) LaO and (0001) LaO$_3$ facets with previous work [24], it can be found that there is a certain degree of error even if we use the same method and model. It is because that there will be a variety of uncertain spin configurations in the calculation process when we ignore the spin states, since there are a lot of local minimums on the energy surface with different spin configurations. In fact, the above result also represents that the spin states extremely affect surface oxygen vacancies. Therefore, considering the spin states of Co$^{3+}$ ions is an essential role in exploring the surface properties of LaCoO$_3$.
Ⅳ. CONCLUSIONWe performed DFT+$U$ calculations to study the effect of Co$^{3+}$ spin states on surface stabilities of several lowindex terminations of the perovskite LaCoO$_3$. And the $E_{\textrm{O}_\textrm{v}}$ of possible exposed facets is calculated in different spin states. The parameter $U$ is employed to correct the onsite Coulomb and the electron interactions for local d orbitals. It is found that the spin states of Co$^{3+}$ ions in hexagonalphase LaCoO$_3$ have an important effect on the region of stability diagrams. For the most cases with different spin states and spin configurations of Co$^{3+}$ ions, the (1102) CoO$_2$ termination can be stabilized in oxygen rich environment, which agrees well with the experimental observation. Moreover, the spin states of Co$^{3+}$ ions can affect surface oxygen vacancy. The higher Co$^{3+}$ spin states give out the lower oxygen vacancy formation energies. Our results reveal that the spin states of Co$^{3+}$ ions should be well considered for studying the properties of LaCoO$_3$ facets.
Supplementary materials: The surface grand potential ($\Omega$) of different terminations in LS, IS (FM, AFM) and HS (FM, AFM) states are displayed, which is expressed as functions of the excess O and Co chemical potentials ($\Delta \mu_{\textrm{Co}}$ and $\Delta \mu_{\textrm{O}}$).
Ⅴ. ACKNOWLEDGMENTSThis work was supported by the National Natural Science Foundation of China (No.U1232118, No.21203099), the National Basic Research Program (No.2014CB932403), the Program of Introducing Talents of Disciplines to China Universities (No.B06006), Research Program for Advanced and Applied Technology of Tianjin (No.13JCYBJC36800), Doctoral Fund of Ministry of Education of China (No.20120031120033), Fok Ying Tung Education Foundation (No.151008), and Special Program for Applied Research on Super Computation of the NSFCGuangdong Joint Fund (the second phase). We appreciate the supports from the National SuperComputing Center at Tianjin and Guangzhou.
Table S1 has listed the grand surface energy (Ω) of lowindex facets in low spin (LS) LaCoO_{3}. Based on the following formulas, we can reveal the relationship between the grand surface energy (Ω) with △μ_{Co} and△μ_{O}. Thus, Figs. 4a and 4b in article can be obtained.
Table S2 has listed the grand surface energy (Ω) of lowindex facets in intermediate spin (IS) LaCoO_{3}. Based on the following formulas, we can reveal the relationship between the grand surface energy (Ω) with △μ_{Co} and△μ_{O}. Thus, Figs. 5a, 5b, 5c and 5d in article can be obtained.
Table S3 has listed the grand surface energy (Ω) of lowindex facets in high spin (HS) LaCoO_{3}. Based on the following formulas, we can reveal the relationship between the grand surface energy (Ω) with △μ_{Co} and △μ_{O}. Thus, Figs. 6a, 6b, 6c and 6d in article can be obtained.
[1]  E. P. Murray, T. Tsai, and S. A. Barnett, Nature 400 , 649 (1999). DOI:10.1038/23220 
[2]  J. Suntivich, K. J. May, H. A. Gasteiger, and J. B. Goodenough, ShaoHorn Y., Science 334 , 1383 (2011). 
[3]  J. Suntivich, H. A. Gasteiger, N. Yabuuchi, H. Nakanishi, and J. B. Goodenough, ShaoHorn Y., Nat. Chem. 3 , 546 (2011). DOI:10.1038/nchem.1069 
[4]  W. Desmond Ng J., Y. Gorlin, T. Hatsukade, and T. F. Jaramillo, Adv. Energy Mater. 3 , 1545 (2013). DOI:10.1002/aenm.v3.12 
[5]  J. Jung, M. Risch, S. Park, M. G. Kim, G. Nam, H. Y. Jeong, Y. ShaoHorn, and J. Cho, Energy Environ. Sci. 9 , 176 (2016). DOI:10.1039/C5EE03124A 
[6]  C. H. Kim, G. S. Qi, K. Dahlberg, and W. Li, Science 327 , 1624 (2010). DOI:10.1126/science.1184087 
[7]  X. G. Li, Y. H. Dong, H. Xian, W. Y. Hernández, M. Meng, H. H. Zou, A. J. Ma, T. Y. Zhang, Z. Jiang, N.Tsubaki, and P. Vernoux, Energy Environ. Sci. 4 , 3351 (2011). DOI:10.1039/c1ee01726h 
[8]  K. S. Song, S. K. Kang, and S. D. Kim, Catal. Lett. 49 , 65 (1997). DOI:10.1023/A:1019072314394 
[9]  G. Saracco, G. Scibilia, A. Iannibello, and G. Baldi, Appl. Catal. B 8 , 229 (1996). DOI:10.1016/09263373(95)000844 
[10]  Y. M. Choi, D. S. Mebane, M. C. Lin, and M. L. Liu, Chem. Mater. 19 , 1690 (2007). DOI:10.1021/cm062616e 
[11]  Y. L. Lee, J. Kleis, J. Rossmeisl, and D. Morgan, Phys. Rev. B 80 , 224101 (2009). DOI:10.1103/PhysRevB.80.224101 
[12]  Y. L. Lee, J. Kleis, J. Rossmeisl, Y. ShaoHorn, and D. Morgan, Energy Environ. Sci. 4 , 3966 (2011). DOI:10.1039/c1ee02032c 
[13]  S. O. Choi, M. Penninger, C. H. Kim, W. F. Schneider, and L. T. Thompson, ACS Catal. 3 , 2719 (2013). DOI:10.1021/cs400522r 
[14]  M. W. Penninger, C. H. Kim, L. T. Thompson, and W. F. Schneider, J. Phys. Chem. C 119 , 20488 (2015). DOI:10.1021/acs.jpcc.5b06351 
[15]  S. Khan, R. J. Oldman, F. Corà, C. R. A. Catlow, S. A. French, and S. A. Axon, Phys. Chem. Chem. Phys. 8 , 5207 (2006). DOI:10.1039/b602753a 
[16]  Z. Z. Chen, C. H. Kim, L. T. Thompson, and W. F. Schneider, Surf. Sci. 619 , 71 (2014). DOI:10.1016/j.susc.2013.09.012 
[17]  S. G. Zhang, N. Han, and X. Y. Tan, RSC Adv. 5 , 760 (2015). DOI:10.1039/C4RA12563K 
[18]  P. G. Radaelli, and S. W. Cheong, Phys. Rev. B 66 , 094408 (2002). DOI:10.1103/PhysRevB.66.094408 
[19]  S. M. Zhou, L. F. He, S. Y. Zhao, Y. Q. Guo, J. Y. Zhao, and L. Shi, J. Phys. Chem. C 113 , 13522 (2009). DOI:10.1021/jp9003032 
[20]  M. Risch, A. Grimaud, K. J. May, K. A. Stoerzinger, T. J. Chen, and A. N. Mansour, ShaoHorn Y., J. Phys. Chem. C 117 , 8628 (2013). DOI:10.1021/jp3126768 
[21]  S. Mukhopadhyay, M. W. Finnis, and N. M. Harrison, Phys. Rev. B 87 , 125132 (2013). DOI:10.1103/PhysRevB.87.125132 
[22]  P. Ravindran, P. A. Korzhavyi, H. Fjellvàg, and A. Kjekshus, Phys. Rev. B 60 , 16423 (1999). DOI:10.1103/PhysRevB.60.16423 
[23]  A. Mineshige, M. Inaba, T. Yao, Z. Ogumi, K. Kikuchi, and M. Kawase, J. Solid State Chem. 121 , 423 (1996). DOI:10.1006/jssc.1996.0058 
[24]  X. Liu, Z. Z. Chen, Y. W. Wen, R. Chen, and B. Shan, Catal. Sci. Technol. 4 , 3687 (2014). DOI:10.1039/C4CY00538D 
[25]  G. Thornton, B. C. Tofield, and A. W. Hewat, J. Solid State Chem. 61 , 301 (1986). DOI:10.1016/00224596(86)900356 
[26]  K. Asai, P. Gehring, H. Chou, and G. Shirane, Phys. Rev. B 40 , 10982 (1989). DOI:10.1103/PhysRevB.40.10982 
[27]  S. Yamaguchi, Y. Okimoto, and Y. Tokura, Phys. Rev. B 55, R8666(1997). 
[28]  I. A. Nekrasov, S. V. Streltsov, M. A. Korotin, and V. I. Anisimov, Phys. Rev. B 68 , 235113 (2003). DOI:10.1103/PhysRevB.68.235113 
[29]  D. P. Kozlenko, N. O. Golosova, Z. Jirk, L. S. Dubrovinsky, B. N. Savenko, M. G. Tucker, Le Godec Y., and V. P. Glazkov, Phys. Rev. B 75 , 064422 (2007). DOI:10.1103/PhysRevB.75.064422 
[30]  W. T. Hong, M. Gadre, Y. L. Lee, M. D. Biegalski, H. M. Christen, and D. Morgan, ShaoHorn Y., J. Phys. Chem. Lett. 4 , 2493 (2013). DOI:10.1021/jz401271m 
[31]  A. M. Ritzmann, M. Pavone, A. B. MuñozGarcıa, J. A. Keith, and E. A. Carter, J. Mater. Chem. A 2 , 8060 (2014). DOI:10.1039/C4TA00801D 
[32]  G. Kresse, and J. Hafner, Phys. Rev. B 47 , 558 (1993). DOI:10.1103/PhysRevB.47.558 
[33]  W. Kohn, and L. J. Sham, Phys. Rev. , A1133 (1965). 
[34]  P. E. Blöchl, Phys. Rev. B 50 , 17953 (1994). DOI:10.1103/PhysRevB.50.17953 
[35]  J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77 , 3865 (1996). DOI:10.1103/PhysRevLett.77.3865 
[36]  H. J. Monkhorst, and J. D. Pack, Phys. Rev. B 13 , 5188 (1976). DOI:10.1103/PhysRevB.13.5188 
[37]  A. Chainani, M. Mathew, and D. D. Sarma, Phys. Rev. B 46 , 9976 (1992). DOI:10.1103/PhysRevB.46.9976 
[38]  V. I. Anisimov, J. Zaanen, and O. K. Andersen, Phys. Rev. B 44 , 943 (1991). DOI:10.1103/PhysRevB.44.943 
[39]  H. Hsu, K. Umemoto, M. Cococcioni, and R. Wentzcovitch, Phys. Rev. B 79 , 125124 (2009). DOI:10.1103/PhysRevB.79.125124 
[40]  L. Wang, T. Maxisch, and G. Ceder, Phys. Rev. B 73 , 195107 (2006). DOI:10.1103/PhysRevB.73.195107 
[41]  C. L. Ma, and J. Cang, Solid State Commun. 150 , 1983 (2010). DOI:10.1016/j.ssc.2010.08.023 
[42]  F. Bottin, F. Finocchi, and C. Noguera, Phys. Rev. B 68 , 035418 (2003). DOI:10.1103/PhysRevB.68.035418 
b. 天津化学化工协同创新中心, 天津 300072;
c. 南开大学物理科学学院, 天津 300071;
d. 天津大学理学院化学系, 天津 300072