Next Article in Journal
First-Principle Investigations on the Electronic and Transport Properties of PbBi2Te2X2 (X = S/Se/Te) Monolayers
Next Article in Special Issue
Quantitative Assessment of Chirality of Protein Secondary Structures and Phenylalanine Peptide Nanotubes
Previous Article in Journal
CD123-Targeted Nano-Curcumin Molecule Enhances Cytotoxic Efficacy in Leukemic Stem Cells
Previous Article in Special Issue
Dual-Spectral Plasmon-Induced Transparent Terahertz Metamaterial with Independently Tunable Amplitude and Frequency
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Iron in Hydroxyapatite: Interstitial or Substitution Sites?

1
Faculty of Physics, Southern Federal University, 5 Zorge St., 344090 Rostov-on-Don, Russia
2
Institute of Mathematical Problems of Biology, Branch of Keldysh Institute of Applied Mathematics, Russian Academy of Sciences, 1 Vitkevicha St., Pushchino, 142290 Moscow, Russia
3
I3N, Department of Physics, Campus Santiago, University of Aveiro, 3810-193 Aveiro, Portugal
4
Université Clermont Auvergne, Clermont Auvergne INP, CNRS, ICCF, F-63000 Clermont-Ferrand, France
*
Author to whom correspondence should be addressed.
Nanomaterials 2021, 11(11), 2978; https://0-doi-org.brum.beds.ac.uk/10.3390/nano11112978
Submission received: 2 October 2021 / Revised: 25 October 2021 / Accepted: 1 November 2021 / Published: 5 November 2021
(This article belongs to the Special Issue Simulation and Modeling of Nanomaterials)

Abstract

:
Iron-doped hydroxyapatite (Fe-HAp) is regarded as a promising magnetic material with innate biocompatibility. Despite the many studies reported in the literature, a detailed theoretical description of Fe inclusions is still missing. There is even no consensual view on what kind of Fe defects take place in Fe-HAp—iron interstitial or calcium substitutions? In order to address these questions, we employ modern first-principles methodologies, including hybrid density functional theory, to find the geometry, electronic, magnetic and thermodynamic properties of iron impurities in Fe-HAp. We consider a total of 26 defect configurations, including substitutional (phosphorus and calcium sites) and interstitial defects. Formation energies are estimated considering the boundaries of chemical potentials in stable hydroxyapatite. We show that the most probable defect configurations are: Fe 3 + and Fe 2 + substitutions of Ca(I) and Ca(II) sites under Ca-poor conditions. Conversely, Fe interstitials near the edge of the hydroxyl channel are favored in Ca-rich material. Substitutional Fe on the P site is also a probable defect, and unlike the other forms of Fe, it adopts a low-spin state. The analysis of Fe K-XANES spectra available in the literature shows that Fe-HAp usually contains iron in different configurations.

Graphical Abstract

1. Introduction

Ceramic materials based on calcium apatites are important substances in medicine, biology, ecology and catalysis [1,2,3,4]. Among the apatites, hydroxyapatite (Ca 10 (PO 4 ) 6 ( OH ) 2 = HAp) is singled out as the main inorganic component of vertebrate bones and teeth. Biogenic HAp differs from the synthetic compound due to the presence of foreign atomic (Mg 2 + , Na + , Cl ) and molecular (CO 3 2 , SiO 4 4 ) ions. Although the role of most impurities is largely unknown, some studies suggest that the presence of a small concentration of iron (below 157 ppm) prevents the loss of calcium from bones and desorption of enamel from teeth [5,6]. Iron-doped HAp (Fe-HAp) has also attracted interest due to its intrinsic magnetic properties [7]. It is a promising material for magneto-resonance imaging agents [8], heat centers for magnetic hyperthermia [4,9], sunscreen filter creams [10], antimicrobial coatings [11], and to be incorporated in drug delivery systems [12,13].
Despite nearly five decades of research on Fe-HAp and related systems [6,14,15,16,17,18,19], our knowledge regarding the atomistic and electronic structure and iron defects is still very limited. The main difficulty arises from the quality of the samples—single-phase products are hard to synthesize, iron metal and iron oxide often segregate at the surface, and the admixture of stable phosphates, in particular tricalcium phosphate (Ca 3 (PO 4 ) 2 = TCP), is usually observed [16,18]. Discrepancies between the results from different authors are often attributed to different preparation conditions, including sintering temperature and atmosphere, the duration and temperature of doping treatments, the source of precursors, etc. A comprehensive review of experimental methods and results can be found in Ref. [6].
The HAp structure contains two symmetrically distinct positions of calcium, Ca(I) and Ca(II) (Figure 1), both of which are candidates for iron replacement. The pioneering study of Fe-HAp using Mössbauer (MB) spectroscopy by Ok [14] indicated the preference of Fe Ca ( I ) over Fe Ca ( II ) , while subsequent MB studies [15,16] concluded on comparable amounts of Fe Ca ( I ) and Fe Ca ( II ) defects. In both cases, electron paramagnetic resonance (EPR) data indicated the presence of isolated high-spin Fe 3 + ions with 4 or 6 neighbors, in coexistence with superparamagnetic iron oxide on the surface. The MB experiments of Boda et al. [6] suggest that Fe Ca ( I ) 3 + defects are more probable when conventional sintering is used during sample preparation, whereas Fe Ca ( II ) 3 + substitutions appear as a result of hot pressing.
Interstitial iron (Fe i ) defects have been reported as well. Analysis of X-ray diffraction (XRD) data shows that iron can occupy 2 b or 12 i Wyckoff positions [5,18,19,20]. The first site corresponds to the formation of a linear O-Fe-O structure in the center of the OH channel (site A at Figure 1), while the latter involves a three-fold coordinated Fe atom displacement towards the wall of the channel (site C at Figure 1). However, from XRD and Raman scattering data, Kato et al. [19] argues on the formation of four-fold coordinated Fe. The absence of Fe Ca in these experiments was postulated based on the lack of CO 3 2 groups [6,16]. These are commonly present in biogenic HAp, and lead to the formation of several defects, including: [CO 3 2 ] [ PO 4 3 ] , [CO 3 2 ] [ OH ] , [O 2 ] [ OH ] and, finally, [Fe 3 + ] Ca 2 + . Here we limit our analysis to iron-only defects, leaving carbonate-related defects outside the scope of this work.
Besides the lack of agreement between different MB studies, XRD data also holds controversy due to a low electronic contrast between Fe 3 + and Ca 2 + ions [18]. Hence, atomic-level simulations are regarded as a source of clarification. Using empiric interatomic potentials, Jiang et al. [15] concluded that the Fe 2 + ion finds its favorite site at the Ca(I) position, while Fe 3 + favors the Ca(II) site. More recently, Zilm et al. [17] investigated this problem using semi-local density functional theory, but the results were limited to the 2+ oxidation state. Despite enlightening, both studies focused on Fe Ca substitutions only, with alternative substitutional and interstitial forms of iron being overlooked. Hence, we may conclude that the starting point of the present work is a largely unknown picture of Fe impurities in HAp regarding their local structure, oxidation and magnetic states, and optoelectronic properties.
Recently we demonstrated [21] that the use of hybrid density functional theory (hybrid-DFT) can describe the electronic properties of HAp with an accuracy that rivals with highly accurate Green’s function (G) method with screened interaction (W), so called G W approximation [22]. The G W and hybrid-DFT methods provided band gaps of HAp close to 7 eV, but generalized gradient functional resulted in a 30 % smaller gap [23,24,25,26,27,28]. Errors of the same magnitude are therefore expected for the calculated optoelectronic properties of defects in HAp, as shown for the case of the OH vacancy [21].
In view of the above, we present a hybrid-DFT study of the most probable Fe-related defects in hydroxyapatite (substitutional and interstitial). The results are interpreted and scrutinized from the perspective of the available experimental data. The paper is organized as follows: Section 2 outlines the methodologies and computational details. Section 3 discusses the results, subdivided in Subsections pertaining the chemical phase diagram of HAp (Section 3.1), the structure and thermodynamics of neutral and charged defects (Section 3.3 and Section 3.4, respectively), and a simulation of the Fe K-XANES spectra of Fe-HAp available in the literature (Section 3.5). Finally, Section 4 summarizes the results.

2. Computation Methods

HAp crystallizes in the form of an molecular ionic crystal of symmetry P 6 3 / m (#176 in the crystallographic tables), whose hexagonal unit cell encloses two Ca 5 ( PO 4 ) 3 OH units [29,30]. Figure 1 illustrates the supercell used in this study, comprising 8 unit cells (16 HAp formula units). HAp primitive cell contains two and four inequivalent calcium and oxygen sites, respectively. The Ca(I) cation columns are surrounded by O(I) and O(II) from PO 4 anion groups, while mirror-symmetric O(III) sites and Ca(II) ions form a hexagonal channel enclosing the hydroxyl anions O(IV)H. Variation of alignment of OH dipoles in these channels lead to different HAp phases. Namely: (i) a hexagonal disordered phase, with random orientations of OH dipoles; (ii) a hexagonal ordered phases, where OH dipoles are all oriented along the same direction, or (iii) a monoclinic ( P 2 1 / b ) phase, made of two cells repeated along a basal lattice vector ( a or b ), with the first possessing a ...-OH-OH-... column, while the second one showing an opposite ...-HO-HO-... ordering [31]. The last phase shows anti-ferroelectric properties [32]. Given that the electronic band structure of the above polytypes is rather similar [27], the impact of OH-flipping on the problems being investigated is expected to be minor. Henceforth, we will consider the introduction of iron defects in the ordered hexagonal phase.

2.1. First-Principles Defect Energetics

The calculations were performed using density functional theory using the QUANTUM ESPRESSO package [33,34]. The many-body electronic potential was evaluated using the hybrid density functional of Heyd, Scuseria and Ernzerhof (HSE) [35,36]. The Kohn-Sham eigenstates obtained using PBE [37] exchange-correlation functional were used to initiate the self-consistent cycle. Due to severe limitations of the semi-local functionals (such as PBE) in describing the electronic structure of HAp in the band gap vicinity [21], the application of hybrid-DFT turns out to be critical in the evaluation of the stability of charged defects.
By opting for a local treatment of the electronic correlation, we expect some error in the description of the d-shell in Fe. Using a GGA + U approach would lead to a narrow gap material, and hence to an artificial hybridization between the d-levels of Fe and those from the conduction and valence band edges of HAp. Unfortunately, it was not possible to combine the non-local HSE functional, with the Hubbard-like local correction to the electronic correlation on the Fe atom.
Core electron states were described by means of optimized norm-conserving Vanderbilt pseudopotentials [38,39], while the Kohn-Sham problem was solved within the plane wave formalism with cut-off energy up to 60 and 240 Ry to expand the wave functions and semi-local potential. The exact exchange operator was evaluated on a grid which corresponded to a plane wave cut-off of 120 Ry. The occupation of states in the vicinity of the Fermi energy were smeared out with a Gaussian function of width 0.002 Ry to improve convergence.
The equilibrium lattice parameters of the unit cell were obtained a = 9.481 Å and c = 6.859 Å using HSE exchange-correlation functional. These values are particularly close to a = 9.417 Å and c = 6.875 Å from an experimental report [30].
Defects were introduced in orthorhombic supercells made up of 8 HAp primitive cells (spacegroup P 6 3 / m ), containing a total of 352 atoms. The smooth dispersion of the band structure of the large HAp supercell allowed us to use a single point ( k = Γ ) for sampling the Brillouin zone (BZ). Convergence tests showed that the total energy of such supercells obtained with Γ only sampling, differs by less than 0.1 eV (0.3 meV/atom) from a calculation where the zone was sampled over a 2 × 2 × 2 grid of k -points.
Defect containing supercells were first relaxed on a PBE level, since claculation of forces on a hybrid-DFT level is prohibitively expensive for large systems. During the relaxations, coordinates of all atoms were varied until the maximum force became less than 0.4 mRy/Bohr 10 meV/Å. On a second step the total energy was calculated on hybrid-DFT level for obtained structures. It was shown [40,41,42] that relative errors in energies obtained within this methodology usually have an order of 10 meV.
The formation energy E f of a defect in a crystalline sample can be expressed as [43]:
E f = E d ( q , R ) E HAp i Δ n i μ i + q E v + E F ,
where the main contribution is the energy difference between defective ( E d ) and pristine ( E HAp ) crystals. R and q denotes the defect configuration and charge state, correspondingly. Last two terms account for stoichiometric and charge differences between the defect and pristine cells. In particular, μ i is the chemical potential of species i which must be added ( Δ n i > 0) or removed ( Δ n i < 0) to or from the ideal crystal to create the defect, respectively. The last term in Equation (1) accounts for the exchange of electrons between the defect with charge q and an electron reservoir with chemical potential μ e = E v + E F . Here E v and E F denotes the energy of the valence band top and Fermi level, respectively. The Fermi energy can vary within the band gap ( 0 E F E g ), where E g is the band gap width), depending on the doping of the crystal.
The variation of chemical potentials μ i is also limited,
i n i ϕ · μ i μ i 0 Δ H f ϕ ,
with the upper limit taking place if HAp becomes unstable with respect to formation of a compound ϕ made of n i ϕ elements of species i and with heat of formation Δ H f ϕ . The chemical potentials of O, H, Ca and P in standard phases, μ i 0 , were calculated from the energy per atom in molecular oxygen, molecular hydrogen, calcium metal and black phosphorous. The molecules (O 2 in the spin-triplet state or H 2 ) were placed in a simulation cubic box of 20 Å size. The chemical potential of iron was found from the bcc-Fe phase with the BZ sampled over a 16 × 16 × 16   k -point mesh. The resulting ground state spin density corresponded to a magnetic moment of 2.80 μ B per cell. The deviation from the experimental value of 2.22 μB results from a known over-localization of the exchange interactions in bulk Fe when using hybrid functionals [44]. An identical result was obtained using projected augmented-wave potentials. Nevertheless, it is noted that hybrid functionals have successfully been used for the study of iron oxides and defects in oxides (see for instance [45,46]).
The heat of formation of the compounds ϕ used in Equation (2) was estimated from their hybrid-DFT total energies ( E ϕ ) with respect to energies of their constituents in their standard phases ( μ i 0 ),
Δ H f ϕ = E ϕ i n i ϕ μ i 0 .
Due to periodic boundary conditions, the calculation of charged defects ( q 0 ) is accompanied by a compensating uniform charge density of opposite sign [47]. The artificial interactions between the periodically repeated charged defects and the background lead to the deviation of calculated periodic total energy E ˜ d from E d , in Equation (1). The energy correction E corr , so that E d = E ˜ d + E corr , is obtain following the recipe of Lany and Zunger [48],
E corr 2 3 α M q 2 ϵ L ,
where α M is the Madelung constant of the HAp supercell with edge length L, and ϵ is the static dielectric constant ( ϵ 11 [21]). For a singly charged defect the correction amounts to about 0.04 eV. However, this figure can grow considerably in the case of multi-ionized defects ( E corr 0.33 eV for q = + 3 ).
The magnetization of the Fe defects was calculated from spin-polarized electron densities n ( r ) and n ( r ) ,
M z = Ω ( n n ) d 3 r ,
while formal oxidation states of iron ions were deduced from the occupations of d-orbitals as proposed by Sit et al. [49].

2.2. X-ray Absorption near Edge Structure (XANES)

In order to compare the defect structures found from the first-principles total energy calculations, with those reported in the literature, we simulated the Fe K-XANES spectra for each structure. The calculations were performed using a full-multiple scattering method and the “muffin-tin” approximation for the interatomic potential as implemented in the FDMNES code [50]. The size of the atomic cluster and spectral convolution parameters were adjusted using the spectra of magnetite Fe 3 O 4 as reference. The radius of the cluster was estimated of 7 Å which gives cluster of Fe 1 Ca 30 P 15 O 70 H 5 , with exact composition dependent on the structure. The used approach cannot reproduce pre-edge features, which is a well known limitation of a description based on dipolar transitions and the single-particle approximation [51,52].

3. Results and Discussion

3.1. Chemical Stability Diagram

The formation of the most probable iron-related defects in HAp depends most notably on the chemical potentials of the several elements involved. Each chemical potential can vary within a certain range, limited by thermodynamic stability conditions of the HAp crystal itself. We estimate the ranges using Equation (2) with respect to a set of boundary phases { ϕ } . The methodology used to find the structure and energy of each phase was identical to that used for the HAp supercell. This includes the exchange-correlation functional and energy cutoffs. The first candidates for the bordering phases are those involved in the HAp synthesis. There is a variety of production routes, but we will only consider reactants whose species are common to those found in HAp, namely CaO [6], Ca(OH) 2 [8], H 3 PO 4 , H 2 O and P 2 O 5 [18]. Next, we considered calcium phosphates: dicalcium phosphate dihydrate CaHPO 4 · 2 H 2 O = DCPD (mineral brushite), anhydrous dicalcium phosphate CaHPO 4 = DCPA (mineral monetite), tricalcium phosphate Ca 3 (PO 4 ) 2 = TCP and tetracalcium phosphate Ca 4 (PO 4 ) 2 O = TTCP. The last two have [Ca]: [P] concentration ratio of 1.5 and 2.0, respectively, while HAp has an intermediate [Ca]: [P] ratio of 5 / 3 1.67 .
Table 1 presents the formation energies Δ H f , calculated according to Equation (3) within HSE, for the materials enumerated above. The overall agreement with the reference data [53,54] (including the energy, cell volume and bulk modulus) is in line with the usual accuracy of hybrid-DFT. The phosphoric acid and water were calculated in a gas phase (single molecules in large periodic box), hence the lack of cell volume and bulk modulus.
From the data of Table 1 and Equation (2), we can estimate the chemical potential ranges within which HAp is thermodynamically stable. We can also reduce the number of independent chemical potentials, for instance to μ H , μ Ca and μ P , by expressing the chemical potential of oxygen as a function of μ i of the remaining elements (c.f. Equation (2)),
μ O μ O 0 = Δ H f HAp i = Ca , P , H n i HAp ( μ i μ i 0 ) / n O HAp .
Let us look first at the range of hydrogen chemical potentials. In this case, phases with stoichiometry (not) satisfying n H ϕ · n O HAp < n O ϕ · n H HAp give lower (upper) bounds for μ H . In particular, the planes corresponding to CaO, TCP and pure O 2 may form lower bounds, while planes for hydrogen-containing phases may be responsible for upper bounds.
According to our calculations, the HAp stability domain is limited by CaO, Ca(OH ) 2 , P, H 2 , DCPD and DCPA phase planes. Figure 2 illustrates the domain of chemical stability of HAp as a convex faceted hull with vertices P i (with i = 1 , , 11 ). These correspond to regions in phase space where three different phases coexists with HAp. Table S1 in Supplementary Materials contains the calculated values of chemical potentials at these intersection points. We note that there is no boundary with Ca metal, with the upper limit for μ Ca μ Ca 0 being 1.8 eV. This means that if μ Ca = μ Ca 0 was assumed, formation energies of defects involving Ca substitutions would be affected by an error of 2 eV. In the analysis regarding substitutional iron defects, we will consider the following chemical phase space conditions:
  • Ca-rich and P-rich, points P 7 and P 8 , μ Ca μ Ca 0 = 1.8 eV, μ P μ P 0 = 0 ;
  • Ca-poor, point P 3 , μ Ca μ Ca 0 = 7.8 eV, μ P μ P 0 = 7.9 eV;
  • P-poor, points P 1 and P 2 , μ Ca μ Ca 0 = 6.0 eV, μ P μ P 0 10.7 eV.

3.2. Notation for Defect Structures

Before describing the structure of the iron defects let us first establish a few notation rules. We denote the defect structure as M S n , where M is a substituting species, most notably iron. However, this can also be a vacancy V or group of atoms (e.g., OH). S denotes the position of the defect, which can be a crystalline atomic site (Ca(I), Ca(II) or P) or an interstitial site (i(A), i(B), etc). The flipping of an hydroxyl unit is denoted as OH HO . The number n is the formal oxidation state of iron, estimated from the occupation of the d-orbitals [49]. Table 2 shows the correspondence between the used notation, Kröger-Vink notation and the net charge of the defective supercell. We found a systematic correspondence between the formal oxidation state of iron atom and the defect charge state. This indicates that the gap states of the Fe defects under investigation (which can trap electrons or holes) belong to the 3d shell.

3.3. Neutral Iron Defects

We start by considering neutral iron defects based on structures previously proposed in the literature. Table 3 summarizes the relevant data collected. The first seven rows correspond to the Fe Ca substitutions, while the next five to iron interstitials. There is no apparent correlation between the defect structure and sample preparation conditions (e.g., Fe concentration or synthesis temperature).
We considered Ca(I), Ca(II), as well as P substitutions, the later being unexplored in the literature. Figure 3a–c illustrate the substitutional defect structures obtained after relaxation, while Figure 3d–g depict the interstitial structures. The relaxation of the structure where iron was initially set up with the Fe i ( A ) configuration (Wyckoff position 2 b ), resulted in another structure, where iron is displaced away from the center of the channel. This new structure is denoted as Fe i ( C ) (Figure 3f) and approximately corresponds to Wyckoff position 12 i . In order to achieve a stable Fe i ( A ) configuration, one of the nearest OH groups has to be flipped (Figure 3d). Another interstitial position with the iron ion displaced from the center of the channel is labeled as Fe i ( B ) (Figure 3e) and corresponds to Wyckoff position 2 c . The Fe i ( C ) structure may be slightly modified by flipping of the neighboring OH group, leading to a minute ( 0.1 eV) benefit in the total energy. Such small difference is lower than the flipping of isolated OH, estimated as 0.22 eV. Our result suggest that OH flipping may be stimulated by the presence of interstitial iron in the channel.
Finally, we considered iron inserted in the region between PO 4 groups (Wyckoff position 6 g ), marked as Fe i ( D ) (Figure 3g). At this location, the negatively charged PO 4 3 groups are expected to screen the positive charge of the iron ion. Although we scanned the relative stability of other structures, those whose formation energy was above 7 eV were discarded and not investigated further (e.g., Fe OH and Fe i in the vicinity of the Ca(I) column).
Table 4 shows the formation energies of the most stable defects, estimated according to Equation (1) using chemical potentials at extreme points of the HAp stability diagram, namely for material grown under Ca- and P-rich ( P 7 and P 8 ), Ca-poor ( P 3 ) and P-poor ( P 1 ) conditions. At Ca- and P-rich conditions, the Fe i ( A ) 0 defect is likely to be the most probable as it shows the lowest formation energy. Other stable defects (within less that 1 eV above the ground state) are Fe i ( B ) 0 and Fe Ca 2 + substitutions. Other interstitial sites ( Fe i ( C ) 0 , Fe i ( D ) 0 ) have higher formation energies. However, they still should be considered since hole or electron trapping may stabilize them.
At Ca-poor conditions Fe Ca 2 + defects are the most stable, with formation energy 5 eV below that of Fe i ( A ) 0 . A low value of μ P at Ca-poor conditions ( P 3 in Figure 2) leads to an easier depletion of phosphorus and to the stabilization of Fe P 5 + as well. Subsequently, at P-poor conditions this effect is further enhanced and Fe P 5 + becomes nearly 4 eV more stable than Fe Ca 2 + , and about 7.5 eV than the most favorable interstitial defect, Fe i ( A ) 0 .
Table 4 reports the number of oxygen atoms neighboring iron (coordination number), their respective Fe-O distances, and the magnetic moment of each defect. The six-fold coordination of Fe Ca ( I ) (Figure 3a) is in line with the results of Jiang et al. [15], Low et al. [16] and Zilm et al. [17]. The five-fold coordinated Fe Ca ( II ) structure differs: six-fold coordinated iron is proposed by Jiang et al, while four-fold – by Zilm et al. However, we may count six neighbors if we include the far oxygen neighbor at distance R Fe - O = 2.55 Å marked with dashed line on Figure 3b. Zilm et al obtained two oxygen atoms at high distances of R Fe - O = 2.73 and 2.97 Å, but excluded them from the neighbors count, resulting in four-fold coordination.
The structure of Fe i ( A ) with two O(IV) neighbors shows Fe-O distances 0.2 Å longer than those observed by XRD [18]. This may result from the XRD analysis, which gives an ordered structure where oxygen atoms are fixed to crystallographic sites, or may reflect the distances of a positively charged state, where the iron cation is closer to the oxygen anions. Fe i ( C ) configurations have four-fold coordinated iron. This is at variance with the XRD study of Gomes et al. [18], where three-fold coordinated iron was found, but is in line with XRD results of Kato et al. [19]. The only three-fold coordination of Fe that we find is in the Fe i ( B ) structure (Figure 4e). However, the geometry differs from that proposed in Ref. [18], where Fe connects to one oxygen atom from a neighboring PO 4 group and two O atoms from OH groups. In general, the calculated Fe-O distances look slightly overestimated. Of course, the picture could improve when considering charged cells.

3.4. Charged Iron Defects

Iron (Fe:4s 2 3d 6 ) is not isoelectronic with respect to the species being replaced (Ca:4s 2 and P:3s 2 3p 3 ). Hence, Fe impurities are expected to create states within the band gap of HAp. These states are localized and may act as hole or electron traps. The description of the Fe-HAp requires the consideration these cases, which can be accounted for by changing the occupation of the highest occupied gap states of the defective supercell. In a real crystal, the capture of a hole (creation of a local positive charge) can be compensated by numerous possibilities: cation vacancies V Ca , V H , foreign anion interstitials like [CO 3 2 ] i , etc., thus leading to a lowering of the Fermi energy (electron chemical potential), We leave the exact mechanisms of charge compensation out of the scope of this work.
An qualitative picture of the charge states allowed for each defect can be found from the respective Kohn-Sham levels in the band gap [61]. Positive charge states (hole trapping) requires the presence of filled states in the gap, while the negatively charged defects (electron trapping) requires the presence of empty states. Figure 4 illustrates the energy of the one-electron states in the band gap obtained from spin-polarized calculations of substitutional and interstitial defects in the neutral charge state.
Substitutional iron on Ca sites, Fe Ca 2 + , have six electrons in the 3d shell, but after the first ionization all filled states move below the valence band top. Hence, further ionization would require an energy equivalent to the band gap width (∼7 eV), which essentially tell us that Fe Ca 2 + defects can be single donors, but not double donors. Unoccupied states of Fe Ca 2 + are rather close to the conduction band, suggesting that they are not acceptors either (can not trap electrons).
Iron on the phosphorous site in the Fe P 5 + state (neutral charge state) has three filled one-electron states and 7 empty states in the gap (3d 3 configuration). We will show below that Fe P 5 + can trap one hole or up to two electrons, becoming Fe P 6 + or Fe P 4 + and Fe P 3 + states, respectively.
Neutral interstitial defects show only filled states in the gap. These can trap up to three holes, thus leading to Fe i + , Fe i 2 + and Fe i 3 + , respectively.
We found that the structure of the defects depend strongly on the local charge. Figure 5 illustrates some of the most remarkable structural changes of the defects induced by the capture of electrons and holes. Additional reconfigurations are depicted in Figure S1 in Supplementary Information. Table 5 presents the number of nearest oxygen neighbors to the Fe ion ( N Fe - O ), the range of Fe-O bond lengths ( R Fe - O ), as well as the magnetic moment of the defect ( M z ). Note, that some states exhibit a high magnetization, up to 5 μ B for Fe Ca ( I ) 3 + and Fe i ( C ) 3 + . These are to be compared to the case of neutral defects, where the magnetic moment was at most 2 μ B .
In general, an increase of the positive charge leads to a decrease of Fe-O bond lengths. This effect is strikingly illustrated in Figure 5h,i, which showcase the Fe P 3 + –Fe P 6 + sequence of defects. This rule is understandably violated when there is a change in the coordination of the Fe ion, and therefore, a significant modification of the local electrostatics. Examples are the increases of R Fe - O in Fe i ( B ) 0 , Fe i ( C ) + , or Fe i ( D ) 2 + upon hole capture,
Fe i ( B ) 0 + h + Fe i ( B ) + (Figure 5f,g),
Fe i ( C ) + + h + Fe i ( C ) 2 + (Figure S1p,q),
Fe i ( D ) 2 + + h + Fe i ( D ) 3 + (Figure S1y,z).
The above processes are reversible, i.e., electron trapping at (or hole emission from) the positively charged defects result in the recovery of the longer bond lengths. As expected, HAp cations repel the iron ion. In the case of Fe i ( C ) 3 + without OH flipping, this effect leads to migration of a proton from OH to a close PO 4 group (see Figure 5d,e). We found that the high-symmetry defect configuration Fe i ( A ) 2 + is not stable and spontaneously transforms to Fe i ( C ) 2 + . Hence, upon hole capture by (or electron emission from) Fe i ( A ) + , the iron moves toward the edge of the OH channel (see Figure 5a,b),
Fe i ( A ) + + h + Fe i ( C ) 2 + + OH HO .
The reverse process involves overcoming a barrier (not calculated), and Fe i ( A ) + is not recovered spontaneously (see Figure 5b,c).
The formation energy of a charged defect is a function of the Fermi energy (c.f. Equation (1)). This dependence is clearly illustrated in Figure 6a,b, which show the results for Ca- and P-rich material and for Ca-poor HAp, respectively. The red shadow area on both diagrams indicates the whole range for the formation energy of Fe P . The upper bound corresponds to P-rich conditions, while the lower bound to P-poor material. The solid red line in Figure 6b shows the formation energy of Fe P under Ca-poor conditions. Thick dashed-dotted lines correspond to formation energies of Fe Ca defects, and thick solid line and other non-solid lines correspond to interstitial defects (see legend).
According to our results, the phosphorous substitutions, especially Fe P 3 + (charge state q = 2 ) is a rather stable species when there is abundance of electrons in the material (n-type HAp). This is even more evidenced in Ca-poor and P-poor conditions, where in p-type material we expect Fe P 4 + , Fe P 5 + and even Fe P 6 + to become more stable and compete with other species, namely Fe i and Fe Ca .
However, the phosphorous substitutions were not previously considered in the literature, since the replacing of phosphorus by iron cation requires the breaking of pretty stable P–O bonds of PO 4 group. Alternatively, the replacement of whole PO 4 group by FeO 4 one will provide the same structural result. The obtained most probable Fe P configuration is tetrahedral coordinated ferric cation ([FeO 4 ] 5 group) is not unusual and can be found in magnetite [52] or in iron-phosphate glass [62].
Regarding substitutional Fe at calcium sites, we find that Fe Ca ( I ) 3 + is the most probable form under Ca-poor and p-type conditions (Figure 6b). We note that despite Fe Ca ( II ) 2 + having lower formation energy than Fe Ca ( I ) 2 + in n-type HAp, the Fe P 3 + species is even more stable, thus making the phosphorus substitution more probable than a replacement of calcium.
In contrast, under Ca- and P-rich conditions (Figure 6a) the interstitial defects are expected to prevail, especially in p-type and intrinsic HAp. Depending on the Fermi level location, the most probable states are Fe i ( C ) 3 + , Fe i ( C ) 2 + and Fe i ( A ) + defects. The two-fold coordinated Fe i ( A ) 0 defect, is more favorable in n-type HAp, yet again, Fe P 3 + is more stable and more likely to form.
We can compare the obtained local atomic structure of iron defects in HAp (see Table 4 and Table 5) with those reported in the literature (Table 3). Iron substitutions with long distances R Fe - O 2.2 Å and 6 oxygen neighbors [15,16] are best described by Fe Ca ( I ) 2 + . The slightly more favorable Fe Ca ( II ) 2 + substitution has 5 oxygen neighbors, which could explain the result of Boda et al. [6]. The iron with only two neighbors in the study of Gomes et al. [18] could be described by Fe i ( A ) + , although the small Fe-O distances of 1.7 Å are only reproduced by Fe P 6 + , which is rather unstable. Iron with three [18] or four [19] oxygen neighbors at distances R Fe - O 1.8 Å can be accounted for by Fe i ( C ) 3 + or Fe P 3 + defects.
In summary, we find that Fe-HAp can contain both substitutional and interstitial defects depending on the preparation conditions. The phosphorus substitutions have iron in low-spin states ( M z μ B ,), making them less useful for many applications envisaged for Fe-HAp. To avoid those defects the synthesis should be performed closer to P-rich and p-type conditions. In that case we expect the formation of high-spin defects with M z = 5 μ B (Fe Ca ( I ) 3 + and Fe i ( C ) 3 + ) and M z = 4 μ B (Fe i ( C ) 2 + ).

3.5. Fe K-XANES of Fe-HAp

The near-edge structure of X-ray absorption (XANES) spectra is particularly sensitive to the details of local atomic structure of the absorbing element [51]. X-ray absorption spectroscopy has been applied to materials without long range order, and that includes the Fe-HAp. We consider the experimental spectra of Fe-HAp published by Gomes et al. [18]. We also keep the notation of the original study regarding the experimental conditions, i.e., 15Fe-500 and 15Fe-1100 which correspond to 15 mol % of Fe per HAp unit cell of samples sintered at 500 °C and 1100 °C, respectively.
The top three (black colored) curves in Figure 7 shows the experimental data of Gomes et al. [18] for Fe-HAp and for magnetite (Fe 3 O 4 ). The latter was used as a reference for the alignment of the theoretical energy scale to the experimental one. Greek letters mark the main spectral features: pre-edge ( α ), the main peak ( β ), its satellite peaks ( β and β ), and a more distant peak ( γ ). The vertical dashed line in Figure 7 provides guidance for the relative positions of the minimum between features β and γ .
In order to determine the types of Fe defects in the 15Fe-500 and 15Fe-1100 samples we simulated the Fe K-XANES spectra for each structure considered in Section 3.3 and Section 3.4. Figure 7 shows the simulated spectra of the most probable defect structures. Figure S2 in Supplementary Information shows all calculated spectra. The comparison of experimental and simulated spectra for magnetite Fe 3 O 4 (dotted curves in Figure 7) reveals the main insufficiencies of the simulation method, including (i) the lack of pre-edge features ( α ) and (ii) a poor reproduction of the β satellite. However, the intensities and energy positions of the main features ( β , β and γ ) of the magnetite spectrum are correctly reproduced. All spectral features of another reference iron oxide, hematite Fe 2 O 3 , could be reproduced (not shown) using exact the same calculation scheme. In this case, the differences between experimental and simulated spectra show qualitatively the same insufficiencies as for magnetite.
The simulated spectra for Fe-HAp show good sensitivity with respect to their atomistic structures and charge states. In most cases the increase of the iron oxidation state shifts the position of the main peak ( β ) to higher energies. This follows from an increase of the local positive charge and binding energy of K-electrons. An exception is seen for the Fe i ( D ) 3 + spectrum (Figure S2d), which is explained by significant changes in its local atomic structure (coordination change from 4 to 6) upon ionization.
The relative intensities I β and I β , corresponding to features β and β from the simulated spectra, are sensitive to the defect type: phosphorus substitutions have I β < I β , calcium substitutions show I β > I β , and most interstitial defects have I β I β . Additionally, the spectra show considerable differences concerning the position and shape of the γ feature.
The spectrum of the 15Fe-1100 sample shows comparable β and β intensities ( I β I β ), suggesting that most Fe defects are of interstitial character. Conversely, the higher intensity of peak β in the spectrum of sample 15Fe-500 may indicate the presence of substitutional Fe on calcium sites or Fe i ( A ) n defect ( n = 0 , 1 + ). The latter scenario is in line with the conclusions of Gomes et al. [18], where the formation of the high-symmetry Fe i ( A ) structure was proposed in the 15Fe-500 and 15Fe-800 samples, whereas Fe i ( C ) -type defects were suggested for the high-temperature treated 15Fe-1100 sample.

4. Conclusions

We performed a hybrid-DFT study of iron defects in hydroxyapatite. The most favorable structures (and respective charge states) were used to simulate Fe K-XANES spectra of Fe-HAp. These were compared to experimental data available in the literature. The results allow us to draw the following conclusions:
  • The chemical stability of HAp (Ca 10 (PO 4 ) 6 ( OH ) 2 ) is limited by a thermodynamic equilibrium with O 2 , CaO, Ca(OH) 2 , H 2 , P, DCPD (CaHPO 4 · 2 H 2 O), DCPA (CaHPO 4 ) and TCP (Ca 3 (PO 4 ) 2 ) phases. Formation of these compounds limit the range of values of Ca and P chemical potentials, and we could identify three extreme cases: Ca- and P-rich, Ca-poor, P-poor.
  • Under P-poor conditions, phosphorous substitutions are the most favorable, resulting in the formation of low-spin Fe defect states. However, p-type HAp may contain high spin iron interstitials Fe i ( C ) 3 + .
  • Under Ca-poor conditions, the high-spin calcium substitution, Fe Ca ( I ) 3 + , is the most favorable species. Fe Ca 2 + defects have higher formation energy comparing to Fe P 3 + and Fe Ca ( I ) 3 + .
  • Under Ca- and P-rich conditions, interstitial iron atoms in the OH channel are the most prominent. Depending on the position of Fermi level the most favorable are Fe i ( C ) 3 + , Fe i ( C ) 2 + and Fe i ( A ) + . When compared to Fe i ( C ) , the last structure involves the flipping of a nearby hydroxyl unit. The OH flipping does not introduce a significant change to the formation energy of Fe i ( C ) 3 + and Fe i ( C ) 2 + defects.
  • High-spin iron defects are Fe Ca ( I ) 3 + and Fe i ( C ) 3 + . These are both expected in p-type HAp. Such configurations are expected to be most useful for materials targeting magnetic hyperthermia or magnetic resonance imaging applications.
  • The comparison of Fe K-XANES spectra of theoretically predicted defect structures with experimental data [18] confirms the interstitial character of iron defects in samples sintered at high (1100 °C) temperature, but does not exclude the substitution defects for samples sintered at lower temperatures.

Supplementary Materials

The following are available online at https://0-www-mdpi-com.brum.beds.ac.uk/article/10.3390/nano11112978/s1, Table S1: Chemical potentials in eV at selected points from the phase stability diagram, Figure S1: Illustration of the local atomic structure of iron atom (≲4 Å) for all considered Fe-HAp configurations, Figure S2: Simulated Fe K-XANES spectra for the considered Fe-HAp structures. Pseudopotential files, example input script and quantum espresso output file. Atomic coordinates of Fe-HAp structures in POSCAR format.

Author Contributions

Conceptualization, J.C. and V.B.; methodology, L.A. and E.P.; software, E.P. and J.C.; validation, S.G., G.R. and J.C.; investigation, L.A. and E.P.; resources, V.B. and J.C.; writing—original draft preparation, L.A.; writing—review and editing, S.G., G.R. and J.C.; visualization, L.A.; funding acquisition, V.B. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by Russian Science Foundation (RSF) grant number 21-12-00251. J.C. thanks the FCT in Portugal for support through Projects UIDB/50025/2020, UIDP/50025/2020.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The atomic coordinates of optimized Fe-HAp structures are available in POSCAR format through Supplementary Materials.

Acknowledgments

Authors are grateful to the I3N at the University of Aveiro and to the Institute of Mathematical Problems of Biology by providing access to the HPC facilities. Authors acknowledge SLS for provision of synchrotron radiation facilities and would like to thank Maarten Nachtegaal for assistance in using beamline SuperXAS.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

The following abbreviations are used in this manuscript:
HApHydroxyapatite
Fe-HApIron-doped HAp
DFTDensity functional theory
HSEHeyd, Scuseria and Ernzerhof / Hybrid Screened Exchange
XANESX-ray absorption near-edge structure

References

  1. Epple, M.; Ganesan, K.; Heumann, R.; Klesing, J.; Kovtun, A.; Neumann, S.; Sokolova, V. Application of calcium phosphate nanoparticles in biomedicine. J. Mater. Chem. 2010, 20, 18–23. [Google Scholar] [CrossRef]
  2. Bystrov, V. Computational studies of the hydroxyapatite nanostructures, peculiarities and properties. Math. Biol. Bioinform. 2017, 12, 14–54. [Google Scholar] [CrossRef]
  3. Duminis, T.; Shahid, S.; Hill, R.G. Apatite Glass-Ceramics: A Review. Front. Mater. 2017, 3, 59. [Google Scholar] [CrossRef] [Green Version]
  4. Degli Esposti, L.; Carella, F.; Adamiano, A.; Tampieri, A.; Iafisco, M. Calcium phosphate-based nanosystems for advanced targeted nanomedicine. Drug Dev. Ind. Pharm. 2018, 44, 1223–1238. [Google Scholar] [CrossRef]
  5. Mayer, I.; Diab, H.; Felner, I. Ferric iron in synthetic carbonate apatites: A Mössbauer effect study. J. Inorg. Biochem. 1992, 45, 129–133. [Google Scholar] [CrossRef]
  6. Boda, S.K.; Basu, B.; Anupama, A.V.; Sahoo, B. Structural and Magnetic Phase Transformations of Hydroxyapatite-Magnetite Composites under Inert and Ambient Sintering Atmospheres. J. Phys. Chem. C 2015, 119, 6539–6555. [Google Scholar] [CrossRef]
  7. Mondal, S.; Manivasagan, P.; Bharathiraja, S.; Moorthy, M.S.; Kim, H.H.; Seo, H.; Lee, K.D.; Oh, J. Magnetic hydroxyapatite: A promising multifunctional platform for nanomedicine application. Int. J. Nanomed. 2017, 12, 8389–8410. [Google Scholar] [CrossRef] [Green Version]
  8. Laranjeira, M.S.; Moço, A.; Ferreira, J.; Coimbra, S.; Costa, E.; Santos-Silva, A.; Ferreira, P.J.; Monteiro, F.J. Different hydroxyapatite magnetic nanoparticles for medical imaging: Its effects on hemostatic, hemolytic activity and cellular cytotoxicity. Colloids Surf. B Biointerfaces 2016, 146, 363–374. [Google Scholar] [CrossRef]
  9. Tampieri, A.; D’Alessandro, T.; Sandri, M.; Sprio, S.; Landi, E.; Bertinetti, L.; Panseri, S.; Pepponi, G.; Goettlicher, J.; Bañobre-López, M.; et al. Intrinsic magnetism and hyperthermia in bioactive Fe-doped hydroxyapatite. Acta Biomater. 2012, 8, 843–851. [Google Scholar] [CrossRef]
  10. Piccirillo, C.; Rocha, C.; Tobaldi, D.M.; Pullar, R.C.; Labrincha, J.A.; Ferreira, M.O.; Castro, P.M.L.; Pintado, M.M.E. A hydroxyapatite–Fe2O3 based material of natural origin as an active sunscreen filter. J. Mater. Chem. B 2014, 2, 5999–6009. [Google Scholar] [CrossRef]
  11. Gambardella, A.; Bianchi, M.; Kaciulis, S.; Mezzi, A.; Brucale, M.; Cavallini, M.; Herrmannsdoerfer, T.; Chanda, G.; Uhlarz, M.; Cellini, A.; et al. Magnetic hydroxyapatite coatings as a new tool in medicine: A scanning probe investigation. Mater. Sci. Eng. C 2016, 62, 444–449. [Google Scholar] [CrossRef] [PubMed]
  12. Lin, K.; Chen, L.; Liu, P.; Zou, Z.; Zhang, M.; Shen, Y.; Qiao, Y.; Liu, X.; Chang, J. Hollow magnetic hydroxyapatite microspheres with hierarchically mesoporous microstructure for pH-responsive drug delivery. CrystEngComm 2013, 15, 2999. [Google Scholar] [CrossRef]
  13. Gu, L.; He, X.; Wu, Z. Mesoporous Fe3O4 /hydroxyapatite composite for targeted drug delivery. Mater. Res. Bull. 2014, 59, 65–68. [Google Scholar] [CrossRef]
  14. Ok, H.N. Mössbauer Studies of Natural Apatite. Phys. Rev. 1969, 185, 477–482. [Google Scholar] [CrossRef]
  15. Jiang, M.; Terra, J.; Rossi, A.M.; Morales, M.A.; Saitovitch, E.M.B.; Ellis, D.E. Fe2+/Fe3+substitution in hydroxyapatite: Theory and experiment. Phys. Rev. B 2002, 66, 224107. [Google Scholar] [CrossRef]
  16. Low, H.R.; Phonthammachai, N.; Maignan, A.; Stewart, G.A.; Bastow, T.J.; Ma, L.L.; White, T.J. The Crystal Chemistry of Ferric Oxyhydroxyapatite. Inorg. Chem. 2008, 47, 11774–11782. [Google Scholar] [CrossRef]
  17. Zilm, M.E.; Chen, L.; Sharma, V.; McDannald, A.; Jain, M.; Ramprasad, R.; Wei, M. Hydroxyapatite substituted by transition metals: Experiment and theory. Phys. Chem. Chem. Phys. 2016, 18, 16457–16465. [Google Scholar] [CrossRef]
  18. Gomes, S.; Kaur, A.; Grenèche, J.M.; Nedelec, J.M.; Renaudin, G. Atomic scale modeling of iron-doped biphasic calcium phosphate bioceramics. Acta Biomater. 2017, 50, 78–88. [Google Scholar] [CrossRef]
  19. Kato, S.; Ikeda, S.; Saito, K.; Ogasawara, M. Fe incorporation into hydroxyapatite channels by Fe loading and post-annealing. J. Solid State Chem. 2018, 265, 411–416. [Google Scholar] [CrossRef]
  20. Renaudin, G.; Gomes, S.; Nedelec, J.M. First-row transition metal doping in calcium phosphate bioceramics: A detailed crystallographic study. Materials 2017, 10, 92. [Google Scholar] [CrossRef] [Green Version]
  21. Avakyan, L.A.; Paramonova, E.V.; Coutinho, J.; Öberg, S.; Bystrov, V.S.; Bugaev, L.A. Optoelectronics and defect levels in hydroxyapatite by first-principles. J. Chem. Phys. 2018, 148, 154706. [Google Scholar] [CrossRef] [Green Version]
  22. Martin, R.; Reining, L.; Ceperley, D. Interacting Electrons; Cambridge University Press: Cambridge, UK, 2016. [Google Scholar]
  23. Calderin, L.; Stott, M.; Rubio, A. Electronic and crystallographic structure of apatites. Phys. Rev. B 2003, 67, 134106. [Google Scholar] [CrossRef] [Green Version]
  24. Rulis, P.; Ouyang, L.; Ching, W.Y. Electronic structure and bonding in calcium apatite crystals: Hydroxyapatite, fluorapatite, chlorapatite, and bromapatite. Phys. Rev. B 2004, 70, 155104. [Google Scholar] [CrossRef]
  25. Rulis, P.; Yao, H.; Ouyang, L.; Ching, W.Y. Electronic structure, bonding, charge distribution, and x-ray absorption spectra of the (001) surfaces of fluorapatite and hydroxyapatite from first principles. Phys. Rev. B 2007, 76, 245410. [Google Scholar] [CrossRef]
  26. Matsunaga, K.; Kuwabara, A. First-principles study of vacancy formation in hydroxyapatite. Phys. Rev. B 2007, 75, 014102. [Google Scholar] [CrossRef] [Green Version]
  27. Slepko, A.; Demkov, A.A. First-principles study of the biomineral hydroxyapatite. Phys. Rev. B 2011, 84, 134108. [Google Scholar] [CrossRef] [Green Version]
  28. Bystrov, V.S.; Coutinho, J.; Bystrova, A.V.; Dekhtyar, Y.D.; Pullar, R.C.; Poronin, A.; Palcevskis, E.; Dindune, A.; Alkan, B.; Durucan, C.; et al. Computational study of hydroxyapatite structures, properties and defects. J. Phys. D Appl. Phys. 2015, 48, 195302. [Google Scholar] [CrossRef]
  29. Kay, M.; Young, R.; Posner, A. Crystal Structure of Hydroxyapatite. Nature 1964, 204, 1050–1052. [Google Scholar] [CrossRef] [PubMed]
  30. Hughes, J.M.; Rakovan, J. The Crystal Structure of Apatite, Ca5(PO4)3(F,OH,Cl). Rev. Mineral. Geochem. 2002, 48, 1–12. [Google Scholar] [CrossRef]
  31. Elliott, J.; Mackie, P.; Young, R. Monoclinic hydroxyapatite. Science 1973, 180, 1055–1057. [Google Scholar] [CrossRef] [PubMed]
  32. Ikoma, T.; Yamazaki, A.; Nakamura, S.; Akao, M. Preparation and dielectric property of sintered monoclinic hydroxyapatite. J. Mater. Sci. Lett. 1999, 18, 1225–1228. [Google Scholar] [CrossRef]
  33. Giannozzi, P.; Baroni, S.; Bonini, N.; Calandra, M.; Car, R.; Cavazzoni, C.; Ceresoli, D.; Chiarotti, G.L.; Cococcioni, M.; Dabo, I.; et al. QUANTUM ESPRESSO: A modular and open-source software project for quantum simulations of materials. J. Phys. Condens. Matter 2009, 21, 395502. [Google Scholar] [CrossRef] [PubMed]
  34. Giannozzi, P.; Andreussi, O.; Brumme, T.; Bunau, O.; Nardelli, M.B.; Calandra, M.; Car, R.; Cavazzoni, C.; Ceresoli, D.; Cococcioni, M.; et al. Advanced capabilities for materials modelling with Quantum ESPRESSO. J. Phys. Condens. Matter 2017, 29, 465901. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  35. Heyd, J.; Scuseria, G.E.; Ernzerhof, M. Hybrid functionals based on a screened Coulomb potential. J. Chem. Phys. 2003, 118, 8207–8215. [Google Scholar] [CrossRef] [Green Version]
  36. Krukau, A.V.; Vydrov, O.A.; Izmaylov, A.F.; Scuseria, G.E. Influence of the exchange screening parameter on the performance of screened hybrid functionals. J. Chem. Phys. 2006, 125, 224106. [Google Scholar] [CrossRef]
  37. Perdew, J.P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865–3868. [Google Scholar] [CrossRef] [Green Version]
  38. Hamann, D.R. Optimized norm-conserving Vanderbilt pseudopotentials. Phys. Rev. B 2013, 88, 085117. [Google Scholar] [CrossRef] [Green Version]
  39. Schlipf, M.; Gygi, F. Optimization algorithm for the generation of ONCV pseudopotentials. Comput. Phys. Commun. 2015, 196, 36–44. [Google Scholar] [CrossRef] [Green Version]
  40. Coutinho, J.; Torres, V.J.; Demmouche, K.; Öberg, S. Theory of the carbon vacancy in 4 H-SiC: Crystal field and pseudo-Jahn-Teller effects. Phys. Rev. B 2017, 96, 174105. [Google Scholar] [CrossRef] [Green Version]
  41. Bystrov, V.S.; Avakyan, L.A.; Paramonova, E.V.; Coutinho, J. Sub-Band Gap Absorption Mechanisms Involving Oxygen Vacancies in Hydroxyapatite. J. Phys. Chem. C 2019, 123, 4856–4865. [Google Scholar] [CrossRef] [Green Version]
  42. Gouveia, J.D.; Coutinho, J. Can we rely on hybrid-DFT energies of solid-state problems with local-DFT geometries? Electron. Struct. 2019, 1, 015008. [Google Scholar] [CrossRef] [Green Version]
  43. Qian, G.X.; Martin, R.M.; Chadi, D.J. First-principles study of the atomic reconstructions and energies of Ga- and As-stabilized GaAs(100) surfaces. Phys. Rev. B 1988, 38, 7649–7663. [Google Scholar] [CrossRef]
  44. Jang, Y.R.; Yu, B.D. Structural, Magnetic, and Electronic Properties of Fe: A Screened Hybrid Functional Study. J. Magn. 2011, 16, 201–205. [Google Scholar] [CrossRef]
  45. Meng, Y.; Liu, X.W.; Huo, C.F.; Guo, W.P.; Cao, D.B.; Peng, Q.; Dearden, A.; Gonze, X.; Yang, Y.; Wang, J.; et al. When Density Functional Approximations Meet Iron Oxides. J. Chem. Theory Comput. 2016, 12, 5132–5144. [Google Scholar] [CrossRef] [PubMed]
  46. Ingebrigtsen, M.E.; Varley, J.B.; Kuznetsov, A.Y.; Svensson, B.G.; Alfieri, G.; Mihaila, A.; Badstübner, U.; Vines, L. Iron and intrinsic deep level states in Ga2O3. Appl. Phys. Lett. 2018, 112, 042104. [Google Scholar] [CrossRef] [Green Version]
  47. Makov, G.; Payne, M.C. Periodic boundary conditions in ab initio calculations. Phys. Rev. B 1995, 51, 4014–4022. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  48. Lany, S.; Zunger, A. Assessment of correction methods for the band-gap problem and for finite-size effects in supercell defect calculations: Case studies for ZnO and GaAs. Phys. Rev. B 2008, 78, 235104. [Google Scholar] [CrossRef] [Green Version]
  49. Sit, P.H.L.; Car, R.; Cohen, M.H.; Selloni, A. Simple, Unambiguous Theoretical Approach to Oxidation State Determination via First-Principles Calculations. Inorg. Chem. 2011, 50, 10259–10267. [Google Scholar] [CrossRef]
  50. Bunau, O.; Joly, Y. Self-consistent Aspects of X-ray Absorption Calculations. J. Phys. Condens. Matter 2009, 21, 345501. [Google Scholar] [CrossRef]
  51. Lamberti, C.; van Bokhoven, J.A. X-ray Absorption and X-ray Emission Spectroscopy: Theory and Applications. In X-ray Absorption and X-ray Emission Spectroscopy; Chapter X-ray Absorption and Emission Spectroscopy for Catalysis; van Bokhoven, J.A., Lamberti, C., Eds.; John Wiley & Sons, Ltd.: Hoboken, NJ, USA, 2016; pp. 351–383. [Google Scholar] [CrossRef]
  52. Soldatov, M.A.; Göttlicher, J.; Kubrin, S.P.; Guda, A.A.; Lastovina, T.A.; Bugaev, A.L.; Rusalev, Y.V.; Soldatov, A.V.; Lamberti, C. Insight from X-ray Absorption Spectroscopy to Octahedral/Tetrahedral Site Distribution in Sm-Doped Iron Oxide Magnetic Nanoparticles. J. Phys. Chem. C 2018, 122, 8543–8552. [Google Scholar] [CrossRef]
  53. Lide, D.R. (Ed.) CRC Handbook of Chemistry and Physics, 88th ed.; CRC Press: Boca Raton, FL, USA, 2005; p. 2661. [Google Scholar]
  54. Gražulis, S.; Daškevič, A.; Merkys, A.; Chateigner, D.; Lutterotti, L.; Quiros, M.; Serebryanaya, N.R.; Moeck, P.; Downs, R.T.; Le Bail, A. Crystallography Open Database (COD): An open-access collection of crystal structures and platform for world-wide collaboration. Nucleic Acids Res. 2012, 40, D420–D427. [Google Scholar] [CrossRef] [PubMed]
  55. Fukui, H.; Ohtaka, O.; Fujisawa, T.; Kunisada, T.; Suzuki, T.; Kikegawa, T. Thermo-elastic property of Ca(OH)2portlandite. High Press. Res. 2003, 23, 55–61. [Google Scholar] [CrossRef]
  56. Brazhkin, V.V.; Gromnitskaya, E.L.; Danilov, I.V.; Katayama, Y.; Lyapin, A.G.; Popova, S.V. High pressure behavior of P2O5crystalline modifications: Compressibility, elastic properties and phase transitions. Mater. Res. Express 2015, 2, 025201. [Google Scholar] [CrossRef]
  57. Zhai, S.; Liu, X.; Shieh, S.R.; Zhang, L.; Ito, E. Equation of state of -tricalcium phosphate, -Ca3(PO4)2, to lower mantle pressures. Am. Mineral. 2009, 94, 1388–1391. [Google Scholar] [CrossRef]
  58. Hughes, J.; Fransolet, A.; Schreyer, W. The atomic arrangement of iron-bearing apatite. In Neues Jahrbuch Fur Mineralogie-Monatshefte; Schweizerbart Science Publishers: Stuttgart, Germany, 1993; pp. 504–510. [Google Scholar]
  59. Carrera, K.; Huerta, V.; Orozco, V.; Matutes, J.; Fernández, P.; Graeve, O.; Herrera, M. Formation of vacancy point-defects in hydroxyapatite nanobelts by selective incorporation of Fe3+ ions in Ca(II) sites. A CL and XPS study. Mater. Sci. Eng. B 2021, 271, 115308. [Google Scholar] [CrossRef]
  60. Bauminger, E.; Ofer, S.; Gedalia, I.; Horowitz, G.; Mayer, I. Iron uptake by teeth and bones: A mossbauer effect study. Calcif. Tissue Int. 1985, 37, 386–389. [Google Scholar] [CrossRef] [PubMed]
  61. Tsuneda, T.; Song, J.W.; Suzuki, S.; Hirao, K. On Koopmans’ theorem in density functional theory. J. Chem. Phys. 2010, 133, 174101. [Google Scholar] [CrossRef]
  62. Stoch, P.; Szczerba, W.; Bodnar, W.; Ciecinska, M.; Stoch, A.; Burkel, E. Structural properties of iron-phosphate glasses: Spectroscopic studies and ab initio simulations. Phys. Chem. Chem. Phys. 2014, 16, 19917–19927. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Figure 1. Top-projection (a) and side-view (b) of an orthorhombic 2 × 2 × 2 HAp supercell with the hexagonal unit cell highlighted. O, Ca and H atoms are shown in red, green and white, respectively. Each PO 4 group contains O(I), O(II) and two O(III) sites. Capital letters A-D indicate the location of the interstitial sites considered.
Figure 1. Top-projection (a) and side-view (b) of an orthorhombic 2 × 2 × 2 HAp supercell with the hexagonal unit cell highlighted. O, Ca and H atoms are shown in red, green and white, respectively. Each PO 4 group contains O(I), O(II) and two O(III) sites. Capital letters A-D indicate the location of the interstitial sites considered.
Nanomaterials 11 02978 g001
Figure 2. Domain of chemical stability of HAp represented by the limits of chemical potentials of Ca, P, H. Symbols P 1 P 11 correspond to points where three boundary phases coexists with the HAp crystal.
Figure 2. Domain of chemical stability of HAp represented by the limits of chemical potentials of Ca, P, H. Symbols P 1 P 11 correspond to points where three boundary phases coexists with the HAp crystal.
Nanomaterials 11 02978 g002
Figure 3. Illustration of neutral defect structures in Fe-HAp. For each defect, we show a top view of its location within supercell and a lateral view of its detailed atomic structure in the vicinity of ∼4 Å around the iron atom.
Figure 3. Illustration of neutral defect structures in Fe-HAp. For each defect, we show a top view of its location within supercell and a lateral view of its detailed atomic structure in the vicinity of ∼4 Å around the iron atom.
Nanomaterials 11 02978 g003
Figure 4. Energy level diagram based on the Khon-Sham states of Fe defects in HAp.
Figure 4. Energy level diagram based on the Khon-Sham states of Fe defects in HAp.
Nanomaterials 11 02978 g004
Figure 5. Illustration of the change in the local atomic structure of selected iron defects (ligands within ≲4 Å from the Fe atom are also depicted), following the capture of electrons or holes.
Figure 5. Illustration of the change in the local atomic structure of selected iron defects (ligands within ≲4 Å from the Fe atom are also depicted), following the capture of electrons or holes.
Nanomaterials 11 02978 g005
Figure 6. Formation energies of Fe-HAp as a function of the Fermi energy at extreme points of HAp stability diagram: (a) Ca-rich and P-rich, (b) Ca-poor. The horizontal segments correspond to neutral simulation cells (Table 4), and slanted – to non-neutral with charge marked by numbers from 2 to + 3 .
Figure 6. Formation energies of Fe-HAp as a function of the Fermi energy at extreme points of HAp stability diagram: (a) Ca-rich and P-rich, (b) Ca-poor. The horizontal segments correspond to neutral simulation cells (Table 4), and slanted – to non-neutral with charge marked by numbers from 2 to + 3 .
Nanomaterials 11 02978 g006
Figure 7. Comparison between experimental Fe K-XANES spectra of Fe-HAp published in Ref. [18] (top three curves) and simulated spectra (all other curves) as described in the text.
Figure 7. Comparison between experimental Fe K-XANES spectra of Fe-HAp published in Ref. [18] (top three curves) and simulated spectra (all other curves) as described in the text.
Nanomaterials 11 02978 g007
Table 1. Comparison of calculated heats of formation Δ H f (eV/formula unit), cell volumes Ω 3 ) and bulk modulus B (GPa) of compounds with available reference data [53,54,55,56,57].
Table 1. Comparison of calculated heats of formation Δ H f (eV/formula unit), cell volumes Ω 3 ) and bulk modulus B (GPa) of compounds with available reference data [53,54,55,56,57].
Compound HSE Ref.
Δ H f Ω B Δ H f Ω B
HAp 64.40 532.187 69.44 528.789
CaO 6.00 110.4114 6.58 110.5116
Ca(OH) 2 9.63 56.4631 10.23 54.7833
P 2 O 5 14.22 330.764 15.60 322.012–40
H 3 PO 4 11.52 13.35
H 2 O 2.61 2.97
TCP 39.2 445.5105 42.8 445.5100
TTCP 45.12 796.079 49.4 797.3
DCPD 23.89 490.351 22.33 497.7
DCPA 17.41 311.975 309.3
Table 2. Correspondence between the defect charge state q (in units of e), Kröger-Vink notation and that used in this work (right-most column).
Table 2. Correspondence between the defect charge state q (in units of e), Kröger-Vink notation and that used in this work (right-most column).
Notation
qKröger-VinkThis Work
0 Fe i × Fe i 0
+ 1 Fe i Fe i +
+ 2 Fe i Fe i 2 +
+ 3 Fe i Fe i 3 +
0 Fe Ca × Fe Ca 2 +
+ 1 Fe Ca Fe Ca 3 +
2 Fe P Fe P 3 +
1 Fe P Fe P 4 +
0 Fe P × Fe P 5 +
+ 1 Fe P Fe P 6 +
Table 3. List of proposed structures for Fe defects in Fe-HAp. Also included (when available) are the iron-oxygen interatomic distances R Fe - O , numbers of oxygen neighbors N Fe - O , number of Fe atoms per HAp unit cell n Fe , temperatures of sample preparation (doping treatment or sintering) and detection methods.
Table 3. List of proposed structures for Fe defects in Fe-HAp. Also included (when available) are the iron-oxygen interatomic distances R Fe - O , numbers of oxygen neighbors N Fe - O , number of Fe atoms per HAp unit cell n Fe , temperatures of sample preparation (doping treatment or sintering) and detection methods.
Structure R Fe - O , Å N Fe - O n Fe T prep . , ° CMethodsRef.
Fe Ca ( I ) 1000MB[14]
Fe Ca ( II ) ≥0.15 XRD[58]
Fe Ca ( I ) 2.2–2.360.590MD, EPR, MB[15]
Fe Ca 2.4–2.56≤0.5600–1000EPR, MB[16]
Fe Ca 0.240XRD, XAS[9]
Fe Ca ( II ) /Fe Ca ( I ) 5 or 6/60.5–2.0950XRD, XPS, MB[6]
Fe Ca ( II ) /Fe Ca ( I ) 2.0–2.2/2.2–2.34 / 60.01225XRD, DFT[17]
Fe Ca ( II ) 0.3–6.0100XPS, XRD[59]
FeOOH biogenic/120MB[5,60]
Fe i ( A ) 1.7 20.15–0.75<1000XRD[18]
Fe i ( C ) 1.80–1.8530.15–0.751100XRD[18]
Fe i ( C ) 1.84–1.94 3 ± 0.7 0.151100EXAFS[18]
Fe i ( C ) 1.8–2.440.1–0.960XRD, IR[19]
Table 4. Formation energies ( E f in eV) under extreme growth conditions of HAp, number of nearest oxygen neighbors ( N Fe - O ), Fe-O bond length ( R Fe - O in Å) and magnetic moment ( M z in μ B units) calculated for several defect configurations in the neutral charge state, q = 0 .
Table 4. Formation energies ( E f in eV) under extreme growth conditions of HAp, number of nearest oxygen neighbors ( N Fe - O ), Fe-O bond length ( R Fe - O in Å) and magnetic moment ( M z in μ B units) calculated for several defect configurations in the neutral charge state, q = 0 .
Defect Structure E f N Fe - O R Fe - O M z
Ca- & P-RichCa-PoorP-Poor
Fe Ca ( I ) 2 + 5.25−0.761.0262.04–2.24 0.0
Fe Ca ( II ) 2 + 5.05−0.960.8251.99–2.21 0.0
Fe P 5 + 6.05−0.09−3.0341.69–1.70 1.0
Fe i ( A ) 0 4.484.484.4821.88–1.89 2.0
Fe i ( B ) 0 5.245.245.2431.86–2.22 2.0
Fe i ( C ) 0 6.056.056.0541.94–2.18 2.0
Fe i ( D ) 0 6.506.506.5042.02–2.06 2.0
Table 5. Charge state, q, Number of oxygen neighbors, N Fe - O , and respective range of distances, R Fe - O (Å), and magnetization M z ( μ B units) calculated for charged iron defects in HAp.
Table 5. Charge state, q, Number of oxygen neighbors, N Fe - O , and respective range of distances, R Fe - O (Å), and magnetization M z ( μ B units) calculated for charged iron defects in HAp.
Defect Structureq N Fe - O R Fe - O M z
Fe Ca ( I ) 3 + + 1 61.96–2.12 5.0
Fe Ca ( II ) 3 + + 1 51.90–2.01 1.0
Fe P 6 + + 1 41.65–1.67 0.0
Fe P 4 + 1 41.74–1.75 0.0
Fe P 3 + 2 41.81–1.84 1.0
Fe i ( A ) + + 1 21.84–1.85 3.0
Fe i ( B ) + + 1 41.89–2.02 3.0
Fe i ( B ) 2 + + 2 41.82–1.95 2.0
Fe i ( B ) 3 + + 3 41.78–1.88 3.0
Fe i ( C ) + + 1 41.89–2.08 3.0
Fe i ( C ) 2 + + 2 51.84–2.16 4.0
Fe i ( C ) 3 + + 3 51.72–2.12 5.0
Fe i ( D ) + + 1 41.98–1.99 3.0
Fe i ( D ) 2 + + 2 41.90–1.93 2.0
Fe i ( D ) 3 + + 3 61.95–2.06 1.0
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Avakyan, L.; Paramonova, E.; Bystrov, V.; Coutinho, J.; Gomes, S.; Renaudin, G. Iron in Hydroxyapatite: Interstitial or Substitution Sites? Nanomaterials 2021, 11, 2978. https://0-doi-org.brum.beds.ac.uk/10.3390/nano11112978

AMA Style

Avakyan L, Paramonova E, Bystrov V, Coutinho J, Gomes S, Renaudin G. Iron in Hydroxyapatite: Interstitial or Substitution Sites? Nanomaterials. 2021; 11(11):2978. https://0-doi-org.brum.beds.ac.uk/10.3390/nano11112978

Chicago/Turabian Style

Avakyan, Leon, Ekaterina Paramonova, Vladimir Bystrov, José Coutinho, Sandrine Gomes, and Guillaume Renaudin. 2021. "Iron in Hydroxyapatite: Interstitial or Substitution Sites?" Nanomaterials 11, no. 11: 2978. https://0-doi-org.brum.beds.ac.uk/10.3390/nano11112978

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop