Evaluation of dielectric and piezoelectric behavior of unpoled and poled barium titanate polycrystals with oxygen vacancies using phase field method

ABSTRACT This article studies the dielectric and piezoelectric behavior of unpoled and poled barium titanate (BaTiO3) polycrystals with oxygen vacancies. A phase field model is employed for BaTiO3 polycrystals, coupled with the time-dependent Ginzburg–Landau theory and the oxygen vacancies diffusion, to demonstrate the interaction between oxygen vacancies and domain evolutions. To generate grain structures, the phase field model for grain growth is also used. The hysteresis loop and butterfly curve are predicted at room and high temperatures. The permittivity, and longitudinal and transverse piezoelectric constants of the BaTiO3 polycrystals are then examined for various grain sizes and oxygen vacancy densities.


Introduction
Piezoelectric ceramics, polymers, and composites demonstrate varying extents of coupling between mechanical and electrical properties, and are recognized for their potential utility in energy-harvesting devices. Lead zirconate titanate (PZT) materials are, in particular, candidates for use as these devices [1,2]. However, due to lead toxicity, the introduction of legislation in Europe to limit the usage of lead in electronic products has led to a worldwide search for lead-free compounds [3]. BaTiO 3 ceramics is one of the widely applied lead-free ferroelectric materials due to high permittivity and low cost [4,5]. Recently, Zheng et al. [6] studied experimentally the effect of grain size on the permittivity and piezoelectric constant of the poled BaTiO 3 ceramics. Huan et al. [7] also prepared the BaTiO 3 ceramics with different grain sizes by two-step sintering method, and investigated the effect of grain size on the permittivity and piezoelectric constant. In the aforementioned two papers, it was shown that the permittivity and piezoelectric constant increase with decreasing average grain size and reach the maximum at about 1 μm. On the other hand, it has been shown that the electrical properties and microstructures of the BaTiO 3 ceramics depend on the oxygen partial pressure [8,9]. Although grain size and oxygen activity have been found to play a significant role in dielectric and piezoelectric behavior of BaTiO 3 ceramics, it is difficult to find the optimal conditions experimentally.
Phase field method is an efficient tool for interpreting the relationship between the effective property and microstructure of ferroelectric materials [10]. Zhang et al. [11] conducted the phase field simulations to study the oxygen-vacancy-induced memory effect in the BaTiO 3 single crystal. However, the investigations for the ferroelectric polycrystals are very few [12]. Recently, Shindo et al. [13] developed the phase field simulation framework for BaTiO 3 polycrystals with oxygen vacancies, and studied the effects of grain size and oxygen vacancy density on the electromechanical response of poled polycrystals.
This article constitutes a continuing study of the previous work [13] on the electromechanical behavior of BaTiO 3 polycrystals with oxygen vacancies. The phase field model of grain growth is used to generate two-dimensional grain structures, and a ferroelectric model coupled with the classical Ginzburg-Landau theory and the oxygen vacancies diffusion is applied. The hysteresis loop (polarization versus electric field) and butterfly curve (strain versus electric field) for unpoled BaTiO 3 polycrystals at room and high temperatures are first calculated, and the dielectric and piezoelectric behavior are evaluated. Then the effects of grain size and oxygen vacancy density on the permittivity, and longitudinal and transverse piezoelectric constants of the poled BaTiO 3 polycrystals are discussed. This information is essential for the establishment of optimal processing conditions, such as temperature, time of annealing, and gas-phase composition, of BaTiO 3 polycrystals.

Phase field model
Consider unpoled and poled BaTiO 3 polycrystals denoted by the global rectangular Cartesian coordinates , where superscript G denotes the global coordinates. The temporal evolution of the polarization vector component P G i and the domain structure is governed by the following time-dependent Ginzburg-Landau (TDGL) equation [14]: where F is the total free energy of the ferroelectric polycrystal, M is the kinetic coefficient related to the domain mobility, and t is the time. In Equation (1), δF=δP G i represents the thermodynamic driving force for the evolution of P G i . The total free energy can be expressed as follows: where f bulk , f grad , f elas , f elec are the local bulk free energy density, gradient energy density which is only nonzero around domain walls and grain boundaries, elastic strain energy density, electrostatic energy density, respectively; and V is the volume of the polycrystal.
The local bulk free energy density f bulk is commonly expressed by a polynomial of the polarization vector component as follows: where superscript L denotes the local coordinates, and α i , α ij , α ijk , and α ijkl are adjusted to the BaTiO 3 single crystal properties. The values of these properties [15] are given by the following equation: where T is the temperature in kelvin. The component of polarization vector in the local coordinate system is related to that in the global coordinate system by the following equation: The transformation matrix is given by the following equation: where θ, φ, and ψ are the Euler angles. The contribution of the domain walls to the total free energy is introduced through the gradient of the global polarization field, and the gradient energy density can be written as follows: where G 11 ; G 12 ; G 44 ; G 0 44 are the gradient energy coefficients. For simplicity, here, G 11 ¼ G 44 and G 12 ¼ G 0 44 ¼ 0 are assumed [12], and G 11 ¼ 0:6 Â 10 À7 m 4 N=C 2 [11] is used. The elastic energy density f elas can be expressed as follows: are the total strain tensor components [16] calculated based on the Khachaturyan's theory [17] 13 are the spontaneous strain tensor components, and c 11 ; c 12 ; c 44 are the elastic stiffnesses. The spontaneous strain components in the global system are related to that in the local system by the following equation: The spontaneous strain components in the local system are given by the following equation: where Q 11 ; Q 12 ; Q 44 are the electrostrictive coefficients. The elastic stiffnesses and electrostrictive coefficients for BaTiO 3 [18] are obtained as follows: The electrostatic energy density f elec is given by the following equation: where E G 1 ; E G 2 ; E G 3 denote the inhomogeneous electric field intensity vector components, E G 01 ; E G 02 ; E G 03 are the applied electric field intensity vector components, and f defect is the vacancy-induced depolarization energy density. A detailed method of solving the inhomogeneous electric field components can be found in Ref. [16] The vacancy-induced depolarization energy density [11] is given by the following equation: are the defect dipolar field vector components. These are given by the following equation: where 0 ¼ 8:85 Â 10 À12 F/m is the permittivity of free space, V 0 % 64 Â 10 À3 nm 3 is the unit cell volume, P defect = 0.515 C/m 2 is the local defect polarization due to the defect dipole, and N d is the oxygen vacancy density of BaTiO 3 at a given time t at each position ðx G 1 ; x G 2 ; x G 3 Þ. The oxygen vacancy density is governed by the following diffusion equation [19]: where β d = 5 Â 10 5 m 2 /sJ is the mobility of the defect, W d is the contribution to the free energy due to defects, e = 1.0621773 Â 10 -19 C is the coulomb charge per electron, Z = 1 is the valency of donors (oxygen vacancies), and ϕ is the electrical potential. The polarization evolution and oxygen diffusion are coupled only from the electrical potential ϕ.
To study the effect of the grain size on the dielectric and piezoelectric properties, a phase field model for grain growth is used and a two-dimensional grain structure is generated. The microstructure of a polycrystal unit cell is specified by a set of grain orientation field variables χ q (q = 1, . . ., Q) defined at a given time t at each position ðx G 1 ; x G 2 ; x G 3 Þwithin the unit cell [20]. The parameter Q denotes the number of distinct grain orientations. The TDGL equations are given by the following equation: where and M q (q = 1, . . ., Q) are kinetic rate coefficients related to the grain boundary mobility and assumed to have the values of M q = 1 μm 3 /Js (q = 1, . . ., Q). Equation (14) can be solved numerically at each position x G i of the unit cell to obtain the time evolution of the microstructure.

Numerical results and discussion
We consider a rectangular BaTiO 3 polycrystal of dimensions 3.3 Â 3.3 μm. In the simulation, 64 Â 64 discrete grid points were employed and periodic boundary conditions were applied along the x G 1 and x G 3 axes. The grid spacing was Δx G 1 =Δx G 3 % 52 nm. The time step for the evolutions in Equations (1), (13), and (14) were Δt/(1/|α 1 |M) = 0.005, Δt/(1/|α 1 |M) = 0.02, and Δt/(1/M q ) = 0.04 (q = 1, . . ., Q), respectively, to obtain reasonable and convergent results, where M = 1 C 2 /Nm 2 s and Q = 20, 30, or 100. Figure 1 shows three examples of simulated grain structures using the grain growth model. Two-dimensional simulations were performed for BaTiO 3 polycrystals without electric field and oxygen vacancy. The average grain radii are d = 0.38, 0.77, and 1.03 μm. We calculate the average polarization and normal strain in the x G 3 direction of unpoled and poled BaTiO 3 polycrystals under uniform electric field E G 03 =E 0 , and discuss the dielectric and piezoelectric behavior for the polycrystals with various grain sizes and oxygen vacancy densities. Figure 2 shows the polarization P G 3 versus electric field E 0 of unpoled BaTiO 3 polycrystals for average grain radius d = 0.89 μm and without oxygen vacancy (N d = 0 ppm) at T = 300 K = 27°C (room temperature) and T = 373 K = 100°C. The hysteresis loops are observed, and it is shown that the coercive electric field decreases as the temperature increases. Figure 3 shows the normal strain ε G 33 versus electric field E 0 of the same BaTiO 3 polycrystals. It is found that the butterfly curves are sensitive to the temperature. The coercive electric field decreases as the temperature increases, and therefore, it is better to perform poling process at a high temperature. This is a well-known fact, but by using the phase field simulation, it has been shown numerically.
It was shown that grain boundaries act as nucleation sites for new domains during switching [12] and a polycrystal with a smaller grain size results in a lower coercive electric field [12,13]. Here, we focus on the effect of oxygen vacancy density on the coercive electric field and polarization/strain versus electric field curve. Figure 4 shows the polarization P G 3 versus electric field E 0 of poled BaTiO 3 polycrystals for average grain radius d = 0.77 μm and oxygen vacancy N d = 0 and 100 ppm at T = 27°C. The coercive electric field decreases due to oxygen vacancy. Figure 5 shows the normal strain ε G 33 versus electric field E 0 of the same BaTiO 3 polycrystals. The slope of the curve under high positive electric fields increases with oxygen vacancy.
To obtain the dielectric and piezoelectric properties at room temperature, the variations of polarization and strain are evaluated under high positive electric fields. The permittivity 33 is defined as 33 ¼ ΔP G 3 =ΔE 0 . The longitudinal and transverse piezoelectric coefficients d 33 and d 31 are also defined as d 33 ¼ Δε G 33 =ΔE 0 and d 31 ¼ Δε G 11 =ΔE 0 , respectively. Figure 6 shows the permittivity 33 versus average grain radius d of poled BaTiO 3 polycrystals for various oxygen vacancy densities N d at T = 27°C. The permittivity is sensitive to the average grain radius and oxygen vacancy density as is expected. The BaTiO 3 polycrystals with N d = 0 and 500 ppm displayed a permittivity peak at  approximately d = 1 μm. This is in agreement with experimental findings [6,7]. On the other hand, for N d = 100 ppm, there is a possibility that in grains with radius less than d = 0.4 μm, the permittivity may display the largest value. In addition, it was found that in the case of N d = 1000 ppm, the value of permittivity was small overall, showing that it is not dependent on the grain radius. Figure 7 shows the longitudinal piezoelectric coefficient d 33 versus average grain radius d of BaTiO 3 polycrystals for various N d at T = 27°C. It was shown that overall, the longitudinal piezoelectric coefficient increases with an increase of the average grain radius and peaks, depending on the oxygen vacancy density. The value of d 33 for commercial poled BaTiO 3 bulk is 140 Â 10 -12 m/V [4], and the longitudinal piezoelectric coefficient of the BaTiO 3 polycrystals with N d = 500 ppm easily exceeds the value of 140 Â 10 -12 m/V. When it is approximately d = 0.9 μm with N d = 500 ppm, the longitudinal piezoelectric coefficient displays the greatest value. By improving the values of the grain radius and oxygen vacancy density, it may be possible to fabricate BaTiO 3 polycrystals that surpass the value of d 33 at 300 Â 10 -12 m/V. Our simulation results point to the possibility of fabricating piezoelectric ceramics that are comparable to PZT ceramics. Figure 8 shows similar results for the transverse piezoelectric coefficient d 31 . Likewise, the absolute value of d 31 peaked at approximately d = 0.9 μm with N d = 500 ppm.

Conclusions
The dielectric and piezoelectric behavior of unpoled and poled BaTiO 3 polycrystals was simulated using the phase field model, which takes into account the oxygen vacancies, and the following conclusions are obtained: (1) The hysteresis and butterfly loops depend on the temperature, and the coercive electric field decreases with increasing the temperature. (2) The permittivity depends on the average grain radius and oxygen vacancy density, and shows the largest value at average grain radius of approximately 1 μm with oxygen vacancy density of 500 ppm. (3) The longitudinal and transverse piezoelectric coefficients depend on the average grain radius and oxygen vacancy density, and at average grain radius of approximately 0.9 μm with oxygen vacancy density of 500 ppm, they show the largest values.
By knowing how the grain boundary and domain interact, we can design a BaTiO 3 polycrystal to improve the overall properties. The results of this study will help to provide useful information for optimizing the BaTiO 3 polycrystal performance by selecting the optimal poling process temperature, grain size, and oxygen vacancy density.

Acknowledgment
This work was supported by JSPS KAKENHI Grant Number 16H04227.

Disclosure statement
No potential conflict of interest was reported by the authors.

Funding
This work was supported by JSPS KAKENHI [Grant Number 16H04227].