Quantum mechanics and molecular mechanics study of the reaction mechanism of quorum quenching enzyme: N-acyl homoserine lactonase with C6-HSL

Abstract

N-Acyl-homoserine lactonase from Ochrobactrum sp. strain (AidH) is a novel AHL (N-acyl-homoserine lactone)-lactonase that hydrolyzes the ester bond of the homoserine lactone ring of AHLs. In this article, on the basis of the high-resolution crystal structure of the mutated AidH in complex with substrate, a combined quantum mechanics/molecular mechanics (QM/MM) approach has been employed to study the detailed catalytic mechanism of AidH using C6-HSL (N-hexanoyl homoserine lactone) as the substrate. The calculation results reveal that the catalytic reaction starts from the abstraction of the hydroxyl hydrogen of Ser102 by His248. Both the formation and cleavage of the covalent bond between Ser102 and the substrate are possible rate-limiting steps, corresponding to energy barriers of 19.2 and 21.7 kcal mol⁻¹, respectively. The ring-opening of the covalent intermediate is calculated to be quite easy. During the catalysis, His248 acts as a dual Lewis acid/base, whereas Glu219 is not directly involved in the chemical process and Tyr160 only plays a role in stabilizing the transition state and orienting the position of the hydrolytic water molecule. In addition, the distant surrounding residues were found to have different influences on the reaction by their electrostatic interaction with the substrate. These results may provide useful information for the novel treatment of plant and animal infections that rely on AHL signaling.

---

1. Introduction

Quorum sensing (QS) is a system of stimuli and responses correlated to population density.¹ In bacteria, the QS system is used to synthesize, secrete and detect small signal molecules in order to perceive the population density and regulate the expression of specific genes in response to the changing environment.^2,3 There are many signal molecules, such as the N-acyl-homoserine lactones (AHLs), p-coumarate, quinolone, oligopeptides in Gram-positive bacteria, and 3OH palmitic acid methyl ester (3OH PAME) in several Gram-negative bacteria.⁴ Since QS processes are widely used by bacteria, the disruption of QS represents a potential strategy to block infection of hosts by pathogens.^5–8 Unlike traditional antibiotics that lead to the emergence of antibiotic-resistant strains of bacteria to limit cell growth or kill the pathogen, a quorum quenching enzyme aims to shut down the virulence expression in pathogenic bacteria.⁹ It should be noted that the enzymatic protection is completely harmless to the environment and its application will lead to the reduction of the use of chemicals.⁹

One of the best characterized types of quorum sensing is the AHL-mediated cell-to-cell communication.^2,10 Structurally, all AHLs share the homoserine lactone ring but differ in length and the C3 substitution of the acyl side chain.¹¹ AHLs are used by various Gram negative bacteria, such as the animal pathogens Burkholderia¹² and Pseudomonas aeruginosa,¹³ plant pathogens Erwinia carotovora,¹⁴ which cause soft rot in a wide range of economically important crops, including potatoes, beets and carrots. In the past few years, several means have been tried for disrupting the QS process, and the enzymatic degradation of QS signal molecules (AHLs) has been proved to be most appreciative and applicable.¹⁵ Up to now, several AHL-degrading enzymes have been discovered in a number of bacteria.¹⁶ Based on their catalytic sites, these enzymes can be divided into three families: AHL-acylases, AHL-oxidases and AHL-lactonases. AHL-acylases inactivate AHLs by hydrolyzing the amide linkage and separate the AHL into homoserine lactone and fatty acid.¹⁷ AHL-oxidases degrade AHLs by catalyzing the oxidation at the ω-1, ω-2, and ω-3 carbons of the acyl chain of AHL.¹⁸ AHL-lactonases degrade AHLs by hydrolyzing the lactone bond.⁷ Recently, AHL-lactonase is a hot focus among these AHL-degrading enzymes. According to the crystal structure, AHL-lactonase can be divided into metallo-lactamase and α/β-hydrolase superfamilies.¹⁶ Structurally, members of metallo-lactamase superfamily are metal-dependent enzymes, which have two zinc ions bound at the conserved HXHXDH motif and both are located in the active site, for example, BTK-AiiA from Bacillus thuringiensis subsp. kurstaki and AiiB from A. tumefaciens.¹⁹ The catalytic reaction is considered to occur via a dual Lewis acid catalysis mechanism.²⁰ Unlike most AHL-lactonases that are members of metallo-lactamase superfamily, AidH from Ochrobactrum sp. strain is a novel AHL-lactonase,¹⁶ which does not contain the conserved HXHXDH motif and Zn²⁺ is not essential for its catalytic activity.

In 2013, the crystal structures of AidH, AidH–C4-HS complex, AidH_S102G–C6-HSL complex, AidH_E219G–C6-HS complex, AidH_S102G and AidH_E219G have been determined by Gao et al.,¹⁶ with PDB codes of 4G5X, 4G9E, 4G8B, 4G8C, 4G8D and 4G9G at resolutions of 1.29, 1.09, 1.30, 1.11, 1.35 and 1.35 Å, respectively. Based on these high-resolution crystal structures and biochemical analysis, a typical acid–base covalent catalytic mechanism has been proposed for AidH.¹⁶ In the catalytic reaction, the lactonase ring is cleaved and a water molecule is incorporated into the product, as shown in Scheme 1. AidH can degrade a number of AHLs with different alkyl substituent (R₁) and functional groups (R₂ = –OH, =O) at the C3 position. Their experimental results revealed that there was no obvious difference in catalytic activity towards different substrates.¹⁶

Scheme 1 Overall reaction catalyzed by AidH.

Although previous experiment has given a rough picture of the catalytic mechanism, open questions still remain, which cannot be solved by experiment alone, for example, the energetics of the whole catalytic reaction, the structures of transition states and intermediates, the detail of each elementary step and the roles of some key residues, such as Tyr160, which was suggested to play an important role in substrate binding and catalysis by the point mutation.¹⁶ These questions inspire us to further investigate the catalytic mechanism of AidH. In this paper, we use a combined quantum mechanical/molecular mechanical (QM/MM) method to study the catalytic mechanism of AidH. This method has been successfully applied in exploring many enzymatic reactions.^{21–25,39–42} To the best of our knowledge, this is the first theoretical report regarding the catalytic reaction of AidH, which may provide important guides for the development of novel therapeutics against quorum sensing.

2. Computational methods

2.1 Setup of the system

The initial structure of the enzyme–substrate complex was prepared on the basis of the recently determined crystal structure of AidH_S102G in complex with C6-HSL (N-hexanoyl homoserine lactone) at 1.30 Å resolution (PDB code: 4G8B).¹⁶ Since AidHS_102G is not a wild type enzyme but a single mutant, Gly102 was firstly mutated to serine with the VMD software.²⁶ The protonation states of titratable residues were determined by H++ program²⁷ together with visual inspections. All histidine residues were singly protonated. For example, His6 and His275 were protonated at the ε position, and His12, His31, His62, His106, His248, His274, His275, His276, His277 and His278 were protonated at the δ position. All missing hydrogen atoms of the crystal structure were added by the HBUILD in CHARMM program package.²⁸ Then the system was solvated into a water sphere of radius 32 Å while all the crystal waters were kept at their original positions. Furthermore, the obtained system was neutralized by adding Na⁺ ions at random positions. After the resulted system was relaxed by a series of energy minimizations, a 20 ns molecular dynamic (MD) simulation was performed with the CHARMM22 force field²⁹ and stochastic boundary condition³⁰ at 298 K and 1 atm. From the calculated root-mean-square deviation (RMSD) of the protein structure, one can see that the system reached a basic plateau after 10 ns (ESI Fig. S1†).

2.2 QM/MM methodology

The QM/MM calculations follow the general procedure as in the review article.⁴² After the enzyme–substrate system was constructed and equilibrated using the CHARMM22 force field,²⁹ the system was divided into two regions: QM and MM. The QM region contains the substrate C6-HSL, two water molecules labeled W1 and W2, the side chains of Ser102, Tyr160, Glu219 and His248, and the backbones of Gly32 and Leu103 (Fig. S2†). All other atoms were defined as MM region. The QM atoms and MM atoms within 15 Å around the carbonyl C1′ atom of the substrate ring were set free, while all the remaining atoms of MM region were kept frozen. All QM/MM calculations were carried out employing ChemShell package,³¹ which combines Turbomole³² for the QM region and DL_POLY³³ for the MM region. An electronic embedding scheme³⁴ was used to include the polarizing effect of the enzymatic environment on the QM region. Hydrogen link atoms with the charge shift scheme were applied to treat the QM/MM boundary.³⁵ B3LYP function with 6-31G (d,p) basis set was applied for QM atoms, while CHARMM22 force field²⁹ was applied for MM atoms. Geometry optimizations were carried out by the hybrid delocalized internal coordinate (HDLC) optimizer.³⁶ The transition state was primarily determined as the highest point on the scanned potential energy surface (PES) along the reaction coordinates, and further optimized by the partitioned rational function optimization (P-RFO) algorithm³⁷ implemented in the HDLC code.³⁶ The last structure from the PES scan was optimized to find the intermediate, which was employed to start another scan to find the consecutive transition state and intermediate on the reaction path. Followed by a high level single point energy calculation at B3LYP function with 6-311++G (2d,2p) basis set, a dispersion correction using DFT-D3 program³⁸ was performed to obtained accurate value of energy. All energies are reported at this level.

3. Results and discussion

3.1 Structure of the reactant complex

Since the crystal structure of the wild-type AidH in complex with its substrate is not available, we used the crystal structure of AidH_S102G–C6-HSL complex to construct the computational model. The final snapshot obtained from the MD simulation was chosen as the initial structure for the following QM/MM calculations.

The active site structure derived from QM/MM optimization is displayed in Fig. 1B. One can see that Tyr160 does not form hydrogen bond with the substrate. Instead, a water molecule (W1) forms a hydrogen bond to the substrate with a distance of 1.88 Å. Other interactions are similar as in the crystal. For example, the backbone NH groups of Gly32 and Leu103 form two hydrogen bonds with the carbonyl O atom of substrate with distances of 1.96 Å and 2.31 Å, respectively. The amino group of Asn33 forms a hydrogen bond to the side chain carbonyl O atom of substrate with a distance of 1.93 Å, and the hydroxyl H atom of Ser102 forms a hydrogen bond with the N^ε of His248 with length of 1.81 Å. This large hydrogen bonding network makes the substrate in a right orientation for the hydrolytic reaction.

Fig. 1 (A) Crystal structure of the active site of mutant S102G in complex with the substrate C6-HSL (PDB code: 4G8B), in which all the hydrogen atoms are hidden; (B) structure of the active site after QM/MM calculations.

3.2 Reaction mechanism

On the basis of our calculations, the overall degradation process of C6-HSL contains five elementary steps, as shown in Fig. 2. Unlike the previous proposal,¹⁶ the ring opening is a one-step reaction (from IM1 to IM2). The reaction details will be discussed in the following. The first step is the covalent binding of substrate with residue Ser102. The optimized structures of transition states and intermediates are displayed in Fig. 3. This reaction is initiated by a His248-mediated nucleophilic attack, in which His248 abstracts a proton from the hydroxyl of Ser102 and simultaneously the generated hydroxyl anion acts as a nucleophile to attack on carbonyl C1′ atom of substrate, forming the intermediate IM1. In TS1, the hydroxyl O–H bond is weakened to 1.35 Å and the distance between the hydroxyl O atom of Ser102 and carbonyl C1′ atom decreases from 2.71 Å to 1.71 Å, which suggests that the proton transfer and the nucleophilic attack are concerted. According to our calculations, this step corresponds to an energy barrier of 19.2 kcal mol⁻¹ (Fig. 4). The relative energy of IM1 is higher than R by 14.1 kcal mol⁻¹, indicating the instability of IM1.

Fig. 2 Catalytic mechanism of AidH based on QM/MM calculations and previous study.¹⁶

Fig. 3 Optimized structures of the species involved in the first two steps of the reaction catalyzed by wild-type enzyme. Distances are given in angstroms.

Fig. 4 Energy profiles of the catalytic reaction catalyzed by the wild-type enzyme and Y160A mutant. Energies are given in kcal mol⁻¹.

The second step is the ring opening of the substrate, which is accompanied by a proton transfer from the doubly protonated His248 to the O3′ atom of the substrate. From IM1 to IM2, the C1′–O3′ bond is cleaved with its length elongated from 1.51 Å to 2.60 Å via 1.58 Å in TS2. By examining the structure of TS2 one can see that TS2 is very early transition state, demonstrating the instability of IM1. We also note that the proton transfer between His248 and substrate is a result of ring-opening process. The energy barrier of this ring-opening step is only 0.4 kcal mol⁻¹, however, IM2 is lower than R by 7.1 kcal mol⁻¹, indicating IM2 to be a stable covalent intermediate.

The next two steps involve the collapse of the covalent intermediate IM2. In the subsequent reaction, IM2 is hydrolyzed by a water molecule (W1) to yield the final product. As shown in Fig. 2, the water molecule (W1) should first adjust its position to facilitate the following hydrolytic reaction. Generally, the movement of a free small molecule is not difficult, we therefore did not perform MD simulate on this process, but manually changed the position of W1, so that W1 was located between the C1′ atom and the His248 N^ε atom. The QM/MM optimized structure is named as IM3, which is shown in Fig. 5. Comparison of the single point energies of IM2 and IM3 suggests that IM3 is only higher than IM2 by 4.3 kcal mol⁻¹, which implies that the slight movement of W1 in the active site is quite easy. In IM3, a strong hydrogen bond is formed between W1 and His248 with a distance of 1.81 Å. In addition, the distance between the C1′ atom of substrate and the O atom W1 is only 3.10 Å. Our calculation results reveal that the activation of W1 and the attack of generated OH⁻ on the carbonyl C1′ are concerted, corresponding to an energy barrier of 21.7 kcal mol⁻¹ (Fig. 4). In TS4, the N–H_w length changes from 1.81 to 1.05 Å, whereas the C1′–O_w distance changes from 3.10 to 1.74 Å, which means the activation of W1 by His248 is earlier than the attack of generated OH⁻.

Fig. 5 Optimized structures of species involved in the last two steps of the reaction catalyzed by wild-type enzyme. Distances are given in angstroms.

The final step corresponds to the cleavage of C′–O bond between the substrate and Ser102, generating the final product N-hexanoyl homoserine (C6-HS) and recovering the enzyme protein. In this process, His248 again plays important role, which acts as a general acid to protonate the leaving group. In TS5, the C′–O bond distance between the substrate and Ser102 increases from 1.49 to 1.68 Å, and the N–H bond length at ε position of His248 changes to 1.33 Å, implying that the C–O bond cleavage and the proton transfer are also concerted. The energy barrier of this step is 3.9 kcal mol⁻¹.

During our calculations, we also note the role of Glu219. Fig. 3 and 5 show that the changes of H–N bond at δ position of His248 are within 0.05 Å, and the distance between the H^δ of His248 and the O atom of Glu219 keeps longer than 1.50 Å (Fig. S3†), which indicate that Glu219 does not act as Lewis acid/base, but plays a role by forming a strong hydrogen bond with His248.

3.3 Y160A mutant

Previous point mutation experiment suggested that Tyr160 plays a role in substrate binding and catalysis.¹⁶ Therefore, we further performed calculations using Y160A mutant. We built the computational model on the basis of the optimized structure of the reactant (R, in Fig. 3) by simply mutating the Tyr160 to Ala160. The optimized geometries are shown in Fig. S4 and S5,† and the corresponding energy profile is shown in Fig. 4. Fig. 4 shows that the energy profile of Y160A mutant is similar to the wild enzyme except a high energy barrier of 29.1 kcal mol⁻¹ in step 4, which is 7.4 kcal mol⁻¹ higher than that of the wild type enzyme. By comparing the structures in Fig. 5 and S5,† we find that, in the wild type enzyme, Tyr160 forms hydrogen bonding interaction to the substrate and hydrolytic water (W1), which plays role in stabilizing the transition state (TS4) and fixing the position of hydrolytic water, whereas in Y160A mutant, there is no hydrogen bonding interaction in this area. Our calculation results are consistent with the previous experimental findings that replacement of Tyr160 has a strong impact on the catalytic efficiency of AidH.

3.4 Residue electrostatic analysis

In our QM/MM calculations, only the essential residues were included in the QM region. The far surrounding residues in MM region may also have significant influence on the catalytic reaction. Considering the relative energies of the three transition states TS1, TS4 and TS5 are higher than other species, we further performed the electrostatic analysis of eleven surrounding residues toward steps 1, 4 and 5. The difference of energy barrier caused by residue i can be described as: ΔE^i–0 = ΔEⁱ − ΔE⁰, where ΔE^i–0 is the change of the energy barrier, ΔEⁱ is the activation energy with charges on the atoms of residue i were set to 0, and ΔE⁰ is the original value of energy barrier. The structures of related species were kept unchanged during the calculations. The negative ΔE^i–0 value means residue i lowers the energy barrier and promotes the catalytic reaction. In contrast, the positive ΔE^i–0 value means residue i stabilizes the reactant or intermediate more than the transition state. The calculated values of ΔE^i–0 are shown in Fig. 6. Among these residues, His106 increases the energy barrier by 7.1 kcal mol⁻¹ in the fourth step. Although Asn33 decreases the energy barrier of the first step about 3.5 kcal mol⁻¹, it increases the energy barrier of the fifth step about 2.7 kcal mol⁻¹. Thus Asn33 has a little effect on the whole catalytic reaction, which is in agreement with the experimental result. The remaining nine residues have minor influence towards these three steps (−2.5 kcal mol⁻¹ < ΔE^i–0 < 2.5 kcal mol⁻¹). Our electrostatic influence analysis suggests His106 to be a disfavored residue for the catalytic reaction, which may be tested by future mutant study.

Fig. 6 Electrostatic influences of eleven residues on the step 1, step 4 and step 5.

4. Conclusion

In this work, we reported the first theoretical study of the degradation mechanism of AidH from Ochrobactrum sp. strain.¹⁶ According to our results, the overall catalytic process contains five elementary steps. The catalytic reaction starts from the abstraction of the hydroxyl hydrogen of Ser102 by His248, accompanies which the deprotonated Ser102 forms covalent intermediate (IM1) with the substrate. This concerted step corresponds to the energy barrier of 19.2 kcal mol⁻¹. The subsequent ring-opening of IM1 (step 2) is calculated to be quite easy with an energy barrier of only 0.4 kcal mol⁻¹. To cleave the covalent bond between Ser102 and substrate, the hydrolytic water (W1) should be firstly activated. Our calculations reveal that His248 performs this task, and the activation of W1 and the attack of generated OH⁻ on the covalent intermediate are also concerted, corresponding to an energy barrier of 21.7 kcal mol⁻¹. During the two steps involving the formation and breaking of covalent bond of Ser102 with the substrate, His248 acts as a dual Lewis base, which twice abstracts the hydroxyl hydrogen of Ser102 and the hydrolytic water (W1). Additionally, it acts as a dual Lewis acid to protonate the leaving groups (–CH₂–O⁻ or –COO⁻). Glu219 is not directly involved in the chemical process but assists the catalysis by forming hydrogen bonding contact with His248. Tyr160 mainly plays role in stabilizing the transition state and orienting the hydrolytic water (W1) in suitable position. Furthermore, analysis of electrostatic effect reveals that the distant His106 may play important role in accelerating the catalytic efficiency of AidH. Our studies may provide useful information for understanding the degradation mechanism of AHLs and developing novel treatments for plant and animal infections that relay on AHL signaling.

Acknowledgements

This work was supported by the National Natural Science Foundation of China (21373125, 21573127).

References

1. V. C. Kalia, Biotechnol. Adv., 2013, 31, 224–245.
2. A. Camilli and B. L. Bassler, Science, 2006, 311, 1113–1116.
3. M. Teplitski, U. Mathesius and K. P. Rumbaugh, Chem. Rev., 2010, 111, 100–116.
4. C. M. Waters and B. L. Bassler, Annu. Rev. Cell Dev. Biol., 2005, 21, 319–346.
5. Y. H. Dong, L. H. Wang, J. L. Xu, H. B. Zhang, X. F. Zhang and L. H. Zhang, Nature, 2001, 411, 813–817.
6. G. Y. Mei, X. X. Yan, A. Turak, Z. Q. Luo and L. Q. Zhang, Appl. Environ. Microbiol., 2010, 76, 4933–4942.
7. Y. H. Dong, J. L. Xu, X. Z. Li and L. H. Zhang, Proc. Natl. Acad. Sci. U. S. A., 2000, 97, 3526–3531.
8. J. Momb, P. W. Thomas, R. M. Breece, D. L. Tierney and W. Fast, Biochemistry, 2006, 45, 13385–13393.
9. F. Chen, Y. X. Gao, X. Y. Chen, Z. M. Yu and X. Z. Li, Int. J. Mol. Sci., 2013, 14, 17477–17500.
10. B. LaSarre and M. J. Federle, Microbiol. Mol. Biol. Rev., 2013, 77, 73–111.
11. S. Schauder and B. L. Bassler, Genes Dev., 2001, 15, 1468–1480.
12. R. L. Ulrich, Appl. Environ. Microbiol., 2004, 70, 6173–6180.
13. L. Passador, J. M. Cook, M. J. Gambello, L. Rust and B. H. Iglewski, Science, 1993, 260, 1127–1130.
14. M. Pirhonen, D. Flego, R. Heikinheimo and E. T. Palva, EMBO J., 1993, 12, 2467–2476.
15. V. C. Kalia and H. J. Purohit, Crit. Rev. Microbiol., 2011, 37, 121–140.
16. A. Gao, G. Y. Mei, S. Liu, P. Wang, Q. Tang, Y. P. Liu, H. Wen, X. M. An, L. Q. Zhang, X. X. Yan and D. C. Liang, Acta Crystallogr., Sect. D: Biol. Crystallogr., 2012, 69, 82–91.
17. J. J. Huang, J. I. Han, L. H. Zhang and J. R. Leadbetter, Appl. Environ. Microbiol., 2003, 69, 5941–5949.
18. P. K. Chowdhary, N. Keshavan, H. Q. Nguyen, J. A. Peterson, J. E. González and D. C. Haines, Biochemistry, 2007, 46, 14429–14437.
19. D. Liu, P. W. Thomas, J. Momb, Q. Q. Hoang, G. A. Petsko, D. Ringe and W. Fast, Biochemistry, 2007, 46, 11789–11799.
20. D. Liu, B. W. Lepore, G. A. Petsko, P. W. Thomas, E. M. Stone, W. Fast and D. Ringe, Proc. Natl. Acad. Sci. U. S. A., 2005, 102, 11882–11887.
21. J. Momb, C. Wang, D. Liu, P. W. Thomas, G. A. Petsko, H. Guo, D. Ringe and W. Fast, Biochemistry, 2008, 47, 7715–7725.
22. R. Z. Liao, J. G. Yu and F. Himo, Inorg. Chem., 2009, 48, 1442–1448.
23. G. C. Ma, W. Y. Zhu, H. Su, N. Cheng and Y. J. Liu, ACS Catal., 2015, 5, 5556–5566.
24. K. E. Ranaghan and A. J. Mulholland, Int. Rev. Phys. Chem., 2010, 29, 65–133.
25. X. Dong, Z. Tian, X. Yang and Y. Xue, Theor. Chem. Acc., 2015, 134, 1–10.
26. W. Humphrey, A. Dalke and K. Schulten, J. Mol. Graphics, 1996, 14, 33–38.
27. J. C. Gordon, J. B. Myers, T. Folta, V. Shoja, L. S. Heath and A. Onufriev, Nucleic Acids Res., 2005, 33, 368–371.
28. B. R. Brooks, R. E. Bruccoleri, B. D. Olafson, D. J. States, S. Swaminathan and M. Karplus, J. Comput. Chem., 1983, 4, 187–217.
29. A. D. MacKerell Jr, D. Bashford, M. Bellott, R. L. Dunbrack Jr, J. D. Evanseck, M. J. Field, S. Fischer, J. Gao, H. Guo, S. Ha, D. Joseph-McCarthy, L. Kuchnir, K. Kuczera, F. T. K. Lau, C. Mattos, S. Michnick, T. Ngo, D. T. Nguyen, B. Prodhom, W. E. Reiher III, B. Roux, M. Schlenkrich, J. C. Smith, R. Stote, J. Straub, M. Watanabe, J. Wiórkiewicz-Kuczera, D. Yin and M. Karplus, J. Phys. Chem. B, 1998, 102, 3586–3616.
30. M. Berkowitz and J. A. McCammon, Chem. Phys. Lett., 1982, 90, 215–217.
31. P. Sherwood, A. H. de Vries, M. F. Guest, G. Schreckenbach, C. R. A. Catlow, S. A. French, A. A. Sokol, S. T. Bromley, W. Thiel, A. J. Turner, S. Billeter, F. Terstegen, S. Thiel, J. Kendrick, S. C. Rogers, J. Casci, M. Watson, F. King, E. Karlsen, M. Sjøvoll, A. Fahmi, A. Schäfer and C. Lennartz, J. Mol. Struct., 2003, 632, 1–28.
32. R. Ahlrichs, M. Bär, M. Häser, H. Horn and C. Kölmel, Chem. Phys. Lett., 1989, 162, 165–169.
33. W. Smith and T. R. Forester, J. Mol. Graphics Modell., 1996, 14, 136–141.
34. S. R. Billeter, A. J. Turner and W. Thiel, Phys. Chem. Chem. Phys., 2000, 2, 2177–2186.
35. A. H. de Vries, P. Sherwood, S. J. Collins, A. M. Rigby, M. Rigutto and G. J. Kramer, J. Phys. Chem. B, 1999, 103, 6133–6141.
36. B. R. Brooks, R. E. Bruccoleri, B. D. Olafson, D. J. States, S. Swaminathan and M. Karplus, J. Comput. Chem., 1983, 4, 187–217.
37. A. Banerjee, N. Adams, J. Simons and R. Shepard, J. Phys. Chem., 1985, 89, 52–57.
38. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, J. Chem. Phys., 2010, 132, 154104.
39. W. Y. Zhu and Y. J. Liu, ACS Catal., 2015, 5, 3953–3965.
40. H. Monajemi, S. M. Zain and W. A. T. W. Abdullah, RSC Adv., 2015, 5, 25489–25503.
41. S. T. Keaveney, J. B. Harper and A. K. Croft, RSC Adv., 2015, 5, 35709–35729.
42. M. G. Quesne, T. Borowski and S. P. de Visser, Chem.–Eur. J., 2016, 22, 2562–2581.

---

Footnote

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

---