Direct and indirect proton–proton coupling in quantum-chemical theory of H-bonded materials

Journal of Molecular Structure 552 (2000) 39±44

www.elsevier.nl/locate/molstruc

Direct and indirect proton±proton coupling in quantum-chemical theory of H-bonded materials A.A. Levin*, S.P. Dolin Kurnakov Institute of General and Inorganic Chemistry, Russian Academy of Sciences, Leninsky pr 31, 117907 Moscow, Russia Received 3 September 1999; revised 16 November 1999; accepted 16 November 1999

Abstract The direct (electrostatic) and indirect mechanisms of proton±proton coupling in H-bonded solids of different nature and dimensionalities are treated. The 3D crystals of the KH2PO4 family, squaric acid (H2C4O4, 2D) as well as the 0D K3H(SO4)2-like compounds are examined as examples. We found out the role of direct and indirect mechanisms in the formation of coupling parameters (Jij) of the Ising model, commonly applied to describe the thermodynamic and dielectric properties of H-bonded materials. It was shown that the evaluated contribution to any Jij due to the direct mechanism does not exceed 30 K for all the Hbonded materials of interest. In the cases of the KH2PO4 family and squaric acid these contributions are minor in comparison with the indirect one, whereas the relative direct contribution in the Ising parameter can be more signi®cant in the case of the K3H(SO4)2-like materials. q 2000 Elsevier Science B.V. All rights reserved. Keywords: H-bonded crystals; Order±disorder ferroelectrics; Ising model; Proton±proton coupling

1. Introduction Solids composed of molecular structural units with `strong' intermolecular H-bonds are of importance from the point of view of both theory and application. This class of compounds includes, in particular, the well-known order±disorder ferroelectrics and antiferroelectrics of the KH2PO4 (KDP) family [1], antiferroelectric H2C4O4 (squaric acid, H2SQ) [2] and the K3H(SO4)2 (TKHS)-like materials with the possible structural phase transitions of the antiferroelectric type [3]. In the last decade the possibilities of application of H-bonded solids as nanomaterials have been discussed [4]. In these solids, of particular interest are those properties conditioned by ordering or reordering of * Corresponding author. E-mail address: [email protected] (A.A. Levin).

the H-bond protons and de®ned by the mechanism of proton±proton coupling. Usually, the latter is not examined in detail, and the phenomenological Ising model, with or without accounting for the effects of proton tunneling, is applied in microscopic theories of H-bonded materials. Such approaches, which are traditionally employed to describe the thermodynamic and dielectric properties of these materials, are based on the Hamiltonian [1]: H ˆ 2V

X i

s ix 2

1 2

X i; j

Jij s iz s jz

…1†

Here, s ix ; s iz are the `transverse' and `longitudinal' components of the ith pseudospin, respectively, corresponding to the proton of the ith H-bond, 2V is the frequency (in atomic units) of proton tunneling between two minima on the potential energy pro®le of a proton …V ˆ 0; if the tunneling effects are

0022-2860/00/$ - see front matter q 2000 Elsevier Science B.V. All rights reserved. PII: S 0022-286 0(00)00457-9

40

A.A. Levin, S.P. Dolin / Journal of Molecular Structure 552 (2000) 39±44 (a)

(b)

sions concerning the relative importance of these interactions are made. 2. Direct proton±proton coupling

Fig. 1. Positions of hydrogen atoms (a) and corresponding `arrows' (b) in a model squaric lattice.

neglected), and J ij denote the Ising coupling parameters describing the interaction of protons of the ith and jth H-bonds. In the frames of the static approximation …V ˆ 0† used below the pseudospin components s iz can be considered as dichotomical variables …s iz ˆ ^1†: Each of these variables corresponds to two possible equilibrium proton positions on H-bond. These variables can be represented by `arrows' or electrical dipoles directed along H-bonds toward the hydrogens, as shown in Fig. 1, where the lattice vertices correspond to non-hydrogen structural units of H-bonded material. In this case the lattice energy (Eq. (1)) can be visualized as a sum of the pair interactions of arrows (only the nearest-neighbor arrows are taken into account as a rule). In the early stages of development of the microscopic theory of H-bonded ferroelectrics it was noted [5] that the Ising parameters involve two contributions of different origins, one caused by direct electrostatic interactions of protons, the other by indirect interactions by virtue of the H-bonded molecular units. The reliable microscopic model of indirect proton±proton interactions for H-bonded ferroelectrics was formulated only recently [6±9], where KDP and its analogues as well as H2SQ were examined as examples. However, the realistic estimations of direct contributions to the Ising parameters were not obtained. In this paper the estimations of contributions in Jij parameters from direct and indirect interactions are considered and compared for several H-bonded solids of different nature and dimensionalities: KDP (3D material), H2SQ layer (quasi-2D crystal) and Rb3H(SeO4)2 (referred to as 0D compounds). Conclu-

The direct, electrostatic contributions in the Jij parameters can be evaluated in a very simple way by utilizing the results of quantum-chemical calculations of the effective charges qH on H atoms. From the diffraction data on H-bonded ferroelectrics it follows that the distance 2so between two minima of the potential energy curve of a proton along H-bond is small in comparison with the distance between different H-bond centers. Then, adopting the point charge approximation for interacting H atoms, one can estimate the maximal value (Jij)direct of the direct contribution in terms of an interaction energy in vacuum of two point dipoles with the dipole moments (qHs0). The dipoles are located in the centers of the ith and jth H-bonds and are directed along these bonds. In order to obtain the values of qH, the results of electronic structure calculations by both semiempirical (MNDO/H, AM1, PM3) and ab initio (HF/ 3-21, 6-31G) methods were exploited. The various cluster models and structural fragments of the above-mentioned crystals (containing up to 17 molecular units) were used in these quantum-chemical examinations. In accordance with semi-empirical calculations, the value of qH ranges from 0.3 to 0.4, and weakly depends on the H-bonded system and the chosen cluster under study as well as on the particular computational procedure. The value of qH obtained from ab initio calculations falls within a more wide range: 0.4±0.5 or 0.5±0.65 as evaluated by using the LoÈwdin or Mulliken procedures, respectively. At the same time the Mulliken charges are usually close to the natural-bond-orbital charges. The evaluated values of qH have a tendency to decrease with an extension of the basis set used (e.g. in passing from 6-31 to 6-31 pp). So, the mean values of qH ˆ 0:35 or 0.60 (found by semi- or non-empirical calculations, respectively) seem to be suitable for our objectives. In other respects the values of Jij for an individual material are due to the geometry of its H-bond network. The crystal framework of the KDP-like compounds (Fig. 2) consists of AO4-tetrahedra (A ˆ P, As), each of which is linked by four H-bonds with its nearest

A.A. Levin, S.P. Dolin / Journal of Molecular Structure 552 (2000) 39±44

41

Fig. 2. Crystal structure of KDP. Fig. 3. Crystal structure of H2SQ (low-temperature phase).

neighbors. In these crystals each H-bond has six nearest-neighboring H-bonds, two of which are located in the same layer and four others in the neighboring layers above and below it, respectively. Therefore, there are only two independent Ising parameters J' and Jk for this class of ferroelectrics. The use of the neutron diffraction data [10] leads to the following estimations of the maximal direct contributions: …J' †direct < 2…Jk †direct < 30…10† K (here and below the ®rst number refers to the non-empirical estimation of qH, the second one in parenthesis to the semiempirical estimation). In the layered crystal H2SQ (Fig. 3), which can be referred to as a 2D ferroelectric and 3D anti-ferroelectric system, we also have two independent Ising parameters Jtr and Jcis which describe the proton interactions in a single molecular layer (H2C4O4)1. They describe the interactions of pseudospins located in the trans- and cis-positions with respect to the molecular fragments [C4O4]. The use of the experimental geometry from the neutron diffraction data [11] gives …Jtr †direct < 2…Jcis †direct < 10…4† K: Finally, the TKHS type materials (Fig. 4) can be mentioned as well, where the anionic subsystem consists of isolated pairs of AO4-tetrahedra (A ˆ S, Se) linked together by one H-bond. For the somewhat idealized structural model of the Rb3H(SeO4)2 crystal, one has, in the nearest-neighbor approximation, u…J†direct u < 15…5† K: The neutron powder diffraction data [12] of this compound were used.

At the same time, the total values of the Ising parameters are known for the solids under consideration. They have been extracted from available experimental (thermodynamic) data by the ®tting procedure [13±15] and for H2SQ there is the estimation obtained by the Monte Carlo method as well [16]. Thus, these `experimental' data of Jij are J' < 2Jk < 350 K for KDP and < 500 K for DKDP [13]; Jtr < 550 K and 2Jcis < 350 K for H2SQ [14] (where the Monte Carlo estimations were taken into account [16]) and uJuz < 95 K for K3D12xHx(SO4)2; where z is the number of the nearest neighbors of a given H-bond [15]. One can see that for the 3D materials of the KDP family and 2D layer of H2SQ the values of (Jij)direct amounts to only a small part of the corresponding values of Jij (mainly because qH essentially differs from the bare proton charge), while for the 0D materials of the TKHS family the value of Jdirect in comparison with the `experimental' J parameter is much more substantial. Nevertheless, the above-mentioned comparison speaks for the principal role of the indirect mechanism in the formation of the Ising parameters in any case.

3. Indirect proton coupling in 3D and 2D H-bond networks The indirect mechanism of proton±proton coupling,

42

A.A. Levin, S.P. Dolin / Journal of Molecular Structure 552 (2000) 39±44

hydrogen atom by a covalent/hydrogen bond. Thus, the Slater parameter energies depend on the distribution of the O(±H) and O(´´ ´H) oxygen atoms within the A[O(±H)]n2k[O(´ ´ ´H)]k units. It is convenient to consider these molecular units as the result of the conversions A[O(--H)]n ) A[O(±H)]n2k[O(´ ´ ´H)]k. Here O(--H) denotes an oxygen atom state in the `transition structure' of crystal, where all protons are located in the centers of H-bonds (O- -H--O). Such conversions result in a shift of the position …Qfn † and the energy (DE) of the minimum point of the lower sheet of the A[O(--H)]n adiabatic potential. The use of the perturbation theory leads to [6±9] X Snm Annm =vnm …2† Qfn ˆ 4…Kn †21 n;m

DE…Qfn † ˆ S0 2 2

Fig. 4. Schematic crystal structure of TKHS.

proposed in Refs. [6±9], has a many-particle nature. It can be conveniently discussed in terms of the Slater parameters associated with the energies of manyproton con®gurations around the non-hydrogen molecular units of a crystal (the Slater con®gurations). It should be marked that the thermodynamic and dielectric properties of H-bonded solids are often described in terms of the Slater parameters [1,2], which have a simple physical meaning as well (see below). However, the nature of the Slater con®guration energies is not discussed in the frames of the conventional approaches because these quantities are ®tted usually using the corresponding experimental data of the treated material. In the quantum-chemical model [6±9] the indirect contributions in the Slater parameters are treated as the (total) energies of the non-hydrogen molecular units `inside' the many-proton Slater con®gurations. These units in the case of typical 3D and 2D H-bonded materials have the form AOn or, more exactly, A[O(±H)]n2k[O(´ ´ ´H)]k, where n ˆ 4 and 0 # k # 4 for KDP and H2SQ. Here, A is the non-oxygen `core' of molecular units (mono- or polyatomic), and O(±H)/O(´ ´´H) denotes an oxygen atom linked with a



X n;m

X n;m

S2nm =vnm 2 8

Snm Annm =vnm

X n

…Kn †21

!2 …3†

Here, Qn are the normal (or symmetry) coordinates of the transition structure unit A[(O- -H)]n; Kn the corresponding force constants; Annm the orbital vibronic constants in the MO basis set; v nm the gaps between occupied (n) and unoccupied (m) MOs. The values Snm in Eqs. (2) and (3) are the matrix elements of the so-called `substitution operator', which describes the substitutions O(±H) à O(- -H) ! O(´´ ´H) by means of the replacement ar ! ar ^ Dar ; where a r is the energy of rth orbital of O(±H) atom, participating in A±O bonding. Eqs. (2) and (3) are valid if the electrons of A-core and those of s and p A±O bonds form a closed-shell system. Therefore, Eq. (3) determines the energies of the Slater con®gurations which are explicitly expressed in terms of electronic and vibronic characteristics. Taking into account the relationships between Slater and Ising parameters [1], the latter can also be expressed in terms of electronic and vibronic characteristics of non-hydrogen units of a material. Let us consider the Slater and Ising parameters for KDP and H2SQ. For two independent Ising parameters in KDP J' and Jk we have [1] J k ˆ w=2 2 e=4 and Jk ˆ e=2 2 w=2: Analogous expressions for the Ising

A.A. Levin, S.P. Dolin / Journal of Molecular Structure 552 (2000) 39±44

parameters of the H2SQ layer are Jtr ˆ w=2 and Jcis ˆ e=4 2 w=2: Here, the quantities e and w are the Slater parameters which, of course, are different for KDP and H2SQ, although their physical meaning is similar for these materials. The Slater parameter e relates to the energy difference between two non-equivalent four-proton groups of the (2/2) type involving two covalently and two weakly bound hydrogen atoms. The parameter w is the energy increase which corresponds to transitions …2=2†0 ! …k=4 2 k†; k ˆ 1 or 3. Here (2/2)0 is the (2/2) proton group of lower energy and …k=4 2 k† denote the proton groups involving k covalently and 4 2 k weakly bound hydrogens. The use of the simple approximations allows us to derive rather visible analytical formulae of the Slater parameters of interest. For KDP these formulae have the form: e < 2…c2a l2t Das2 =v' † £‰…1 1 4A2' =K' v' † 2 …1 1 4A2k =Kk vk †…v' =vk †Š …4† w < …c2a l2t Das2 =2v' † £‰3…1 1 4A2' =K' v' † 2 2…1 1 4A2k =Kk vk †…v' =vk †Š

…5† and for the H2SQ layer have a similar form as follows: e < 2…c2a c2e Da2p =vcis £‰…1 1 4A2cis =Kcis vcis † 2 …c2b =c2e †…1 1 4A2tr =Ktr vtr †…vcis =vtr †Š

…6† w < …c2a c2e Da2p =2vcis † £‰2…1 1 4A2cis =Kcis vcis † 2 …c2b =c2e †…1 1 4A2tr =Ktr vtr †…vcis =vtr †Š

…7† Eqs. (4)±(7) were derived by using the frontier orbitals approximation. They are e, b and a p s MOs in the case of PO4-tetrahedron with the local S4 point symmetry (which are related to t2 and ap1 s MOs of regular tetrahedron) and apu ; epg ; bpu p MOs in the case of the non-hydrogen unit [C4O4]. In Eqs. (4)±(7) the quantities K, A and v with different indexes are the force constants, orbital vibronic constants and energy

43

gaps between MOs, respectively. The products ca lt Das in Eqs. (4) and (5) and ca ce Dap and ca cb Dap in Eqs. (6) and (7) are related to matrix elements Snm in the MO basis set. 4. Ising parameters for indirect mechanism All the quantities required to receive the quantitative estimations of the Slater parameters are extracted from the above-mentioned quantum-chemical calculations. The values of Da s and Da p for the O(±H) and O(´´ ´H) atoms were calculated by using the values of the intramolecular charge transfers accompanied the proton reordering between minima of their potential energy pro®les. In this way the numerical estimations of J', Jk for KDP and Jtr, Jcis for the H2SQ layer were obtained. (in fact, in the last case we applied more complex expressions in which all of the p MOs are taken into consideration). Our quantum-chemical estimations of the Ising parameters for the KDP-like materials show that J' and Jk are close to each other in magnitude and lie in the range 450±900 K [6±8]. The corresponding estimations of the Ising parameters for H2SQ layer are equal to J tr < 450±600 K and 2J cis < 100±150 K [9]. These quantum-chemical estimations of Jij must be compared with those, which have been referred to as `experimental' (and Monte Carlo), as well as with the direct contributions (Jij)direct. It is evident that the reasonable description of the indirect mechanism can be achieved for the KDP family and H2SQ layer in the frame of model, proposed in Refs. [6±9], in which the reorganization of the electronic and geometrical structures of non-hydrogen framework of a material were taken into account. Besides, the correct explanation of the geometry changes of non-hydrogen units in going from the paraelectric phase to the ferroelectric one [17], given by Eq. (2), also speaks for the proposed model. In addition, it should be noted that the proposed analytical model [6±9] permits us to elucidate the molecular mechanism of ferroelectric behavior formation for the treated and related solids. Two remarks of general character should be made. The developed model [6±9] refers to materials with continuous H-bond networks, for example, to 3D or quasi-2D systems. In the case of 0D systems (TKHS

44

A.A. Levin, S.P. Dolin / Journal of Molecular Structure 552 (2000) 39±44

and related materials) the electrostatic interactions of the (neighboring) isolated dimers [O3AO± H´ ´´OAO3] 32 play, perhaps, an important role in the Jij formation, where the charge redistribution within these dimers is caused by proton transfer [15] as well as by transfer-induced reorganization of electronic density. Of course, the proposed analytical model is not a unique way of obtaining the numerical information of interest. For instance the estimates of the Slater and Ising parameters for H-bonded materials can be extracted directly from quantum-chemical calculations of the total energies of their suitable structural fragments [18] without the use of the explicit formulae like Eqs. (4)±(7). The Monte Carlo method [16], as well as the modern techniques, such as the ab initio molecular dynamics simulation [19] and the path integral Car±Parrinello molecular dynamics [20±23], are also applicable for H-bonded solids. 5. Conclusion Thus, the quantum-chemical comparative study of the direct and indirect mechanisms of proton±proton coupling in H-bonded solids demonstrates an important and sometimes a principal role of the latter mechanism. According to our ab initio calculations the relative contributions of the indirect mechanism in Jij formation are not less than 50±70% for the TKHS-like materials and more than 90% for the KDP family. For the treated solids of higher dimensionality (the KDP family, H2SQ) the indirect mechanism of proton±proton coupling is reasonably described by the supposed quantum-chemical model. Acknowledgements This work is supported by the Russian Fundamental Science Foundation, project 99-03-33234.

References [1] M.E. Lines, A.M. Glass, Principles and Application of Ferroelectrics and Related Materials, Clarendon Press, Oxford, 1977. [2] R. Blinc, in: H. Ratajczak, W.J. Orville-Thomas (Eds.), Molecular Interactions, vol. 2, Wiley, Chichester, 1981 (chap. 5). [3] M. Ichikawa, T. Gustafsson, L. Olovsson, Ferroelectrics 159 (1994) 73. [4] V.A. Russell, M.D. Ward, Chem. Mater. 8 (1996) 1654. [5] P.G. de Gennes, Solid State Commun. 1 (6) (1963) 132. [6] A.A. Levin, S.P. Dolin, Dokl. Phys. Chem. 341 (1995) 99. [7] A.A. Levin, S.P. Dolin, Proc. Est. Acad. Sci. Phys. Math. 44 (2/3) (1995) 144. [8] A.A. Levin, S.P. Dolin, J. Phys. Chem. 100 (1996) 6258. [9] A.A. Levin, S.P. Dolin, V.L. Lebedev, Mikhailova T.Yu, Chem. Phys. Reps. (2000) (in press). [10] R.J. Nelmes, Z. Tun, W.E. Kuhs, Ferroelectrics 71 (1987) 125. [11] M.I. McMahon, R.J. Nelmes, W.F. Kuhs, D. Semmingsen, Z. Kristallogr. 195 (1991) 231. [12] R. Melzer, R. Sonntag, K.S. Knight, Acta Crystallogr. C 52 (1996) 1061. [13] V.G. Vaks, Introduction to Microscopic Theory of Ferroelectrics, Nauka, Moscow, 1973 (in Russian). [14] V.G. Vaks, V.I. Zinenko, V.E. Schneider, Usp. Fiz. Nauk. 141 (1983) 629 (in Russian). [15] Y. Moritomo, Y. Tokura, N. Nagaosa, T. Suzuki, K. Kumagai, Phys. Rev. Lett. 71 (1993) 2833. [16] H.E. Muser, Ferroelectrics 55 (1984) 275. [17] A.A. Levin, S.P. Dolin, Russian J. Coord. Chem. 24 (1998) 270. [18] S.P. Dolin, V.L. Lebedev, A.A. Levin, Dokl. Phys. Chem. 341 (1995) 109. [19] R. Car, M. Parrinello, Phys. Rev. Lett. 55 (1985) 2471. [20] M.E. Tuckerman, D. Marx, M.L. Klein, M. Parrinello, J. Chem. Phys. 104 (1996) 5579. [21] M.E. Tuckerman, D. Marx, M.L. Klein, M. Parrinello, Science 275 (1997) 817. [22] S. Miura, M.E. Tuckerman, M.L. Klein, J. Chem. Phys. 109 (1998) 5290. [23] M.E. Tuckerman, M.L. Klein, Chem. Phys. Lett. 283 (1998) 147.