Complex interplay of interatomic bonding in a multi-component pyrophosphate crystal: K2Mg (H2P2O7)2·2H2O

The electronic structure and interatomic bonding of pyrophosphate crystal K2Mg (H2P2O7)2·2H2O are investigated for the first time showing complex interplay of different types of bindings. The existing structure from single-crystal X-ray diffraction is not sufficiently refined, resulting in unrealistic short O─H bonds which is rectified by high-precision density functional theory (DFT) calculation. K2Mg (H2P2O7)2·2H2O has a direct gap of 5.22 eV and a small electron effective mass of 0.14 me. Detailed bond analysis between every pair of atoms reveals the complexity of various covalent, ionic, hydrogen bonding and bridging bonding and their sensitive dependence on structural differences. The K--O bonds are much weaker than Mg--O bonds and contributions from the hydrogen bonds are non-negligible. Quantitative analysis of internal cohesion in terms of total bond order density and partial bond order density divulges the relative importance of different types of bonding. The calculated optical absorptions show multiple peaks and a sharp Plasmon peak at 23 eV and a refractive index of 1.44. The elastic and mechanical properties show features unique to this low-symmetry crystal. Phonon calculation gives vibrational frequencies in agreement with reported Raman spectrum. These results provide new insights indicating that acidic pyrophosphates could have a variety of unrealized applications in advanced technology.


Introduction
Phosphates is a vast and complex subject with a long history intersecting multiple disciplines of physics, chemistry and biology. However, it is still not fully understood compared to another type of well-known inorganic crystals, the silicates [1]. The first research on phosphates can be traced to two centuries ago with the preparation of Na 4 P 2 O 7 [2]. The vastly different classifications and overlapping nomenclatures of phosphates without consensus across different communities cause great difficulties and confusion among researchers. It is generally agreed that it originated with phosphoric acid H 3 PO 4, and the P atom is pentavalently bonded to four O atoms in forming a tetrahedral [PO 4 ] −3 unit. Within this context, condensed phosphates can be described as compounds formed with cations consisting of alkali or other metallic elements and phosphoric anions denoted by (P n O 3n+1 ) −(n+2) (n < 20) or their various combinations with or without the presence of water molecules [1].
Biologists believe that Nature chooses phosphates to be part of the genetic material such as DNA, RNA and ATP because phosphoric acid can meet the stringent condition required for the formation of linkage via hydrolytically stable ester bonds [3]. This leads to the assertion of evolutionary centrality of phosphates [3], as pointed out by Todd that 'Where there is life, there is phosphorus' [4]. Organic phosphates and phosphate-based glasses are also an important part of biomaterials which is another extremely active arena of research [5,6]. There are essentially two types of bio-glasses: (i) a network glass such as ternary phosphate glass (CaO) x (Na 2 O) 1−x−y (P 2 O 5 )y in the form of a continuous random network similar to silicate glasses and (ii) aqueous solution of orthophosphate species H n (PO 4 ) 3−n (n = 0-3) important for biological reactions. Thus phosphates play an important role in biomedical and life sciences.
Applications of phosphates are ubiquitous because of their diverse compositions and unique complex structures. There are many biomedical and biochemical applications such as special drugs and their deliveries [7], bone repair, tissue engineering [8] and periodontal ligament fibroblasts for cell regeneration [9], etc. In energy science, they are known as predominating materials for Li or Na ion batteries [10,11]. They are also used in a variety of sensors [12,13] and in water oxidation catalysts [14]. Some phosphates such as KTP are known for their exceptionally large second harmonic generation for special electro-optical applications [15,16]. Hence, phosphates are critical to the development of modern technology.
In modern materials research, computational study always accompanies experimentation in any class of materials, and phosphates are no exception. A limited amount of computational works were done in tri-and diphosphates by Hansia et al. [17], interaction of pyrophosphate with hydroxyapatite by Rivas et al. [18] and cluster model studies on vanadyl pyrophosphates [19]. However, most of the studies in the past concentrate on the structures and geometry with little insights on the fundamental interactions at the atomic level and seldom touch on their properties which are important in their existing and emerging applications. This actually does not come as a surprise, given the confusing status of their definitions and many unsolved issues such as the accuracy of the measured structure, and the appropriate methods for different and diverse structures of phosphate. Most of the computational studies focus on molecular clusters rather than crystalline band structures. They also concentrate more on the vibrational spectra rather than the electronic structure. Quantitative information on different types of bonding in multicomponent phosphates are missing. It turns out that we actually have several studies on different types of phosphates within the past decade with the earlier version of the methods we developed. They cover inorganic phosphates [20][21][22], apatite crystals and their spectroscopic or mechanical properties [23][24][25], and some major phosphor-olivines Li x MPO 4 [26,27]. More recently, calculation has been extended to complex biomolecular systems containing phosphorus [28][29][30][31][32].
In this study, we focus on a particular acidic pyrophosphate crystal K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O because of its well-characterized crystal structure and detailed vibrational analysis on Raman and infrared  [33]. To the best of our knowledge, there have not been any electronic structure information for this crystalate or similar pyrophosphates. Here, we report the first such calculation with detailed interatomic bonding together with its optical and mechanical properties. Moreover, we show that, for such crystals, the experimentally determined structure may not be sufficiently accurate due to the difficulty in locating exact positions of H. We first describe its complex structure followed by a brief description of the computational methods used. The calculated electronic structure and properties are presented and discussed in §4. We end with a brief conclusion and discussion on the future prospect of research and potential applications of pyrophosphates in general.  (2), Z = 2 with a total of 31 atoms in the unit cell as shown in figure 1a. The original synthesis of this pyrophosphate was done by Harcharras et al. [33] in 2003 and the single-crystal structure was obtained using X-ray diffractometry (XRD). The unit cell refinement used 1055 of the 1446 reflections observed by omitting those with intensity I less than 2σ /(I obs ). The final structure was resolved using the direct method of SIR97 [34] and refined by a full-matrix least square technique based on F2-SHELXL-97 [35] with residue R = 0.0368 and supposed to be extremely accurate. Four years later, one of us (R.K.) used a new synthesis method, the co-precipitation method according to the reaction scheme: 2 K 4 P 2 O 7 + MgCl 2 + 4HCl to give K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O + 6KCl (dissolved). The recovered single crystal was then analysed in the same laboratory using the same high-precision procedure. The structure differs only slightly from the one published earlier by Harcharras et al. [33] (see electronic supplementary material for comparison). This structure is now referred to as the experiment XRD or unrelaxed structure in the present discussion. Figure 1b shows the polyhedral picture of this crystal in a 2 × 2 × 2 supercell showing layers of MgO 6 octahedra and PO 4 tetrahedra, and the locations of the K ions. However, it was observed that there are two very short O−H bonds of 0.636 Å, 0.650 Å in the experimentally measured data which appear to be unrealistic. We proceed to investigate this problem using density functional theory (DFT) by firstly optimizing the measured structure theoretically using  [36]. We used PAW (projector augmented wave) with PBE potential [37] for the exchange and correlation functional within the generalized gradient approximation of the DFT. To achieve high accuracy for the structure optimization, a high cut-off energy of 600 eV is used. The electronic convergence and ionic force convergence are set at of 10 -9 eV and 10 -7 eV Å −1 , respectively. A 3 × 3 × 3 k-point sampling using the standard Monkhorst scheme [38] as implemented in VASP is adopted. Test calculations by using other exchange correlation functionals such as hybrid functionals PBE0, HSE03 and also applying van der Waals correction show only very minor change in the relaxed structure.

Crystal
In table 1, we compare the experimental structure (XRD) and the final optimized structure (DFT) for K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O. It can be seen that there are some notable changes in the lattice constants with a slight increase in cell volume by 3 figure 7 in blue circles), respectively. Another pair of HBs with BLs of 1.897 Å and 1.906 Å are also decreased to 1.807 Å and 1.768 Å, respectively. These changes have significant consequence in the subsequent analysis of interatomic bonding in this crystal (see §4.4). The total energy of the relaxed crystal is 19.12 eV lower than in the unrelaxed structure, or 59.51 KJ mol −1 per atom, indicating that the relaxed structure is more stable than the unrelaxed one. All different bond lengths for all atomic pairs are listed in table 1 for direct one-to-one comparison with those O−H, O···H and O−H−O bridging bonds with significant changes mentioned above shown in bold. Figure 2 shows a two-dimensional sketch of the relaxed structure with different bond angles and bond lengths clearly marked. In the MgO 6 octahedron, Mg forms ionic bonds with four O atoms from the four PO 4 tetrahedra and two O atoms from the two water molecules. In the two PO 4 tetrahedra within the (H 2 P 2 O 7 ) 2− group, one of its O atoms (O7) is bonded with another P atom in forming the strong P−O−P bond, which has been the centre of discussion in many pyrophosphate crystals and molecules. The

Method of electronic structure and properties calculation
For the electronic structure, interatomic bonding and optical properties of the titled crystal, we use the Orthogonalized linear combination of atomic orbital (OLCAO) method [39], which is also a DFT-based       method. The VASP-optimized structure described above is used as input for the OLCAO calculation. The OLCAO method uses atomic orbitals for the basis expansion and is very efficient as well as economical, especially for large complex systems. It has been used previously to study the electronic structure and spectroscopic properties of many phosphate-related inorganic crystals and their surfaces [21][22][23][24][25][26][27]. The use of the OLCAO method with a VASP-optimized structure has resulted in many successful calculations in electronic and optical properties for large complex materials including inorganic glasses [40][41][42][43][44] and biomolecules [30][31][32].
In most OLCAO calculations, a more localized minimal basis (MB) is used for the calculation of effective charges Q * α and the bond-order (BO) values ρ αβ between a pair of atoms (α, β) using Mulliken population analysis [45,46]. and In equations (3.1) and (3.2), S iα,jβ are the overlap integrals between the ith orbital in the αth atom and the jth orbital in the βth atom; and C m jβ are the eigenvector coefficients of the mth band, jth orbital in the βth atom. We have calculated the partial charge Q or the charge transfer which is a deviation from neutral charge (Q o ) of the effective charge (Q*) on the same atom i.e. Q = Q o -Q*. The MB set used in the present study has the atomic configurations of P: 1s 2 2s 2 2p 6 3s 2 3p 3 ; Mg: 1s 2 2s 2 2p 6 3s 2 ; K: 1s 2 2s 2 2p 6 3s 2 3p 6 3d 4s 1 ; O: 1s 2 2s 2 2p 4 and H: 1s 1 . A full-basis (FB) set which has one more shell of atomic orbitals than the MB is used for calculation of self-consistent potential, the total density of states (TDOS), partial density of states (PDOS), band structure and optical properties.
The bond order (BO) ρ αβ in equation (3.2) provides a quantitative measure for the strength of the bond. The summation of all BOs normalized by the cell volume gives us the total bond order density (TBOD), which is a single metric to assess the internal cohesion in the crystal [47]. TBOD can be resolved into partial components, or the partial BO density (PBOD) from different types of atomic pairs or groups of atoms in a structural unit.
For the interband optical properties in the form of frequency-dependent complex dielectric function (hω) = ε 1 (hω) + iε 2 (hω), the imaginary part ε 2 (hω) is calculated first according to  where l and n are for the occupied and unoccupied states, respectively, ψ n (k, r) are the ab initio Bloch functions from OLCAO calculation using FB and a large k-point sampling, and f 1 (k) and f n (k) are the Fermi distribution functions. The real part ε 1 (hω) is obtained from the imaginary part ε 2 (hω) through Kramers-Kronig transformation [48]. The energy loss function, F(ω) is obtained from the imaginary part of (1/ε).
. (3.4) The elastic and mechanical properties of the pyrophosphate crystal are calculated using the stress versus strain approach to obtain the complete set of elastic coefficient C ij and compliance tensor. The fully relaxed crystal structure from VASP is used to calculate the second-order elastic tensors for the crystal using an efficient scheme [49,50] by applying a strain of ±0.5% to the cell to obtain the stress data σ j . The elastic coefficients C ij are obtained by solving the set of linear equation C ij ε j (i, j = 1, 2, 3, 4, 5, 6). (3.5) The usual mechanical properties K (bulk modulus), G (shear modulus), E (Young's modulus) and η (Poisson's ratio) are obtained based upon the Voigt-Reuss-Hill approximation for polycrystals [51,52].This approach for mechanical properties has been successfully used by us in many different inorganic crystals and glasses [44,[53][54][55][56][57]. 4. Results and discussion 4

Partial charge
The partial charges (PCs) obtained from the effective charges according to equation           figure 2 where the positions of the labelled atoms in the sequential order are marked. Figure 5 shows that there are two data points each for O1 to O8 1.124 Å that has been emphasized in §2 above in the structure change due to optimization. O6 is the atom bonded only to P2, whereas O7 is forming the bridging bond between two phosphorus atoms P1 and P2. Finally, O8 is the atom in the two water molecules but they also participate in the formation of the octahedra centred at Mg. The above discussion on the PC distribution of all atoms within the unit cell of triclinic K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O crystal provided the most detailed and quantitative analysis of the charge transfer in the crystal in relation to their geometric arrangements.

Interatomic bonding
The TBOD is the sum of total bond order in the system normalized by volume. We use the TBOD to assess the cohesion in the system. Table 3 compares the TBO and TBOD and their components, the PBOD for   [33] is shown in electronic supplementary material, table S1. It is quite obvious that the optimized structure has the largest TBO and TBOD indicating a stronger crystal cohesion.

Optical properties and refractive index
The optical properties of K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O are calculated in the form of complex dielectric function based on interband optical transitions. Figure 8a shows the calculated real (ε 1 ) and imaginary (ε 2 ) parts of frequency-dependent dielectric function in black and red lines.    Figure 8b shows the energy loss function (ELF) obtained from the complex dielectric function. The main peak in ELF is identified as the plasma frequency (ω p ) and is 22.98 eV. This is the frequency for the collective excitation of the electrons in a solid and is experimentally an easily measurable quantity. Thus, calculation of optical properties of different pyrophosphates can be effectively used as a predictor for the possible phases in an unknown sample.

Mechanical properties
As pyrophosphates may have many diversified uses, it is essential to know their mechanical properties like bulk modulus (K), shear modulus (G), Young's modulus (E) and Poisson's ratio (η). Table 4 lists the calculated elastic tensor C ij from equation (3.5). K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O has a triclinic structure of very low symmetry, so all the 22 second-order tensor elements are non-zero. The largest one is C 11 followed by C 22 and C 33, which are very close. It was noted that the off-diagonal element C 13 is actually negative but small, which could be related to the rather complex structure of K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O, with the presence of water molecules discussed in §2. The calculated K and G are 22.88 GPa and 15.63 GPa, respectively. The low K and G show K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O with a G/K ratio of 0.691, which shows the materials to be brittle and not soft. This seems to be counterintuitive to the fact that K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O contains water molecules and a larger number of weak ionic bonds. The calculated E is 38.19 GPa and the Poisson's ratio η is 0.2218. This value is smaller than similar calculations on HAP (0.262) [58] which contains the unlinked PO 4 tetrahedrons. The low Poisson's ratio η and high E for K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O also suggest that the material is less ductile and stiff in nature. This could be linked to the rigidity of the pyrophosphate units in contrast to the isolated PO 4 unit. This is the first time that the mechanical properties of a pyrophosphate crystal have been calculated and much deeper analysis of their intriguing behaviour needs further investigation. We are not aware of any experimental measurements on the mechanical properties of pyrophosphate crystals.

Connection to experimental vibrational spectrum
It is not clear at this moment if the above new findings may influence the interpretation of the vibrational properties because Raman and infrared measurements were conducted on actual samples and did not explicitly depend on the reported crystal parameters [33]. It would be of great interest to establish a    . . . . connection between the electronic structure and bonding to the discussion of the vibrational spectra, which have many interesting features at the low-frequency region. To this end, we have calculated the phonon spectrum of the K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O crystal using the phonopy package [59], with VASP for the force calculation. This requires even higher accuracy because it involves a higher vibration frequency of O−H modes and the HBs. The convergence criteria for electronic and ionic force are now set at 10 −9 eV and 10 −7 eV Å −1 , respectively. The calculated phonon spectrum and phonon DOS are shown in electronic supplementary material, figure S1a,b, respectively. The vibrational frequencies at the zone centre Γ gives the information on Raman and infrared active modes. The vibrational frequency range can be divided into three regions I (0-37.72 THz), II (38.34-57 THz) and III (97.67-105.79 THz) (electronic supplementary material, figure S1b). This is in very good agreement with the measured Raman spectrum as shown in fig. 5 of Harcharras et al. [33], both showing multiple peak structures in region I, very pronounced double peak structure in region II and H-related vibration modes in region III separated from region II with a large gap.

Conclusion
In this article, we have presented a detailed analysis of the electronic structure, interatomic bonding, optical and mechanical properties calculation for the acidic diphosphate K 2 Mg (H 2 P 2 O 7 ) 2 ·2H 2 O. This is the first such study in this pyrophosphate crystal of very low symmetry. It is concluded that the existing crystal structure data from XRD is not sufficiently refined, which results in unrealistically The most significant part of the above finding is that conventional XRD data and its structural analysis may be questionable for crystals involving light H atoms because of their weak intensity signals in the measurement. This could affect many of the published data on the structures of pyrophosphates. High-level DFT calculation can ameliorate the situation and should be part of future crystallographic analysis. The combination of an optical spectrum with large absorption at much higher energy beyond the ultraviolet region and the interesting mechanical properties unique to this low-symmetry crystal, and their possible connection to the low-frequency vibrational spectra offer new insights and opportunities for novel applications of the pyrophosphate materials.
We end this section with comments on possible extension of the present work and what specific aspects of the properties calculation can be related to real applications. Firstly, there are many similar diphosphate crystals such as K 2 [63], where the presence of transition metals could significantly affect its electronic structure and properties, leading to some unexpected applications. Secondly, the calculation can be extended to triphosphates or even to non-crystalline pyrophosphate glasses where the application in biological or medical areas are abundant. The OLCAO method is fully capable of treating such complex systems as demonstrated in many recent applications to large biomolecular systems [28,29,31,32] and inorganic glasses [41][42][43]56].