On Estimation of Adhesive Strength of Implants Bioactive Coating with Titanium by Density Functional Theory and Molecular Dynamics Simulations

Abstract

One of the ways to improve and accelerate osseointegration of a surgical implant with bone is application of biocompatible coatings, in particular, hydroxyapatite (HAp). Since the cases of delamination of the coating take place in dental practice, it is very important to estimate the adhesive strength of HAp with the implant. A measure of the coating-to-substrate bond strength is the energy of this bond. In this research, quantum chemistry is used to calculate the binding energy of functional groups (anions) of hydroxyapatite and titanium 2+, which is a standard implant material. First, using Density Functional Theory with Becke three-parameter Lee-Yang-Parr hybrid exchange-correlation functional, the lowest potential energy surface is calculated. Then, by ab initio molecular dynamics, the reaction path, the reaction products, frequencies of oscillations, the activation energies and binding energies between various combinations of component anions of HAp and Ti(II) are calculated.

Keywords:
hydroxyapatite coating; bond strength; MD simulation

1. Introduction

The key requirement of dental implantation is strong osseointegration (fusion of an implant with a bone - Figure 1)^1,2. To achieve this, implants are often covered with bioactive coatings^3,4,5. One of the most popular materials is hydroxyapatite (HAp), with the molecular formula Ca₁₀ (PO₄)₆ (OH)₂ (Figure 2)^6,7,8. The composition of HAp is similar to the mineral constituent of bone. The nanocrystalline structure of HAp gives a microrelief of the surface, favorable for osseointegration^9-13.

Figure 1
General view of a dental implant screwed into the jaw.

Figure 2
The unit cell of hydroxyapatite(image courtesy University of Liverpool)¹⁴.

However, such coatings showed a number of drawbacks: cases of peeling of the coating from the titanium base^15,16, a moderate rate of osseointegration, and exposure to the environment. An important task therefore, among others, is to study the bond strength of the HAp coating with titanium^17,18. The strength characteristics of the bond-substrate interaction is the energy of this bond¹⁹. The aim of the present work is to determine the binding energies between hydroxyapatite functional groups (anions) and Ti(II) (titanium cation with plus 2 charge) - the standard material for implants. This is accomplished with computational methods. In the future, these components will be used to calculate the total binding energy of the unit cell HAp and the cation Ti(II).

2. Method: Calculation of Binding Energy by Density Functional Theory and Molecular Dynamics Simulations

2.1 General computational method

In computational chemistry, the Density Functional Theory (DFT) is used to describe the state of many-particle systems and to determine their geometric and chemical properties^20,21. In this paper it is used to determine the ground state energies of polyatomic complexes in the Ti (II) -hydroxyapatite system. The ultimate goal is the theoretical calculation of the binding energy of the HAp coating and Ti(II).

Using the formulas obtained with the help of DFT, we calculate the ground-state energies in the computational chemistry software suite, Gaussian 09 (Revision C.01 was used²²). The suite allows us to perform geometric optimizations of structures of interest, which result in obtaining incredibly specific nuclei coordinates in three-dimensional space. These geometry optimizations give the necessary positions of nuclei that deliver a global minimum on the potential energy surface. The program calculates and interprets the partial differential equations that deliver the solution with minimum energy to the Schrodinger equation. All calculations for the polyatomic complexes in the Ti (II) - hydroxyapatite system were performed using the DFT functional B3LYP (Becke three-parameter Lee-Yang-Parr functional) in Gaussian 09 system.

2.2 Mathematical and physical model: DFT

2.2.1 The time independent Schrödinger equation

The time independent (general) Schrödinger equation describes the interactions of N particles. It states that when the Hamiltonian, Ĥ, acts on a wave function, Ψ, it is proportional to the stationary state of the same wave function with a proportionality constant E. For our HAp calculations, we determine the lowest possible energy that the systems can have:

(1) $\hat{\mathit{H}}\mathit{\psi }\left({\mathit{x}}_1,{\mathit{x}}_2\mathrm{,...}{\mathit{x}}_{\mathit{N}}\right)=\mathit{E}\mathit{\psi }\left({\mathit{x}}_1,{\mathit{x}}_2\mathrm{,...}{\mathit{x}}_{\mathit{N}}\right)$

where E is the total energy for the system, and Ĥ is a Hamilton operator for a molecular system containing N nuclei and electrons:

(2) $\hat{\mathit{H}}=\mathit{T}+\mathit{V}$

where T is a kinetic energy and V is the potential energy, or

$\begin{array}{c}\hat{\mathit{H}}={\mathit{T}}^{\mathit{elec}}\left(\mathit{r}\right)+{\mathit{T}}^{\mathit{nuc}}\left(\mathit{R}\right)+{\mathit{V}}^{\mathit{nuc}-\mathit{elec}}\left(\mathit{R},\mathit{r}\right)+\hfill \\ {\mathit{V}}^{\mathit{nuc}-\mathit{nuc}}\left(\mathit{R}\right)+{\mathit{V}}^{\mathit{elec}-\mathit{elec}}\left(\mathit{r}\right)\hfill \end{array}$

where the two first terms are kinetic energies for electrons and nuclei, and the last three terms are potential energies of nuclei-electron attractions, nuclei-nuclei repulsions and electron-electron repulsions correspondingly.

2.2.2 Density functional theory and Born-Oppenheimer approximation

In Born-Oppenheimer (BO) approximation, which is an important aspect of Density Functional Theory (DFT), we assume that the kinetic energy of electrons is much greater than the kinetic energy of nuclei, and because of that the nuclei are believed to be fixed²¹. Their kinetic energy is zero, and the potential energy due to nucleus-nucleus repulsion is a constant. From this concept, for each arrangement of nuclei, the ground-state energy, E, (or full energy of a many-particle system) can be found on the base of DFT. Thus, at each moment in time the system is settled at the ground state energy, or at absolute zero, and therefore the velocity of nuclei are close to zero and the kinetic energy of nuclei is negligible small ( T^nuc = 0). Then the expression (3) becomes

(4) $\hat{\mathit{H}}={\mathit{T}}^{\mathit{elec}}\left(\mathit{r}\right)+{\mathit{V}}^{\mathit{nuc}-\mathit{elec}}\left(\mathit{R},\mathit{r}\right)+{\mathit{V}}^{\mathit{nuc}-\mathit{nuc}}\left(\mathit{R}\right)+{\mathit{V}}^{\mathit{elec}-\mathit{elec}}\left(\mathit{r}\right)$

or

$\begin{array}{c}\hat{\mathit{H}}=\sum _{\mathit{i}}\left(\left(\frac{{\partial }^2}{{\partial }^2\mathit{x}}\right)+\left(\frac{{\partial }^2}{{\partial }^2\mathit{y}}\right)+\left(\frac{{\partial }^2}{{\partial }^2\mathit{z}}\right)\right)-\sum _{\mathit{i}}\sum _{\mathit{J}}\left(\frac{{\mathit{Z}}_{\mathit{j}}}{\left|{\mathit{R}}_{\mathit{j}}-{\mathit{r}}_{\mathit{i}}\right|}\right)+\hfill \\ \sum _{\mathit{i}}\sum _{\mathit{j}<\mathit{i}}\left(\frac{1}{\left|{\mathit{r}}_{\mathit{i}}-{\mathit{r}}_{\mathit{j}}\right|}\right)+\sum _{\mathit{I}}\sum _{\mathit{J}<\mathit{I}}\left(\frac{{\mathit{Z}}_{\mathit{i}}{\mathit{Z}}_{\mathit{j}}}{\left|{\mathit{R}}_{\mathit{i}}-{\mathit{R}}_{\mathit{j}}\right|}\right)\hfill \end{array}$

is called the Electronic Hamiltonian. It is easy to see that, in this case, in the Schrödinger equation

(5) ${\hat{\mathit{H}}}^{\mathit{elec}}{\mathit{\psi }}^{\mathit{elec}}\left(\overrightarrow{\mathit{r}},\overrightarrow{\mathbf{R}}\right)={\mathit{E}}^{\mathit{eff}}\left(\overrightarrow{\mathbf{R}}\right){\mathit{\psi }}^{\mathit{elec}}\left(\overrightarrow{\mathit{r}},\overrightarrow{\mathbf{R}}\right)$

the eigenvalues E^eff will give the ground state energy of the nuclei, or ab-initio lowest energy E _min, and that is a unique function of $\overrightarrow{\mathit{R}}$, the coordinates of the nuclei. The ab initio lowest total energies can be evaluated for many configurations, which form the minimum energy surface, or ab-initio potential-energy surface (see Figure 3), for a given reaction^23-28.

According to DFT, the reactants and products are both at the positions of the absolute minima on the potential energy surface. From energy principals, the path between reactants and products goes in the direction of minimum spending energy, which is viewed as a path through a saddle point between two global minima on the potential energy surface (see Figure 3). The difference between the saddle point energy and the reactant energy will give the activation energy barrier, E _act, or the energy that is needed to be overcome to form a bond, and the difference in reactants energy and the product energy, ΔG, will give the binding energy, or the energy stored in the product compared with the reactants. In this paper, the subject of our interest is the binding energy only. In the future however, we will be analyzing our results for the activation energies as well, which are related with the manufacturing technologies and processes for spattering nano-coatings on the implant.

Figure 3
Potential-energy surface (PES), activation energy and binding energy.

2.2.3 Density functional theory basics (electronic optimization)

Multiplying Eg. (5) by ${\stackrel{\leftharpoonup }{\mathit{\psi }}}^{\mathit{elec}}$

${\stackrel{\leftharpoonup }{\mathit{\psi }}}^{\mathit{elec}}{\hat{\mathit{H}}}^{\mathit{elec}}{\mathit{\psi }}^{\mathit{elec}}\left(\overrightarrow{\mathit{r}},\overrightarrow{\mathbf{R}}\right)={\mathit{E}}^{\mathit{eff}}\left(\overrightarrow{\mathbf{R}}\right){\stackrel{\leftharpoonup }{\mathit{\psi }}}^{\mathit{elec}}{\mathit{\psi }}^{\mathit{elec}}$

and taking into account the normalization of wave functions, this expression can be rewritten

(6) ${\stackrel{\leftharpoonup }{\mathit{\psi }}}^{\mathit{elec}}{\hat{\mathit{H}}}^{\mathit{elec}}{\mathit{\psi }}^{\mathit{elec}}\left(\overrightarrow{\mathit{r}},\overrightarrow{\mathbf{R}}\right)={\mathit{E}}^{\mathit{eff}}$

which means the minimization of quadratic functional with respect to Ψ^elec will give us the true energy minimum, E^eff. In Density Functional Theory Energy ρ (r), the electron density is introduced as

(7) $\mathit{\rho }\left(\mathit{r}\right)=2\sum _{\mathit{i}}^{\raisebox{1ex}{$\mathit{N}$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}{\mathit{\psi }}_{\mathit{i}}\left(\mathit{r}\right){\mathit{\psi }}_{\mathit{i}}\left(\mathit{r}\right),$

Then the quadratic functional taking into account expression for the Hamiltonian (4) can be written as

(8) $\begin{array}{c}\mathit{E}\left[\mathit{\rho }\right]=\mathit{T}\left[\mathit{\rho }\right]+\int {\mathit{V}}_{\mathit{ext}}\left(\mathit{r}\right)\mathit{\rho }\left(\mathit{r}\right)\mathit{dr}\hfill \\ +\frac{1}{2}\int \mathit{dr}\int \mathit{d}{\mathit{r}}^{\prime }\frac{\mathit{\rho }\left(\mathit{r}\right)\mathit{\rho }\left({\mathit{r}}^{\prime }\right)}{\left|\mathit{r}-{\mathit{r}}^{\prime }\right|}+{\mathit{E}}_{\mathit{xc}}\left[\mathit{\rho }\right]\hfill \end{array}$

Hohenberg-Kohn Theorems²¹:

HK1E[ ρ ] is a unique functional of ρ (r).

HK2 Functional minimized with respect to ρ (r) will give true ground state energy (which means E_min), or E[ ρ ] reaches its minimum value for the real ground state density, ρ (r).

In the expression (7)Ψ _i is the Kohn-Sham orbital satisfying the following orbital equation

(9) $\mathit{H}{\mathit{\psi }}_{\mathit{i}}\left(\mathit{r}\right)=\left(-\frac{1}{2}{▽}^2+{\mathit{V}}_{\mathit{eff}}\left(\mathit{r}\right)\right){\mathit{\psi }}_{\mathit{i}}\left(\mathit{r}\right)={\epsilon }_{\mathit{i}}{\mathit{\psi }}_{\mathit{i}}\left(\mathit{r}\right)$

with V_eff as the Kohn-Sham effective potential that depends in turn on the electron density. The expression (8) for total energy can be expressed as

(10) $\mathit{E}\left[\mathit{\rho }\right]=2\sum _{\mathit{i}}^{\raisebox{1ex}{$\mathit{N}$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}{\epsilon }_{\mathit{i}}+\mathit{Q}\left[\mathit{\rho }\right]+\sum _{\mathit{A},\mathit{B}}\raisebox{1ex}{${\mathit{Z}}_{\mathit{A}}{\mathit{Z}}_{\mathit{B}}$}\!\left/ \!\raisebox{-1ex}{${\mathit{R}}_{\mathit{AB}}$}\right.$

where

(11) $\mathit{Q}\left[\mathit{\rho }\right]=\int \mathit{\rho }\left[\raisebox{1ex}{$-\Phi \left(\mathit{r}\right)$}\!\left/ \!\raisebox{-1ex}{$2-{\mathit{V}}_{\mathit{xc}}\left(\mathit{r}\right)$}\right.\right]\mathit{dr}+{\mathit{E}}_{\mathit{xc}}\left[\mathit{\rho }\right]$

and

(12) $\Phi \left(\mathit{r}\right)=\int \frac{\mathit{\rho }\left({\mathit{r}}^{\prime }\right)\mathit{d}{\mathit{r}}^{\prime }}{\left|\mathit{r}-{\mathit{r}}^{\prime }\right|}$

E_xc is the exchange-correlation functional and

(13) ${\mathit{V}}_{\mathit{xc}}\left(\mathit{r}\right)=\raisebox{1ex}{$\mathit{\delta }{\mathit{E}}_{\mathit{xc}}\left[\mathit{\rho }\right]$}\!\left/ \!\raisebox{-1ex}{$\mathit{\delta }\mathit{\rho }\left(\mathit{r}\right)$}\right.$

3. Results: Vibrational Analysis of Stationary Structures for Different Configurations of HAp Component ions

3.1 Calculation method

This paper represents the first step in solving the general problem of determining the binding energy of HAp with titanium. We determine the binding energies between individual fragments (anions) of hydroxyapatite and titanium. In general, a Ti²⁺ atom was paired with various anions from HAp (see Figure 2). In Gaussian, it requires the initial “guess” set up of positions of the atoms. The more complex the final structure is, the more important the initial guess. This is due to the fact that as the geometry gets more complicated, the charge distribution among the atoms in the structure plays an important role in the structure of the final optimized geometry. The charge distribution may cause atoms in the final geometry to repel and attract one another, resulting in particular geometries (e.g. angles and bond lengths).

First, by using the B3LYP hybrid-exchange-correlation functional, we calculated ground state energies for different configurations of hydroxyapatite, namely for OH^-, PO₄ ^3- (which are anions), and then for Ti²⁺ (which we took with charge 2+ to have attraction between ions) (Figure 2). Then, we calculated various combinations of Ti²⁺, OH^-, PO₄ ^3- (that is when the initial guesses become more important). Once these energies were found, we calculated the binding energies (BE) by subtracting the ground state energies of products from the ground state energies of reactants.

To achieve convergence during calculations, we used fine integration grid and basis set 3-21G during optimization. Also, frequency calculations were performed before geometric optimization. If frequencies are negative, the initial configurations are not correct and ground state energy is not found and the optimized position of nuclei is not at global minimum. If that is the case, we need to figure out the new initial guess configurations and repeat the optimization process from the beginning. For all the combinations considered, the main energy levels and binding energies are calculated, with particular emphasis on the equilibrium bond lengths and the angles of Ti-O bonds. In the next stage of work, this data will be used in calculating the total binding energy of the anion Ti (II) and the whole unit cell HAp. Depending of the complexity of the final configuration, the calculation time varied from several minutes to several hours.

3.2 Results of calculations

The results are shown in Fig. 4-8. For example, let us briefly summarize the general results of one of the calculations. According to calculations of the complex [TiOH]+ (Figure 4), the charge of the Ti (II) ion decreases. In terms of the Density Functional Theory (DFT), this means that Ti (II) experiences an increase in the electron density due to its attraction with negatively charged hydroxide ions, OH^-. The length of the reference bond is about 1.8 Å. For the vibrational frequency, which indicates as interaction between Ti (II) and the oxygen atom, the calculation gives a reasonable value of about 970 cm^-1. The Hessian index contains all positive eigenvalues ​​of a constant-force matrix and corresponds to real vibrational frequencies of bonds within molecules, or in this case, interacting ions.

Figure 4
Stable interaction configuration of Ti (II) with oxygen in the complex [TiOH] ⁺ (Ti - larger grey atom, O - red atom, H - smaller grey atom) - linear structure.

3.3 Analysis of results: differences of electron densities of titanium for component ions of HA are capital

3.3.1 Titanium(II) and the hydroxide ion: [TiOH]⁺

In the first structure, the Titanium (II) ion experiences a decrease of almost half of its initial charge and demonstrates an expected charge distribution. The initial charge of Ti(II) is 2+, however, after its reaction with the negatively charged hydroxide (the oxygen atom in particular), the Titanium ion loses approximately .876 charge (Table 1 and Figure 4). This information provides several implications about the function of titanium in the [TiOH]⁺ structure and its centralized relationship with the oxygen atom. The orbital space around titanium is now host to more electrons, which suggests an increase in electron density around the titanium atom²⁹. The [TiOH]⁺ structure is a useful reference as each configuration grows in complexity, especially to make inferences about the maximum amount of phosphate ions that can form a bond with a Ti (II) ion and the overall prediction of the binding energy for HAp. This structure demonstrates the lowest binding energy with Ti (II) at .734 atomic units (Table 1).

Table 1
Binding Energy Calculations of Ti (II) With Constituent Ions of Hydroxyapatite.

	Ti2+	OH-	[TiOH]+	Ti(OH)2	PO43-	[TiPO4]-	[HO5PTi]-2	[H2O6PTi]3-
Picture
Ti-O Bond Lengths (Å)	N/A	N/A	1.8	1.73	N/A	1.94 per bond to Oxygen	1.82	1.87
Ti-O Bond Angles (degrees)	N/A	180	180	126	N/A	Approx. 74 per bond to Oxygen	81.4	78.0
Ground Energy of Structure (a.u.)	-844.42	-75.246	-920.4	-996.2	-638.2	-1484.5	-1559.9	-1635.2
Charge of Titanuim	2+	N/A	1.124+	0.486+	N/A	0.124+	-0.215	-0.888
Binding Energy (a.u.)	N/A	N/A	.734	1.28	N/A	1.88	2.03	2.09

3.3.2 Titanium (II) and two Hydroxide ions: Ti(OH)₂

A neutrally charged compound, Ti(OH)₂ introduces electrons from another oxygen atom to the orbital space of the titanium atom. The point charge of Ti (II) in the Ti(OH)₂ structure is approximately 0.486 compared to the 1.12 charge of Ti (II) in the [TiOH]⁺ structure (Table 1 and Figure 5). As mentioned earlier, the bond length between oxygen and titanium in [TiOH]⁺ is approximated at 1.8 Ǻ. The bond length between oxygen and titanium in Ti(OH)₂ is approximated at 1.73 Ǻ, which suggests a stronger bonding strength between titanium and oxygen in the structure. This is additionally supported by a stronger binding energy than [TiOH]⁺ at 1.28 atomic units. The binding energy of the structure nearly doubles and this can primarily be attributed to the strong electronic differences between the oxygen and titanium atoms.

Figure 5
Stationary geometry of the Ti (II) interaction with oxygen in the [Ti (OH) ₂ ] complex (Ti - larger grey atom, O - red atom, H - smaller grey atom) (a) symmetric stretch (b) antisymmetric stretch; slightly larger than trigonal planar geometry due to O-O repulsion.

3.3.3 Titanium (II) and the phosphate ion: [TiPO₄]^-

In this structure, the titanium atom loses a significant amount of charge due to electrons from three oxygen atoms occupying its orbital space. The oxygen atoms lose much of their electrons and the titanium atom results in a charge of 0.124 during the interaction (Figure 6). The effect that this change in electron density has on the bond lengths between the oxygen atoms, titanium, and phosphorus is such: three of the oxygen atoms that contain most of the negative charge in the phosphate ion are pulled towards the titanium atom and, in consequence, the remaining oxygen atom develops a stronger bond to the phosphorus atom in [TiPO₄]^2-. As a reference for later binding energies between phosphate, titanium, and hydroxide, it is important to note that the binding energy for [TiPO4]^2- is calculated at approximately 1.88 atomic units.

Figure 6
Stationary geometry of the Ti²⁺ interaction with Oxygen in the complex [TiPO ₄ ] ^- (Ti - larger grey atom, O - red atom, P - orange atom, H - smaller grey atom) (a) three O atoms oscillate (b) Ti, P and a single O oscillate - structure becomes even more stable due to the presence of additional negatively charged O atoms compared to the previous configurations.

3.3.4 Titanium (II), hydroxide, and phosphate: [Ti(OH)PO₄]^2-

In this configuration, this is the first instance where the titanium atom will have a negative charge (Figure 7). The trend in electron density of titanium is such that titanium consistently gains electron density at a much more dramatic rate than the phosphorus atom. The phosphorus atom is also consistently gaining electron density from the oxygen atoms, but at much smaller amounts than the titanium atom.

Figure 7
Stationary geometry of Ti (II) interaction with oxygen in the complex [Ti(OH)PO ₄ ]^2- (Ti - larger grey atom, O - red atom, P - orange atom, H - smaller grey atom) (a) antisymmetric centralized oscillation (b) symmetric centralized oscillation; quasi-chair structure, which is a favorable structure.

3.3.5 Titanium (II), two hydroxides, and One phosphate: [Ti(OH)₂PO₄]^3-

In the final and most complex structure, the oxygen atoms are contributing even more electron density to the orbital space of the primarily titanium. This results in a centrally negative structure, with all of the positive charge belonging to the hydrogen and phosphorus atoms (Figure 8). Most of the positive charge belongs to phosphorus with a charge of 1.01. The charge of phosphorus in the [TiPO₄]^- was 1.19, which means phosphorus will experience a maximum charge loss of only 0.180. This is considerably smaller than the maximum charge loss of titanium in the [TiOH]⁺ structure, with a maximum charge decrease of approximately 2.00. All of this seems to suggest very favorable trends and interactions with the titanium (II) ion and the component ions of Hydroxyapatite. This particular structure may also demonstrate the maximum amount of phosphate that can form a bond with Ti (II), due to an insufficient positive charge and energy of Ti (II).

Figure 8
Stationary geometry of Ti(II) with oxygen in [Ti(OH) ₂ PO ₄ ]^3- (Ti - larger grey atom, O - red atom, P - orange atom, H - smaller grey atom) (a) longitudinal oscillation centralized between Ti and P (b) transverse oscillation between Ti and Phosphate (c) longitudinal oscillation centralized between Ti and P with strong Hydrogen-Hydrogen repulsion - boat structure, which is the most stable configuration, since Titanium and P interact because of Ti having now a negative charge.

4. Discussion: Differences in Geometries and Vibrational Behavior of Titanium (II) Bonds With Component Ions of HAp

4.1 Titanium(II) and the hydroxide ion: [TiOH]⁺

Once again [TiOH]⁺ is a convenient initial structure to act as a reference for the following structures. [TiOH]⁺ is a linear structure with a bond length of 1.8 Ǻ between the oxygen and the titanium atoms (Figure 4). The repulsion between the positively charged atoms titanium and hydrogen are likely the cause for this linearity. The calculated frequency interaction between titanium and oxygen is approximately 970 cm^-1, a reasonable vibration³⁰.

4.2 Titanium (II) and two hydroxide ions: Ti(OH)₂

In the Ti(OH)₂ structure, higher frequencies are calculated between the titanium and oxygen atoms. The geometry of the structure changes to that of being slightly larger than the trigonal planar geometry (Figure 5). The geometry is slightly larger due to the oxygen-oxygen repulsion of Ti(OH)₂. Additionally, it is important to recognize that the bond length between the oxygen atoms and the titanium atom has decreased to 1.73 Ǻ when compared to the [TiOH]⁺ structure, whose Titanium-oxygen bond length is 1.80 Ǻ. A larger interaction exists between oxygen and titanium in the Ti(OH)₂ structure, which suggests a favorably higher binding energy to the titanium atom.

4.3 Titanium (II) and the phosphate ion: [TiPO₄]^-

[TiPO₄]^- (Figure 6) is useful as a reference for the other phosphate configurations in terms of their geometries. Initially, the bond lengths between titanium and oxygen are 1.94 Ǻ per bond. Additionally, the bond angles between titanium and oxygen are 74 degrees per bond. There exists one oxygen atom that is not paired to the titanium atom and instead develops a stronger and shorter bond to the phosphorus atom after reaction (Figure 6).

4.4 Titanium (II), hydroxide, and phosphate: [Ti(OH)PO₄]^2-

This structure is a special geometry of interest, as it exhibits a quasi-chair structure (Figure 7). This is a favorable geometry due to a more mechanically reliable structure³¹. When compared to the standard vibrational behavior of phosphate, the increase in frequency between the titanium and oxygen interaction suggests a higher energy in the structure. Additionally, the bond lengths between titanium and oxygen have shortened, also suggesting a more favorable binding energy.

4.5 Titanium (II), two hydroxides, and one phosphate: [Ti(OH)₂PO₄]^3-

Now let us consider the case where another hydroxide reacts with [Ti(OH)PO ₄]^2- . A tetrahedral geometry is formed and the titanium atom continues to pull oxygen atoms away from the phosphorus atom. However, since the titanium atom also has electrons from two hydroxide ions occupying its orbital space, the bond lengths between the titanium and oxygen atoms are slightly larger than that of the [Ti(OH)PO ₄]^2- structure.

5. Conclusions

5.1 Trends of electron density

Regarding the information from Table 1, the change in charge of Ti(II) appears to occur at a constant rate. As the structure gains more oxygen atoms, the electron density increases much more significantly than that of the phosphorous atom’s charge. This seems to suggest that the most influential interaction regarding the binding energy of Ti(II) and HAp is the interaction between titanium and oxygen. As the electron density of titanium increases at a constant rate for each structure, there is a direct relationship with other aspects of the structures such as stronger bonds, higher interaction frequencies, and more favorable geometries.

5.2 Geometries and binding energies for each structure.

Additionally, the geometries become more favorable as the complexity of the structures grow. With the [Ti(OH)PO₄]^2- (Figure 7) the quasi-chair geometry contains a four-member ring with low bond angles and high energy. There is an open space on top of the structure near the hydroxide, which would be an area available and energetically favorable for additional interactions. This most likely explains why the [Ti(OH)₂PO₄]^3- structure (refer to Figure 8) contains the final hydroxide on the opposite side of the titanium atom, as the titanium atom contains a more positive charge than either of the two oxygen atoms at the top of the structure. This is the closest spot the hydroxide can bond to the lower energy part of [Ti(OH)₂PO₄]^3- (the top space that contains the open oxygen and hydroxide). The [Ti(OH)₂PO₄]^3- structure displays the most complex geometry as well as the lowest ground state energy, which according to DFT suggests the highest binding energy compared to the other structures.

5.3 There seems to be a constant relationship between vibrational frequencies, bond lengths, geometries, and electron densities with the overall binding energy for Ti(II) for each compound.

As stated in 5.1, as the electron density of titanium increases with each structure, the bonding strength increases between titanium and its neighbor oxygen atoms. Titanium is predominantly influential with its drastic change in charge compared to the phosphorous atom, which has very little change in charge per structure. As the structures increase in complexity, the titanium atom’s charge becomes increasingly more negative, which results in a more centrally negative structure and a lower ground state energy. Since the orbital space of titanium is increasingly occupied with electrons from neighbor oxygen atoms, there is a limited number of phosphate ions that can form a bond with titanium. This suggests that the phosphate-titanium (II) interactions may act as a parameter for this model, which will constrain its scope.

5.4 What can the conclusions suggest about the overall strength of HAp and what is the “right” bonding strength for the interaction with bone, so that it doesn’t interfere with necessary bodily healing functions.

In biological and medical practices, it is typical to consider if the binding energy of chemicals that attach to compounds and tissues found in the human body will interfere with necessary biological functions of the body^32,33. Determining an accurate maximization optimization of the binding energy of HAp to Ti(II) will be useful for practicalities such as ensuring that the bonding strength of HAp and Ti(II) doesn’t interfere with the natural healing process. However, considering the grafting and healing capabilities that the HAp has with bone, it may be of best interest to maximize the binding energy of HAp with Ti(II) for medical procedures.

6. Acknowledgements

The work was done with partial support of RFBR grants No. 17-08-01579 and No. 17-08-01312 and of state grant #AAAA-A17-117021310386-3. This work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number ACI-1548562. We thank consultants for their assistance with porting, optimization and visualization, which was made possible through the XSEDE Extended Collaborative Support Service (ECSS) program ^34,35.

7. References

• ¹
  Misch CE, Perel ML, Wang HL, Sammartino G, Galindo-Moreno P, Trisi P, et al. Implant success, survival, and failure: the International Congress of Oral Implantologists (ICOI) Pisa Consensus Conference. Implant Dentistry 2008;17(1):5-15.
• ²
  Albrektsson T, Johansson C. Osteoinduction, osteoconduction and osseointegration. European Spine Journal 2001;Suppl 2:S96-101.
• ³
  Xuereb M, Camilleri J, Attard NJ. Systematic review of current dental implant coating materials and novel coating techniques. International Journal of Prosthodontics 2015;28(1):51-59.
• ⁴
  Menicucci G, Mussano F, Schierano G, Rizzati A, Aimetti M, Gassino G, et al. Healing properties of implants inserted concomitantly with anorganic bovine bone. A histomorphometric human study. Australian Dental Journal 2013;58(1):57-66.
• ⁵
  Smeets R, Stadlinger B, Schwarz F, Beck-Broichsitter B, Jung O, Precht C, et al. Impact of Dental Implant Surface Modifications on Osseointegration. BioMed Research International 2016;2016:6285620.
• ⁶
  Eliaz N, Metoki N. Calcium Phosphate Bioceramics: A Review of Their History, Structure, Properties, Technologies and Biomedical Applications. Materials 2017;10(4):334.
• ⁷
  Dorozhkin SV. Nanodimensional and Nanocrystalline Hydroxyapatite and Other Calcium Orthophosphates. In: Gshalaev VS, Demirchan AC, eds. Hydroxyapatite: Synthesis, Properties, and Applications New York: Nova Science Publishers; 2012. p. 1-90.
• ⁸
  Makarova E, Rubstein A, Trakhtenberg I, Kisilyova N, Bliznets D. Osteointegration of Titanium Porous Implants with Carbon Nanocoating and Hydroxyapatite Particles into the Pores. In: Gshalaev VS, Demirchan AC, eds. Hydroxyapatite: Synthesis, Properties, and Applications New York: Nova Science Publishers; 2012. p. 431-446.
• ⁹
  Mendes V, Davies JE. Discrete Calcium Phosphate Nanocrystalline Deposition Enhances Osteoconduction on Titanium-based Implant Surfaces In: Canadian Biomaterials Society, 25th Annual Meeting; 2006 May 26-28; Calgary, AB, Canada.
• ¹⁰
  Mendes V, Davies JE. Discrete Calcium Phosphate Nanocrystals Render Titanium Surfaces Bone Bonding. International Journal of Oral & Maxillofacial Implants 2007;22:484
• ¹¹
  Kenealy J, Pharm D, Berckmans B, Stach R. Nanometer-scale CaP Enhances Early Implant-bone Fixation in an Animal Model In: European Association for Osseointegration, 15th Annual Meeting; 2006 Oct. 5-7; Zurich, Switzerland.
• ¹²
  Suttin Z, Gubbi P. Surface Area Increase Due to Discrete-Crystalline-Deposition of Nanometer-Scale CaP Crystals In: European Association for Osseointegration, 15th Annual Meeting; 2006 Oct. 5-7; Zurich, Switzerland.
• ¹³
  Bigi A, Bracci B, Cuisinier F, Elkaim R, Fini M, Mayer I, et al. Human osteoblast response to pulsed laser deposited calcium phosphate coatings. Biomaterials 2005;26(15):2381-2389.
• ¹⁴
  University of Liverpool. ChemTube3D Available from: <http://www.chemtube3d.com>. Access in: 03/2018.
  » http://www.chemtube3d.com
• ¹⁵
  Roy M, Bandyopadhyay A, Bose S. Induction plasma sprayed nano hydroxyapatite coatings on titanium for orthopaedic and dental implants. Surface and Coatings Technology 2011;205(8-9):2785-2792.
• ¹⁶
  Pezeshki P, Lugowski S, Davies J. Dissolution of Discrete Calcium Phosphate Crystals From Candidate Ti-based Implant Surfaces In: Canadian Biomaterials Society, 25th Annual Meeting; 2006 May 26-28; Calgary, AB, Canada.
• ¹⁷
  Shekhawat G, Suttin Z. Determination of CaP Crystal Shear Strength Using contact-Mode Atomic Force Microscopy (AFM) Palm Beach Gardens: BIOMET 3i(tm).
• ¹⁸
  Suttin Z, Gubbi P. Adhesion Shear Strength of Nanometer-scale CaP Crystals Applied by Discrete-Crystalline-Deposition In: European Association for Osseointegration, 15th Annual Meeting; 2006 Oct. 5-7; Zurich, Switzerland.
• ¹⁹
  Balueva AV, Dashevskiy IN. Atomistic Modeling and Strength Calculations of Ni-Al Submicron and Nanosized Composites In: Proceedings of the 7th International Conference on Mechanics and Materials in Design; 2017 Jun 11-15; Albuferia, Portugal. p. 1501-1502.
• ²⁰
  Laird BB, Ross RB, Ziegler T. Chemical Applications of Density Functional Theory Washington: American Chemical Society; 1996.
• ²¹
  Koch W, Holthausen M. A Chemist's Guide to Density Functional Theory Berlin, Heidelberg: Wiley-VCH; 2001. 295 p.
• ²²
  Frisch MJ, Trucks GW, Schlegel HB, Scuseria GE, Robb MA, Cheeseman JR, et al. Gaussian 09, Revision A02, Wallingford: CT: Gaussian, Inc.; 2016. Available from: http://gaussian.com/glossary/g09/
  » http://gaussian.com/glossary/g09/
• ²³
  Gross A, Wilke S, Scheffer M. Six-dimensional quantum dynamics of adsorption and desorption of H2 at Pd(100): no need for a molecular precursor adsorption state. Surface Science 1996;357-358:614-618.
• ²⁴
  Gross A, Wei CM, Scheffler M. Poisoning of hydrogen dissociation at Pd (100) by adsorbed sulfur studied by ab-initio quantum dynamics and ab-initio molecular dynamics. Surface Science 1998;416(1-2):L1095-L1100.
• ²⁵
  Jónsson H, Mills G, Jacobsen KW. Nudged elastic Band Method for finding minimum energy paths of transitions. In: Berne BJ, Ciccotti G, and D.F. Coker DF, eds. Classical and Quantum Dynamics in Condensed Phase Simulations Chapter 6. Hackensack, NJ: World Scientific; 1998.
• ²⁶
  Elber R, Karplus M. Enhanced sampling in molecular dynamics: use of the time-dependent Hartree approximation for a simulation of carbon monoxide diffusion through myoglobin. Journal of the American Chemical Society 1990;112(25):9161-9175.
• ²⁷
  Carr JM, Trygubenko SA, Wales DJ. Finding pathways between distant local minima. Journal of Chemical Physics 2005;122(23):234903.
• ²⁸
  Trygubenko SA, Wales DJ. A Doubly Nudged Elastic Band Method for Finding Transition States. Journal of Chemical Physics 2004;120(5):2082-2094.
• ²⁹
  Balueva AV, Dashevskiy IN. Mathematical Modelling of Technologies for the Design of Innovative Functional Bio-coatings for Dental Implants In: Proceedings IRF2018: 6th International Conference Integrity-Reliability-Failure; 2018 Jul 22-26; Lisbon, Portugal.
• ³⁰
  Kwon KD, Kubicki JD. Molecular Orbital Theory Study on Surface Complex Structures of Phosphates to Iron Hydroxides: Calculation of Vibrational Frequencies and Adsorption Energies. Langmuir 2004;20(21):9249-9254.
• ³¹
  Dashevskiy I, Balueva A, Todebush P. Calculation of the Binding Energy of Implant Coating Components with a Substrate by Quantum Chemistry Methods In: 3rd International Symposium on Physics, Engineering and Technologies for Biomedicine; 2018 Oct 15-17; Moscow, Russia. p. 94-96. Available from:http://plasma.mephi.ru/ru/uploads/files/conferences/2018PBS/20181017_PBS2018_BookOfAbstracts.pdf
  » http://plasma.mephi.ru/ru/uploads/files/conferences/2018PBS/20181017_PBS2018_BookOfAbstracts.pdf
• ³²
  Fernandes GVO, Alves G, Linhares ABR, da Silva MHP, Granjeiro JM. Evaluation of Cytocompatibility of Bioglass-Niobium Granules with Human Primary Osteoblasts: A Multiparametric Approach Key Engineering Materials. 2012;493-494:37-42.
• ³³
  Dashevskiy IN, Gribov DA. Patient-Specific Biomechanical Analysis in Computer Planning of Dentition Restoration with the Use of Dental Implants In: Second International Conference of Artificial Intelligence, Medical Engineering, Education (AIMEE 2018); 2018 Oct 6-8; Moscow, Russia.
• ³⁴
  Towns J, Timothy C, Dahan M, Foster I, Gaither K, Grimshaw A, et al. XSEDE: Accelerating Scientific Discovery. Computing in Science & Engineering 2014;16(5):62-74. DOI: 10.1109/MCSE.2014.80
  » https://doi.org/10.1109/MCSE.2014.80
• ³⁵
  Wilkins-Diehr N, Sanielevici S, Alameda J, Cazes J, Crosby L, Pierce M, et al. An Overview of the XSEDE Extended Collaborative Support Program. In: Gitler I, Klapp J, eds. High Performance Computer Applications 6th International Conference, ISUM 2015, Mexico City, Mexico, March 9-13, 2015, Revised Selected Papers Cham: Springer International Publishing; 2016. DOI: 10.1007/978-3-319-32243-8
  » https://doi.org/10.1007/978-3-319-32243-8

Publication Dates

• Publication in this collection
  30 May 2019
• Date of issue
  2019

History

• Received
  16 Jan 2019
• Reviewed
  22 Mar 2019
• Accepted
  08 Apr 2019