Communication: Different behavior of Young's modulus and fracture strength of CeO2: Density functional theory calculations Ryota Sakanoi, Tomomi Shimazaki, Jingxiang Xu, Yuji Higuchi, Nobuki Ozawa, Kazuhisa Sato, Toshiyuki Hashida, and Momoji Kubo Citation: The Journal of Chemical Physics 140, 121102 (2014); doi: 10.1063/1.4869515 View online: http://dx.doi.org/10.1063/1.4869515 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/140/12?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Communication: Improving the density functional theory+U description of CeO2 by including the contribution of the O 2p electrons J. Chem. Phys. 136, 041101 (2012); 10.1063/1.3678309 Residual stress-dependent electric conductivity of sputtered co-doped CeO2 thin-film electrolyte J. Appl. Phys. 109, 084321 (2011); 10.1063/1.3573669 Critical tensile and compressive strains for cracking of Al2O3 films grown by atomic layer deposition J. Appl. Phys. 109, 084305 (2011); 10.1063/1.3567912 On the correlation of Young’s modulus and the fracture strength of metallic glasses J. Appl. Phys. 109, 033515 (2011); 10.1063/1.3544202 Temperature and strain-rate dependent fracture strength of graphene J. Appl. Phys. 108, 064321 (2010); 10.1063/1.3488620
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
THE JOURNAL OF CHEMICAL PHYSICS 140, 121102 (2014)
Communication: Different behavior of Young’s modulus and fracture strength of CeO2 : Density functional theory calculations Ryota Sakanoi, Tomomi Shimazaki, Jingxiang Xu, Yuji Higuchi, Nobuki Ozawa, Kazuhisa Sato, Toshiyuki Hashida, and Momoji Kuboa) Fracture and Reliability Research Institute, Graduate School of Engineering, Tohoku University, 6-6-11 Aoba, Aramaki, Aoba-ku, Sendai 980-8579, Japan
(Received 28 January 2014; accepted 13 March 2014; published online 24 March 2014) In this Communication, we use density functional theory (DFT) to examine the fracture properties of ceria (CeO2 ), which is a promising electrolyte material for lowering the working temperature of solid oxide fuel cells. We estimate the stress-strain curve by fitting the energy density calculated by DFT. The calculated Young’s modulus of 221.8 GPa is of the same order as the experimental value, whereas the fracture strength of 22.7 GPa is two orders of magnitude larger than the experimental value. Next, we combine DFT and Griffith theory to estimate the fracture strength as a function of a crack length. This method produces an estimated fracture strength of 0.467 GPa, which is of the same order as the experimental value. Therefore, the fracture strength is very sensitive to the crack length, whereas the Young’s modulus is not. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4869515] I. INTRODUCTION
Solid oxide fuel cells (SOFCs), which generate electrical energy from the reaction of H2 with O2 , are regarded as a key industrial technology for clean energy conversion systems that do not emit pollutants and greenhouse gases such as NOx and CO2 .1, 2 However, the high working temperatures of SOFCs have so far prevented their widespread use. Standard zirconia electrolytes for SOFCs operate at 600 ◦ C–800 ◦ C.3 This means expensive heat-proof components are required, and it also increases the number of breakdowns and failures from repetitive start-stop processes causing large temperature differences.4, 5 Doped ceria (CeO2 ) is a promising electrolyte material for SOFCs with low working temperatures (500 ◦ C– 650 ◦ C).6 Reducing the working temperatures would make SOFCs a more competitive technology. There are many experimental studies on the physical properties of CeO2 -based materials, such as the ionic conductivity.7–9 The mechanical properties of CeO2 -based electrolytes are also important for the reliability of SOFCs.10, 11 It is essential to understand fracture properties at the atomic scale, including the bond stretching, bending, and breaking where fractures occur. Theoretical studies are required for analyzing fracture properties on an atomic scale because experimental studies alone are not sufficient. First-principles calculations are a powerful computational tool for understanding atomic-scale phenomena in fields such as chemistry, physics, and materials science.12–16 We have also used first-principles calculations to investigate metal-oxide systems such as ZnO,17 PtO2 ,18 MgO,19 SnO2 ,20 Al2 O3 ,21 and CeO2 .22 First-principles calculations based on density functional theory (DFT) techniques have also been used to investigate fracture properties such as Young’s modulus and fracture strength.23–27 Deyirmenjian et al. pioneered DFT-based tensile tests for solid-state materials by publishing a) [email protected]
0021-9606/2014/140(12)/121102/4/$30.00
the stress-strain curve for defective aluminum.28 Kohyama studied cubic SiC to explain the grain boundary effects.29 Li and Tang also investigated the elasticity and ideal strength of SiC with a plane-wave DFT-based tensile test.30 These theoretical tensile tests show the atomic-level structural changes and predict fracture properties under uniform loads. However, the actual stresses would be concentrated around cracks and there would be much higher stresses around meso-sized cracks. It is difficult to measure local stresses experimentally because only an averaged stress is measured against the applied strain. Therefore, it is necessary to predict and compare the fracture properties of an ideal material with no cracks and a realistic material with cracks. In this Communication, we develop a theoretical method for predicting the fracture properties of CeO2 containing a crack. First, we calculate the theoretical stress-strain curve for CeO2 without cracks, and use DFT to discuss the Young’s modulus and the fracture strength. Because first-principles methods usually only consider an ideal model with no defects such as cracks, they are inadequate for discussing the fracture properties of CeO2 . Thus, we model the effect of the cracks and estimate the real fracture strength by using a combination of DFT and Griffith theory. We explain our theoretical methods and calculation conditions in Sec. II, present and discuss the calculation results in Sec. III, and summarize our work in Sec. IV.
II. CALCULATION METHODS
We employ the Perdew-Wang 91 (PW91) functional for the DFT calculations in order to obtain the optimized geometrical structures of the CeO2 bulk model and their total energies. The electronic structures for core electrons are described by the effective core potential. We employed DMol3 software package for these DFT calculations, which has been used for various investigations of CeO2 such as surface defects,31
140, 121102-1
© 2014 AIP Publishing LLC
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
121102-2
Sakanoi et al.
J. Chem. Phys. 140, 121102 (2014)
water adsorption,32 and catalytic reduction.33 The validity of DMol3 was already confirmed by the previous studies. We used the double-numeric polarized basis set with a globalspace cutoff of 5.8 Å and periodic boundary conditions with a Monkhorst-Pack grid of 4 × 4 × 4. Carter and co-workers34, 35 showed the influence of tensile direction on mechanical properties of chromia and iron oxides by DFT calculations. In CeO2 , the structure is cubic system and symmetry. Thus, we do not need considering the effect of tensile direction. To calculate the Young’s modulus and fracture strength, we assume that the energy-density, Ecell /Vcell , can be described in the following polynomial equation as a function of the strain, , 1 2 (1) Y + A 3 + B 4 . 2 Here, Ecell and Vcell are the total energy and the volume of the unit cell, respectively. Y is the Young’s modulus, and A and B are the coefficients of the third- and fourth-order terms, respectively. In the equation, the first-order term disappears, because of the equilibrium structure at = 0, which is obtained by relaxing the atom positions and the lattice constants of the CeO2 model. The stress, σ , can be obtained from the first derivative of Eq. (1) Ecell /Vcell =
σ = Y + 3A 2 + 4B 3 .
(2)
This method is also used in other studies for the second- or third-order elastic modulus to estimate Young’s modulus.36–38 To determine the yield point and fracture strength, we take into account the polynomial terms up to the fourth-order term. The yield point can be determined from the condition δσ /δ = 0, and the yield point, yield . The fracture strength, σ yield , is obtained analytically by using the coefficients of the polynomial equation as follows: √ −3A − 9A2 − 12BY yield = , (3) 12B 2 3 + 4Byield , σyield = Y yield + 3Ayield
=
2 2 Y yield + Ayield . 3
(4) (5)
III. RESULTS AND DISCUSSIONS
First, we use a first-principles calculation technique to obtain the Young’s modulus and fracture strength of crackfree CeO2 from the variation of the total energy with the applied strain. In this calculation, the CeO2 bulk model consists of 4 Ce4 + and 8 O2 − atoms, as shown in Figure 1. The lattice constant of CeO2 is 5.4707 Å, which is optimized value at the strain of 0.00 by DFT calculation. The lattice constant depends on basis set and a Monkhorst-Pack grid value. The lattice constant value is 5.47 Å in PAW-PBE,39 FPLAPWPBE,39 and GTO-PBE,40 5.40 Å in PAW-HSE,39 5.49 Å in PAW-PBE+U,39 5.39 Å in PAW-PBE0,39 5.45 Å in GTOTPSS,40 and 5.41 Å in GTO-HSE.40 The experimental value is 5.41 Å.41 According to these results, our value of 5.47 Å is in agreement with the previous results. On the other hand, the
O
2-
Ce
4+
FIG. 1. CeO2 calculation model.
lattice constant does not depend on a Monkhorst-Pack grid value much because we calculated with a Monkhorst-Pack grid of 4 × 4 × 4 and 8 × 8 × 8 and the change of the lattice constant is about 0.0002%. Figure 2(a) shows the change in the energy density with the strain on the [001] axis, and each point is determined from the structural relaxation with a fixed lattice constant. Strains from 0.00 to 0.13 at 0.01 increments are applied. The Young’s modulus and coefficients of polynomial equation (2) are estimated by the least squares method, and are summarized in Table I. In this fitting procedure, we use strains from 0.00 to 0.13 in one direction. The Young’s modulus obtained from these calculations is 221.8 GPa on the [001] axis and the experimental Young’s modulus is 175 GPa.42 These values are similar. In this model, a crystalline material with no faults, such as grain boundaries or cracks, is used. Experimental materials usually contain many defects and the Young’s modulus is thought to decrease as the number of defects increases.43 Our calculated Young’s modulus is similar to the experimental Young’s modulus, suggesting that the effect of defects is small. Jiang et al.44 showed pioneered work on the effect of on-site Coulomb repulsion term U on the electronic structure and the bond length of CeO2 . Then, we also discussed the effect of on-site Coulomb repulsion term U on the fracture properties. We used CASTEP code for the calculation because DMol3 code, which was used for the calculation of fracture properties in this Communication, cannot consider on-site Coulomb repulsion term U. We employ CA-PZ functional for the DFT calculations by CASTEP code. We calculated the energy-densities at the strain of 0.00 and 0.09 with the standard DFT and DFT+U. Here, the strain of 0.00 and 0.09 corresponds to the origin and yield point of the stress-strain curve in Fig. 2, respectively. The reasonable U-J value is set to 6.3 eV and applied to 4f electrons of Ce atoms following to the earlier study.44 The difference TABLE I. Young’s modulus, Y, and coefficients of Eq. (2) for the energy density vs. strain curve of CeO2 based on the DFT method. Y (GPa)
A (GPa)
B (GPa)
221.8
1198.5
−9052.6
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
Sakanoi et al.
121102-3
J. Chem. Phys. 140, 121102 (2014)
Fracture strength [GPa]
(a) Ecell/Vcell [GPa]
3.0 Calculation result
2.0 1.0 0.0 0.00
Stress [GPa]
(b)
Fitting
0.05 0.10 Strain
0.15
3.0 2.5 2.0 1.5 1.0 0.5 0.0 0.0
1.0 2.0 Crack length [µm]
3.0
FIG. 3. Fracture strength, σ f , as a function of crack length, a.
25.0 20.0
(c)
15.0 10.0
(d) 2-
5.0 0.0 0.00
O 4+ Ce
0.05 0.10 Strain
0.15
FIG. 2. (a) Changes in energy density vs. strain. Black circles indicate the DFT calculation results, and the black solid line is the polynomial equation (2) fitted by the least squares method. (b) Stress-strain curve obtained from the derivative of Eq. (2). Figures (c) and (d) show CeO2 structures before and after the yield point, respectively.
of change in the energy density at the strain of 0.09 is 0.08 GPa. Therefore, the effect of Coulomb repulsion term U on the fracture property is not critical. Next, we calculate the fracture strength of crack-free CeO2 and compare it with the experimental fracture strength. The stress-strain curve is obtained from the first derivative of Eq. (2) as shown in Figure 2(b), where the fitting is calculated from strains of 0.00–0.13. The curve shows a yield point of 0.09 and a fracture strength of 22.7 GPa. Interestingly, the stress value is not zero after the curve reaches the yield point. At the yield point, the cerium atom, which has an eightfold oxygen coordination number in CeO2 , changes to a sixfold one with the strain (Figs. 2(c) and 2(d)). Thus, a phase transition from a fluorite CeO2 crystal structure is observed, which is in good agreement with the classical molecular dynamics simulation reported by Sayle et al.45 However, there are large differences between the calculated and experimental fracture strengths. A fracture strength of 22.7 GPa is obtained from the DFT calculation, whereas the experimental value is 0.250 GPa.42 The mechanical properties generally depend on factors such as the grain size, defects, flaws, cracks, voids, and grain boundaries. These factors are affected by experimental preparation techniques, thus it is difficult to directly compare calculated and experimental results. The theoretical tensile tests yield ideal values for the fracture strength, because simulation models do not possess any macro- or meso-sized cracks and voids. Wong et al. prepared a micrometer-scale SiC whisker,
which had no macro- or meso-sized cracks and no voids.46 They used atomic force microscopy to determine the experimental fracture strength as 53.4 GPa. Li and Wang30 estimated fracture strength of 50.8 GPa by using a DFT calculation. These results strongly suggest that DFT-based tensile tests are reliable, and provided that the samples do not contain any macro- or meso-sized cracks and voids, experiments can provide the ideal fracture properties that are calculated by DFT. However, to estimate the fracture properties of CeO2 materials, we must take into account the effect of defects. Next, we discuss the discrepancy between the calculation results and the experimentally observed fracture strength caused by cracks or other flaws. Stress is concentrated around the end of cracks contained in experimental materials, and the fracture strength is estimated by the following equation, based on Griffith theory,47 2γs Y σf = , (6) π a2 where σ f is the value of the stress at which fracture occurs, a is the length of the crack, γ s is the surface energy, and Y is the Young’s modulus. This equation indicates that fracture occurs when the incremental release of stored elastic strain energy caused by the increase in crack length becomes larger than the increase in surface energy from the creation of new surface area. We used a surface energy, γ s , of 1.545 J/m2 and a Young’s modulus, Y, of 221.8 GPa, which were obtained by DFT, in Griffith theory calculations. The surface energy is calculated by DFT method: γ s = (Esurface − Ebulk )/(2 × S), where Esurface is energy of CeO2 (001) with vacuum phase, Ebulk is energy of bulk model, and S is surface area. We employed the constant value against crack length. Figure 3 shows the dependence of fracture strength on crack length, a given by this equation. The fracture strength decreases sharply when the length of the crack is small. This indicates that the fracture strength is very sensitive to the crack length. A scanning electron microscopy (SEM) image of gadolinia-doped ceria published by Tsoga et al.48 showed that the defect size is about 1.0 μm. Using SEM, Leng et al.49 showed that the gadoliniadoped ceria defect was about 2.0 μm. According to these experimental data, we estimate the fracture strength as follows by using Eq. (6). When the crack lengths are 1.0 and 2.0 μm, the calculated fracture strengths are 0.660 and 0.467 GPa, respectively, which are close to the experimental value of
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
121102-4
Sakanoi et al.
0.250 GPa. When the fracture strength is 0.250 GPa, we estimate that the crack length is about 7.0 μm. Therefore, our method for estimating the fracture strength from the crack size and the crack size from the fracture strength is valid. Finally, we conclude that the fracture strength is very sensitive to crack length, although the Young’s modulus is not sensitive to defects; thus, it is important to take account of cracks in the theoretical estimation of the fracture strength. IV. SUMMARY
We have investigated the difference in the behavior of the Young’s modulus and fracture strength of CeO2 with DFT. The calculated Young’s modulus of 221.8 GPa was of the same order as the experimental value of 175 GPa. However, the fracture strength of 22.7 GPa was two orders of magnitude larger than the experimental value of 0.250 GPa. We combined DFT and Griffith theory to estimate the fracture strength as a function of crack length. For a crack length of 2.0 μm, the fracture strength estimated by our method was 0.467 GPa, which was close to the experimental value of 0.250 GPa. These results indicate that the Young’s modulus is not sensitive to defects whereas the fracture strength is very sensitive to defects. Thus, it is important to take account of defects around which a large stress concentration arises when calculating fracture properties with DFT. Our method is effective for predicting both the fracture strength from crack length and the crack length from fracture strength. 1 B.
Steele and A. Heinzel, Nature 414, 345 (2001). J. Jacobson, Chem. Mater. 22, 660 (2010). 3 R. Nedelec, S. Uhlenbruck, D. Sebold, V. Haanappel, H.-P. Buchkremer, and D. Stoever, J. Power Sources 205, 157 (2012). 4 Y. Du, N. Sammes, G. Tompsett, D. Zhang, J. Swan, and M. Bowden, J. Electrochem. Soc. 150, A74 (2003). 5 B. Steele, Solid State Ionics 86–88, 1223 (1996). 6 H. Inaba and H. Tagawa, Solid State Ionics 83, 1 (1996). 7 M. Mogensen, N. Sammes, and G. Tompsett, Solid State Ionics 129, 63 (2000). 8 H. Yoshida, H. Deguchi, K. Miura, M. Horiuchi, and T. Inagaki, Solid State Ionics 140, 191 (2001). 9 W. Zajac and J. Molenda, Solid State Ionics 179, 154 (2008). 10 C. Sevik and T. Cagin, Phys. Rev. B 80, 014108 (2009). 11 T. Zhang, J. Ma, L. Kong, P. Hing, and J. Kilner, Solid State Ionics 167, 191 (2004). 12 T. Gurel and R. Eryigit, Phys. Rev. B 74, 014302 (2006). 13 L. Petit, A. Svane, Z. Szotek, and W. Temmerman, Phys. Rev. B 72, 205118 (2005). 14 R. Pornprasertsuk, P. Ramanarayanan, C. Musgrave, and F. Prinz, J. Appl. Phys. 98, 103513 (2005). 15 H. Yao, L. Ouyang, and W.-Y. Ching, J. Am. Ceram. Soc. 90, 3194 (2007). 2 A.
J. Chem. Phys. 140, 121102 (2014) 16 M.
Nakayama and M. Martin, Phys. Chem. Chem. Phys. 11, 3241 (2009).
17 S. Ito, T. Shimazaki, M. Kubo, H. Koinuma, and M. Sumiya, J. Chem. Phys.
135, 241103 (2011). Shimazaki, T. Suzuki, and M. Kubo, Theor. Chem. Acc. 130, 1031 (2011). 19 K. Serizawa, H. Onuma, H. Kikuchi, K. Suesada, M. Kitagaki, I. Yamashita, R. Miura, A. Suzuki, H. Tsuboi, N. Hatakeyama, A. Endou, H. Takaba, M. Kubo, H. Kajiyama, and A. Miyamoto, Jpn. J. Appl. Phys. 49, 04DJ14 (2010). 20 Z. Zhu, R. C. Deka, A. Chutia, R. Sahnoun, H. Tsuboi, M. Koyama, N. Hatakeyama, A. Endou, H. Takaba, C. A. Del Carpio, M. Kubo, and A. Miyamoto, J. Phys. Chem. Solids 70, 1248 (2009). 21 F. Ahmed, M. K. Alam, A. Suzuki, M. Koyama, H. Tsuboi, N. Hatakeyama, A. Endou, H. Takaba, C. A. Del Carpio, M. Kubo, and A. Miyamoto, J. Phys. Chem. C 113, 15676 (2009). 22 M. K. Alam, F. Ahmed, K. Nakamura, A. Suzuki, R. Sahnoun, H. Tsuboi, M. Koyama, N. Hatakeyama, A. Endou, H. Takaba, C. A. Del Carpio, M. Kubo, and A. Miyamoto, J. Phys. Chem. C 113, 7723 (2009). 23 G. Lu, S. Deng, T. Wang, M. Kohyama, and R. Yamamoto, Phys. Rev. B 69, 134106 (2004). 24 O. Nielsen and R. Martin, Phys. Rev. B 32, 3792 (1985). 25 J. Chen, Y. Xu, P. Rulis, L. Ouyang, and W. Ching, Acta Mater. 53, 403 (2005). 26 J. Chen, L. Ouyang, P. Rulis, A. Misra, and W. Ching, Phys. Rev. Lett. 95, 256103 (2005). 27 W. Y. Ching, P. Rulis, and A. Misra, Acta Biomater. 5, 3067 (2009). 28 V. Deyirmenjian, V. Heine, M. Payne, V. Milman, R. LyndenBell, and M. Finnis, Phys. Rev. B 52, 15191 (1995). 29 M. Kohyama, Phys. Rev. B 65, 184107 (2002). 30 W. Li and T. Wang, Phys. Rev. B 59, 3993 (1999). 31 M. Fronzi, A. Soon, B. Delley, E. Traversa, and C. Stampfl, J. Chem. Phys. 131, 104701 (2009). 32 M. Fronzi, S. Piccinin, B. Delley, E. Traversa, and C. Stampfl, Phys. Chem. Chem. Phys. 11, 9188 (2009). 33 Y. Peng, J. Li, L. Chen, J. Chen, J. Han, H. Zhang, and W. Han, Environ. Sci. Technol. 46, 2864 (2012). 34 N. J. Mosey and E. A. Carter, J. Mech. Phys. Solids 57, 287 (2009). 35 P. Liao and E. A. Carter, J. Mater. Chem. 20, 6703 (2010). 36 J. Diao, K. Gall, and M. Dunn, J. Mech. Phys. Solids 52, 1935 (2004). 37 C. Lee, X. Wei, J. W. Kysar, and J. Hone, Science 321, 385 (2008). 38 H. Bu, Y. Chen, M. Zou, H. Yi, K. Bi, and Z. Ni, Phys. Lett. A 373, 3359 (2009). 39 P. J. Hay, R. L. Martin, J. Uddin, and G. E. Scuseria, J. Chem. Phys. 125, 034712 (2006). 40 J. L. F. Da Silva, M. V. Ganduglia-Pirovano, J. Sauer, V. Bayer, and G. Kresse, Phys. Rev. B 75, 045121 (2007). 41 V. Esposito and E. Traversa, J. Am. Ceram. Soc. 91, 1037 (2008). 42 K. Sato, H. Yugami, and T. Hashida, J. Mater. Sci. 39, 5765 (2004). 43 H. Kim and M. Bush, Nanostruct. Mater. 11, 361 (1999). 44 Y. Jiang, J. B. Adams, and M. van Schilfgaarde, J. Chem. Phys. 123, 064701 (2005). 45 T. X. T. Sayle, B. J. Inkson, A. Karakoti, A. Kumar, M. Molinari, G. Moebus, S. C. Parker, S. Seal, and D. C. Sayle, Nanoscale 3, 1823 (2011). 46 E. Wong, P. Sheehan, and C. Lieber, Science 277, 1971 (1997). 47 D. P. Miannay, Fracture Mechanics, Mechanical Engineering Series (Springer, 1998). 48 A. Tsoga, A. Gupta, A. Naoumidis, and P. Nikolopoulos, Acta Mater. 48, 4709 (2000). 49 Y. Leng, S. Chan, S. Jiang, and K. Khor, Solid State Ionics 170, 9 (2004). 18 T.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
THE JOURNAL OF CHEMICAL PHYSICS 140, 121102 (2014)
Communication: Different behavior of Young’s modulus and fracture strength of CeO2 : Density functional theory calculations Ryota Sakanoi, Tomomi Shimazaki, Jingxiang Xu, Yuji Higuchi, Nobuki Ozawa, Kazuhisa Sato, Toshiyuki Hashida, and Momoji Kuboa) Fracture and Reliability Research Institute, Graduate School of Engineering, Tohoku University, 6-6-11 Aoba, Aramaki, Aoba-ku, Sendai 980-8579, Japan
(Received 28 January 2014; accepted 13 March 2014; published online 24 March 2014) In this Communication, we use density functional theory (DFT) to examine the fracture properties of ceria (CeO2 ), which is a promising electrolyte material for lowering the working temperature of solid oxide fuel cells. We estimate the stress-strain curve by fitting the energy density calculated by DFT. The calculated Young’s modulus of 221.8 GPa is of the same order as the experimental value, whereas the fracture strength of 22.7 GPa is two orders of magnitude larger than the experimental value. Next, we combine DFT and Griffith theory to estimate the fracture strength as a function of a crack length. This method produces an estimated fracture strength of 0.467 GPa, which is of the same order as the experimental value. Therefore, the fracture strength is very sensitive to the crack length, whereas the Young’s modulus is not. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4869515] I. INTRODUCTION
Solid oxide fuel cells (SOFCs), which generate electrical energy from the reaction of H2 with O2 , are regarded as a key industrial technology for clean energy conversion systems that do not emit pollutants and greenhouse gases such as NOx and CO2 .1, 2 However, the high working temperatures of SOFCs have so far prevented their widespread use. Standard zirconia electrolytes for SOFCs operate at 600 ◦ C–800 ◦ C.3 This means expensive heat-proof components are required, and it also increases the number of breakdowns and failures from repetitive start-stop processes causing large temperature differences.4, 5 Doped ceria (CeO2 ) is a promising electrolyte material for SOFCs with low working temperatures (500 ◦ C– 650 ◦ C).6 Reducing the working temperatures would make SOFCs a more competitive technology. There are many experimental studies on the physical properties of CeO2 -based materials, such as the ionic conductivity.7–9 The mechanical properties of CeO2 -based electrolytes are also important for the reliability of SOFCs.10, 11 It is essential to understand fracture properties at the atomic scale, including the bond stretching, bending, and breaking where fractures occur. Theoretical studies are required for analyzing fracture properties on an atomic scale because experimental studies alone are not sufficient. First-principles calculations are a powerful computational tool for understanding atomic-scale phenomena in fields such as chemistry, physics, and materials science.12–16 We have also used first-principles calculations to investigate metal-oxide systems such as ZnO,17 PtO2 ,18 MgO,19 SnO2 ,20 Al2 O3 ,21 and CeO2 .22 First-principles calculations based on density functional theory (DFT) techniques have also been used to investigate fracture properties such as Young’s modulus and fracture strength.23–27 Deyirmenjian et al. pioneered DFT-based tensile tests for solid-state materials by publishing a) [email protected]
0021-9606/2014/140(12)/121102/4/$30.00
the stress-strain curve for defective aluminum.28 Kohyama studied cubic SiC to explain the grain boundary effects.29 Li and Tang also investigated the elasticity and ideal strength of SiC with a plane-wave DFT-based tensile test.30 These theoretical tensile tests show the atomic-level structural changes and predict fracture properties under uniform loads. However, the actual stresses would be concentrated around cracks and there would be much higher stresses around meso-sized cracks. It is difficult to measure local stresses experimentally because only an averaged stress is measured against the applied strain. Therefore, it is necessary to predict and compare the fracture properties of an ideal material with no cracks and a realistic material with cracks. In this Communication, we develop a theoretical method for predicting the fracture properties of CeO2 containing a crack. First, we calculate the theoretical stress-strain curve for CeO2 without cracks, and use DFT to discuss the Young’s modulus and the fracture strength. Because first-principles methods usually only consider an ideal model with no defects such as cracks, they are inadequate for discussing the fracture properties of CeO2 . Thus, we model the effect of the cracks and estimate the real fracture strength by using a combination of DFT and Griffith theory. We explain our theoretical methods and calculation conditions in Sec. II, present and discuss the calculation results in Sec. III, and summarize our work in Sec. IV.
II. CALCULATION METHODS
We employ the Perdew-Wang 91 (PW91) functional for the DFT calculations in order to obtain the optimized geometrical structures of the CeO2 bulk model and their total energies. The electronic structures for core electrons are described by the effective core potential. We employed DMol3 software package for these DFT calculations, which has been used for various investigations of CeO2 such as surface defects,31
140, 121102-1
© 2014 AIP Publishing LLC
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
121102-2
Sakanoi et al.
J. Chem. Phys. 140, 121102 (2014)
water adsorption,32 and catalytic reduction.33 The validity of DMol3 was already confirmed by the previous studies. We used the double-numeric polarized basis set with a globalspace cutoff of 5.8 Å and periodic boundary conditions with a Monkhorst-Pack grid of 4 × 4 × 4. Carter and co-workers34, 35 showed the influence of tensile direction on mechanical properties of chromia and iron oxides by DFT calculations. In CeO2 , the structure is cubic system and symmetry. Thus, we do not need considering the effect of tensile direction. To calculate the Young’s modulus and fracture strength, we assume that the energy-density, Ecell /Vcell , can be described in the following polynomial equation as a function of the strain, , 1 2 (1) Y + A 3 + B 4 . 2 Here, Ecell and Vcell are the total energy and the volume of the unit cell, respectively. Y is the Young’s modulus, and A and B are the coefficients of the third- and fourth-order terms, respectively. In the equation, the first-order term disappears, because of the equilibrium structure at = 0, which is obtained by relaxing the atom positions and the lattice constants of the CeO2 model. The stress, σ , can be obtained from the first derivative of Eq. (1) Ecell /Vcell =
σ = Y + 3A 2 + 4B 3 .
(2)
This method is also used in other studies for the second- or third-order elastic modulus to estimate Young’s modulus.36–38 To determine the yield point and fracture strength, we take into account the polynomial terms up to the fourth-order term. The yield point can be determined from the condition δσ /δ = 0, and the yield point, yield . The fracture strength, σ yield , is obtained analytically by using the coefficients of the polynomial equation as follows: √ −3A − 9A2 − 12BY yield = , (3) 12B 2 3 + 4Byield , σyield = Y yield + 3Ayield
=
2 2 Y yield + Ayield . 3
(4) (5)
III. RESULTS AND DISCUSSIONS
First, we use a first-principles calculation technique to obtain the Young’s modulus and fracture strength of crackfree CeO2 from the variation of the total energy with the applied strain. In this calculation, the CeO2 bulk model consists of 4 Ce4 + and 8 O2 − atoms, as shown in Figure 1. The lattice constant of CeO2 is 5.4707 Å, which is optimized value at the strain of 0.00 by DFT calculation. The lattice constant depends on basis set and a Monkhorst-Pack grid value. The lattice constant value is 5.47 Å in PAW-PBE,39 FPLAPWPBE,39 and GTO-PBE,40 5.40 Å in PAW-HSE,39 5.49 Å in PAW-PBE+U,39 5.39 Å in PAW-PBE0,39 5.45 Å in GTOTPSS,40 and 5.41 Å in GTO-HSE.40 The experimental value is 5.41 Å.41 According to these results, our value of 5.47 Å is in agreement with the previous results. On the other hand, the
O
2-
Ce
4+
FIG. 1. CeO2 calculation model.
lattice constant does not depend on a Monkhorst-Pack grid value much because we calculated with a Monkhorst-Pack grid of 4 × 4 × 4 and 8 × 8 × 8 and the change of the lattice constant is about 0.0002%. Figure 2(a) shows the change in the energy density with the strain on the [001] axis, and each point is determined from the structural relaxation with a fixed lattice constant. Strains from 0.00 to 0.13 at 0.01 increments are applied. The Young’s modulus and coefficients of polynomial equation (2) are estimated by the least squares method, and are summarized in Table I. In this fitting procedure, we use strains from 0.00 to 0.13 in one direction. The Young’s modulus obtained from these calculations is 221.8 GPa on the [001] axis and the experimental Young’s modulus is 175 GPa.42 These values are similar. In this model, a crystalline material with no faults, such as grain boundaries or cracks, is used. Experimental materials usually contain many defects and the Young’s modulus is thought to decrease as the number of defects increases.43 Our calculated Young’s modulus is similar to the experimental Young’s modulus, suggesting that the effect of defects is small. Jiang et al.44 showed pioneered work on the effect of on-site Coulomb repulsion term U on the electronic structure and the bond length of CeO2 . Then, we also discussed the effect of on-site Coulomb repulsion term U on the fracture properties. We used CASTEP code for the calculation because DMol3 code, which was used for the calculation of fracture properties in this Communication, cannot consider on-site Coulomb repulsion term U. We employ CA-PZ functional for the DFT calculations by CASTEP code. We calculated the energy-densities at the strain of 0.00 and 0.09 with the standard DFT and DFT+U. Here, the strain of 0.00 and 0.09 corresponds to the origin and yield point of the stress-strain curve in Fig. 2, respectively. The reasonable U-J value is set to 6.3 eV and applied to 4f electrons of Ce atoms following to the earlier study.44 The difference TABLE I. Young’s modulus, Y, and coefficients of Eq. (2) for the energy density vs. strain curve of CeO2 based on the DFT method. Y (GPa)
A (GPa)
B (GPa)
221.8
1198.5
−9052.6
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
Sakanoi et al.
121102-3
J. Chem. Phys. 140, 121102 (2014)
Fracture strength [GPa]
(a) Ecell/Vcell [GPa]
3.0 Calculation result
2.0 1.0 0.0 0.00
Stress [GPa]
(b)
Fitting
0.05 0.10 Strain
0.15
3.0 2.5 2.0 1.5 1.0 0.5 0.0 0.0
1.0 2.0 Crack length [µm]
3.0
FIG. 3. Fracture strength, σ f , as a function of crack length, a.
25.0 20.0
(c)
15.0 10.0
(d) 2-
5.0 0.0 0.00
O 4+ Ce
0.05 0.10 Strain
0.15
FIG. 2. (a) Changes in energy density vs. strain. Black circles indicate the DFT calculation results, and the black solid line is the polynomial equation (2) fitted by the least squares method. (b) Stress-strain curve obtained from the derivative of Eq. (2). Figures (c) and (d) show CeO2 structures before and after the yield point, respectively.
of change in the energy density at the strain of 0.09 is 0.08 GPa. Therefore, the effect of Coulomb repulsion term U on the fracture property is not critical. Next, we calculate the fracture strength of crack-free CeO2 and compare it with the experimental fracture strength. The stress-strain curve is obtained from the first derivative of Eq. (2) as shown in Figure 2(b), where the fitting is calculated from strains of 0.00–0.13. The curve shows a yield point of 0.09 and a fracture strength of 22.7 GPa. Interestingly, the stress value is not zero after the curve reaches the yield point. At the yield point, the cerium atom, which has an eightfold oxygen coordination number in CeO2 , changes to a sixfold one with the strain (Figs. 2(c) and 2(d)). Thus, a phase transition from a fluorite CeO2 crystal structure is observed, which is in good agreement with the classical molecular dynamics simulation reported by Sayle et al.45 However, there are large differences between the calculated and experimental fracture strengths. A fracture strength of 22.7 GPa is obtained from the DFT calculation, whereas the experimental value is 0.250 GPa.42 The mechanical properties generally depend on factors such as the grain size, defects, flaws, cracks, voids, and grain boundaries. These factors are affected by experimental preparation techniques, thus it is difficult to directly compare calculated and experimental results. The theoretical tensile tests yield ideal values for the fracture strength, because simulation models do not possess any macro- or meso-sized cracks and voids. Wong et al. prepared a micrometer-scale SiC whisker,
which had no macro- or meso-sized cracks and no voids.46 They used atomic force microscopy to determine the experimental fracture strength as 53.4 GPa. Li and Wang30 estimated fracture strength of 50.8 GPa by using a DFT calculation. These results strongly suggest that DFT-based tensile tests are reliable, and provided that the samples do not contain any macro- or meso-sized cracks and voids, experiments can provide the ideal fracture properties that are calculated by DFT. However, to estimate the fracture properties of CeO2 materials, we must take into account the effect of defects. Next, we discuss the discrepancy between the calculation results and the experimentally observed fracture strength caused by cracks or other flaws. Stress is concentrated around the end of cracks contained in experimental materials, and the fracture strength is estimated by the following equation, based on Griffith theory,47 2γs Y σf = , (6) π a2 where σ f is the value of the stress at which fracture occurs, a is the length of the crack, γ s is the surface energy, and Y is the Young’s modulus. This equation indicates that fracture occurs when the incremental release of stored elastic strain energy caused by the increase in crack length becomes larger than the increase in surface energy from the creation of new surface area. We used a surface energy, γ s , of 1.545 J/m2 and a Young’s modulus, Y, of 221.8 GPa, which were obtained by DFT, in Griffith theory calculations. The surface energy is calculated by DFT method: γ s = (Esurface − Ebulk )/(2 × S), where Esurface is energy of CeO2 (001) with vacuum phase, Ebulk is energy of bulk model, and S is surface area. We employed the constant value against crack length. Figure 3 shows the dependence of fracture strength on crack length, a given by this equation. The fracture strength decreases sharply when the length of the crack is small. This indicates that the fracture strength is very sensitive to the crack length. A scanning electron microscopy (SEM) image of gadolinia-doped ceria published by Tsoga et al.48 showed that the defect size is about 1.0 μm. Using SEM, Leng et al.49 showed that the gadoliniadoped ceria defect was about 2.0 μm. According to these experimental data, we estimate the fracture strength as follows by using Eq. (6). When the crack lengths are 1.0 and 2.0 μm, the calculated fracture strengths are 0.660 and 0.467 GPa, respectively, which are close to the experimental value of
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
121102-4
Sakanoi et al.
0.250 GPa. When the fracture strength is 0.250 GPa, we estimate that the crack length is about 7.0 μm. Therefore, our method for estimating the fracture strength from the crack size and the crack size from the fracture strength is valid. Finally, we conclude that the fracture strength is very sensitive to crack length, although the Young’s modulus is not sensitive to defects; thus, it is important to take account of cracks in the theoretical estimation of the fracture strength. IV. SUMMARY
We have investigated the difference in the behavior of the Young’s modulus and fracture strength of CeO2 with DFT. The calculated Young’s modulus of 221.8 GPa was of the same order as the experimental value of 175 GPa. However, the fracture strength of 22.7 GPa was two orders of magnitude larger than the experimental value of 0.250 GPa. We combined DFT and Griffith theory to estimate the fracture strength as a function of crack length. For a crack length of 2.0 μm, the fracture strength estimated by our method was 0.467 GPa, which was close to the experimental value of 0.250 GPa. These results indicate that the Young’s modulus is not sensitive to defects whereas the fracture strength is very sensitive to defects. Thus, it is important to take account of defects around which a large stress concentration arises when calculating fracture properties with DFT. Our method is effective for predicting both the fracture strength from crack length and the crack length from fracture strength. 1 B.
Steele and A. Heinzel, Nature 414, 345 (2001). J. Jacobson, Chem. Mater. 22, 660 (2010). 3 R. Nedelec, S. Uhlenbruck, D. Sebold, V. Haanappel, H.-P. Buchkremer, and D. Stoever, J. Power Sources 205, 157 (2012). 4 Y. Du, N. Sammes, G. Tompsett, D. Zhang, J. Swan, and M. Bowden, J. Electrochem. Soc. 150, A74 (2003). 5 B. Steele, Solid State Ionics 86–88, 1223 (1996). 6 H. Inaba and H. Tagawa, Solid State Ionics 83, 1 (1996). 7 M. Mogensen, N. Sammes, and G. Tompsett, Solid State Ionics 129, 63 (2000). 8 H. Yoshida, H. Deguchi, K. Miura, M. Horiuchi, and T. Inagaki, Solid State Ionics 140, 191 (2001). 9 W. Zajac and J. Molenda, Solid State Ionics 179, 154 (2008). 10 C. Sevik and T. Cagin, Phys. Rev. B 80, 014108 (2009). 11 T. Zhang, J. Ma, L. Kong, P. Hing, and J. Kilner, Solid State Ionics 167, 191 (2004). 12 T. Gurel and R. Eryigit, Phys. Rev. B 74, 014302 (2006). 13 L. Petit, A. Svane, Z. Szotek, and W. Temmerman, Phys. Rev. B 72, 205118 (2005). 14 R. Pornprasertsuk, P. Ramanarayanan, C. Musgrave, and F. Prinz, J. Appl. Phys. 98, 103513 (2005). 15 H. Yao, L. Ouyang, and W.-Y. Ching, J. Am. Ceram. Soc. 90, 3194 (2007). 2 A.
J. Chem. Phys. 140, 121102 (2014) 16 M.
Nakayama and M. Martin, Phys. Chem. Chem. Phys. 11, 3241 (2009).
17 S. Ito, T. Shimazaki, M. Kubo, H. Koinuma, and M. Sumiya, J. Chem. Phys.
135, 241103 (2011). Shimazaki, T. Suzuki, and M. Kubo, Theor. Chem. Acc. 130, 1031 (2011). 19 K. Serizawa, H. Onuma, H. Kikuchi, K. Suesada, M. Kitagaki, I. Yamashita, R. Miura, A. Suzuki, H. Tsuboi, N. Hatakeyama, A. Endou, H. Takaba, M. Kubo, H. Kajiyama, and A. Miyamoto, Jpn. J. Appl. Phys. 49, 04DJ14 (2010). 20 Z. Zhu, R. C. Deka, A. Chutia, R. Sahnoun, H. Tsuboi, M. Koyama, N. Hatakeyama, A. Endou, H. Takaba, C. A. Del Carpio, M. Kubo, and A. Miyamoto, J. Phys. Chem. Solids 70, 1248 (2009). 21 F. Ahmed, M. K. Alam, A. Suzuki, M. Koyama, H. Tsuboi, N. Hatakeyama, A. Endou, H. Takaba, C. A. Del Carpio, M. Kubo, and A. Miyamoto, J. Phys. Chem. C 113, 15676 (2009). 22 M. K. Alam, F. Ahmed, K. Nakamura, A. Suzuki, R. Sahnoun, H. Tsuboi, M. Koyama, N. Hatakeyama, A. Endou, H. Takaba, C. A. Del Carpio, M. Kubo, and A. Miyamoto, J. Phys. Chem. C 113, 7723 (2009). 23 G. Lu, S. Deng, T. Wang, M. Kohyama, and R. Yamamoto, Phys. Rev. B 69, 134106 (2004). 24 O. Nielsen and R. Martin, Phys. Rev. B 32, 3792 (1985). 25 J. Chen, Y. Xu, P. Rulis, L. Ouyang, and W. Ching, Acta Mater. 53, 403 (2005). 26 J. Chen, L. Ouyang, P. Rulis, A. Misra, and W. Ching, Phys. Rev. Lett. 95, 256103 (2005). 27 W. Y. Ching, P. Rulis, and A. Misra, Acta Biomater. 5, 3067 (2009). 28 V. Deyirmenjian, V. Heine, M. Payne, V. Milman, R. LyndenBell, and M. Finnis, Phys. Rev. B 52, 15191 (1995). 29 M. Kohyama, Phys. Rev. B 65, 184107 (2002). 30 W. Li and T. Wang, Phys. Rev. B 59, 3993 (1999). 31 M. Fronzi, A. Soon, B. Delley, E. Traversa, and C. Stampfl, J. Chem. Phys. 131, 104701 (2009). 32 M. Fronzi, S. Piccinin, B. Delley, E. Traversa, and C. Stampfl, Phys. Chem. Chem. Phys. 11, 9188 (2009). 33 Y. Peng, J. Li, L. Chen, J. Chen, J. Han, H. Zhang, and W. Han, Environ. Sci. Technol. 46, 2864 (2012). 34 N. J. Mosey and E. A. Carter, J. Mech. Phys. Solids 57, 287 (2009). 35 P. Liao and E. A. Carter, J. Mater. Chem. 20, 6703 (2010). 36 J. Diao, K. Gall, and M. Dunn, J. Mech. Phys. Solids 52, 1935 (2004). 37 C. Lee, X. Wei, J. W. Kysar, and J. Hone, Science 321, 385 (2008). 38 H. Bu, Y. Chen, M. Zou, H. Yi, K. Bi, and Z. Ni, Phys. Lett. A 373, 3359 (2009). 39 P. J. Hay, R. L. Martin, J. Uddin, and G. E. Scuseria, J. Chem. Phys. 125, 034712 (2006). 40 J. L. F. Da Silva, M. V. Ganduglia-Pirovano, J. Sauer, V. Bayer, and G. Kresse, Phys. Rev. B 75, 045121 (2007). 41 V. Esposito and E. Traversa, J. Am. Ceram. Soc. 91, 1037 (2008). 42 K. Sato, H. Yugami, and T. Hashida, J. Mater. Sci. 39, 5765 (2004). 43 H. Kim and M. Bush, Nanostruct. Mater. 11, 361 (1999). 44 Y. Jiang, J. B. Adams, and M. van Schilfgaarde, J. Chem. Phys. 123, 064701 (2005). 45 T. X. T. Sayle, B. J. Inkson, A. Karakoti, A. Kumar, M. Molinari, G. Moebus, S. C. Parker, S. Seal, and D. C. Sayle, Nanoscale 3, 1823 (2011). 46 E. Wong, P. Sheehan, and C. Lieber, Science 277, 1971 (1997). 47 D. P. Miannay, Fracture Mechanics, Mechanical Engineering Series (Springer, 1998). 48 A. Tsoga, A. Gupta, A. Naoumidis, and P. Nikolopoulos, Acta Mater. 48, 4709 (2000). 49 Y. Leng, S. Chan, S. Jiang, and K. Khor, Solid State Ionics 170, 9 (2004). 18 T.
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 89.25.153.189 On: Fri, 04 Apr 2014 08:19:42
