Influence of metal ions (Zn²⁺, Cu²⁺, Ca²⁺, Mg²⁺ and Na⁺) on the water coordinated neutral and zwitterionic L-histidine dimer

Abstract

The interaction of metal cations with mono hydrated, neutral and zwitterionic histidine dimer complexes are studied using the density functional theory (DFT-B3LYP) method in both gas and liquid phase. It is found that the metal interaction with histidine is governed by two types of bonds (I) N–Mⁿ⁺ (Mⁿ⁺ = Zn²⁺, Cu²⁺, Ca²⁺, Mg²⁺ and Na⁺), O–Mⁿ⁺ bonds and (II) N⋯H–O and N–H⋯O hydrogen bonds. For both the liquid and gas phase, the coordination of metals with the nitrogen atom is strong in neutral systems, whereas the coordination with the oxygen atom is strong in the zwitterions systems. Among all the metal ions considered, the Cu²⁺ ion strongly coordinates with both the neutral and the zwitterionic forms of histidine dimer. AIM analysis indicates that the N–Mⁿ⁺, O–Mⁿ⁺ and hydrogen bonds are partially covalent and electrostatic in nature. The enthalpy of the reaction indicates that the reaction is energetically favoured and exothermic in nature. For all the metal complexes, pK_a values decrease with the increase in hydrogen bond strength. Redox potential is high for the zwitterionic complex because of the transfer of electrons from the carboxylic group to amine moiety in all the metal complexes. From the MD calculations, we find that the backbone Mg²⁺ substituted complexes are stabilized around 0.9 nm, while for the other metal ions it is more than 1.0 nm. Kirkwood's potential indicates the presence of antiparallel dipole–dipole interactions.

---

Introduction

Metal cations are crucial for the basic physical and chemical operations in living organisms. Metal cations have attracted considerable attention in recent years due to their significance in biomedicine. Most of these metals belong to transition elements, alkali and alkaline earth metals, namely Zn²⁺, Cu²⁺, Ca²⁺, Mg²⁺, Mn²⁺, K⁺, Na²⁺, and possess catalytic and structural properties.^1–4 Several divalent metal based compounds are being used as potent antibacterial, antifungal and anticancer drugs.^5,6 Moreover, metals play a vital role in the bio synthesis and metabolism of certain bio active peptides and are relevant to a variety of critical diseases, including arthritis and cancer.^7,8

Amino acids are good metal-complexing agents forming a chelate ring via N, O and S atoms present in the side chains, amino and carboxylate groups of the amino acids. Backbone peptide groups and the side chains of several amino acid residues, such as aspartate (Asp), glutamate (Glu), cysteine (Cys), histidine (His), asparagine (Asn), glutamine (Gln), serine (Ser) and threonine (Thr), frequently coordinate with the metal cation(s) in metalloproteins.^9,10 The interaction of metal cations with amino acid residues can largely modify its pK_a value, inducing important changes in the structure and reactivity of bio molecules.^11,12 Further, redox-inactive metal ions also play a role in non-biological electron transfer.^13–18

Histidine is one of the major ligands in metalloproteins involved in several coordination modes. Histidine exhibits two isomeric forms, D and L, and has two protonation sites, namely N_ε and N_δ. The D-histidine is not bio-active, whereas the L-histidine is bioactive and can be converted to histamine, which is a neurotransmitter in the brain.¹⁹ In addition, the abovementioned zwitterionic form of histidine with the anionic carboxyl group and a protonated amino group finds strong applications in solvent phase. Both the neutral and zwitterionic forms of histidine are shown in Fig. 1 along with their two protonation sites. According to X-ray crystallographic images of some metalloproteins, it is observed that the histidine side chain is favourable for metal coordination.^20–24 Furthermore, Hunsicker et al.²⁵ have suggested that the nitrogen atom of the imidazole ring in the side chain of histidine coordinates with the transition metals. Among the various metal ions, the coordination of zinc is essential due to its vital role in metabolism and gene expression.²⁶ In general, the coordination of zinc within proteins is with histidine and cysteine. For instance, in Cys₂His₂, zif268 and TATA binding protein, the Zn²⁺ ion makes four coordinations, two each with cysteine and histidine, respectively.^7,27 Furthermore, in order to bind DNA, a water molecule is coordinated with zinc.²⁸

Fig. 1 Optimized structure of neutral and zwitterionic histidine.

From the abovementioned paragraphs, we see that zinc coordination is common with histidine in zinc finger protein. Therefore the search for an alternate metal coordination is essential for inhibiting biological processes in the absence of zinc. Hence, this study focuses on the interaction of L-histidine dimer in both neutral and zwitterionic forms with metal cations (Zn²⁺, Cu²⁺, Ca²⁺, Mg²⁺ and Na⁺) and its coordination with single water in both gas and aqueous phase. We believe the results of this study will provide a better understanding regarding the influence of metal cations on the properties pertaining to the chemistry and biological aspects of peptides, which are implicated in several biochemical and physiological processes in a living cell.

Computational details

Quantum mechanical studies

The density functional theory methods (DFT), particularly B3LYP functional method, has been widely adopted as a cost effective method and provides feasible results for studying metal cation-ligands.^29–32 The geometries of the L-histidine dimer hydrated with one water molecule interacting with metals, such as Zn²⁺, Cu²⁺, Ca²⁺, Mg²⁺ and Na⁺ systems, were optimized using the hybrid density functional method such as B3LYP^33–36 with 6-311G** basis set in both gas and aqueous phase. The 6-311G** basis set is a well established Pople basis set and is widely used for the pK_a calculation and redox potential.³⁷ The formation of metal–histidine complexes can be described by the reaction:

Mⁿ⁺ (g) + 2His (g) + H₂O → 2His·Mⁿ⁺ (g)·H₂O (1)

Mⁿ⁺ = Zn²⁺, Cu²⁺, Ca²⁺, Mg²⁺ and Na⁺

The gas-phase interaction enthalpy ΔH (cation binding affinity) for the abovementioned reaction is defined by the following equations:

ΔH = ΔE + pΔV (2)

ΔH²⁹⁸ = [E²⁹⁸_M+2His+H2O − (E²⁹⁸_M + E²⁹⁸_His + E²⁹⁸_H2O)] + ΔpV (3)

ΔH²⁹⁸ = [E²⁹⁸_M+2His+H2O − (E²⁹⁸_M − E²⁹⁸_His − E²⁹⁸_H2O)] − RT (4)

where E²⁹⁸_M+2His+H2O is the energy of the complex, E²⁹⁸_M is the energy of the respective metal cation and E²⁹⁸_His is the energy of the two histidines at T = 298.15 K. Cation basicity is defined as the increment in Gibbs energy for the reaction at T = 298.15 K. Therefore, Gibbs free energies in the gaseous and aqueous phases were calculated based on the estimated standard enthalpies and entropies. Natural population analysis (NPA) has been carried out to describe the charge transfer within the complexes. The pK_a and redox potential for the metallated histidine dimer have been calculated using quantum chemical methods. All the calculations are carried out using the G09 package.³⁸

Molecular dynamics simulation

All the molecular dynamics were carried out by GROMACS 4.6.2 ³⁹ package using the OPLS-AA force field.^40,41 OPLS-AA treats every atom explicitly including the hydrogen atoms in the system. Water molecules were represented using a simple point charge (SPC) model. Appropriate amounts of Cl⁻ counter ions were added by replacing water particles to form the overall charge neutrality of the false organization. Periodic boundary conditions were used in all the three directions of space. At first, an energy minimization procedure was carried out using the steepest descent method. After the energy minimization process, a position restraint procedure was performed in association with NVT and NPT ensembles. The first phase of equilibration involved simulating for 200 ps under the canonical ensemble (NVT), and the temperature was maintained at 300 K with the Nose–Hoover thermostat.^42–44 In the second phase of equilibration, a 200 ps NPT equilibration was performed using the Parrinello–Rahman Barostat⁴⁵ to maintain the pressure at 1.0 bar. In molecular dynamics simulation, the particle mesh Ewald (PME) method interaction was used,^46,47 and the Lincs algorithm was applied for covalent bond constraints.⁴⁸ Cut off distances for computing the electrostatic and van der Waals interaction were 0.9 and 1.0 mm, respectively.

Results and discussion

Coordination mode

L-Histidine molecule has two tautomers (N_ε–H and N_δ–H) by the presence of a hydrogen atom at the N_ε and N_δ positions of the imidazole part of the histidine. Previous literature⁴⁹ showed that histidine dimer binds well with metal at the N_δ position, and hence the metal cations interacted at this site alone in the present study. The labelling of the bonded atoms of the system is shown in Fig. 2. However, the metal cations can simultaneously interact with the carbonyl oxygen (C=O) part of the carboxylic acid, being an electron rich site, in addition to the nitrogen (N_δ) atom of the imidazole. Further, to tie down the DNA in the transcriptional protein, a water molecule is positioned at the tip of the metal cation. In a similar manner, metal binding in the zwitterionic form is also between the imidazole ring and carbonyl oxygen.

Fig. 2 Pictorial representation of the histidine dimer with atom labelling and numbering.

As a steric requirement, the histidine complexes were coordinated by placing a water molecule to arrange in a square pyramidal geometry. The nitrogen and oxygen atoms of the imidazole and carboxylic parts in histidine bonded with the metal atom to form the corners of the square on the basal plane of the structure, while the water molecule bonded to the metal above the plane making the pyramidal structure. This type of geometry is common for complexes where the metal has d⁰ configurations. In the gas phase, for both neutral and zwitterionic forms, all the complexes are coordinated with metals in a pendate manner. The square pyramidal structure of the neutral and zwitterionic complexes are shown in Fig. 3 and 4.

Fig. 3 Optimized structure of the neutral histidine dimer with metal complexes coordinated with one water molecule in gas phase.

Structural parameters

Neutral complex. The optimized neutral complexes are shown in Fig. 3. Bond distances between the metal cations and binding sites of the structures are given in Table 1. The overall bond length values for the metal coordination are in the range of 1.985 Å to 2.477 Å in the gas phase. From the bond distance values, we infer that the metal cations are strongly bound with the N_δ position (Mⁿ⁺–N) of the imidazole. The strongest bond is observed for the N25–Mⁿ⁺ in Cu²⁺ complex, and weakest bond is for N25–Mⁿ⁺ in the Ca⁺ complex. Cu²⁺ strongly binds with nitrogen, whose bond length is 1.985 Å in the imidazole part of the histidine, 2.013 Å in carboxylate oxygen and 2.341 Å with the water molecule. Zinc prefers softer ligands for coordination and can adopt five or six coordinate geometries. The Zn–ligand bond distances are 2.029 Å for Zn–N (His) and 2.136 Å for Zn–O (water), and our results agree well with the values reported in earlier studies.⁵⁰ In contrast, calcium prefers to bind to a hard base such as oxygen containing ligands. The coordination number of calcium in the protein varies from 6 to 8 depending upon the availability of ligands. Ca²⁺ strongly binds with the oxygen atom of carboxylate and water molecule with the bond lengths of 2.266, 2.290 and 2.377 Å, respectively. Similar to calcium, Mg²⁺, being a hard ion, prefers hard ligands, such as oxygen, for the coordination. Here, the bond between Mg²⁺ and oxygen is (2.019 Å) stronger compared with the nitrogen atom in the imidazole (2.169 Å). In the Na⁺ complex, there is a slight deformation in the structure, and the metal binding bonds are weaker than the abovementioned complexes. The intra molecular interaction, such as hydrogen bonding, plays an important role in the stabilization of the dimers. In this study, a strong intramolecular hydrogen bond (N⋯H–O) is formed between the amine (NH₂) and O–H of the carboxylic terminal of the histidine. The hydrogen bond length in all the complexes is in the range of 1.747–1.864 Å (shown in Fig. 3), particularly in the Ca⁺ complex, and the hydrogen bond is strong with a bond length value of 1.747 Å.

Table 1 Optimized bond length (Å) between metals and ligand in both neutral and zwitterionic forms calculated using the B3LYP level of theory with 6-311G** basis set

Bond	Gas phase					Solvent phase
	M^n+ = Zn^2+	M^n+ = Cu^2+	M^n+ = Ca^2+	M^n+ = Mg^2+	M^n+ = Na^+	M^n+ = Zn^2+	M^n+ = Cu^2+	M^n+ = Ca^2+	M^n+ = Mg^2+	M^n+ = Na^+
a The value in parentheses indicates the zwitterionic complex.b Zwitterionic structure not converged.
N25–M^n+	2.029 (2.095)^a	1.985 (2.016)	2.477 (2.521)	2.169 (2.195)	2.430^b	2.012 (2.065)^a	1.969 (1.987)	2.511 (2.552)	2.166 (2.172)	2.446^b
N14–M^n+	2.029 (2.095)	1.985 (2.016)	2.452 (2.484)	2.169 (2.195)	2.429	2.012 (2.065)	1.969 (1.987)	2.493 (2.482)	2.166 (2.172)	2.446
O6–M^n+	2.122 (2.012)	2.013 (1.955)	2.266 (2.211)	2.019 (1.964)	2.326	2.177 (2.06)	2.040 (1.982)	2.347 (2.278)	2.056 (1.992)	2.365
O38–M^n+	2.122 (2.012)	2.013 (1.955)	2.290 (2.284)	2.019 (1.964)	2.326	2.177 (2.06)	2.040 (1.982)	2.376 (2.552)	2.056 (1.992)	2.367
O21–M^n+	2.136 (2.138)	2.341 (2.389)	2.377 (2.347)	2.084 (2.095)	2.321	2.140 (2.177)	2.285 (2.378)	2.401 (2.411)	2.073 (2.085)	2.347

---

To further understand the nature of complexes in aqueous medium, the complexes are optimized in the aqueous phase using the PCM model.⁵¹ The bond distance values of N–Mⁿ⁺, O–Mⁿ⁺ in the aqueous phase are shown in Table 1. In aqueous phase, neutral complexes have the stronger N–Mⁿ⁺ bond for the Zn²⁺, Cu²⁺ and Mg²⁺ complexes than in the gas phase. In contrast, for the Ca²⁺ and Na⁺ complexes, the N–Mⁿ⁺ bond are weaker in comparison to the gas phase. Likewise, the O–Mⁿ⁺ bond of all the complexes in the aqueous phase are weaker than the gas phase. The O–Mⁿ⁺ bond between the oxygen atom of the water molecule and metal complexes is weaker in the Zn²⁺, Ca²⁺ and Na²⁺ complexes and is stronger in the Cu²⁺ and Mg²⁺ complexes. The binding order of the complex with the metals is Cu²⁺ > Zn²⁺ > Mg²⁺ > Ca²⁺ > Na⁺. Both the gas and solvent phase follow the same binding order. The hydrogen bond present in the aqueous phase falls in the range of 1.73 to 1.90 Å, which is weaker than the gas phase.

Further, comparing the bond length values for C=O in the carbonyl group, N–H in amine, C=C, C–N_δ, and C–N_ε in the imidazole ring between the metal histidine complex and its isolated form, we observe significant dissidences in bond length during complexation. In all the complexes, C=O of the carbonyl bond becomes weaker by 0.02 Å, and the bond (C–N_δ) strengthens by 0.04–0.05 Å. The contraction of the C–N_ε bond may be because of the presence of the Mⁿ⁺ atom, which polarizes the carboxylic bond in such a manner that the oxygen atom transfers some of its electron density to the metal cation. Similarly, in the case of the C–N bond, the strengthening is in the order of about 0.04–0.06 Å. The C–N bond contraction subsequently elongated the other C–N bond length of the imidazole by 0.03–0.08 Å from the isolated form. Similarly, because of the interaction of metal, O–H bonds in the carboxylic and water are found to be elongated compared with their corresponding isolated phase.

Zwitterionic complexes. Zwitterionic forms of histidine interact with metals as neutral complexes and are optimized in both gas and aqueous phase using the same level of theory. The structures of the zwitterionic complexes are shown in Fig. 4. The bond lengths between the metal and ligand are weaker than the bond lengths observed in the neutral complex (1.955 to 2.521 Å). The metal and histidine bonds indicate that Cu²⁺ strongly binds with histidine than other metals, as observed in the neutral form. The Zn²⁺ complex has a slightly weaker bond than the Cu²⁺ complex. The Ca²⁺ substituted complex possesses the weak N25–Ca²⁺ bond length of 2.521 Å. However, Zn²⁺, Cu²⁺, Ca²⁺ and Mg²⁺ have higher affinity toward the oxygen atom (being hard base) in the ligand and eventually have a stronger coordination through the Mⁿ⁺–O bond than the Mⁿ⁺–N bond. In contrast, the metal cations bind slightly weaker with an oxygen atom in water than the Mⁿ⁺–N bond of the complex. The binding order of the complex with the metal cation in the gas phase is as follows Cu²⁺ > Zn²⁺ > Mg²⁺ > Ca²⁺. The presence of a strong intramolecular interaction (shown in Fig. 4) between the amine and carboxylic group leads to the formation of the N–H⋯O hydrogen bond. In all the complexes (gas phase), the N–H⋯O hydrogen bonds formed are in the range of 1.674 to 1.848 Å. A strong hydrogen bond is observed for the Mg²⁺ complex, whereas it is weak for Ca²⁺.

Fig. 4 Optimized structure of the zwitterionic histidine dimer with metal complexes coordinated with one water molecule in gas phase.

In solvent phase, the metal interacted complexes, such as Zn²⁺, Cu²⁺, Ca²⁺, have weaker metal–nitrogen bonds (N25–Mⁿ⁺, N14–Mⁿ⁺) than the gas phase. Nevertheless, in the case of O6–Mⁿ⁺, O38–Mⁿ⁺ bonds, the bond lengths are stronger. Similar to gas phase, all the metals in solvent phase have greater affinity toward the oxygen atom (being hard base) in the ligands except for the Ca²⁺ cation, which prefers nitrogen (also hard base) rather than the oxygen atom of the ligands. On comparing the zwitterionic metal interacted complexes with their isolated form, the carbonyl C=O is elongated by ∼0.03 Å, while the C–N_δ bond becomes strengthened. The other C–N_ε bond is weakened through 0.02 Å and 0.06 Å, respectively. This trend is a common feature in both gas and solvent phases. Finally the order of the M–X [X = N, O] distances for Mg²⁺, Ca²⁺ and Zn²⁺ can be approximately correlated with the ionic radii for the three charged metals.

Topological analysis

In addition to the structural parameters, a more insightful appeal comes from the consideration of electron density based on topological parameters, such as the electron density (ρ(r)) and its Laplacian (∇²ρ (r)) at the bond critical point (BCP), which provides valuable insights into the nature of bonding. The bonding nature between the metals and the ligands are analyzed in terms of electron density and its derivatives using Bader's topological atom in molecular analysis (AIM).⁵² The Laplacian of electron density (∇²ρ(r)) indicates whether electron density is locally concentrated (∇²ρ(r) < 0) or depleted (∇²ρ(r) > 0). Therefore, a negative value of ∇²ρ(r) at a BCP is related to the covalent character of a bond, where ∇²ρ(r) > 0 implies “closed-shell-type” interaction as found in ionic bonds, hydrogen bonds, and van der Waals molecules. Further, the Laplacian of electron density can be used to predict reactivity and gain insights into various chemical properties, such as lone pairs, unpaired electrons, electrophilic and nucleophilic reactive sites, and hence provides interesting information about chemical bonds.⁵³ If the Laplacian is negative, the potential energy is dominant, negative charge is concentrated, kinetic energy dominates positively, and negative charge depletes. The electronic energy density (H(r)) at the BCP is defined as H(r) = G(r) + V(r), where G(r) and V(r) correspond to kinetic and potential energy densities, respectively.^54–56 The sign of H(r) determines whether the accumulation of charge at a given point r is stabilizing [H(r) < 0] or destabilizing [H(r) > 0]. The orbital model of the electronic structure, ellipticity provides a quantitative measure of π-bond character and delocalization electronic charge. In addition, ellipticity also measures bond stability; high ellipticity values indicate the instability of the bond.^52,57–60

The values of the topological indices were calculated only for optimized geometries in the gas phase and are presented in the ESI, Table S1.† The values of electron densities obtained for all the complexes are within a range of 0.01922 to 0.08907 a.u. (see Table S1†), and their Laplacian values are positive, indicating the typical closed-shell type of interaction in the complexes. Chen et al. and Sastry et al.^61,62 reported that the electron density of alkali and alkaline earth metals are in the range of 0.010–0.074 a.u. and is subsequently observed for Ca²⁺, Mg²⁺ and Na⁺ complexes in the present study. For transition metals, such as Zn²⁺ and Cu²⁺, N–Mⁿ⁺ bonds for electron density are slightly higher than the abovementioned range in neutral and zwitterionic complexes. The N–Mⁿ⁺ bond in zwitterionic complex behaves like an alkali earth metal ion binding. On comparing various metal cations, Cu²⁺ complex has a greater electron density value, and also has a strong (N26–Mⁿ⁺ and N15–Mⁿ⁺) bond in both the forms. Na⁺ substituted complex has the lowest ρ value due to a weaker bond (N25–Na⁺, N14–Na⁺) between nitrogen atoms and the metal. In general, stronger (N26–Mⁿ⁺ and N15–Mⁿ⁺) bonds are found to have large electron density values.

The positive values of ∇²ρ(r) indicate that the interaction between the metal and the histidine complex is electrostatic in nature. The values of ∇²ρ(r) are more positive in the entire neutral complex than the zwitterionic complex for the N26–Mⁿ⁺ and N15–Mⁿ⁺ bond. The Cu²⁺ complex has the largest Laplacian of the electron density in both the forms, followed by the Zn²⁺ complex. Similar to electron density, stronger bonds have large Laplacian of the electron density. Further, the calculated values of H(r), reported in Table S1,† are negative for N26–Mⁿ⁺ and N15–Mⁿ⁺ bonds in the zinc and copper coordinated complexes, indicating the presence of covalent bonds. However, the positive value of ∇²ρ(r) and the negative value of H(r) conclude that N26–Mⁿ⁺ and N15–Mⁿ⁺ bonds in zinc and copper coordinated complexes are partially covalent and electrostatic in nature. Despite the abovementioned result, the Laplacian and H(r) are positive for the N26–Mⁿ⁺ and N15–Mⁿ⁺ bonds in Mg²⁺, Ca²⁺ and Na⁺ complexes, which suggest electrostatic nature. Ellipticity values, which indicate bond stability for all the complexes and metal nitrogen bond, are in the range of 0.0227–0.0971 a.u. The Cu²⁺ substituted complex has the lowest ellipticity value and indicates that the N26–Cu²⁺ and N15–Cu²⁺ bonds are more stable.

However, the bond between the metal and the oxygen atom of the carboxylic group in amino acids, such as O6–Mⁿ⁺ and O38–Mⁿ⁺, shows different behaviour unlike the N–Mⁿ⁺ bond discussed above. The electron density ρ(r) is higher for the O–Mⁿ⁺ bond in zwitterionic complexes compared to neutral form. The Laplacian of the electron density of O–Mⁿ⁺ is positive for both neutral and zwitterionic forms and does not have the same trend as the N–Mⁿ⁺ bond. Nevertheless, the zwitterionic complex has the higher Laplacian of the electron density for O–Mⁿ⁺ bonds for all the substituted metals. However, the exceptions are for the bonds O38–Ca²⁺, O6–Mg²⁺ and O38–Mg²⁺, which have lower Laplacian of electron density values than the neutral complex. In the Zn²⁺ and Cu²⁺ coordinated complex, the ∇²ρ(r) and H(r) are positive and negative, respectively, which indicates that the bond is partially covalent and electrostatic in nature. On the contrary, for the complexes with Ca²⁺, Mg²⁺ and Na⁺ cations, the Laplacian and H(r) is positive, which indicates that the interaction between the metal and histidine is electrostatic in nature. The ellipticity values of O6–Mⁿ⁺ and O38–Mⁿ⁺ for all the complexes are in the range of 0.0192–0.1393 a.u. and does not support the bond stability criteria.

Further, the topological parameter based on atoms in molecular theory (AIM) has been applied to characterize hydrogen bond strength in various molecular systems.^58,63–72 For the hydrogen bond, ρ and ∇²ρ(r) values are in the range of 0.002–0.34 and 0.016–0.13 a.u., respectively. In the present study, there are two conventional hydrogen bonds present in each complex, i.e. the N⋯H–O and N–H⋯O hydrogen bonds, formed between the carboxylate oxygen and amine of histidine in neutral and zwitterionic complex, respectively. The electron density and Laplacian of the electron density is positive for all the complexes, which shows that the interaction between the carboxylic and amine terminal is electrostatic in nature. For all the hydrogen bonds present in the complex, the Laplacian of the electron density is higher in the zwitterionic complex than the neutral complex. Overall, the ellipticity values are in the range of 0.0796–0.1555 a.u. for the hydrogen bond present in the complexes. The ellipticity value of the complexes indicates that the bond formed in most of the zwitterionic complexes is more stable than the neutral complexes due to low ellipticity values. In addition, the correlation between the hydrogen bond distances, electron density and the Laplacian of electron density are shown in Fig. 5 and 6, which indicate that the bond length and electron density are inverse to each other, i.e., an increase in hydrogen bond length corresponds to a decrease in electron density because an increase in distance results in reduced orbital overlap and hence low electron density. The correlation coefficient values obtained are mostly equal to unity except for the Laplacian of the electron density in neutral complexes. The H(r) is negative for almost all the complexes except for the single hydrogen bond (N2⋯H5–O8) in the calcium coordinated complex. The negative value of H(r) indicates that the bond formed between the carboxylate oxygen atom and NH of the amine group is electrostatically dominant. In the Ca²⁺ zwitterionic complex, among the two hydrogen bonds, electron density of one hydrogen bond (N41⋯H44–O40) is higher than the neutral one. The N41⋯H44–O40 hydrogen bond also has positive Laplacian and H(r) values, which indicates that the bond is purely electrostatic in nature. Thus, from the abovementioned results, the positive Laplacian and negative H(r) show that the hydrogen bonds formed in the complexes possess partially covalent and partially electrostatic nature.

Fig. 5 Correlation between electron density at bond critical point and hydrogen bond distance for neutral and zwitterionic complexes in gas phase.

Fig. 6 Correlation between the Laplacian of electron density at bond critical point and hydrogen bond distance for neutral and zwitterionic complexes in gas phase.

In general, from the abovementioned AIM analysis in both the forms, the metal binds with imidazole nitrogen partially covalently and electrostatically. The electron density value indicates strong N–Mⁿ⁺ and O–Mⁿ⁺ bonds in neutral complexes and zwitterionic complexes, respectively. The Laplacian of the electron density and H(r) values conclude that hydrogen bond is partially covalent and electrostatic in nature.

Proton dissociation enthalpies and Gibbs energies

To describe the metallization process of an amino acid, it is important to study the affinities, enthalpies and Gibbs free energies of the metal cation interacting with amino acids. In general, metals can strongly interact with heteroatoms in histidine and affect electron distribution. These interactions can change the basicity and acidity of the groups involved in hydrogen bonding by affecting the hydrogen bond strength. The Gibbs free energy, enthalpy and entropy difference (i.e., ΔG, ΔH and ΔS) for the histidine dimer interacting with metals, which is hydrated with one water molecule, are shown in Table 2.

Table 2 Computed gas and solvent phase enthalpies (ΔH), entropies (ΔS) and Gibbs free energies (ΔG) (in kcal mol⁻¹) of the metal histidine hydrated with one water molecule

Complex	Gas phase			Solvent phase
	ΔH^298	ΔS^298	ΔG^298	ΔH^298	ΔS^298	ΔG^298
a The value in parentheses indicates the zwitterionic complex.
Zn^2+(His)2·H2O	−404.40 (−379.25)^a	−56.70 (−55.07)	−370.10 (−346.626)	−83.33 (−93.72)	−57.00 (−55.81)	−49.56 (−60.65)
Cu^2+(His)2·H2O	−545.14 (−521.45)	−57.77 (−55.60)	−510.91 (−488.51)	−194.20 (−208.82)	−57.35 (−54.23)	−160.22 (−176.69)
Ca^2+(His)2·H2O	−311.85 (−292.67)	−51.97 (−56.10)	−281.06 (−259.43)	−48.90 (−63.87)	−52.51 (−54.20)	−17.79 (−31.75)
Mg^2+(His)2·H2O	−380.52 (−358.79)	−55.52 (−55.10)	−347.63 (−326.14)	−84.53 (−99.31)	−55.14 (−54.57)	−51.85 (−22.24)
Na^+(His)2·H2O	−117.01	−44.58	−90.598	−29.65	−43.72	−3.74

---

Standard enthalpy values shows that all the complexes formed are energetically favoured exothermically. In both the neutral and zwitterionic forms, the copper bonded complex forms a strong base in both the phases with the proton affinity values of −545.14 kcal mol⁻¹ (neutral, gas), −521.45 kcal mol⁻¹ (zwitterionic, gas), −194.20 kcal mol⁻¹ (neutral, solvent), and −208.826 kcal mol⁻¹ (zwitterionic, solvent). Cu²⁺ complex has strong Cu²⁺–X (X = N, O) bonds with histidine and water, and hence its proton affinity is high. Next to the copper binding complex, the Zn²⁺ complex has high proton affinity, which is around ∼140 kcal mol⁻¹ lesser than the Cu²⁺ complex. The Na⁺ complex in neutral form has low proton affinity in all the phases. This implies that proton affinity of the complexes depends upon the interaction affinity of the metal and ligand. In both the neutral and zwitterionic forms, the tendency of binding enthalpies in gas phase has the same trend as follows:

Cu²⁺ > Zn²⁺ > Mg²⁺ > Ca²⁺ > Na⁺ and Cu²⁺ > Zn²⁺ > Mg²⁺ > Ca²⁺

On comparing the neutral and zwitterionic enthalpies in the gas phase, zwitterionic forms have slightly lower values than the neutral ones. The differences between the neutral and zwitterionic forms are around ∼20–25 kcal mol⁻¹. However, the scenario changes in the solvent phase, where the order of the enthalpies is Cu²⁺ > Mg²⁺ > Zn²⁺ > Ca²⁺ > Na⁺ (neutral) and Cu²⁺ > Mg²⁺ > Zn²⁺ > Ca²⁺ (zwitterionic) in both the forms. Thus, Cu²⁺ has more binding affinity toward the histidine dimer.

For the metal interacted complexes, the absolute value of entropy (ΔS) is lower than the absolute value of enthalpy (ΔH). During complexation, the entropy of the system decreases from its isolated phase. The entropy of the Cu²⁺ complex is maximum in both the neutral and zwitterionic forms, while for Na⁺ it is minimum in neutral form for gas and solvent phases. The ΔG value of all the complexes is in the range from −3.74 to −510.91 kcal mol⁻¹. In both the gas and solvent phases of neutral complexes, the Na⁺ complex has the minimum Gibbs free energy, while the Cu²⁺ complex has the maximum Gibbs free energy. The calculated Gibbs free energies are negative, indicating that all the complex formations in both the phases considered here are chemically feasible. Based on the Gibbs free energy, the stability and order of both the forms in gas phase is Na⁺ > Ca²⁺ > Mg²⁺ > Zn²⁺> Cu²⁺. However, in the solvent phase, the stability and order behave differently for both the forms. The stability order is Na⁺ > Ca²⁺ > Zn²⁺ > Mg²⁺ > Cu²⁺ for neutral form and Mg²⁺ > Ca²⁺ > Zn²⁺ > Cu²⁺ for zwitterionic forms.

pK_a values

In proteins, polar amino acid residues form hydrogen bonds, which contribute to the stabilization of entire proteins. When a hydrogen atom is removed from one such residue to estimate its pK_a value, two negative atoms in the amino acid residues repel each other, leading to the destabilization of the protein structure. However, small differences in the pK_a scale may entail large changes in the degree of ionization and related properties such as enzyme inhibition and environmental fate. Therefore, the knowledge of the pK_a for each compound is of fundamental importance to obtain a better comprehension of the compounds. The interaction of metal cations with the amino acid can largely modify its pK_a value, inducing important changes in the structure and reactivity of biomolecules. Thus, the knowledge of how metal cation coordination influences the acidity/basicity of relevant residues is of significant importance. For the calculation of pK_a value of the metal interaction complex, we use the thermodynamic cycle shown in Scheme 1. Scheme 1 is commonly used for calculating pK_a values:

ΔG_aq,1 = G_g(AH) + G_g(H⁺) − G_g(A⁻) + ΔG_s(AH) + ΔG_s(H⁺) − ΔG_s(A⁻) + RTIn(24.46) (5)

Scheme 1 Scheme for calculating pK_a and redox potential.

For G_gas(H⁺), we employed the value of −6.28 kcal mol⁻¹, as obtained from the Sackur–Tetrode equation.^73,74 G_gas(H⁺) uses a reference state of 1 atm:

The thermodynamic cycle for calculating the value of ΔG_aq follows Scheme 1, which is subsequently used for calculating pK_a values using eqn (6).

The calculated values of ΔG_aq for the dissociation and pK_a values are presented in Table 3. The deprotonation of aqueous Gibbs free energy generally increases when the acidity of the compound decreases. For both the neutral and zwitterionic forms, the calculated pK_a values are in the range of 45.39 to 47.79. The calculated pK_a values for the neutral complex are in the range of 45.96 to 47.79. In the neutral form, Mg²⁺ and Na⁺ substituted complexes have high (47.79) and low (45.96) pK_a values, respectively. In the case of the zwitterionic complex, Cu²⁺ complex has the maximum pK_a value of about 47.42, while Zn²⁺ complex shows the minimum value of about 45.39. The value also varies due to different metals in the complexes. Earlier studies^75–77 have reported the influence of the hydrogen bond on the pK_a values. In the complexes considered here, there are two hydrogen bonds (N–H⋯O and N⋯H–O) whose bond length could affect the pK_a value. Based on the hydrogen length (shown in Fig. 2 and 3), we find that the neutral complexes possess slightly weaker hydrogen bonds than the zwitterionic complexes. Thus, the pK_a value for neutral complexes are high compared to the zwitterionic complexes. This agrees well with the observation that the stronger hydrogen bond leads to decreasing pK_a value.⁷⁵ The correlation between the ΔG_aq and pK_a is linear (Fig. 7) with a correlation coefficient of 0.9999.

Table 3 Calculated pK_a and redox potential of the neutral and zwitterionic complexes at the B3LYP/6-311G** level of theory

Complex	ΔGaq		pKa		Redox potential
	Neutral	Zwitterionic	Neutral	Zwitterionic	Neutral	Zwitterionic
a Not available.
Zn^2+(His)2·H2O	179.706	170.825	47.75	45.39	−0.338	−0.322
Cu^2+(His)2·H2O	178.681	178.483	47.48	47.42	−0.337	−0.336
Ca^2+(His)2·H2O	177.533	173.091	47.17	45.99	−0.334	−0.326
Mg^2+(His)2·H2O	179.875	177.890	47.79	47.27	−0.339	−0.335
Na^+(His)2·H2O	172.968	^a	45.96	^a	−0.326	^a

---

Fig. 7 Correlation between ΔG_aq values with the calculated pK_a values in both neutral and zwitterionic complexes.

Redox potential

The redox potential of a protein is interesting for many reasons. The redox potential of proteins influence the interaction of the protein with other molecules, enzymatic reactions, charge transfer reaction and protein stability.^78–81 It provides a quantitative measure for determining the feasibility and direction of the reduction and oxidation reactions, where the positive charge is transferred spontaneously from an oxidizing agent having a higher reduction potential to a reducing agent having a lower reduction potential.⁸² The calculated redox potential of the complexes is listed in Table 3, and the scheme for calculating the redox potential is shown in Scheme 1. In accordance with the Nernst's law, the redox potential (E_redox) can be calculated as:where F is the Faraday constant, ΔG_redox and n are the Gibbs reaction energy and the number of electrons participating in a given redox reaction, respectively. The redox potential values range from −0.326 to −0.339 for neutral complexes, while for zwitterionic complexes it range from −0.322 to −0.336. Redox potential for the zwitterionic complexes is higher than the neutral complexes because of electron transfer from the carboxylic acid group to the amine terminal of the histidine. A negative redox potential value indicates that the complexes are strong reducing agents and have the tendency to lose electrons to new species. Thus, the redox potential values indicate that the reaction between the metal and ligand is spontaneous and reverse in direction. Further, the Cu²⁺ complex has more negative redox potential in the zwitterionic complex because copper is a strong reducing agent. In the case of Zn²⁺ and Mg²⁺ complex, redox value is more negative in the neutral complex than the zwitterionic complex. The redox potential order for the neutral and zwitterionic form is as follows:

Na⁺ > Ca²⁺ > Cu²⁺ > Zn²⁺ > Mg²⁺ and Zn²⁺ > Ca²⁺ > Mg²⁺ > Cu²⁺

Ligand to metal charge transfer

Charge transfer (CT) is one of the key components that stabilize histidine and metal complexes. In this study, the natural population analysis (NPA) was used to discuss the qualitative behaviour of charge transfer during the complexation. NPA charges only for selected atoms of metal coordinated complexes are presented in Table S2.† The amount of CT between metal and ligands are easily determined as the difference between the charges of the metal in the corresponding complexes. Comparing the metal with and without the ligand, the metal cation coordination shifted the electron density towards the metals, which is indicated by the decrease in positive charges in the metal cation. In all the complexes, the electron charge transfer occurs from the Lewis base (ligand) to the Lewis acid (metal cation), which was also reported in an earlier study.⁸³ Overall, the electronic charge of metals decreases in the order Mg²⁺ > Cu²⁺ > Ca²⁺ > Zn²⁺ > Na⁺ in neutral complexes and Mg²⁺ > Cu²⁺ > Ca²⁺ > Zn²⁺ in zwitterionic complexes.

Further, the extent of charge transfer in coordination complexes can be significant^84–89 because there is a small energy difference between the (n − 1)d and ns valence orbitals and the Rydberg orbitals (such as np and nd) of the metal ion.^90,91 These orbitals participate in hyper conjugative interactions with the ligand, facilitating charge transfer between the metal and the ligand.^92,93 From Table 4, the NPA result indicates that both the valence orbitals (4s and 3d) and the near-valent Rydberg orbitals (4p and 4d) are all partially occupied for Zn²⁺, Cu²⁺, Mg²⁺, and Ca²⁺. This indicates the transfer of electron charge between the metal, Lewis acid, and the ligand, Lewis base, upon the formation of the complex. However, in the Na⁺ complex, the 3s and 3p orbitals are occupied, which indicates that some amount of charge is transferred between the ligand and metal.

Table 4 Natural occupancy of orbitals (in e) of Mⁿ⁺ ions in histidine dimer complexes from natural population analysis (NPA) using the B3LYP/6-311G** level of theory

M^n+	4s	3d	4p	4d	Total
Zn^2+	0.33	9.96	0.01	0.001	10.30
Cu^2+	0.30	9.19	0.004	0.001	9.49
Ca^2+	0.12	0.091	0.18	0.000	0.391
Mg^2+	0.002	0.004	0.01	0.0	0.034
	3s	3p	—	—	—
Na^+	0.10	0.16			0.26

---

Structural stability

Molecular dynamics (MD) simulations of metal interacted histidine complexes were carried out to characterize structural fluctuations due to the introduction of metal cations. Root mean square deviation (RMSD) calculations indicate that the change in the metal cation altered the structure and flexibility of the histidine dimer system. To confirm the stability of the system, RMSD of the backbone atoms were analyzed. Fig. 8 shows the mean-square fluctuations plot of alpha carbons obtained from the ps trajectory. Interestingly, during the simulation, all the metals attain the equilibrated trajectory around 300–500 ps, after which the system loses its equilibrium trajectory. The RMSD of the zinc substituted metal complexes fluctuate randomly at about 1.2 Å, 1.4 Å and 1.6 Å for the 400, 700 and 950 ps, respectively. The copper substituted complex exhibits 2.0 Å around 800 ps and is stabilized around 1.5 Å. The calcium substituted complex also behaves similar to that of the copper substituted complex, thereby attaining equilibrium at around 1.6 Å. The magnesium substituted complex retains its structural stability around 1.0 Å and maintains this throughout the simulation. However, the trajectory of the sodium substituted complex steadily increases up to 0.8 Å, fluctuates up to 1.8 Å, and finally falls down at 1.2 Å. This indicates that Na⁺ ion is unable to form a stable complex with the histidine complex.

Fig. 8 Plot for root mean square deviation of various metal substituted complexes.

Kirkwood's potential

The Kirkwood–Froehlich correlation factor of dielectric polarization (g) is an important parameter, which provides insight into short range intermolecular forces that arise due to dipole–dipole interactions and hence provides information regarding the local order of the polar solute. Thermodynamic properties, such as kinetic energy, potential energy, total energy, and dielectric properties, such as finite Kirkwood factor G_k and infinite Kirkwood factor g_k and dielectric constant (ε₀), are given in Table 5. The finite system Kirkwood factor G_k is determined by the formula:where M is the dipole moment of the total simulation system, N is the total number of molecules and μ is the dipole moment of a single molecule. The infinite Kirkwood factor g_k is determined by using the formula:⁹⁴

Table 5 Finite and infinite system Kirkwood's factors (G_k and g_k) and dielectric constant (ε) calculated from MD calculations

Complex	Gk	gk	ε
Zn^2+(His)2·H2O	0.100	0.095	1.168
Cu^2+(His)2·H2O	0.132	0.126	1.177
Ca^2+(His)2·H2O	0.266	0.246	1.286
Mg^2+(His)2·H2O	0.090	0.086	1.153
Na^+(His)2·H2O	1.205	0.987	2.185

---

Dielectric constants ε(0) were calculated from the fluctuations in the total dipole moment 〈M²〉 of the system⁹⁵ using a Clausius–Mossotti type equation for reaction fields:

where ε₀ is the vacuum permittivity, V is the volume, K_B is the Boltzmann's constant, and T is the temperature. To estimate the uncertainty ε(0), a graph of the convergence of ε(0) with respect to the simulation time is shown in Fig. 9. The divalent cation substituted complexes show that the dielectric constant is stable throughout the simulation; however, for the monovalent (Na⁺) cation, the value increases linearly and attains stability around 800 ps with a dielectric constant of about 2.15. From Table 5, it is also evident that the ε values of the metal coordinated complexes are in the range of 1.153 to 2.185. The minimum value of ε is observed in Mg²⁺ substituted complexes, which is shown in the graph with a blue line. The ε value of the Zn²⁺ complex increases slightly around 180 ps and stabilizes at 380 ps. The ε value of the Ca²⁺ complex is 1.286 and stabilizes at 300 ps. For the Cu²⁺ complex, the ε value stabilizes over 500–900 ps and thereafter increases to 1.4. The Na⁺ cation behaves slightly different from other divalent cations in such a manner that the ε value stabilizes after 800 ps.

Fig. 9 Changes in ε graph for metal interacted histidine complexes with respect to time.

The values of g_k relate directly to the correlation between the directions of a particular molecule under consideration with the neighbouring molecules in solution. If the dipole–dipole interactions lead to a positive correlation and parallel orientation and give g > 1, they are called α-multimers. If the dipole–dipole interactions lead to a negative correlation, anti-parallel orientation, and give g < 1, they are called β-multimers^96,97 molecules in solutions. The g_k values are shown in the Table 5, which indicates that all the values are below unity, which indicates the presence of anti-parallel dipole–dipole interactions in the system.

Conclusion

From this DFT study of the interaction between metal cations and histidine dimer, we have obtained the following major conclusions:

1. Metal coordinated histidine dimer with one water hydrated complex exhibits square pyramidal geometry, and its coordination is in pendate manner.

2. The N–Mⁿ⁺ and O–Mⁿ⁺ bonds are stronger in neutral and zwitterionic complexes possessing partially covalent and electrostatic character.

3. Coordination of metal cations with the nitrogen atom (N–Mⁿ⁺) is strong in neutral complexes, however; the binding of zwitterionic system with the oxygen atom (O–Mⁿ⁺) is strong.

4. Cu²⁺ ion coordinates strongly with histidine dimer and is followed by Zn²⁺ ion in both the gas and aqueous phases. The sodium monovalent cation exhibits weak coordination.

5. C–N_δ and C–N_ε bonds in histidine are strengthened and weakened, respectively, due to metal cation substitutions.

6. From AIM analysis, the N–Mⁿ⁺, O–Mⁿ⁺ bonds are partially covalent while the hydrogen bonds are electrostatic in nature.

7. The negative Gibbs free energy for all the complexes in both gas and liquid phase indicate that the reaction is chemically feasible.

8. The pK_a value of the complex depends upon the hydrogen bond strength.

9. Negative redox potential indicates that the redox reaction is spontaneous and reverse in direction.

10. The charge transfer occurs from Lewis base (ligand) to the Lewis acid (metal) and is large for Mg²⁺ ion.

11. MD calculations show that the backbone stability obtained for Mg²⁺ substituted complexes around is 0.9 nm, while it is more than 1.0 nm for other metal ions.

12.The values of Kirkwood's potential suggest the presence of antiparallel dipole–dipole interactions in zwitterionic form.

Biological implication

This work would provide some insights into the metal binding nature of histidine in some metallo proteins. The mimicry of metal cations into amino acids may offer opportunities for developing novel probes to reveal the functions of proteins, design sensors, and developing new drugs for specific diseases. The presence of various coordination numbers and geometry leads to structural diversity and could be exploited for designing pharmacophores.

Acknowledgements

Authors Dr L. Senthilkumar and P. Umadevi gratefully acknowledge the DST-SERB, New Delhi, India, for granting the project and fellowship under the Fast-Track Scheme.

References

1. H. Sigel, A. Sigel, Probing of proteins by metal ions and their low- molecular weight complexes, 2000, vol. 38.
2. N. S. Poonia and A. V. Bajaj, Chem. Rev., 1979, 79, 389–445.
3. P. G. Daniele, C. Foti, A. Gianguzza, E. Prenesti and S. Sammartano, Coord. Chem. Rev., 2008, 252, 1093–1107.
4. Advances in Protein Chemistry, ed. J. P. Glusker, C. B. Anfinsen and S. E. David, 1991, vol. 42, pp. 1–76.
5. C. Orvig and M. J. Abrams, Chem. Rev., 1999, 99, 2201–2204.
6. S. Jurisson, D. Berning, W. Jia and D. Ma, Chem. Rev., 1993, 93, 1137–1156.
7. W. N. Lipscomb and N. Sträter, Chem. Rev., 1996, 96, 2375–2434.
8. C. E. Brinckerhoff and L. M. Matrisian, Nat. Rev. Mol. Cell Biol., 2002, 3, 207–214.
9. Y. M. Lee and C. Lim, J. Mol. Biol., 2008, 379, 545–553.
10. S. Campbell, M. T. Rodgers, E. M. Marzluff and J. L. Beauchamp, J. Am. Chem. Soc., 1995, 117, 12840–12854.
11. P. D. Schnier, W. D. Price, R. A. Jockusch and E. R. Williams, J. Am. Chem. Soc., 1996, 118, 7178–7189.
12. S. G. Summerfield, A. Whiting and S. J. Gaskell, Int. J. Mass Spectrom. Ion Processes, 1997, 162, 149–161.
13. A. L. Feig and S. J. Liard, Chem. Rev., 1994, 94, 759–805.
14. T. Inoue, H. Sugawara, S. Hamanaka, H. Tsukui, E. Suzuki, T. Kohzuma and Y. Kai, Biochemistry, 1999, 38, 6063–6069.
15. M. J. Ellis, F. E. Dodd, R. W. Strange, M. Prudencio, G. Sawers, R. R. Eady and S. S. Hasnain, Acta Crystallogr., Sect. D: Biol. Crystallogr., 2000, 57, 1110–1118.
16. M. Namazian and S. Halvani, J. Chem. Thermodyn., 2006, 38, 1495–1502.
17. M. Risch, K. Klingan, F. Ringleb, P. Chernev, I. Zaharieva, A. Fischer and H. Dau, ChemSusChem, 2012, 5, 542–549.
18. M. W. Kanan, J. Yano, Y. Surendranath, M. Dincă, V. K. Yachandra and D. G. Nocera, J. Am. Chem. Soc., 2010, 132, 13692–13701.
19. A. M. Wade and H. N. Tucker, J. Nutr. Biochem., 1998, 9, 308–315.
20. S. Gubba, V. Cipriano and J. M. Musser, Infect. Immun., 2000, 68, 3716–3719.
21. C. C. Tseng, M. Miyamoto, K. Ramalingam, K. C. Hemavathy, M. J. Levine and N. Ramasubbu, Arch. Oral Biol., 1999, 44, 119–127.
22. N. Wolff, C. Deniau, S. Létoffé, C. Simenel, V. Kumar, I. Stojiljkovic, C. Wandersman, M. Delepierre and A. Lecroisey, Protein Sci., 2002, 11, 757–765.
23. T. Soulimane, G. Buse, G. P. Bourenkov, H. D. Bartunik, R. Huber and M. E. Than, EMBO J., 2000, 19, 1766–1776.
24. E. Frieden, S. Osaki and H. Kobayashi, J. Gen. Physiol., 1965, 49, 213–252.
25. L. M. Hunsicker-Wang, A. Heine, Y. Chen, E. P. Luna, T. Todaro, Y. M. Zhang, P. A. Williams, D. E. McRee, J. Hirst, C. D. Stout and J. A. Fee, Biochemistry, 2003, 42, 7303–7317.
26. B. L. Vallee and D. S. Auld, Biochemistry, 1990, 29, 5647–5659.
27. A. Belker and H. Azeddoug, Microbiol. Res., 2013, 168, 99–105.
28. J. E. Coleman, Curr. Opin. Chem. Biol., 1998, 2, 222–234.
29. A. Ignaczak and J. A. N. F. Gomes, J. Electroanal. Chem., 1997, 420, 209–218.
30. A. Ignaczak and J. N. F. Gomes, Chem. Phys. Lett., 1996, 257, 609–615.
31. J. Jover, R. Bosque and J. Sales, Dalton Trans., 2008, 45, 6441–6453.
32. D. M. Popovic and A. A. Stuchebrukhov, Photochem. Photobiol. Sci., 2006, 5, 611–620.
33. C. Lee, W. Yang and R. G. Parr, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785–789.
34. B. Miehlich, A. Savin, H. Stoll and H. Preuss, Chem. Phys. Lett., 1989, 157, 200–206.
35. A. D. Becke, Phys. Rev. A, 1988, 38, 3098–3100.
36. A. D. Becke, J. Chem. Phys., 1993, 98, 5648–5652.
37. J.-y. Hasegawa and K. Matsuda, Theor. Chem. Acc., 2011, 130, 175–185.
38. 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. Hratchian, 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. J. Bearpark, J. Heyd, E. N. Brothers, K. N. Kudin, V. N. Staroverov, R. Kobayashi, J. Normand, K. Raghavachari, A. P. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, N. Rega, N. J. 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 and D. J. Fox, Gaussian, Inc., Wallingford, CT, USA, 2009.
39. B. Hess, C. Kutzner, D. van der Spoel and E. Lindahl, J. Chem. Theory Comput., 2008, 4, 435–447.
40. W. L. Jorgensen, D. S. Maxwell and J. Tirado-Rives, J. Am. Chem. Soc., 1996, 118, 11225–11236.
41. G. A. Kaminski, R. A. Friesner, J. Tirado-Rives and W. L. Jorgensen, J. Phys. Chem. B, 2000, 105, 6474–6487.
42. W. G. Hoover, Phys. Rev. A, 1985, 31, 1695–1697.
43. S. Nosé, Mol. Phys., 1984, 52, 255–268.
44. S. Nosé and M. L. Klein, Mol. Phys., 1983, 50, 1055–1076.
45. M. Parrinello and A. Rahman, J. Appl. Phys., 1981, 52, 7182–7190.
46. T. Darden, D. York and L. Pedersen, J. Chem. Phys., 1993, 98, 10089–10092.
47. U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen, J. Chem. Phys., 1995, 103, 8577–8593.
48. B. Hess, H. Bekker, H. J. C. Berendsen and J. G. E. M. Fraaije, J. Comput. Chem., 1997, 18, 1463–1472.
49. L. Zhou, S. Li, Y. Su, X. Yi, A. Zheng and F. Deng, J. Phys. Chem. B, 2013, 11, 78954–78965.
50. T. Dudev and C. Lim, Chem. Rev., 2003, 103, 773–788.
51. J. Tomasi, B. Mennucci and R. Cammi, Chem. Rev., 2005, 105, 2999–3094.
52. R. F. W. Bader, Chem. Rev., 1991, 91, 893–928.
53. S. Miertuš, E. Scrocco and J. Tomasi, Chem. Phys., 1981, 55, 117–129.
54. D. Cremer and E. Kraka, Angew. Chem., 1984, 23, 627–628.
55. K. B. Wiberg, Tetrahedron, 1968, 24, 1083–1096.
56. S. Jenkins and I. Morrison, Chem. Phys. Lett., 2000, 317, 97–102.
57. P. Popelier, Atoms in molecules, an introduction, Pearson Education, London, 2000.
58. S. Mahanta, B. K. Paul, R. Balia Singh and N. Guchhait, J. Comput. Chem., 2011, 32, 1–14.
59. R. F. W. Bader, Atoms in Molecules, a Quantum Theory, Oxford, 1990.
60. B. K. Paul, S. Mahanta, R. B. Singh and N. Guchhait, J. Phys. Chem. A, 2010, 114, 2618–2627.
61. W.-Z. Xu, F.-D. Ren, J. Ren, S.-N. Liu, Y. Yue, W.-L. Wang and S.-S. Chen, J. Mol. Model., 2010, 16, 615–627.
62. J. R. Premkumar, D. Vijay and G. N. Sastry, Dalton Trans., 2012, 41, 4965–4975.
63. R. Bianchi, G. Gervasio and D. Marabello, Inorg. Chem., 2000, 39, 2360–2366.
64. J. N. Woodford, J. Phys. Chem. A, 2007, 111, 8519–8530.
65. U. Koch and P. Popelier, J. Phys. Chem., 1995, 99, 9747–9754.
66. P. Popelier, J. Phys. Chem., 1998, 102, 1873–1878.
67. L. Senthilkumar, T. K. Ghanty and S. K. Ghosh, J. Phys. Chem. A, 2005, 109, 7575–7582.
68. L. Senthilkumar, T. K. Ghanty, S. K. Ghosh and P. Kolandaivel, J. Phys. Chem. A, 2006, 110, 12623–12628.
69. L. Senthilkumar, P. Umadevi, K. Nithya and P. Kolandaivel, J. Mol. Model., 2013, 19, 3411–3425.
70. R. Shankar, P. Kolandaivel and L. Senthilkumar, J. Phys. Org. Chem., 2011, 24, 553–567.
71. V. Umadevi, L. Senthilkumar and P. Kolandaivel, Mol. Simul., 2013, 39, 908–921.
72. A. Mano Priya, L. Senthilkumar and P. Kolandaivel, Struct. Chem., 2014, 25, 139–157.
73. M. D. Liptak and G. C. Shields, J. Am. Chem. Soc., 2000, 123, 7314–7319.
74. F. Blanco, I. Alkorta and J. Elguero, Magn. Reson. Chem., 2007, 45, 797–800.
75. M. D. Joshi, G. Sidhu, J. E. Nielsen, G. D. Brayer, S. G. Withers and L. P. McIntosh, Biochemistry, 2000, 40, 10115–10139.
76. P.-F. Wang, M. J. McLeish, M. M. Kneen, G. Lee and G. L. Kenyon, Biochemistry, 2000, 40, 11698–11705.
77. P. A. Tishmack, D. Bashford, E. Harms and R. L. Van Etten, Biochemistry, 1997, 36, 11984–11994.
78. J. M. Antosiewicz and D. Shugar, Mol. BioSyst., 2011, 7, 2923–2949.
79. D. M. Popović and A. A. Stuchebrukhov, Biochim. Biophys. Acta, Bioenerg., 2012, 1817, 506–517.
80. M. Ubbink, Photosynth. Res., 2004, 81, 277–287.
81. G. M. Ullmann, J. Phys. Chem. B, 2000, 104, 6293–6301.
82. N. Yamamoto and K. Kuwata, J. Biol. Inorg. Chem., 2009, 14, 1209–1218.
83. P. Lipkowski, S. J. Grabowski and J. Leszczynski, J. Phys. Chem. A, 2006, 110, 10296–10302.
84. P. R. Varadwaj and H. M. Marques, Phys. Chem. Chem. Phys., 2010, 12, 2126–2138.
85. P. R. Varadwaj, A. Varadwaj and H. M. Marques, J. Phys. Chem. A, 2011, 115, 5592–5601.
86. P. Varadwaj and H. Marques, Theor. Chem. Acc., 2010, 127, 711–725.
87. P. R. Varadwaj, I. Cukrowski and H. M. Marques, J. Phys. Chem. A, 2008, 112, 10657–10666.
88. A. E. Reed, L. A. Curtiss and F. Weinhold, Chem. Rev., 1988, 88, 899–926.
89. P. R. Varadwaj, A. Varadwaj, G. H. Peslherbe and H. M. Marques, J. Phys. Chem. A, 2011, 115, 13180–13190.
90. S. G. Wang and W. H. E. Schwarz, Angew. Chem., 2009, 48, 3404–3415.
91. W. H. E. Schwarz, J. Chem. Educ., 2010, 87, 444–448.
92. M. M. Francl, W. J. Pietro, R. F. Hout and W. J. Hehre, Organometallics, 1983, 2, 281–286.
93. I. V. Alabugin, K. M. Gilmore and P. W. Peterson, Comput. Mol. Biosci., 2011, 1, 109–141.
94. M. Neumann, Mol. Phys., 1986, 57, 97–121.
95. M. Neumann, Mol. Phys., 1983, 50, 841–858.
96. N. E. Hill, W. E. Vaughan, A. H. Price and M. Davies, Dielectric Properties and Molecular Behaviour (Physical Chemistry), Van Nostrand Reinhold Co., London, 1969.
97. N. V. P. Rao, Int. J. Macromol. Sci., 2011, 1, 23–25.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/c4ra08155b

---