HAp can be poled into an electret by applying a DC bias at elevated temperature and then quenching to room temperature. The electret formation, i.e., the freezing of the dipole, has been explained by the displacement of protons under bias at high temperature and subsequent quenching of its motion [2]. However, the microscopic polarization and ion conduction mechanisms remain unclear, hindering the optimization of this material for various applications mentioned above. In this work, we undertook a computational approach to this problem [3]. We first performed ab-initio molecular dynamics simulations in HAp with/without point defects to identify major polarization and ion conduction mechanisms. Then, the minimum energy paths and activation energies of those mechanisms were studied in detail using the climbing image nudged elastic band (CI-NEB) method [4]. VASP [5] code was used in combination with VTSTtools [6] for the calculations.
The Figure shows the hexagonal structure of hydroxyapatite often observed at elevated temperature. In our simulations, we used a supercell where this unit cell was multiplied by three times in the c-axis direction. Based on previous works on defect formation energetics [7, 8], we considered three types of defects, a proton vacancy, OH- vacancy, and proton interstitial, which are all located within the column of OH- ions along the c-axis. We then performed molecular dynamics simulations at 1000-1500 K for > 20 ps, where we can expect elementary processes with activation energies less than 1 eV to occur several times during the simulation according to transition state theory. From these simulations, we found that the dominant polarization/diffusion mechanisms mainly occur in the column of OH- ions aligned along the c-axis, and that the main processes are the flipping of the direction of OH- ions, exchange of protons between oxide ions, and the hopping of OH- vacancies. When proton interstitials are present as in the case of Ca-deficient HAp, interstitial proton hopping occurs mainly on the OH- column, although diffusion in the ab plane is also viable under bias. The activation energies for these processes were calculated using the CI-NEB method and are tabulated in the Table. The variations in the values originate from different environments around the local hopping/flipping site, which we estimate to be less than few tenths of an eV from the calculations that we have performed. The calculated activation energies are consistent with processes that have been detected in various temperature stimulated depolarization current experiments and conductivity measurements in the literature. This allows us to make a correspondence between measured activation energies and the dominant polarization/diffusion mechanisms at various experimental conditions. With such insight, materials processing conditions can be tailored for utilizing specific microscopic mechanisms in various applications.
[1] S. B. Lang, Phase Transitions 89, 678 (2016).
[2] Y. Tanaka et al., J. Appl. Phys. 107, 014107 (2010).
[3] S. Kasamatsu and O. Sugino, Phys. Chem. Chem. Phys. 20, 8744 (2018).
[4] G. Henkelman et al., J. Chem. Phys. 113, 9901 (2000).
[5] G. Kresse and J. Furthmueller, Phys. Rev. B 54, 11169 (1996).
[6] http://theory.cm.utexas.edu/vtsttools/index.html