340 Acta Chim. Slov. 2014, 61, 340-349 Scientific paper Molecular Modeling of Organic Corrosion Inhibitors: Why Bare Metal Cations are Not Appropriate Models of Oxidized Metal Surfaces and Solvated Metal Cations Anton Kokalj* Department of Physical and Organic Chemistry, Jozef Stefan Institute, Jamova 39, SI-1000 Ljubljana, Slovenia * Corresponding author: E-mail: tone.kokalj@ijs.si, URL: http://www.ijs.si/ijsw/K3-en/Kokalj Tel: +386-1-477-35-23; Fax: +386-1-477-38-22, Received: 19-11-2013 Paper based on a presentation at the 4th RSE-SEE 2013 Symposium on Electrochemistry in Ljubljana, Slovenia Abstract The applicability of various models of oxidized metal surfaces - bare metal cations, clusters of various size, and extended (periodic) slabs - that are used in the field of quantum-chemical modeling of corrosion inhibitors is examined and discussed. As representative model systems imidazole inhibitor, MgO surface, and solvated Mg2+ ion are considered by means of density-functional-theory calculations. Although the results of cluster models are prone to cluster size and shape effects, the clusters of moderate size seem useful at least for qualitative purposes. In contrast, the bare metal cations are useless not only as models of oxidized surfaces but also as models of solvated cations, because they bind molecules several times stronger than the more appropriate models. In particular, bare Mg2+ binds imidazole by 5.9 eV, while the slab model of Mg0(001) by only 0.35 eV. Such binding is even stronger for 3+ cations, e.g., bare Al3+ binds imidazole by 17.9 eV. The reasons for these fantastically strong binding energies are discussed and it is shown that the strong bonding is predominantly due to electron charge transfer from molecule to metal cation, which stems from differences between molecular and metal ionization potentials. Keywords: corrosion inhibition, oxidized metal surfaces, molecular modeling, DFT, adsorption 1. Introduction Quantum chemical calculations of corrosion inhibitors have become very fashionable in corrosion inhibition studies (e.g., see the Review 1). The purpose of such calculations is aimed at either better understanding of atomic-scale details, interpretation of experimental findings, or screening of new corrosion inhibitors utilizing kind of QSAR/QSPR (quantitative-structure-acti-vity-relationship/quantitative-structure-property-rela-tionship) approach. In large majority of cases such screening consists of calculating several electronic parameters of inhibitor molecules - either in gas- or in aqueous-phase with the solvent treated implicitly by continuum models - and associating them to the experimentally determined corrosion inhibition efficiency via some correlation analysis. Such approach was critically examined in Ref. 2, where it was pointed out that, in general, molecular electronic properties cannot be directly related to corrosion inhibition efficiency, because the actual rela- tion is far more involved; anticipations based solely on molecular properties become especially unreliable when bond-breaking and bond-making take place, for example, during the adsorption of corrosion inhibitor molecu-les.3'4 To go beyond this approach several researchers attempted to explicitly model the inhibitor-surface interaction, where the materials were either represented by metal ion,5-7 atom,7-13 cluster,14-19 or slab3 4 20-37 models. Among these models, the bare metal ion models stand out, because fantastically strong binding energies (Eb) were reported for interaction between organic molecules and metal ions with the Eb magnitudes reaching a few tens of eV.56 38 These magnitudes are indeed incredibly large if compared to the bond strength of the N=N triple bond (9.8 eV),39 which is among the strongest known molecular chemical bonds.40 Several technologically important metals are of interest in corrosion inhibition studies, e.g., iron, copper, zinc, aluminum, magnesium, and their alloys. In addition to metallic state (i.e., formal oxidation state of zero), metals Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: Acta Chim. Slov. 2014, 61, 340-349 341 in oxidized states are also relevant in this context, e.g., Fe2+, Fe3+, Cu+, Cu2+, Zn2+, Al3+, and Mg2+, where the superscripts indicate the formal oxidation states (not to be confused with the actual charge as deduced from, say, quantum-mechanical calculations; see the Results section). The listed examples of metal cations hence span the oxidation states from 1+ to 3+, yet for the sake of definite-ness the presentation will focus on oxidation state of 2+ exemplified for Magnesium (one example of Al3+ will be also presented). The Mg was chosen, because it is a simple s-type metal thus having much simpler electronic structure than transition metals. As a model of oxidized surface the Mg0(001) is considered, even though the bulk Mg(0H)2 is more stable than the bulk MgO.41 The reason for this choice is that bulk Mg(0H)2 is a layered material consisting of HO-Mg-OH five-layers, where the Mg ions extend along the (001) plane with the perpendicular OH groups fully decorating them from below and above;42,43 the OH groups therefore completely cover the layer of Mg ions and hence sterically prevent the interaction of foreign molecular species with them. Figure 1. Ball-and-stick model of imidazole and its skeleton structure with pyrrol- and pyridine-type N atoms indicated explicitly (top). The coloring of atoms as used in this work is also shown (bottom). Imidazole will be used as a representative organic corrosion inhibitor. Its derivatives - e.g., benzimidazole and merkapto-benzimidazole - are known as efficient inhibitors for versatile metallic materials,44-46 but because the main point of this paper is conceptual, a simpler imi-dazole was chosen instead. Its structure is shown in Figure 1. The purpose of this publication is to examine and discuss the usability of various models that are used to represent oxidized metal surfaces in the field of corrosion inhibitors and to point out why the isolated metal ions are inadequate models of oxidized metals; the origin of the above mentioned fantastically strong binding energies will be also explained. The solvated metal cations and their interaction with corrosion inhibitor molecules are also of interest in corrosion inhibition context, because it is known that some organic inhibitors form complexes with metal cations, benzotriazole being a notable example.47 For this reason, the modeling of interaction between corrosion inhibitors and solvated metal cations is also addressed. 2. Technical Details Calculations were performed within the framework of Density-Functional-Theory (DFT) using the hybrid B3LYP functional48 and the local Gaussian-type-orbital basis set. Molecular calculations in the gas- and aqueous-phases were performed using the Gaussian program,49 while calculations with periodic-boundary-conditions were performed using the Crystal program.50 Molecular graphics were produced by Xcrysden graphical package.51 Unless explicitly stated otherwise, a 6-311++G(d,p) basis set was used for molecular calculations, which is a triple-zeta basis set augmented with polarization and diffuse functions. This basis set contains too diffuse functions for periodic calculations, hence a smaller basis set was used correspondingly (see below). structures of molecular and smaller cluster systems were fully relaxed, whereas for larger cluster models only the imidazole molecule and the atoms close to the adsorption site were relaxed (more details are specified below). some molecular calculations were also performed in aqueous-phase with the solvent described implicitly by the SMD continuum solvation model52 using the integral equation formalism.53,54 Chemical species in gas- and aqueous-phase are labeled as X(g) and X(aq), while the physical properties corresponding to these phases are labeled as Y(g) and Y(aq), respectively (currently, the label (aq) implies the implicit description of aqueous-phase by the continuum sMD model). For brevity reasons the label (g) will be often omitted, yet the label (aq) will be always stated explicitly. Binding energy (Eb) between two fragments that constitute a given system in the gas-phase is calculated as: Eb = ea-b ea eb, (1) where EA B is the total energy of the A-B system, while EA and Eb are the total energies of the isolated A and B fragments, respectively (for example, A = imidazole and B = Mg0(001) surface). The Eb will be also designated as E® to explicitly indicate its reference to gas-phase. All the reported binding energies are corrected for basis-set superposition-error (BSSE) using the Boys-Bernardi counterpoise correction.55 For the 6-311++G(d,p) calculations a typical BSSE correction of Eb is about 0.05 eV, whereas for larger clusters and slabs, where smaller basis sets were used, it is about 0.2 eV. The following models representing the Mg in 2+ oxidation state - labeled also as Mg(II) - will be considered: Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: 342 Acta Chim. Slov. 2014, 61, 340-349 II V Isolated Mg2+ ion: it can be seen either as a simple model of an oxidized magnesium surface (e.g., MgO or Mg(OH)2) or as solvated Mg2+ ion. It will be shown that isolated Mg2+ ion is a very poor model for both of them. [Mg(H2O)n]2+ complexes, n e [1, 6]: they can be seen as simple models of solvated Mg2+ ion. It is known that metallic cations are very strongly solvated in water, i.e., they are chemically bonded with the first shell of water molecules and form complexes,56 e.g., [Mg(H2O)6]2+,57 III IV [Al(H2O)6] , and [Fe(H2O)6]3+.59 The formation energy of the first-shell complexes represents a significant part of the whole solvation energy.60 Five-atom Mg(OH)2 cluster: it can be seen as a very small cluster representation of either the Mg(OH)2 or the MgO surfaces, where for the latter the cluster's oxygen dangling bonds are saturated with hydrogens (such saturation is often used in cluster model representation of surfaces, albeit mainly for covalent ma-terials).61 MgO(001)[9,i], MgO(001)[99], and MgO^)^^ clusters: these are two- and three-layer cluster model representations of MgO(001) surface, where the subscripts indicate the number of atoms in each layer. These clusters are cut from the optimized structure of bulk MgO as calculated by the Crystal program (see below). For the smaller [9,1] and [9,9] clusters only the central Mg ion in the first layer was relaxed and other ions were fixed, whereas for the larger [21,21,9] cluster six atoms were allowed to relax, i.e., the central Mg ion in the first-layer and its five nearest neighbor O atoms. These six atoms were described with 6-311++G(d) basis set (for all the clusters), whereas other atoms were described by a smaller 6-31G(d) basis set. MgO(001) slab: this is a periodic slab model of MgO(001) surface. The slab consists of four (001) layers. The in-plane lattice spacing was fixed at calculated equilibrium MgO bulk lattice parameter of 4.23 A, which is in good agreement with the experimental value of 4.21 A.39 The bottom layer of the slab was kept frozen to the bulk positions, while all other degrees of freedom were relaxed. Adsorption calculations were modeled with (3 x 3), (4 x 4), and (5 x 5) supercells; the Eb was then extrapolated to zero coverage, because the lateral dipole-dipole interactions between adsorbed imidazoles are very long ranged due to imidazole's large dipole moment of 3.8 D.27,62 These calculations were performed with the following basis sets: 311G(p) of Bredow et al.63 for H; 6-31G(d) of Gatti et al.64 for C and N; 8-411G(d) of Bredow et al.65 for O; and 8-511G(d) of Valenzano et al.66 for Mg. 3. Results and Discussion Optimized structures of imidazole molecule interacting via its pyridine-type N atom with various models of Mg(II) are shown in Figure 2. The corresponding binding energies (Eb) and N-Mg bond-lengths are also stated. The following two observations are evident: (i) the charged models of Mg(II) bind imidazole considerably stronger than the neutral models. The strong binding is in particularly exaggerated for the bare Mg2+ ion model. (ii) The Eb of the neutral models of Mg(II) does not converge smoothly to the value given by the extended slab model as the cluster size increases, but shows oscillations instead, i.e., the smallest five-atom Mg(OH)2 cluster binds imidazole considerably stronger, while the [9,1] and [9,9] cluster-models of Mg0(001) bind it weaker than the extended slab model, yet the largest three-layer [21,21,9] cluster-model binds again stronger. Also the N-Mg distance displays analogous oscillations versus cluster-size, because it correlates well with the Eb, i.e., the stronger is the Eb the shorter is the distance. Such oscillations of binding (adsorption) energy with cluster size (and shape) are well known and are perhaps the most pronounced for ionic and metallic surfaces.67-71 Several methods have been devised in nineties to diminish them, e.g., bond-preparation rule,72 the method of Russier-Mijoule,69 multipole-expansion embedding -which is particularly suitable for ionic surfaces, because it is able to properly account for the Madelung poten-tial70,73 - and embedding schemes based on Green function formalism.67,68,74 Figure 2. B3LYP optimized structures of imidazole interacting via its pyridine-type N atom with various models of Mg(II). The binding energies (Eb) and the imidazole-Mg distances (i.e., N-Mg bond lengths) are also stated. I Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: ... Acta Chim. Slov. 2014, 61, 340-349_ 343 Despite the cluster-size oscillations the clusters of moderate size, treated without any special preparation, seem useful at least for qualitative purposes,75 while the usability of very small clusters appears questionable. In contrast, the bare metal cations (Mg2+ in the current case) are useless as models of oxidized surfaces. Contrary to the finite-size cluster models, the periodic slab model circumvents the problem of the cluster-size oscillations, while its adsorption and other properties converge rapidly with the slab thickness.76,77 A convenient feature of the slab model is that the coverage dependence can be easily accounted for by either changing the size of the surface supercell or by changing the number of adsorbed molecules in a sufficiently large supercell. As for disadvantages, let us mention that a special care must be taken when treating charged systems, because long range Coulomb interactions dictate the use of charge neutrality in periodic calculations. It is also not possible to directly perform calculations in the limit of zero coverage and calculations at very low coverage are computationally expensive. This is relevant for modeling the adsorption of polar molecules with large dipole moment, because the resulting lateral dipole-dipole interactions can be very long ranged.25,27,62 In the current case, the calculated adsorption energies of imidazole on Mg0(001) are -0.26, -0.31, and -0.33 eV at 1/9, 1/16, and 1/25 monolayer coverage, res- pectively, while the value extrapolated to zero coverage is 0.35 eV; the coverage is defined as number of imidazole molecules per surface Mg ion. Having roughly established the usability of various models, it is of some interest to compare the computational costs of these models. Such considerations can be helpful when planning extensive studies on a number of different corrosion inhibitors on different materials. Figure 3 displays the actual time (wall time) taken by a modern desktop computer to perform a converged single-point SCF (self-consistent-field) calculation of systems shown in Figure 2; the corresponding number of electrons of these systems is also plotted. It can be seen that the wall time of small cluster models is on the order of 1 minute, for larger [21,21,9] cluster model it is on the order of 1 hour, while for slab models it is in the range from several hours to 1 day, depending on the size of the slab. 3. 1. Bare Mg2+ as a Model of Oxidized Mg Surface Why the bare Mg2+ binds so much stronger than other models? The reason is mainly due to instability of isolated Mg2+ ion. The first ionization potential (IP) of Mg is 7.7 eV, while the second IP is 15.0 eV.39 Hence, it costs 22.7 eV to create Mg2+ ion from Mg, i.e.: Mg ^ Mg2+ + 2e~ AE = IP1St + IP2nd : 22.7 eV (2) where the labels IP1st and IP2nd stand for first and second IP, respectively. The B3LYP calculated IP1st and IP2nd are 7.7 and 15.5 eV, respectively, resulting in the AE value (cf. Eq (2)) of 23.2 eV. Calculated values are thus in fair agreement with experimental values. The AE value of 23 eV implies that isolated Mg2+ ion is highly electronegative. It will thus seize some fraction of electrons from a mo- Figure 3. Top panel: wall clock times (on 16 CPU cores of dual CPU AMD 0pteron(tm) 6128) of converged single-point SCF calculations for systems shown in Figure 2. Bottom panel: number of electrons in the considered systems. 0rdinate axes are in logarithmic scale with base 4. Figure 4. Energies of isolated Mg (left, red) and imidazole (right, blue) as functions of their charge; the AE = 0 corresponds to neutral species. The energy change of each species due to transfer of x = 0.62 electrons from imidazole to Mg2+, reaction (3a), is schematically indicated by arrows. Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: 344 Acta Chim. Slov. 2014, 61, 340-349 lecular species it interacts with; note that the IP1st + IP2nd value of 23 eV is in between the IP1st of He (24.6 eV) and Ne (21.6 eV), which have the largest IP among all the elements. To foresee that the strong binding between molecule and Mg2+ ion is mainly due to electron charge transfer from molecule to Mg2+ ion, let us decompose the formation of the [imidazole-Mg]2+ complex into two conceptual steps: (i) transfer of electron charge from molecule to Mg2+ (keeping the two species completely separated) and (ii) formation of the actual bond, i.e.: charge transfer: imidazole + Mg2+ ^ imidazole^ + Mg(2-x)+, (3a) bond formation: imidazole^ + Mg(2-x)+ ^ [imidazole-Mg]2+, (3b) o.o >—1 0 1 2 3 15 3D 45 d- o^ (A) d- dv (A) Figure 5. Energy of [imidazole-Mg]2+ system in singlet (red) and triplet (blue) states versus the N-Mg distance; d0 stands for equilibrium N-Mg distance of [imidazole-Mg]2+ complex, while EM is the energy of well separated (d ^ imidazole+ and Mg+ cations. The curve labeled 1/d represents the Coulomb electrostatic repulsion between two unit charges (see text); d corresponds to the distance between Mg ion and barycenter of imidazole. where the x+ and (2-x)+ are the charges of imidazole and Mg, respectively, at relaxed geometry of [imidazo-le-Mg]2+ complex. The Hirshfeld population analysis78 was utilized to estimate the respective charges and it gives the value of x = 0.62. The binding energy between imida-zole and Mg2+ can be correspondingly written as: F = F + F ^b et ^ ^bf' (4) where subscripts ct and bf stand for charge-transfer and bond-formation, respectively, while the Ect and Ebf indicate the energy change of reactions (3a) and (3b), respectively. Because for isolated systems the total energy is a piecewise linear function of the electron number (N) having the dE/dN discontinuities at integer values of N,79 the energy change due to electron charge transfer of reaction (3a) can be written as: Fet = ^(IPimidazole ^ ) f°T X £ [0,1]. (5) This concept is schematically illustrated in Figure 4. Because the IP2nd of Mg is considerably larger than the IP1st of imidazole, the energy gain due to the electron charge transfer is substantial; the B3LYP calculated adia-batic IP1st of imidazole is 8.77 eV, hence Ect = 0.62(8.77 -15.46) eV = -4.15 eV. This implies that the Ect contributes 70% to the binding energy between imidazole and Mg2+, while further 30% (corresponding to -1.75 eV) comes from the subsequent bond formation, Eq (3b). It should be noted that situation is radically different when Mg(II) ion is embedded in either crystal lattice or water solvent, because there it is stabilized considerably by the surrounding species (i.e., by oxygen counter-ions in the case of bulk MgO and by water molecules in the case of solvated Mg2+). Consequently, the actual charge of Mg(II) in such environments is not +2, but can be significantly smaller; e.g., the calculated charge of Mg(II) is +1.75 and +0.5 in bulk MgO and [Mg(H2O)6]2+ complex, respectively. Let us point out an important issue with respect to stability of [imidazole-Mg]2+ complex. According to Eq (5), the Ect contribution is optimal at x = 1 (this corresponds to Mg+ and imidazole+ cations), where Ect = IP1midazole - IPMngd = -6.7 eV.* This value is larger in magnitude than the Eb of [imidazole-Mg]2+ complex, hence the complex is meta-stable because the well separated Mg+ and imidazole+ cations are by 0.8 eV more stable. This is shown graphically in Figure 5, which plots the energy as a function of the N-Mg distance (d) for [imida-zole-Mg]2+ system. A minimum can be seen around the equilibrium N-Mg distance (d0) for the singlet state (red curve). However, for distances d - d0 > 1 A the singlet state becomes less stable than the triplet state (blue curve); this latter state leads to standalone Mg+ and imidazo-le+ cations as d ^ ^ (note that both cations have a single unpaired electron). Indeed, at large distances, say d - d0 > 10 A, the triplet state closely follows the c ■ 1/d curve, which represents Coulomb electrostatic interaction between two unit charges, where d corresponds to the distance between the Mg ion and the barycenter of imidat-zole molecule, and c is a conversion factor such that d can be specified in A (c = a0Eh = 14.39 eV A, where a0 is Bohr radius, 0.529 A, and Eh is Hartree atomic unit of energy, 27.21 eV). r Although the 1/d approximation is not valid at short N-Mg distances as can be apprehended from Figure 5, it is conceptually useful to understand why for the [imidazole-Mg]2+ complex the x <1 (cf. Eq (3)); note that the x = 1 is favored by charge ttransfer contribution as discussed above. According to 1/d approximation, the Coulomb repulsion between the two fragments of [imi- * Note that Eq (5) is valid only for x e [0, 1], but because the IP2nd of imidazole (16.1 eV) is larger than the IP1st of Mg (7.7 eV) there will be no charge transfer beyond the x = 1. Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: ... Acta Chim. Slov. 2014, 61, 340-349 345 dazole-Mg]2+ complex is proportional to x(2-x)/d. Hence, for x e [0, 1] it is minimal at x = 0 and maximal at x = 1, thus showing the opposite trend to charge transfer contribution. The actual value of x is therefore a compromise (Hirshfeld population analysis gives x = 0.62 as mentioned above). The bottom line of all these arguments is that the bare Mg2+ is rather useless as a model of oxidized metal surface. To further support this claim, let us mention that the nominal binding between molecule and bare metal cations is even far larger in magnitude for 3+ cations. Take the Al3+ as an example. The B3LYP calculated binding energy between imidazole and Al3+ within the [imidazole-Al]3+ complex is the incredible -17.9 eV. Again, this is largely due to charge transfer from imida-zole to Al3+, yet despite the enormous magnitude of Eb the [imidazole-Al]3+ complex is nevertheless meta-stable (for similar reason as the [imidazole-Mg]2+ complex), because the well separated imidazole2+ and Al+ cations are by 5.2 eV more stable; namely, the first three IPs of Al are 6.0, 18.9, and 29.0 eV, while the first two adiabatic IPs of imidazole are 8.7 eV and 16.1 eV (B3LYP values). The IP2nd of Al is thus larger than the IP2nd of imidazole, which implies that the most convenient charge redistribution is the following: imidazole + Al3+ ^ imidazole2+ + Al+. The corresponding change of energy equals ff1^ + ^azol - IPAld - IPA? = -23.1 eV, which is by 5.2 eV more negative than the binding energy between imidazole and Al(III) within the [imidazole-Al]3+ complex. 3. 2. Bare Mg2+ and [Mg(H2O)„]2+ as Models of Solvated Ion Bare metal cations are thus useless as models of oxidized metal surfaces, but are they at least useful as models of solvated metal ions? The answer is again no and the ar- guments go along the lines discussed above. Namely, bare metal cations are so reactive that they will chemically interact with whatever molecular species they encounter, the water molecules being no exception. The B3LYP calculated adiabatic IP1st of water molecule is 12.7 eV, hence being smaller than the IP2nd of Mg. This implies a non-negligible charge transfer and a strong chemical interaction between Mg2+ and water molecules; e.g., the calculated Eb between a single water molecule and Mg2+ ion is -3.5 eV. Indeed, it is known that metallic cations are so strongly solvated in water as to be chemically bonded with the first shell of water molecules thus forming com-plexes,56,80 in particular case the [Mg(H2O)6]2+.57 The interplay between the chemical bonding and charge deloca-lization within the [Mg(H2O)n]2+ complexes is shown in Figure 6, which plots the Hirshfeld charge of Mg(II) and cumulative binding energy between water molecules and Mg2+ ion as functions of the number of coordinated water molecules; the cumulative binding energy is calculated as: E cumulative _ e _ E i _ nE ^b _ complex Mg2+ 't£,H2O' (6) where E^^, EMg2+, and £h2o are total energies of [Mg(H2O)n]2+ complex, Mg2+ ion, and water molecule, respectively. It can be seen that the magnitude of cumulative binding energy increases and the Hirshfeld charge of Mg(II) decreases as the number of coordinated water molecules increases. The Mg2+ is strongly bound by 13.6 eV within the [Mg(H2O)6]2+ complex, while the positive charge is delocalized over the whole complex: according to Hirshfeld population analysis the charge of Mg(II) is about +0.5, while the charge on each water molecule is +0.25. The [Mg(H2O)6]2+ complex would be therefore much better model of solvated Mg2+ than the bare Mg2+. Let us now calculate the interaction of imidazole with sol-vated Mg2+ ion using the [Mg(H2O)6]2+ model. Within this Figure 6. Left: B3LYP calculated cumulative binding energy between water molecules and Mg ion (red curve) and the Hirshfeld charge of Mg(II) ion (blue curve) versus the number of water molecules for [Mg(H2O)n]2+ complexes, n e [0, 6]. Right: optimized structures of [Mg(H2O)n]2+ complexes. Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: 346 Acta Chim. Slov. 2014, 61, 340-349 complex the Mg(II) is highly coordinated with water molecules, hence imidazole will have to substitute one water molecule resulting in [imidazole-Mg(H2O)5]2+ complex A relevant question at this point is how well the purely implicit solvent model describes the solvated Mg2+ ion in comparison with the hybrid implicit-explicit (shown in Figure 2). The binding energy between the imi- [Mg(H2O)l2+ model. It is known that adding the explicit dazole and [Mg(H2O)5] is -2.1 eV. This value is considerably smaller in magnitude than the Eb of the plain [imi-dazole-Mg]2+, yet it does not take into account the work needed to substitute one water molecule. The corresponding net energy of the interaction can be therefore estimated by considering the substitution reaction: imidazole + [Mg(H2O)6]2+ ^ ^ [imidazole-Mg(H2O)5]2+ + H2O, AE = AE = [E - E ] sr L products reactantsJ' (7) where Eproducts (Ereactants) corresponds to sum of tOfel energies of product (reactant) species; the energy change of reaction is labeled as AEsr, where the subscript sr stands for substitution-reaction. The calculated AEsr is -0.9 eV. Although, this value is considerably smaller in magnitude than the binding energy between imidazole and [Mg(H2O)5]2+, its exothermicity is nevertheless still significantly overestimated (see below). Further improvement can be achieved by treating all the four species involved in reaction (7) as solvated in water, i.e.: imidazole^ + [Mg^O)/^ ^ [imidazole-Mg(H2O)5]2 (aq^ ^(aq)' + HO,. (8) where the solvation effects are estimated implicitly by continuum description of solvent (the SMD model is currently utilized); the subscript (aq) hence stands for implicit description of aqueous-phase. The so calculated AEsr of reaction (8) is only -0.08 eV. The magnitude of this value is too small to claim that imidazole would form complex with Mg2+ ion in aqueous solution, taking into account the uncertainties due to approximations of the computational model and the neglect of entropic effects that may relatively disfavor the reaction (8), i.e., the entropy likely reduces during this reaction. first-shell water molecules considerably improves the accuracy of the calculated ion solvation free energies,81 but in the current context the relevant quantity is the net interaction between imidazole and solvated Mg2+ ion. Within the purely implicit solvent model, the reaction (8) reads: imidazole(aq) + Mg2(a+q) [imidazole- Mg]2+ (aq) (9) Although this reaction does not look as substitution reaction, it should be noted that the substitution is implicitly modeled, i.e., both imidazole and Mg2+ have to partially loose its solvation shell before they can interact with each other. The resulting calculated value of this reaction is -1.2 eV, which is 30% larger in magnitude than the AEsr of reaction (7) and an order of magnitude more exothermic that the value of -0.08 eV given by the hybrid implicit-explicit [Mg(H2O)J2+ model, reaction (8). The reason for the huge difference between purely implicit and hybrid implicit-explicit models is due to strong interaction of metal cation with first-shell water molecules (cf. Figure 6), i.e., continuum solvent models cannot adequately describe chemical bonds that form between the Mg(II) cation and first-shell water molecules. 3. 3. Implicit Solvent Description of Solvated Cluster Models of Oxidized Mg(II) Surfaces In contrast to bare Mg2+, other models of oxidized Mg(II) surfaces - i.e., cluster and slab models of MgO -interact considerably weaker with imidazole (cf. Figure 2) and consequently also with water; the calculated Eb between a single water molecule and the MgO(001)[21219] cluster is -0.4 eV (to be compared to -3.5 eV for Mg2+). Consequently, the pure implicit description of solvent should be more adequate in this case. For this reason, the Table 1. B3LYP calculated interaction energies (E® and Einaq)) and N-Mg distances (d®Mg and dNq:Mg) of imidazole interacting with various models of Mg(II) in gas-phase and aqueous-phase with the solvent treated implicitly by the SMD model. b int N-Mg N-Mg (eV) (eV) (A) (A) [imidazole-Mg]2+ -5.9 -1.2 -0.3 (-0.08^) 1.97 2.04 [imidazole-Mg(H2O)5]2+ -2.1 (-0.9a) 2.14 2.18 imidazole-Mg(OH)2 -1.3 -0.7 2.11 2.11 imidazole-MgO(0 0 1)[91] -0.1 does not bind 2.50 / imidazole-MgO(0 0 1)[9 9] -0.25 +0.02 2.32 2.42 imidazole-MgO(0 0 1)[21219] -0.8 -0.6 2.23 2.39 imidazole-MgO(0 0 1)slab ' -0.35 / 2.24 / "Calculated according to reaction (7). 'Calculated according to reaction (8). Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: ... Acta Chim. Slov. 2014, 61, 340-349 347 interaction of imidazole with cluster models of oxidized Mg surface was also modeled in the implicitly described aqueous-phase using the continuum SMD solvent model. The corresponding net interaction energies (E(ianqt)) were calculated according to the reaction: imidazole(aq) + Mg(II), ^ imidazole-Mg(II) hence: (aq) '(aq)' P(aq) _ E(ag) — E(ag) — E(ag) int imidazole- Mg(II) ^imidazole Mg(II)' (10a) (10b) where Mg(II) stands for the given model of oxidized Mg surface, while various E(aq) terms are the total energies of corresponding species as obtained from the implicit solvent calculations. Although the E(ianqt) is the aqueous-phase analogue to gas-phase Eg), it is subscripted by "int" to indicate that it measures the net energy of interaction and not the gross binding energy, i.e., both fragments have to partially loose its solvation shell before they can interact and the corresponding energy cost is included in the E(ianqt). Table 1 tabulates the imidazole-Mg(II) interaction energies and N-Mg distances as calculated in gas- and aqueous-phases; for the extended slab model only the gasphase calculations were performed. The comparison between the gas- and aqueous-phase values reveals that (i) E(ianqt) displays smaller magnitudes than Eb and (ii) the N-Mg distances are elongated as passing from gas- to aqueous-phase. The reason is that the adsorption at solid/water interface is influenced by competitive interplay between molecule-surface, molecule-water, and water-surface interactions.424 4. Conclusions In this paper various models of oxidized metal surfaces that are used in the field of corrosion inhibitors were examined and discussed. It was shown that it matters how the model of the surface is chosen, because not all the models provide adequate description. The usability of very small cluster models appears questionable, while the clusters of moderate size, although susceptible to cluster-size effects, seem useful at least for qualitative purposes. In contrast, the bare metal cations - in particular the 2+ and 3+ cations - are useless as models of oxidized surfaces. It was demonstrated that bare Mg2+ and Al3+ cations give exceedingly strong bonding with imidazole molecule, which is predominantly due to electron charge transfer contribution that stems from differences between molecular and metal ionization potentials; the situation is radically different when metal cations are embedded in either crystal lattice or water solvent, because there they are stabilized considerably by surrounding species. Furthermore, the [molecule-Mg]2+ and [molecule-Al]3+ complexes have little physical significance because they are less stable than the well separated cations: imidazole+ /Mg+ and imi-dazole2+/Al+. In addition to these unrealistic binding energies also note that metal cation is a very small spherical (or roughly so) object, while real surfaces are extended planar objects (or nearly so if atomic-scale surface defects are considered). A metal cation or a very small cluster model may therefore not provide adequate steric environment in particular for the parallel adsorption modes. 5. Acknowledgments This work has been supported by the Slovenian Research Agency (Grant No. P2-0148). The author thanks Prof. Ingrid Milošev for bringing the subject of corrosion inhibitors to his attention and Dunja Peca for her careful reading of the manuscript. 6. References 1. G. Gece, Corros. Sci., 2008, 50, 2981-2992. 2. A. Kokalj, Electrochim. Acta, 2010, 56, 745-755. 3. A. Kokalj, Corros. Sci., 2013, 68, 195-203. 4. N. Kovacevic, A. Kokalj, Corros. Sci., 2013, 73, 7-17. 5. M. C. S. S. Macedo, O. E. Barcia, E. C. Da Silva, J. D. O. Mendes, O. R. Mattos, J. Electrochem. Soc., 2012, 159, C160-C169. 6. D. Seifzadeh, S. Hamzedoust-Hasankiadeh, A. Shamkhali, Prot. Met. Phys. Chem. Surf., 2013, 49, 229-239. 7. C. Li, L. Li, C. Wang, Y. Zhu, W. Zhang, Corros. Sci., 2014, 80, 511-516. 8. G. Lendvay-Gyorik, G. Meszaros, B. Lengyel, G. Lendvay, Corros. Sci., 2003, 45, 1685-1702. 9. R. M. Issa, M. K. Awad, F. M. Atlam, Appl. Surf. Sci., 2008, 255, 2433-2441. 10. M. Lashgari, M. R. Arshadi, G. A. Parsafar, Corrosion, 2005, 61, 778-783. 11. A. Lesar, I. Milošev, Chem. Phys. Lett., 2009, 483, 198-203. 12. M. S. Masoud, M. K. Awad, M. A. Shaker, M. M. T. El-Ta-hawy, Corros. Sci., 2010, 52, 2387-2396. 13. M. Sahin, G. Gece, F. Karc, S. Bilgic, J. App. Electrochem., 2008, 38, 809-815. 14. M. Arshadi, M. Lashgari, G. Parsafar, Mater. Chem. Phys., 2004, 86, 311-314. 15. E. Jamalizadeh, S. Hosseini, A. Jafari, Corros. Sci., 2009, 51, 1428-1435. 16. M. Lashgari, A. M. Malek, Electrochim. Acta, 2010, 55, 5253-5257. 17. M. Lashkari, M. R. Arshadi, Chem. Phys., 2004, 299, 131137. 18. M. Mousavi, M. Mohammadalizadeh, A. Khosravan, Corros. Sci., 2011, 53, 3086-3091. 19. M. Mousavi, H. Safarizadeh, A. Khosravan, Corros. Sci., 2012, 65, 249-258. Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: 348 Acta Chim. Slov. 2014, 61, 340-349 20. Y. Jiang, J. B. Adams, Surf. Sci., 2003, 529, 428-442. 21. Y. Jiang, J. B. Adams, D. Sun, J. Phys. Chem. B, 2004, 108, 12851-12857. 22. O. Blajiev, A. Hubin, Electrochim. Acta, 2004, 49, 27612770. 23. A. Kokalj, S. Peljhan, Langmuir, 2010, 26, 14582-14593. 24. A. Kokalj, S. Peljhan, M. Finšgar, I. Milošev, J. Am. Chem. Soc., 2010, 132, 16657-16668. 25. S. Peljhan, A. Kokalj, Phys. Chem. Chem. Phys., 2011, 13, 20408-20417. 26. A. Kokalj, N. Kovacevic, S. Peljhan, M. Finšgar, A. Lesar, I. Milošev, ChemPhysChem, 2011, 12, 3547-3555. 27. N. Kovacevic, A. Kokalj, J. Phys. Chem. C, 2011, 115, 24189-24197. 28. N. Kovacevic, A. Kokalj, Mater. Chem. Phys., 2012, 137, 331-339. 29. S. Sun, Y. Geng, L. Tian, S. Chen, Y.Yan, S. Hu, Corros. Sci., 2012, 63, 140-147. 30. X. Chen, H. Häkkinen, J. Phys. Chem. C, 2012, 116, 2234622349. 31. S. Hong, W. Chen, Y. Zhang, H. Q. Luo, M. Li, N. B. Li, Corros. Sci., 2013, 66, 308-314. 32. W. Wang, Z. Li, Q. Sun, A. Du, Y. Li, J. Wang, S. Bi, P. Li, Corros. Sci., 2012, 61, 101-110. 33. F. Grillo, D. W. Tee, S. M. Francis, H. Fruchtl, N. V. Richardson, Nanoscale, 2013, 5, 5269-5273. 34. J. Radilla, G. E. Negron-Silva, M. Palomar-Pardave, M. Romero-Romo, M. Galvan, Electrochim. Acta, 2013, 112, 577586. 35. S. Peljhan, J. Koller, A. Kokalj, J. Phys. Chem. C, 2014, 118, 933-943. 36. A. Kokalj, S. Peljhan, J. Koller, J. Phys. Chem. C, 2014, 118, 944-954. 37. M. Özcan, D. Toffoli, H. Üstünel, Ilyas Dehri, Corros. Sci., 2014, 80, 482-486. 38. S. Singh, J. Singh, S. Gulia, R. Kakkar, J. Theor. Chem., 2013, 2013, 9. 39. D. R. Lide (Ed.), CRC Handbook of Chemistry and Physics, CRC Press, Boca Raton, Florida USA, 85th ed., 2005. 40. W. Reusch, Virtual textbook of organic chemistry, sec. Bond energy, http://www2.chemistry.msu.edu/faculty/reusch/Virt-TxtJml/react2.htm#rx6, 2013. 41. M. Pourbaix, Atlas of Electrochemical Equilibria in Aqueous Solutions, NACE, Cebelcor, Houston, Texas, 2nd ed., 1974. 42. P. D'Arco, M. Causà, C. Roetti, B. Silvi, Phys. Rev. B, 1993, 47, 3522-3529. 43. M. Catti, G. Ferraris, S. Hull, A. Pavese, Phys. Chem. Miner., 1995, 22, 200-206. 44. Y. I. Kuznetsov, L. P. Kazansky, Russ. Chem. Rev, 2008, 77, 219-232. 45. K. F. Khaled, Electrochim. Acta, 2003, 48, 2493-2503. 46. D.-Q. Zhang, L.-X. Gao, G.-D. Zhou, Corros. Sci., 2004, 46, 3031-3040. 47. M. Finšgar, I. Milošev, Corros. Sci., 2010, 52, 2737-2749. 48. A. D. Becke, J. Chem. Phys., 1993, 98, 5648-5652. 49. M. J. Frisch, G. W. Trucks, H.B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, B. Mennucci, G. A. Petersson, H. Nakatsuji, M. Caricato, X. Li, H. P. Hratc-hian, A. F. Izmaylov, J. Bloino, G. Zheng, J. L. Sonnenberg, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, J. A. Montgomery, Jr., J. E. Peralta, F. Ogliaro, M. Bearpark, J. J. Heyd, E. Brothers, K. N. Kudin, V. N. Staro-verov, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, N. Rega, J. M. Millam, M. Klene, J. E. Knox, J. B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R. E. Stratmann, O. Yazyev, A. J. Austin, R. Cammi, C. Pomelli, J. W. Ochterski, R. L. Martin, K. Morokuma, V. G. Zakrzewski, G. A. Voth, P. Salvador, J. J. Dannenberg, S. Dapprich, A. D. Daniels, Ö. Farkas, J. B. Foresman, J. V. Ortiz, J. Cioslowski, D. J. Fox, Gaussian 09 Revision A.1, 2009, Gaussian Inc. Wallingford CT. 50. R. Dovesi, V. R. Saunders, C. Roetti, R. Orlando, C. M. Zi-covich-Wilson, F. Pascale, B. Civalleri, K. Doll, N. M. Harrison, I. Bush, P. D'Arco, M. Llunell, Crystal06 Users Manual, University of Torino, Torino, http://www.crystal.uni-to.it/, 2006. 51. A. Kokalj, J. Mol. Graphics Modelling, 1999, 17, 176-179, code available from http://www.xcrysden.org/. 52. A. V. Marenich, C. J. Cramer, D. G. Truhlar, J. Phys. Chem. B, 2009, 113, 6378-6396. 53. E. Cancès, B. Mennucci, J. Tomasi, J. Chem. Phys., 1997, 107, 3032-3041. 54. B. Mennucci, E. Cances, J. Tomasi, J. Phys. Chem. B, 1997, 101, 10506-10517. 55. S. F. Boys, F. Bernard, Mol. Phys., 1970, 19, 553. 56. H. Ohtaki, T. Radnai, Chem. Rev., 1993, 93, 1157-1204. 57. S. E. Rodriguez-Cruz, R. A. Jockusch, E. R. Williams, J. Am. Chem. Soc., 1999, 121, 1986-1987. 58. M. B. Hay, S. C. B. Myneni, J. Phys. Chem. A, 2008, 112, 10595-10603. 59. R. Caminiti, M. Magini, Chem. Phys. Lett., 1979, 61, 40-44. 60. J. O. Bockris, A. K. N. Reddy, M. Gamboa-Aldeco, Modern Electrochemistry, vol. 2A, Kluwer Academic/Plenum Publishers, New York, Boston, Dordrecht, London, Moscow, 2nd ed., 2000. 61. J. Sauer, Chem. Rev., 1989, 89, 199-255. 62. A. Kokalj, Phys. Rev. B, 2011, 84, 045418. 63. M. F. Peintinger, D. V. Oliveira, T. Bredow, J. Comput. Chem., 2013, 34, 451-459. 64. C. Gatti, V. R. Saunders, C. Roetti, J. Chem. Phys., 1994, 101, 10686-10696. 65. T. Bredow, K. Jug, R. A. Evarestov, Phys. Status Solidi B, 2006, 243, R10-R12. 66. L. Valenzano, Y. Nol, R. Orlando, C. Zicovich-Wilson, M. Ferrero, R. Dovesi, Theor. Chem. Acc., 2007, 117, 9911000. 67. J. L. Whitten, H. Yang, Surf. Sci. Rep., 1996, 24, 55-124. 68. G. P. Brivio, M. I. Trioni, Rev. Mod. Phys., 1999, 71, 231265. 69. V. Russier, C. Mijoule, J. Phys.: Condens. Matter., 1991, 3, 3193-3198. Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: ... Acta Chim. Slov. 2014, 61, 340- 349 349 70. G. Pacchioni, A. M. Ferrari, A. M. Márquez, F. Illas, J. Com-put. Chem., 1997, 18, 617-628. 71. A. Kokalj, A theoretical study of chemisorption and reactions on transition metal surfaces, Ph.D. thesis, University of Ljubljana, Ljubljana, Slovenia, 2000. 72. I. Panas, J. Schüle, P. Siegbahn, U. Wahlgren, Chem. Phys. Lett., 1988, 149, 265-272. 73. P. Reinhardt, M. Causa, C. M. Marian, B. A. Heß, Phys. Rev. B, 1996, 54, 14812-14821. 74. T. B. Grimley, C. Pisani, J. Phys. C: Solid State, 1974, 7, 2831-2848. 75. A. M. Ferrari, G. Pacchioni, Int. J. Quantum Chem., 1996, 58, 241-250. 76. R. A. Evarestov, A. V. Bandura, Int. J. Quantum Chem., 2004, 100, 452-459. 77. A. Kokalj, Corros. Sci., 2014, 79, 215-220. 78. F. Hirshfeld, Theor. Chim. Acta, 1977, 44, 129-138. 79. J. P. Perdew, R. G. Parr, M. Levy, J. L. Balduz, Phys. Rev. Lett., 1982, 49, 1691-1694. 80. L.-A. Naslund, M. Cavalleri, H. Ogasawara, A. Nilsson, L. G. M. Pettersson, P. Wernet, D. C. Edwards, M. Sandstrom, S. Myneni, J. Phys. Chem. A, 2003, 107, 6869-6876. 81. V. S. Bryantsev, M. S. Diallo, W. A. Goddard III, J. Phys. Chem. B, 2008, 112, 9709-9719. Povzetek Proučili smo uporabnost različnih modelov oksidirane kovinske površine (tj. gole kovinske katione, gruče različnih velikosti in periodične modele plošč), ki se uporabljajo na področju kvantnokemijskega modeliranja inhibitorjev korozije. Za vzorčni modelni sistem smo na podlagi teorije gostotnega funkcionala obravnavali imidazol kot inhibitor korozije, površino MgO, in solvatiran kation Mg2+. Čeprav so rezultati gručnega modela odvisni od velikosti in oblike gruče, se zdi, da so gruče srednje velikosti uporabne vsaj za kvalitativne namene. Nasprotno pa so goli kovinski kationi neuporabni ne le kot model oksidiranih površin, ampak tudi kot model solvatiranih kationov, saj vežejo molekule nekajkrat močnejše kot bolj ustrezni modeli. Na primer, goli Mg2+ se veže na imidazol z jakostjo 5,9 eV, medtem ko model plošče Mg0(001) veže imidazol samo z jakostjo 0,35 eV. Takšna vezava je še močnejša za 3+ katione; tj. goli Al3+ se veže na imidazol z jakostjo 17,9 eV. Pokazali smo, da je razlog za te fantastične jakosti vezi predvsem v prenosu znatne količine elektronskega naboja od molekule do kovinskega kationa, ki je posledica razlik med ionizacijskimi potenciali molekule in kovine. Kokalj: Molecular Modeling Of Organic Corrosion Inhibitors: