First-Principles Study of Structure, Vibrational, and Elastic Properties

Apr 29, 2014 - The Members of CaPO 4 Family. 2016,17-58. Calcium orthophosphates (CaPO4): occurrence and properties. Sergey V. Dorozhkin. Progress in ...
0 downloads 0 Views 4MB Size
Article pubs.acs.org/crystal

First-Principles Study of Structure, Vibrational, and Elastic Properties of Stoichiometric and Calcium-Deficient Hydroxyapatite Soumya S. Bhat,† Umesh V. Waghmare,‡ and Upadrasta Ramamurty*,†,§ †

Department of Materials Engineering, Indian Institute of Science, Bangalore 560012, India Theoretical Science Unit, Jawaharlal Nehru Centre for Advanced Scientific Research, Jakkur, Bangalore 560064, India § Center of Excellence for Advanced Materials Research, King Abdulaziz University, Jeddah 21589, Saudi Arabia ‡

ABSTRACT: Hydroxyapatite (HAp), a primary constituent of human bone, is usually nonstoichiometric with varying Ca/P molar ratios, with the well-known fact that Ca deficiency can cause marked reductions in its mechanical properties. To gain insights into the mechanism of this degradation, we employ first-principles calculations based on density functional theory and determine the effects of Ca deficiency on structure, vibrational, and elastic properties of HAp. Our simulation results confirm a considerable reduction in the elastic constants of HAp due to Ca deficiency, which was experimentally reported earlier. Stress-induced transformation of the Cadeficient defected structure into a metastable state upon the application of stress could be a reason for this. Local structural stability of HAp and Ca-deficient HAp structures is assessed with full phonon dispersion studies. Further, specific signatures in the computed vibrational spectra for Ca deficiency in HAp can be utilized in experimental characterization of different types of defected HAp. missing Ca2+ is compensated by two protons. Winand et al.24 proposed an alternative formula Ca10−x(HPO4)x(PO4)6‑x(OH)2‑x, which involves the charge compensation of Ca2+ vacancies by OH− vacancies and protons. Thermal analysis and infrared spectroscopic measurements on CDHAps with Ca/P molar ratios ranging between 1.5 and 1.67 by Berry et al.16 confirm the results of Winand et al.24 A similar chemical formula is reported by Kuhl and Nebergall.19 More recently, Rietveld refinements and spectroscopic studies on CDHAp by Wilson et al.26 confirmed the presence of (HPO 4) 2− ions in the structure. Using first-principles calculations, Matsunaga27,28 predicted that the proton-related defects exhibit much smaller formation energies than the other defects. The signatures in vibrational spectra can be effective in characterization of the defects and their structure and a number of experimental and theoretical studies on vibrational properties of SHAp crystal are conducted. Fowler et al.,29 Rehman and Bonfield,30 Cusco et al.31 and Markovic et al.32 have experimentally determined the infrared (IR) and Raman active vibrational modes. Theoretical studies of vibrational spectra of HAp are demanding due to the complicated crystal structure with large number of atoms per cell and are limited. Calderin et al.33 developed a shell model and studied the lattice dynamics of HAp and predicted IR spectra and thermal factors. Phonon frequencies at Γ-point (vanishing wave vector or long wavelength) for SHAp are reported by Corno et al.34 and Slepko et al.35 using ab initio methods and by Pedone et al.36

1. INTRODUCTION Hydroxyapatite (HAp), with the chemical formula Ca10(PO4)6(OH)2, is the basic mineral constituent of human bone and tooth enamel.1 In stoichiometric HAp (SHAp here afterward) the Ca/P molar ratio is 1.67. However, HAp’s constitution in human bones is usually nonstoichiometric and Ca deficient with Ca/P ratios varying from 1.5 to 1.67.2−4 This deficiency is commonly attributed to the low driving force available for growth under biomineralization conditions.5 Since biological and mechanical properties of HAp are associated with Ca deficiency, fundamental understanding of the defect chemistry and mechanical properties of Ca-deficient HAp (CDHAp here afterward) is essential for not only understanding how nature works but also for designing synthetic biomaterials. In view of its importance, many theoretical and experimental efforts were devoted to study the mechanical properties of SHAp.6−13 However, theoretical analysis of the mechanical behavior of CDHAp remains scarce. Viswanath et al.14 investigated the effect of Ca deficiency on mechanical properties of single crystals of HAp using nanoindentation. They found a striking 80% drop in its elastic modulus and hardness, and a marked 75% reduction in toughness, when the Ca/P ratio changes from 1.67 to about 1.5. More recently, Sun et al.15 reported a reduction in the elastic constants and moduli of CDHAp due to vacancies using first-principles calculations. However, their work is limited in the sense that it did not explore the defect structures that are directly relevant to experiments. A number of mechanisms for point defects formation in CDHAp have been reported16−25 and are summarized in Table 3.2 of ref 2, Posner et al.21−23 suggested a general formula Ca10−xH2(PO4)6(OH)2 for CDHAp, where the charge due to a © XXXX American Chemical Society

Received: March 27, 2014 Revised: April 27, 2014

A

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

cell, with ε ranging from −0.02 to 0.02 such that ε is within the elastic limit. The second derivative of Etotal with respect to ε is obtained by fitting second order polynomial, from which the elastic constants are evaluated. The polycrystalline bulk modulus (B) and shear modulus (G) are estimated using the calculated elastic constants of the single crystal, using both the Voigt45 and the Reuss46 approximations, given by the following equations:

using classical interatomic force field model. However, similar studies on CDHAp have not yet been reported in the literature. As discussed above, the structure, vibrational, and mechanical properties of CDHAp are keys for understanding of its biological properties, which we investigate in this paper theoretically using first-principles calculations based on density functional theory (DFT). This approach allows for a systematic study of the precise role played by defects. We consider two types of defect structures, which are experimentally relevant. Through comparison of lattice structure, electronic structures of stoichiometric and nonstoichiometric HAps and defect formation energies, we analyze their consequences to properties. Determining complete phonon spectra of SHAp and CDHAp structures, we assess their local structural stability. Elastic constants of single crystal and polycrystalline forms are estimated, and elastic anisotropy is examined for SHAp and CDHAps.

2. THEORY AND COMPUTATIONAL DETAILS 2.1. Ab Initio Total Energy Calculations. Our calculations are based on first-principles DFT as implemented in the Quantum Espresso code,37 with ultrasoft pseudopotentials38 to represent interaction between ionic cores and valence electrons. The exchange-correlation energy of electrons is treated within a generalized gradient approximated functional (GGA) of the Perdew−Burke−Ernzerhof parametrized form.39 Kohn−Sham wave functions are represented with a plane-wave basis truncated at an energy cutoff of 30 Ry and a charge density with a cutoff of 180 Ry. Brillouin zone integration is sampled on a Monkhorst−Pack mesh40 of 3 × 3 × 6 k-points. Structural optimization is carried out by relaxing positions of all the atoms in the unit cell to minimum energy, until forces are less than 0.001 Ry/Bohr. DFT linear response is used to determine dynamical matrices on a 2 × 2 × 2 mesh of wave vectors, which are Fourier interpolated to obtain full phonon dispersion. 2.2. Vacancy Formation Energy. The formation energy of an isolated vacancy in HAp is calculated using the definition: δH = Etotal(CDHAp‐VX ) + μ X − Etotal(SHAp)

2 1 ∂ Etotal V0 ∂εi ∂εj

1 [2(C11 + C12) + C33 + 4C13] 9

(3)

GV =

1 [7C11 − 5C12 + 12C44 + 2C33 − 4C13] 30

(4)

BR =

(C11 + C12)C33 − 2C132 C11 + C12 + 2C33 − 4C13

(5)

GR =

5 2 {⎡⎣[(C11 + C12)C33 − 2C13 ]C44C66⎤⎦ 2 ⎡⎣3B V C44C66 + [(C11 + C12)C33 − 2C132] (C44 + C66)⎤⎦}

(6)

Since these approximations are true only for isotropic crystals and give theoretical upper and lower bounds on the moduli respectively, the Hill’s averages47 are used to evaluate B and G of the polycrystalline aggregates, which are compared with the experimental values. Hill’s average is given by 1 1 (BR + B V ) and G = (G R + G V ) 2 2

B=

(7)

For an isotropic material, the Young’s modulus, E, and Poisson’s ratio, ν, are given by the following formulas: E=

(1)

9BG 3B − 2G and ν = 3B + G 2(3B + G)

(8)

The elastic anisotropy is estimated by using the parameters proposed by Chung and Buessem,48 which take a value between 0 and 1. They are

where Etotal(SHAp), Etotal(CDHAp-VX), and μX, respectively, represent the total energy of SHAp, CDHAp, and chemical potential of atomic species X which is removed from SHAp to get a vacancy in the structure. Chemical potentials of O and Ca are calculated using the relations: μO = 1/2 μO2 and μCaO = μCa + μO.27 The chemical potential μX equals the total energy, Etotal(X) at 0 K and 0 Pa.41 The defect formation energies of the considered systems are listed in Table 2 and discussed in Section 3.1.2. 2.3. Elastic Properties. There are two methods to compute single crystal elastic constants from first-principles calculations: the energy-strain approach and the stress−strain approach.42,43 In the present work, we use the former. The elastic stiffness constants are given by the relation Cij =

BV =

AB =

(B V − B R ) (G V − G R ) and A G = (B V + B R ) (G V + G R )

(9)

For isotropic materials, AB = AG = 0. Further the values of AB and AG from zero, greater is the anisotropy. We also report the universal elastic anisotropy index AU suggested by Ranganathan and Ostioja-Starzewski49 given by the relation AU =

5G V B + V −6 GR BR

(10)

Using the equation for E in a specific direction and the relation between the stiffness and compliance constants for hexagonal crystals,44,50 E of basal and prism plane in terms of stiffness constants can be derived as

(2)

where Etotal is the total energy, V0 is the volume of unit cell and ε is the strain. Crystals with the hexagonal symmetry have five independent elastic constants C11, C12, C13, C33, C44 = C55, and an additional dependent constant C66 = (C11 − C12)/2.44 To determine the first five, Etotal − ε curves are obtained by applying different types of ε on the equilibrium hexagonal unit

E basal

⎞−1 ⎛ C33 C 1 and Eprism = 2⎜ = + ⎟ C11 + C12 C11−C12 ⎠ ⎝ C (11)

C12)−2C213.

where C = C33(C11 + prism plane are given by B

Similarly, G and ν of basal and

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design G basal = C44 and Gprism = νbasal = =

Article

2C44(C11 − C12) 2C44 + C11 − C12

within the typical errors associated with the calculations made using GGA. Further, the estimated values are in good agreement with earlier theoretical estimates by Menendez et al.9 (a = 9.58 Å, c = 6.88 Å), Mostafa et al.10 (a = 9.41 Å, c = 6.85 Å), Zhu et al.54 (a = 9.52 Å, c = 6.87 Å), Leeuw et al.7 (a = 9.56 Å, c = 6.83 Å), and Sun et al.15 (a = 9.52 Å, c = 6.89 Å). 3.1.2. CDHAPs. Two types of CDHAp structures are considered. (a) One Ca atom from the Ca2 site is removed, as Ca2 atoms are most likely to be deficient than Ca1 atoms.25,26 This will be referred to as CDHAp-VCa here afterward. (b) Early works of Berry et al.16 and Winand et al.24 give a clear evidence for the presence of OH− vacancies and protons attaching to (PO4)3− groups as a mechanism for charge compensation of Ca2+ vacancies, which results in the general compositional formula of Ca10−x(HPO4)x(PO4)6−x(OH)2−x, 0 ≤ x ≤ 1, for CDHAp with the Ca/P molar ratio ranging from 1.5 to 1.67. In the present work, this is modeled by taking away the Ca atom from Ca2 site and OH group along c-axis and adding an H atom, as suggested by Zahn et al.25 The charge neutrality is assumed to be maintained as locally as possible, in order to minimize the energetically unfavorable Coulomb interactions. That is, the charge compensating hydroxide defect is incorporated in close neighborhood of deficient Ca2 site and H is added to an adjacent phosphate ion. This structure is referred to as CDHAp-VCaO here afterward. For geometrical optimization, the atomic positions of CDHAp-VCa and CDHAp-VCaO structures are fully relaxed while maintaining the c/a ratio equal to that of the SHAp. The assumption that the c/a ratio of CDHAp structures will not vary significantly from that of SHAp is justified as the force between the atoms are less than 0.001 Ry/Bohr and pressure on the unit cell is less than −4 kbar in the optimized CDHAp structures. Figure 2a shows the variation of Etotal with a for

(12)

C13 and νprism C11 + C12 (C − C12)(C13 − C33) 1 + 11 2 C + C33(C11 − C12)

(13)

3. RESULTS AND DISCUSSION 3.1. Structural Analysis. 3.1.1. SHAp. As an initial guess, we use the crystallographic structure of SHAp reported by Posner et al.51 and Kay et al.,52 with a hexagonal primitive cell of 44 atoms that has two formula units of Ca5(PO4)3OH in a unit cell with space group P63/m (No. 176 in the International X-ray Tables) and 50% partial occupancy of OH sites along the c axis. The crystal consists of tightly bonded PO4 tetrahedral units, two types of Ca atoms, and OH groups. Ca atoms are distributed among two nonequivalent crystallographic sites Ca1 and Ca2, which are columnar Ca and axial Ca, respectively. Four Ca1 atoms are symmetrically located about the 3-fold axes and are coordinated with three oxygen ions forming triangles. Six Ca2 atoms are symmetrically positioned along the 6-fold screw axes and each is coordinated with six oxygen and one OH group. In the present calculation, we use the modification suggested by Mostafa and Brown10 by removing two of the OH groups and placing the other two along the c axis, without mirror symmetry, which changes the space group to P63. The OH group orientation along the c axis was set as -OH−OHand not as -OH−HO-, as the former is energetically more favorable than later27,53 (see Figure 1).

Figure 1. Hexagonal unit cell of SHAp.

The lattice constants of the SHAp calculated using GGA are listed in Table 1 together with experimental data reported in literature. It shows that our estimate of lattice constant a is larger by ∼1.6% than the experimental value, whereas B is smaller by ∼11% (see Table 3). These differences are well Figure 2. Energy as a function of lattice parameter a for (a) CDHApVCaO (b) CDHAp-VCaO-1 and CDHAp-VCaO-2.

Table 1. Calculated Cell Constants (in Å) and Volume (in Å3) of the Unit Cell for SHAp, with P63 Symmetry, Using DFT Calculations Based on GGA with Experimental Values a (= b) c c/a volume a

our calculation

experimentsa

9.58 6.93 0.724 551.2

9.418 6.875 0.7299 528.1

CDHAp-VCaO. A discontinuity is noted between a = 9.525 and 9.575 Å. To examine the reasons for this in detail, the relaxed atomic positions of the structure corresponding to a = 9.525 Å are taken as inputs and full relaxation of the structure was carried out for a up to 9.75 Å. Likewise, the relaxed structure corresponding to a = 9.575 Å is taken as input and full relaxation of atomic positions was carried out for a down to

Ref 55. C

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

Figure 3. (a) Displacement vectors (green arrows) of CDHAp-VCaO-2 structure with respect to CDHAp-VCaO-1 and (b) optimized structure of CDHAp-VCaO-2.

9.45 Å. The resulting Etotal vs a curves, displayed in Figure 2b, correspond to the optimized structures, which are denoted as CDHAp-VCaO-1 and CDHAp-VCaO-2 respectively. These results show the existence of a metastable structure under applied pressure, whose local structural stability is discussed in section 3.2.2. The structural distortion that takes CDHAp-VCaO-1 to CDHAp-VCaO-2 structure is shown in Figure 3a. Significant distortion is observed around the Ca vacancy. H of HPO4 group and its nearest PO4 group along with nearest Ca atom are displaced to form CDHAp-VCaO-2 structure. The zero temperature transition pressure between CDHAp-VCaO-1 and CDHAp-VCaO-2 structures of CDHAp-VCaO is estimated to be ∼7.2 GPa. We also note that this is relatively high pressure. Couple of points is noteworthy in this context. (a) As the temperature is increased, one can expect the pressure required for transformation to decrease [see for example ref 56]. (b) This is purely hydrostatic pressure whereas the stress state underneath the indenter is such that significant shear components will also be present [see for example ref 57]. This reduces the stress required for phase transformation [see for example ref 58]. As a consequence of these two factors, it is reasonable to expect that the pressure required for inducing phase transformation during indentation experiments conducted at room temperature is much lower than 7.2 GPa, estimated at 0 K and under pure hydrostatic conditions. Our geometrical analysis of the local atomic structures of CDHAps reveals the presence of hydrogen bonding in all the defect structures considered here, consistent with the work of Posner,21 who reported the experimental evidence for the presence of hydrogen bonds in CDHAp samples. It is important to note here that the existence of weak hydrogen bonds could significantly alter the elastic moduli. Table 2 lists the lattice constants and vacancy formation energies for the three CDHAp structures together with that of SHAp. We find that a is larger by 0.02, 0.04, and 0.03 Å in

CDHAp-VCa, CDHAp-VCaO-1, and CDHAp-VCaO-2 structures respectively as compared to that of SHAp. Young et al.59 studied a series of samples of CDHAps and reported that the lattice parameters ranged from a = 9.4157 and c = 6.8777 Å for the material with no (HPO4)−2 ions detectable by IR to a = 9.4389 and c = 6.8865 Å for those with easily detectable (HPO4)−2. Similar results were obtained by Trautz et al.,60 who reported an increase in a by 0.01−0.02 Å for CDHAp from that of SHAp whereas Rietveld structure refinement studies of CDHAp by Wilson et al.61 showed an increase in a by ∼0.04 Å. Our results are indeed in good agreement with these experimental findings and trends. The above results also show that CDHAp-VCaO defect, which involves protonation of the PO4 group, has a much smaller vacancy formation energy, (δH) than CDHAp-VCa. The δH for CDHAp-VCa estimated by Sun et al.15 using GGA is 3.72 eV, which in close agreement with our result. Our estimates of δH are all positive, which is partly due to the higher vacancy concentrations in our model (1.8 × 1021 and 3.58 × 1021 cm−3 for CDHAp-VCa and CDHAp-VCaO respectively). These are nearly equal to the highest equilibrium concentration of defects in HAp which is 5 × 1021 cm−3.28 Also, these formation energies are calculated for 0 K temperature and it has been reported that δH decreases with the increase in temperature and it is positive below 1000 K for HAp.27 Figure 4 shows total density of electronic states (TDOS) along with the partial density of states (PDOS) of SHAp and CDHAp-VCaO-1 structures. The PDOS of CDHAp-VCa and CDHAp-VCaO-2 are not shown as they are similar to that of SHAp and CDHAp-VCaO-1 respectively. The top of the valence band of SHAp shows four significant peaks (1−4), which are primarily associated with states of PO4 group, in agreement with the results of Matsunaga and Kuwabara.27 Peak (1) arises from P 3s, O 2s and 2p states, while peak (2) is due to P 3p and O 2p states with small contribution from O 2s. Peaks (3) and (4) are dominated by O 2p states with only small Ca 4s and 3d contribution. The bottom of the conduction band is dominated by 3d states of Ca. In comparison of CDHAp-VCaO-1 structure with SHAp, the PDOS of Ca and OH remain unchanged whereas PDOS of PO4 changed significantly. Along with four peaks (1−4), three extra peaks (A, B and C) are observed in the valence band of CDHAp-VCaO-1 (see Figure 4b). These three peaks arise mainly from O 2p, with a small contribution of H 1s to peaks A and B. From the total electronic density of states (DOS) (see Figure 5), the energy band gap of SHAp was estimated as 5.32 eV, which is in good agreement with the literature data of 5.3 and 5.23 eV estimated using GGA27,35 and 5.40 eV using LDA.33 In

Table 2. Lattice Constants (in Å) and Volume (in Å3) of the Unit Cell for HAp with and without Vacancy and Vacancy Formation Energy, δH (in eV) Calculated Using eq 1, for CDHAp-VCa, CDHAp-VCaO-1, and CDHAp-VCaO-2

a (= b) c volume δH

SHAp

CDHAp-VCa

CDHAp-VCaO-1

CDHAp-VCaO-2

9.58 6.93 551.2

9.60 6.95 554.3 3.48

9.62 6.96 558.8 2.66

9.61 6.95 556.2 2.87 D

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

Figure 4. Combined partial density of states (PDOS) for (a) SHAp and (b) CDHAp-VCaO-1. The Fermi level is set at zero energy and marked by a vertical line.

132 degrees of freedom, of which 3 are translational acoustic modes and remaining 129 are optical modes. The classification of phonon states based on the symmetry group representations is nontrivial due to defects and a physical picture is not simple. However, the availability of the atomic displacements associated with the eigenvectors provide information on the vibration and characteristic frequencies (related to bond stiffness) of group of atoms or molecules in the system. In this study, the analysis of the computed frequencies has been carried out through graphical visualization of atomic displacements associated with the eigenvector of each mode. For SHAp, the bands below 413 cm−1, which have been classified as lattice modes, are mixed vibrational modes involving Ca atoms and PO4 groups. The modes ranging from 420 to 620 cm−1 and 900 to 1080 cm−1 of the spectrum involve mainly the PO4 vibrations. The modes at 684 and 3662 cm−1 correspond to libration and stretching modes of OH, respectively. We find a reasonable agreement with the experimental results:29 the frequencies of the phosphate group vibrations (1028 and 1078 cm−1) are underestimated by our calculations by less than 2% (experimental values 1040 and 1092 cm−1 respectively). On the other hand, theoretical frequencies of the vibrations of OH groups at 684 and 3662 cm−1 are overestimated by 8% and 2% respectively, compared with experimental results (630 and 3572 cm−1 respectively). Overall, the agreement between our phonon spectra and experimental results is quite reasonable. A similar trend is reported by Alexander et al.,35 and our results are in good agreement with earlier reported theoretical results.33−36 A broad discussion of the band assignments for HAp has already been reported elsewhere.29,34 Therefore, we limit the discussion below mainly to the differences in the vibrational spectra of SHAp and CDHAps. 3.2.2. CDHAps. First, we establish that the CDHAp structures that are considered in this work are stable, as there are no imaginary frequencies (see Figure 7). There are some

Figure 5. Calculated electronic density of states (DOS) for SHAp and CDHAp structures. The Fermi level is set at zero energy and marked by a vertical line.

contrast, experimentally estimated values of band gaps range from 3.95 eV62 to more than 6 eV.63 This wide variance in experimental values does indeed indicate that defects must be common in HAp, as they can alter the band gaps considerably. In the case of CDHAp-VCa, extra electronic levels appear just above the valence band maximum, which can be considered as acceptor-like levels induced by vacancies. Indeed, our estimate of the band gap for CDHAp-VCa (5.67 eV) is higher, whereas for CDHAp-VCaO-1 (5.04 eV) and CDHAp-VCaO-2 (4.81 eV), it is lower than that of SHAp (5.32 eV). 3.2. Lattice Dynamics. 3.2.1. SHAp. The phonon density of states and full phonon dispersion curves for SHAp and CDHAp structures are illustrated in Figures 6 and 7 respectively. As SHAp has 44 atoms in the unit cell, there are E

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

Figure 6. Density of phonon states for SHAp and CDHAp structures for (a) frequency range of 0−1200 cm−1 and (b) frequency range of 1200− 3800 cm−1.

Figure 7. Full phonon dispersion for (a) SHAp, (b) CDHAp-VCa, (c) CDHAp-VCaO-1, and (d) CDHAp-VCaO-2. Only low frequency modes are shown for clarity.

3662 cm−1 of SHAp. Second, the modes due to OH are observed at 650, 814, and 842 cm−1 in CDHAp-VCa, which correspond to those of 613 and 684 cm−1 for SHAp. It is evident from figure. 6 that the mode at 1078 cm−1 due to PO4

significant differences in the phonon spectrum of SHAp and CDHAps. Phonon density of states (Figure 6) clearly shows that there are modes at 3307 and 3683 cm−1 corresponding to OH group of CDHAp-VCa, in contrast to the narrow band at F

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

Figure 8. Density of states of IR active modes for SHAp and CDHAp structures for (a) frequency range of 0−1200 cm−1 and (b) frequency range of 1200−3800 cm−1. Distinct peaks of HPO4 are marked as *.

calculated around 1160 cm−1 corresponding to δOH mode of HPO4, is experimentally observed at frequency 1210 cm−1 for Ca deficient samples. Underestimation of these frequencies here is partly due to higher concentration of defects in our analysis. Second, temperature dependence of these modes has not been accounted for. Anharmonicity results in temperature dependence of frequencies of vibrational modes. Such effects of anharmonicity are stronger for lower frequency (soft) modes. Since modes involving H atoms are relatively hard (higher frequency, >800 cm−1), we expect these effects on their frequencies to be weak, which would be reflected in slight decrease in frequencies with increasing temperature. While calculated frequencies are slightly underestimated, the changes in frequencies due to Ca-deficiency are consistent with the experimental results, and should be quite useful in characterization of Ca defects in HAp. Finally, our estimated IR spectra shows extra peaks with small intensities around 1350 and 1920 cm−1 (refer Figure 8b), corresponding to δOH mode of HPO4 for both CDHAp-VCaO-1 and CDHAp-VCaO-2 structures, which have not been reported yet in experiments. 3.3. Elastic Properties. 3.3.1. SHAp. The five independent elastic constants are estimated from the energy-strain curves as described in section 2.3. Table 3 lists all the five independent Cij’s and B together with earlier experimental and theoretical estimates. All the Cij’s except C13 are slightly smaller than that reported by Katz et al.,8 but in agreement with the ab initio results of Menendez et al.9 Recently reported results of Tofail et al.,12 where the elastic constants are determined for a textured polycrystalline sample of SHAp that has transverse isotropy, is also listed for comparison. Note that the disparity between these experimental results is quite large, which can be attributed to the different experimental techniques used to extract the elastic data, and also likely due to different defect structure. The present B of SHAp agrees well with earlier theoretical and some of the experimental results. 3.3.2. CDHAps. Accurate estimation of all the elastic constants of the CDHAp crystals is important for under-

vibration in SHAp is missing in the case of CDHAp-VCa. These constitute the vibrational signatures, and are the basis for spectroscopic characterization of Ca-vacant HAp. The major difference in the phonon spectra of SHAp and CDHAp-VCaO-1 (also CDHAp-VCaO-2) is in their HPO4 vibrations. The modes corresponding to HPO4 vibrations are observed at 801, 879, 1084, 1162, 1351, and 1923 cm−1 in CDHAp-VCaO-1. The OH libration mode at 684 cm−1 of SHAp is significantly altered in the case of CDHAp-VCaO-1. The band near 3662 cm−1 corresponding to OH stretching mode in SHAp hardens to 3696 cm−1 in CDHAp-VCaO-1. Also, extra modes at 700 and 888 cm−1, associated with OH and PO4 vibrations respectively, are found in CDHAp-VCaO-1. Modes near 590, 900, and 1070 cm−1 are mixed vibrational modes due to HPO4 and PO4 groups of CDHAp-VCaO-1. The phonon dispersion of CDHAp-VCaO-2 is only slightly different from that of CDHAp-VCaO-1. The OH stretching mode at frequency at 3696 cm−1 for CDHAp-VCaO-1 is at 3676 cm−1 for CDHAp-VCaO-2. Modes corresponding to HPO4 vibrations of CDHAp-VCaO-2 are shifted to 830, 862, 1128, and 1327 cm−1 (which are at 801, 879, 1162, and 1351 cm−1 respectively for CDHAp-VCaO-1). These results should form the vibrational signatures and basis for spectroscopic analysis of CDHAp structures, in which a Ca vacancy is accompanied by OH vacancy and protonation of PO4 group. To compare vibrational features of SHAp and CDHAp structures with the experimental observations, the density of states of IR active modes for all the four structures are calculated and illustrated in Figure 8. The IR spectra of CDHAp-VCaO-1 and CDHAp-VCaO-2 are consistent with the experimentally reported spectra for CDHAp.14,17,26,61 The P− OH stretching mode (that is ν5 mode) of HPO4 is observed at 800 and 830 cm−1 respectively for CDHAp-VCaO-1 and CDHAp-VCaO-2, which is experimentally reported at 870 cm−1. The ν3 mode of HPO4 which is seen at 1130 cm−1 in experiments is estimated at 1070 and 1075 cm−1 respectively for CDHAp-VCaO-1 and CDHAp-VCaO-2. The IR active mode G

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

3.3.3. Anisotropy. The anisotropy indices of B and G, AB, and AG respectively, are estimated as explained in section 2.3. CDHAp-VCaO-1 exhibits relatively small values of AB and AG compared to other structures. From the universal anisotropy elastic index AU, we note that SHAp shows a larger value of AU than of HAp with vacancies. Thus, the vacancies reduce the anisotropy in elastic moduli, as expected intuitively.

Table 3. Elastic Constants Cij’s and Bulk Modulus B (in GPa) of SHAp with Experimental and Theoretical Results of Others present work C11 C12 C13 C33 C44 BH a

experiment 1a experiment 2b

118.3 31.6 63.7 156.8 33.5 76.4

137 42.5 54.9 172 39.6 82.6

other theoretical calculationc

137.2 53.0 55.1 123.2 42.2

117.1−157.5 26.2−48.9 55.6−68.5 147.3−231.8 38.5−56.4 76−90.7

Table 5. Calculated Anisotropy Ratios AB, AG, AU, and E, B (in GPa), ν of Basal and Prism Planes for SHAp, CDHApVCa, CDHAp-VCaO-1, and CDHAp-VCaO-2

Reference 8. bReference 12. cReference 6, 9, 11, and 15. AB AG AU Ebasal Eprism Gbasal Gprism νbasal νprism

standing their anisotropic elastic behavior. Computed Cij’s for CDHAp-VCa, CDHAp-VCaO-1 and CDHAp-VCaO-2 crystals are listed in Table 4. Since experimental or theoretical data are not Table 4. Calculated Single Crystal Elastic Constants Cij’s, Bulk Moduli B, and Shear Moduli G (in GPa) in Voigt− Reuss−Hill Approaches, Young’s ModuliE (in GPa), and Poisson’s Ratios ν in Hill Approach for SHAp, CDHAp-VCa, CDHAp-VCaO-1, and CDHAp-VCaO-2 elastic constants

SHAp

C11 C12 C13 C33 C44 C66 BV BR BH GV GR GH BH/GH E ν

118.3 31.6 63.7 156.8 33.5 43.4 79.0 73.8 76.4 37.7 36.6 37.1 2.1 95.9 0.29

CDHAp-VCa CDHAp-VCaO-1 87.3 30.6 50.3 128.7 28.5 28.3 62.9 58.1 60.5 28.5 28.1 28.3 2.1 73.5 0.30

104.5 36.0 39.7 110.9 29.6 34.3 61.2 61.1 61.2 32.3 32.2 32.3 1.9 82.3 0.28

SHAp

CDHAp-VCa

CDHAp-VCaO‑1

CDHAp-VCaO-2

0.035 0.015 0.223 102.7 92.1 33.5 37.8 0.42 0.22

0.040 0.008 0.163 85.8 65.8 28.5 28.4 0.43 0.24

0.001 0.003 0.028 88.5 85.1 29.6 31.8 0.28 0.26

0.010 0.014 0.163 95.1 85.2 26.2 30.6 0.32 0.22

Estimated E of basal plane for all the structures considered here are larger than that of prism plane. Higher Ebasal than Eprism by ∼10% for SHAP is consistent with the nanoindentation results reported earlier.13,65 Due to close-packed nature of basal planes, these planes are stiffer than other planes. Anisotropy in E is less pronounced in case of CDHAp-VCaO-1 (∼4%) and more for CDHAp-VCa (∼23%). Calculated values of ν for basal and prism also exhibit similar trends. On the other hand, estimated Gbasal is lower than Gprism for all the structures except for CDHAp-VCa, for which they are almost same.

CDHAp-VCaO-2 101.9 28.5 42.0 122.1 26.2 36.7 61.2 60.0 60.6 32.0 31.2 31.6 1.9 80.8 0.28

4. DISCUSSION The above results illustrate that the stiffness of HAp decreases considerably due to Ca vacancies. This is consistent with the results of first-principles calculations, showing the reduction in the modulus due to point defects in ionic/covalent solids.66−69 Recently, first-principle calculations showed a drop over 60% in E and G of HAp with Ca ionic vacancy.15 Nanoindentation study on the effect of Ca deficiency showed an 80% reduction in elastic modulus and hardness in the plate-shaped CDHAp crystals.14 Our calculations show that the CDHAp-VCaO structure, which is more relevant to such experiments, undergoes a structural phase transformation from CDHApVCaO-1 to CDHAp-VCaO-2 phase under applied pressure. This kind of structural phase transformation under pressure could lead to drastic reduction in the mechanical properties and could be a possible reason for the experimentally reported drastic reduction in the moduli and hardness of CDHAp. In this context it is worth noting that Varughese et al.70 using nanoindentation, observed stress induced transformation of aspirin from metastable polymorph II to the stable polymorph I and reported that the metastable form is considerably softer than the stable form. The CDHAp-VCaO structure, which has a Ca/P ratio 1.5, includes one Ca (out of 10) and one OH (out of two) vacancies per unit cell, and one (out of six) PO4 group substituted by HPO4 group. The presence of Ca vacancies could lead to replacement of strong ionic Ca−O bond, connecting the covalently bonded PO4 tetrahedra and OH groups with weak hydrogen bonding. This could also contribute to the decrease in the strength of CDHAp crystals.

available for elastic constants of precisely the CDHAp structures considered in the present work, we cannot make a comparison. From Table 4, it is clear that the difference between the elastic constants of SHAp and those of CDHAps is significant. In general, all the Cij’s of CDHAps are smaller than that of SHAp except C12 of CDHAp-VCaO-1. The percentage decrease in Cij’s varies from 3% to 35%, and the most significant decrease is seen for C11, C33 and C66. Elastic properties such as B, G, E, and ν for the SHAp and CDHAPs are determined using VRH scheme as described in section 2.3 (see Table 4). All the elastic moduli for CDHAp are considerably smaller than that of defect free HAp. The bulk modulus is ∼20% smaller than that of SHAp. A reduction of ∼23% and ∼15% is found for the E, whereas G reduces by ∼24% and ∼15% due to defect structure of CDHAp-VCa and CDHAp-VCaO respectively. In order to estimate the ductile and brittle behavior, the B/G ratio is also listed in Table 4. The high value of B/G ratio corresponds to ductility and low value to brittleness, with the critical value being 1.75.64 Note that the difference in the B/G ratio between SHAp and CDHAps is marginal, which indicates that the effect of Ca deficiency on the brittle nature of HAp is negligible. H

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

Notes

The reduction in elastic modulus presented in this work is less than that reported in experiment,14 which could be accounted partly by a large scattering in the values of modulus obtained for different samples with varying Ca/P ratios and strong anisotropy of the crystals. Second, the difference in our results and experiments can also arise from the temperature effects and metastability. Due to soft and hard vibrations in the defective HAp structure, vibrational free energy nonlinearity could account for some of this. Nanoindentation measurements on human cortical bone have disclosed anisotropy in its elastic behavior.71−76 Also Viswanath et al.13 and Saber-Samandari and Gross65 reported the elastic anisotropy in single crystals of SHAp. However, no experimental data are available to check the anisotropy behavior of CDHAp. Such data are crucial for explaining the orientational dependence of the CDHAp crystals in collagen fibrils which is important to the biomedical industry. In this regard, results presented in this paper will be useful to develop a theoretical understanding of the anisotropic behavior of bone. One of the major objectives of this study is to establish a realistic calculation to understand the properties of bone, which will be instrumental in designing the artificial bone material, since these materials have their properties very different from that of pure HAp crystals. The full phonon spectra and IR spectra of CDHAps, which is presented for the first time, will be useful in the characterization of Ca deficiency in HAp. The structural details and also electronic structure studies will be effective in the analysis of CDHAps. Using the results presented in this paper for SHAp and CDHAps, a deeper insight into the Ca deficiency of bone can be anticipated.

The authors declare no competing financial interest.



5. SUMMARY We have determined the effects of calcium deficiency on structure, vibrational, and elastic properties of hydroxyapatite using first-principles calculations based on DFT. Estimated elastic constants of SHAp agree well with experimental and theoretical results reported earlier. Calculated elastic constants reveal 30% reduction in the elastic constants and moduli of HAp due Ca deficiency, and more such reduction is expected at higher strains and temperatures. Our calculations reveal the presence of a metastable state, for one of the defect structures considered, under applied pressure. We propose this as a cause for experimentally observed drastic reduction in the modulus and hardness for CDHAp, seen in indentation experiment. The presence of hydrogen bonds in CDHAp structures could also account for some of these. Elastic anisotropic behavior of SHAp and CDHAp are analyzed, which will be useful in understanding the effect of crystal orientation in designing synthetic bone. The defect formation energies are estimated for various defect structures under consideration and their electronic structures are compared with that of SHAp. Local structural stability of the defect structures assessed with full phonon spectra confirms the stability of all these structures including the reported metastable state. The comparative analysis of phonon frequencies of SHAp and of CDHAp structures will form the basis of experimental spectroscopic characterization of CDHAp.



REFERENCES

(1) Young, R. A.; Brown, W. E. Structures of biological minerals. In Biological Mineralization and Demineralization; Springer: New York, 1982; pp 101−141. (2) Elliott, J. C. Structure and Chemistry of the Apatites and Other Calcium Orthophosphates; Elsevier: Amsterdam, 1994; Vol. 4. (3) Eppell, S. J.; Tong, W.; Katz, J. L.; Kuhn, L.; Glimcher, M. J. Shape and size of isolated bone mineralites measured using atomic force microscopy. J. Orthop. Res. 2001, 19 (6), 1027−1034. (4) Weiner, S.; Wagner, H. D. The material bone: structuremechanical function relations. Annu. Rev. Mater. Sci. 1998, 28 (1), 271−298. (5) Viswanath, B.; Ravishankar, N. Controlled synthesis of plateshaped hydroxyapatite and implications for the morphology of the apatite phase in bone. Biomaterials 2008, 29 (36), 4855−4863. (6) Ching, W. Y.; Rulis, P.; Misra, A. Ab initio elastic properties and tensile strength of crystalline hydroxyapatite. Acta Biomaterialia 2009, 5 (8), 3067−3075. (7) de Leeuw, N. H.; Bowe, J. R.; Rabone, J. A. L. A computational investigation of stoichiometric and calcium-deficient oxy-and hydroxyapatites. Faraday Discuss. 2007, 134, 195−214. (8) Katz, J. L.; Ukraincik, K. On the anisotropic elastic properties of hydroxyapatite. J. Biomech. 1971, 4 (3), 221−227. (9) Menéndez-Proupin, E.; Cervantes-Rodríguez, S.; Osorio-Pulgar, R.; Franco-Cisterna, M.; Camacho-Montes, H.; Fuentes, M. E. Computer simulation of elastic constants of hydroxyapatite and fluorapatite. J. Mech. Behav. Biomed. Mater, 2011, 4 (7), 1011−1020. (10) Mostafa, N. Y.; Brown, P. W. Computer simulation of stoichiometric hydroxyapatite: structure and substitutions. J. Phys. Chem. Solids 2007, 68 (3), 431−437. (11) Snyders, R.; Music, D.; Sigumonrong, D.; Schelnberger, B.; Jensen, J.; Schneider, J. M. Experimental and ab initio study of the mechanical properties of hydroxyapatite. Appl. Phys. Lett. 2007, 90 (19), 193902−193902. (12) Tofail, S. A. M.; Haverty, D.; Cox, F.; Erhart, J.; Hána, P.; Ryzhenko, V. Direct and ultrasonic measurements of macroscopic piezoelectricity in sintered hydroxyapatite. J. Appl. Phys. 2009, 105 (6), 064103−064103. (13) Viswanath, B.; Raghavan, R.; Ramamurty, U.; Ravishankar, N. Mechanical properties and anisotropy in hydroxyapatite single crystals. Scr. Mater. 2007, 57 (4), 361−364. (14) Viswanath, B.; Shastry, V. V.; Ramamurty, U.; Ravishankar, N. Effect of calcium deficiency on the mechanical properties of hydroxyapatite crystals. Acta Mater. 2010, 58 (14), 4841−4848. (15) Sun, J. P.; Song, Y.; Wen, G. W.; Wang, Y.; Yang, R. Softening of hydroxyapatite by vacancies: A first principles investigation. Mater. Sci. Eng.: C 2013, 33 (3), 1109−1115. (16) Berry, E. E. The structure and composition of some calciumdeficient apatites. J. Inorg. Nuclear Chem. 1967, 29 (2), 317−327. (17) Berry, E. E. The structure and composition of some calciumdeficient apatitesII. J. Inorg. Nuclear Chem. 1967, 29 (7), 1585− 1590. (18) Calderin, L.; Stott, M. J.; Rubio, A. Electronic and crystallographic structure of apatites. Phys. Rev. B 2003, 67 (13), 134106. (19) Kühl, G.; Nebergall, W. H. Hydrogenphosphat-und Carbonatapatite. Z. Anorg. Allg. Chem. 1963, 324 (5−6), 313−320. (20) Miquel, J. L.; Facchini, L.; Legrand, A. P.; Rey, C.; Lemaitre, J. Solid state NMR to study calcium phosphate ceramics. Colloids Surfaces 1990, 45, 427−433. (21) Posner, A. S. Hydrogen-bonding in calcium-deficient hydroxyapatites. Nature 1960, 188, 486−487. (22) Posner, A. S.; Perloff, A. Apatites deficient in divalent cations. J. Res. Natl. Bur. Stand. 1957, 58, 279. (23) Stutman, J. M.; Posner, A. S.; Lippincott, E. R. Hydrogen bonding in the calcium phosphates. Nature 1960, 188, 486−487.

AUTHOR INFORMATION

Corresponding Author

*Tel: +91-080-22933241. Fax: +91-080-23600472. E-mail: [email protected]. I

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

(24) Winand, L.; Dallemagne, M. J. Hydrogen bonding in the calcium phosphates. Nature 1962, 193 (4813), 369−370. (25) Zahn, D.; Hochrein, O. On the composition and atomic arrangement of calcium-deficient hydroxyapatite: an ab-initio analysis. J. Solid State Chem. 2008, 181 (8), 1712−1716. (26) Wilson, R. M.; Elliott, J. C.; Dowker, S. E. P.; RodriguezLorenzo, L. M. Rietveld refinements and spectroscopic studies of the structure of Ca-deficient apatite. Biomaterials 2005, 26 (11), 1317− 1327. (27) Matsunaga, K.; Kuwabara, A. First-principles study of vacancy formation in hydroxyapatite. Phys. Rev. B 2007, 75 (1), 014102. (28) Matsunaga, K. Theoretical investigation of the defect formation mechanism relevant to nonstoichiometry in hydroxyapatite. Phys. Rev. B 2008, 77 (10), 104106. (29) Fowler, B. O. Infrared studies of apatites. I. Vibrational assignments for calcium, strontium, and barium hydroxyapatites utilizing isotopic substitution. Inorg. Chem. 1974, 13 (1), 194−207. (30) Rehman, I.; Bonfield, W. Characterization of hydroxyapatite and carbonated apatite by photo acoustic FTIR spectroscopy. J. Mater. Sci.: Mater. Med. 1997, 8 (1), 1−4. (31) Cusco, R.; Guitian, F.; Aza, S. d.; Artús, L. Differentiation between hydroxyapatite and β-tricalcium phosphate by means of μRaman spectroscopy. J. Eur. Ceram. Soc. 1998, 18 (9), 1301−1305. (32) Markovic, M.; Fowler, B. O.; Tung, M. S. Preparation and comprehensive characterization of a calcium hydroxyapatite reference material. J. Res. Natl. Inst. Stand. Technol. 2004, 109 (6), 553−568. (33) Calderin, L.; Dunfield, D.; Stott, M. J. Shell-model study of the lattice dynamics of hydroxyapatite. Phys. Rev. B 2005, 72 (22), 224304. (34) Corno, M.; Busco, C.; Civalleri, B.; Ugliengo, P. Periodic ab initio study of structural and vibrational features of hexagonal hydroxyapatite Ca10(PO4)6(OH)2. Phys. Chem. Chem. Phys. 2006, 8 (21), 2464−2472. (35) Slepko, A.; Demkov, A. A. First-principles study of the biomineral hydroxyapatite. Phys. Rev. B 2011, 84 (13), 134108. (36) Pedone, A.; Corno, M.; Civalleri, B.; Malavasi, G.; Menziani, M. C.; Segre, U.; Ugliengo, P. An ab initio parameterized interatomic force field for hydroxyapatite. J. Mater. Chem. 2007, 17 (20), 2061− 2068. (37) Giannozzi, P.; Baroni, S.; Bonini, N.; Calandra, M.; Car, R.; Cavazzoni, C.; Ceresoli, D.; Chiarotti, G. L.; Cococcioni, M.; Dabo, I. QUANTUM ESPRESSO: a modular and open-source software project for quantum simulations of materials. J. Phys.: Condens. Matter 2009, 21 (39), 395502. (38) Vanderbilt, D. Soft self-consistent pseudopotentials in a generalized eigenvalue formalism. Phys. Rev. B 1990, 41 (11), 7892. (39) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized gradient approximation made simple. Phys. Rev. Lett. 1996, 77 (18), 3865. (40) Monkhorst, H. J.; Pack, J. D. Special points for Brillouin-zone integrations. Phys. Rev. B 1976, 13 (12), 5188−5192. (41) Astala, R.; Stott, M. J. First-principles study of hydroxyapatite surfaces and water adsorption. Phys. Rev. B 2008, 78 (7), 075427. (42) Le Page, Y.; Saxe, P. Symmetry-general least-squares extraction of elastic coefficients from ab initio total energy calculations. Phys. Rev. B 2001, 63 (17), 174103. (43) Le Page, Y.; Saxe, P. Symmetry-general least-squares extraction of elastic data for strained materials from ab initio calculations of stress. Phys. Rev. B 2002, 65 (10), 104104. (44) Nye, F. Physical Properties of Crystals; Clarendon Press: Oxford, 1964. (45) Voigt, W. Lehrb. Kristallphys.; Teubner, B. G., Ed.; The University of Michigan: Ann Arbor, MI, 1928; Vol. 34. (46) Reuss, A. Calculation of the flow limits of mixed crystals on the basis of the plasticity of monocrystals. Z. Angew. Math. Mech. 1929, 9, 49−58. (47) Hill, R. The elastic behaviour of a crystalline aggregate. Proc. Phys. Soc., Sect. A 1952, 65 (5), 349. (48) Chung, D. H.; Buessem, W. R.; Vahldiek, F. W.; Mersol, S. A. Anisotropy in Single Crystal Refractory Compounds; Plenum Press, New York, 1968; p 217.

(49) Ranganathan, S. I.; Ostoja-Starzewski, M. Universal elastic anisotropy index. Phys. Rev. Lett. 2008, 101 (5), 055504. (50) Tromans, D. Elastic anisotropy of hcp metal crystals and polycrystals. Int. J. Res. Rev. Appl. Sci. 2011, 6, (4). (51) Posner, A. S.; Perloff, A.; Diorio, A. F. Refinement of the hydroxyapatite structure. Acta Crystallogr. 1958, 11 (4), 308−309. (52) Mi, K. A. Y.; Ra, Y. Crystal structure of hydroxyapatite. Nature 1964, 204, 1050−1052. (53) 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 (15), 155104. (54) Zhu, W.; Wu, P. Surface energetics of hydroxyapatite: a DFT study. Chem. Phys. Lett. 2004, 396 (1), 38−42. (55) Hughes, J. M.; Cameron, M.; Crowley, K. D. Structural variations in natural F, OH, and Cl apatites. Am. Mineral. 1989, 74 (7− 8), 870−876. (56) Xia, H.; Ruoff, A. L.; Vohra, Y. K. Temperature dependence of the ω-bcc phase transition in zirconium metal. Phys. Rev. B 1991, 44 (18), 10374. (57) Patnaik, M.; Narasimhan, R.; Ramamurty, U. Spherical indentation response of metallic glasses. Acta Mater. 2004, 52 (11), 3335−3345. (58) Bhat, S. S.; Waghmare, U. V.; Ramamurty, U. Pressure induced structural phase transformation in TiN: A first-principles study. J. Appl. Phys. 2013, 113 (13), 133507. (59) Young, R. A.; Holcomb, D. W. Variability of hydroxyapatite preparations. Calcif. Tissue Int. 1981, 34, S17−32. (60) Trautz, O. R. X-ray diffraction of biological and synthetic apatites. Ann. N. Y. Acad. Sci. 1955, 60 (5), 696−712. (61) Wilson, R. M.; Elliott, J. C.; Dowker, S. E. P. Formate incorporation in the structure of Ca-deficient apatite: Rietveld structure refinement. J. Solid State Chem. 2003, 174 (1), 132−140. (62) Rosenman, G.; Aronov, D.; Oster, L.; Haddad, J.; Mezinskis, G.; Pavlovska, I.; Chaikina, M.; Karlov, A. Photoluminescence and surface photovoltage spectroscopy studies of hydroxyapatite nano-bioceramics. J. Lumin. 2007, 122, 936−938. (63) Tsukada, M.; Wakamura, M.; Yoshida, N.; Watanabe, T. Band gap and photocatalytic properties of Ti-substituted hydroxyapatite: Comparison with anatase-TiO2. J. Mol. Catal. A: Chem. 2011, 338 (1), 18−23. (64) Pugh, S. F. Philos. Mag. 1954, 45, 823. (65) Saber-Samandari, S.; Gross, K. A. Micromechanical properties of single crystal hydroxyapatite by nanoindentation. Acta Biomater. 2009, 5 (6), 2206−2212. (66) Hu, Q.-M.; Kádas, K.; Hogmark, S.; Yang, R.; Johansson, B.; Vitos, L. Hardness and elastic properties of covalent/ionic solid solutions from first-principles theory. J. Appl. Phys. 2008, 103 (8), 083505. (67) Jhi, S.-H.; Ihm, J.; Louie, S. G.; Cohen, M. L. Electronic mechanism of hardness enhancement in transition-metal carbonitrides. Nature 1999, 399 (6732), 132−134. (68) Jhi, S.-H.; Louie, S. G.; Cohen, M. L.; Ihm, J. Vacancy hardening and softening in transition metal carbides and nitrides. Phys. Rev. Lett. 2001, 86 (15), 3348. (69) Lucas, M.; Wang, Z. L.; Riedo, E. Combined polarized Raman and atomic force microscopy: In situ study of point defects and mechanical properties in individual ZnO nanobelts. Appl. Phys. Lett. 2009, 95 (5), 051904. (70) Varughese, S.; Kiran, M.; Solanko, K. A.; Bond, A. D.; Ramamurty, U.; Desiraju, G. R. Interaction anisotropy and shear instability of aspirin polymorphs established by nanoindentation. Chem. Sci. 2011, 2 (11), 2236−2242. (71) Carnelli, D.; Lucchini, R.; Ponzoni, M.; Contro, R.; Vena, P. Nanoindentation testing and finite element simulations of cortical bone allowing for anisotropic elastic and inelastic mechanical response. J. Biomech. 2011, 44 (10), 1852−1858. (72) Fan, Z.; Swadener, J.; Rho, J.; Roy, M.; Pharr, G. Anisotropic properties of human tibial cortical bone as measured by nanoindentation. J. Orthop. Res. 2002, 20 (4), 806−810. J

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX

Crystal Growth & Design

Article

(73) Franzoso, G.; Zysset, P. K. Elastic anisotropy of human cortical bone secondary osteons measured by nanoindentation. J. Biomech. Eng. 2009, 131 (2), 021001. (74) Rho, J. Y.; Currey, J. D.; Zioupos, P.; Pharr, G. M. The anisotropic Young’s modulus of equine secondary osteones and interstitial bone determined by nanoindentation. J. Exp. Biol. 2001, 204 (10), 1775−1781. (75) Rho, J.-Y.; Tsui, T. Y.; Pharr, G. M. Elastic properties of human cortical and trabecular lamellar bone measured by nanoindentation. Biomaterials 1997, 18 (20), 1325−1330. (76) Swadener, J.; Rho, J. Y.; Pharr, G. Effects of anisotropy on elastic moduli measured by nanoindentation in human tibial cortical bone. J. Biomed. Mater. Res. 2001, 57 (1), 108−112.

K

dx.doi.org/10.1021/cg5004269 | Cryst. Growth Des. XXXX, XXX, XXX−XXX