Ab-Initio Investigation of Polarization and Ion Conduction Mechanisms in the Bone Mineral and Electret Material Hydroxyapatite
Authors: Kasamatsu Shusuke , Sugino Osamu 
Published: February 2020

Hydroxyapatite [Ca10(PO4)6(OH)2] (HAp) comprises up to 50% of human bone and is a subject of intense study in the context of biocompatibility for, e.g., bone transplants. It is also known to exhibit ionic conductivity, as well as piezo-, pyro-, and ferroelectricity, and this fact has stimulated research towards its possible application in various in vivo and ex vivo energy harvesting and nanoelectronics devices [1]. On a related note, it has also been shown that HAp can be turned into remarkably good “electrets”, or materials with a semi-permanent dipole. Such electrets are now of interest due to its potential application in vibration energy harvesting for Internet of Things (IoT), in addition to biological tissue engineering such as in vivo tissue culture under an external electric field.
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

Figure 1