Chemical shifts of p-methylstyrene (δ in ppm) (cf. Figure 2).
Applications of density functional theory (DFT) calculations to organic chemistry are shown, beginning with geometry optimization and the calculation of vibrational frequencies, infrared (IR) intensities, and thermodynamic properties. The isotropic chemical shielding values and anisotropies relevant to nuclear magnetic resonance (NMR) can be calculated using gauge-invariant atomic orbitals (GIAOs); the calculation of spin-spin couplings is possible but time-consuming. For free radicals, hyperfine couplings and g tensors pertaining to EPR can be obtained. Regarding UV/vis spectra, wavelengths and oscillator strengths can be calculated by using a time-dependent Hamiltonian. In addition to gas-phase acidities, approximate pKa values can be obtained, provided that solvation is taken into account. Several sets of substituent parameters have been calculated: Hammett σ and σ+ parameters and inductive and mesomeric effects. Regarding reaction mechanisms, geometries and energies of intermediates and transition structures have been calculated for pericyclic reactions, nucleophilic aliphatic substitutions, electrophilic aromatic substitutions, additions, and eliminations.
- density functional theory
- magnetic resonance
- Hammett parameters
- reaction mechanisms
- pericyclic reactions
Focusing on density functional theory (DFT) calculations with Gaussian 09  and the B3LYP/6-311G(d,p) method, several applications to organic chemistry will be shown. After geometry optimization, which yields the total energy, a frequency calculation can be done, yielding the infrared spectrum (wave numbers and intensities) and, if requested, the Raman intensities and the thermodynamic properties (enthalpy, entropy, and Gibbs free energy).
Using a time-dependent Hamiltonian, UV/vis spectra can be calculated (wave lengths and oscillator strengths). Nuclear magnetic resonance (NMR) spectra can be calculated, providing isotropic shielding values as well as tensor data (anisotropies) of all magnetic nuclei, using gauge-invariant atomic orbitals (GIAOs). The calculation of spin-spin coupling constants is also possible but requires much more computational time. For free radicals, EPR data can be calculated: isotropic hyperfine coupling constants, hyperfine tensors, and g tensors; in this case, the restricted B3LYP method has to be replaced by the unrestricted UB3LYP method.
Substituent effects such as the σ parameters in the Hammett equation can also be estimated. Although a calculation of changes in the charge distribution might seem to be a promising method for that purpose, it was found that calculated 19F shielding values (“virtual NMR experiments”) yielded much more convincing results. The calculation of gas-phase acidities or basicities is straightforward, and the calculated data show a good correlation with experimental data. However, the correlation with pKa values, which refer to aqueous solutions, is very poor. A reasonable correlation was obtained by taking a few water molecules explicitly into account, in addition to the bulk solvent properties of water.
Regarding organic reaction mechanisms, pericyclic reactions are particularly well amenable to DFT calculations. Usually, the transition structure can be obtained which is characterized by a single imaginary frequency, which belongs to the reaction coordinate. For many other reaction types (substitutions, additions, eliminations, and rearrangements), at least an approximation to the transition structure can be calculated. Moreover, starting with such a structure and performing an optimization, the approximate dynamics of the reaction can be followed.
2. Geometries, energies, and thermodynamic data
2.1. Geometry optimization
As a starting point, a reasonable approximation to the geometry of the target molecule is required. Preferably, the coordinate file should be given as Z matrix, and standard bond lengths and angles may be used. A convenient tool for the generation of Z matrices is molden : in the Z-mat editor, start with methane, substitute by phenyl and finally by vinyl, and save as Z matrix (GAMESS). Next, the input file for the quantum-chemical calculation has to be created by supplementing the Z matrix file with the necessary parameters (see Appendix A).
After a successful calculation, the log file contains the energy (in Hartree) and the coordinates of the optimized structure. Again, it is advantageous to use a tool such as molden for analyzing the log file.
2.2. Calculation of thermodynamic properties
For a determination of the thermodynamic properties, it is necessary to calculate the (vibrational) frequencies. In the Gaussian input file, the preliminary coordinates have to be replaced by the optimized ones and the task “Opt” by “Freq”. (Actually, the request for “Freq Prop Pop = Full” additionally provides useful information such as charges and dipole moment. By default, the calculation is done for 298 K and 1.000 atm, but a different temperature or pressure may be specified.)
For the example molecule (p-methylstyrene), the salient results are as follows:
|SCF Done: E(RB3LYP)||−349.054882687 AU|
|Zero-point correction||0.159885 (Hartree/particle)|
|Thermal correction to energy||0.167687|
|Thermal correction to enthalpy||0.168631|
|Thermal correction to Gibbs free energy||0.127033|
|Sum of electronic and zero-point energies||−348.894997|
|Sum of electronic and thermal energies||−348.887196|
|Sum of electronic and thermal enthalpies||−348.886252|
|Sum of electronic and thermal free energies||−348.927850|
|Total E (thermal)||105.225 kcal/mol|
|Total CV||30.240 cal/mol-K|
|Total S (entropy)||87.551 cal/mol-K|
Most data are given in Hartree (see Appendix A), they refer to the formation from atomic nuclei and electrons. It is fairly easy to calculate the energy of formation from the atoms by subtracting the energies obtained for respective calculations of free atoms. In order to obtain approximate values for standard enthalpies of formation, bond energies and possibly enthalpies of phase changes (to the gas phase) have to be taken into account. It should be mentioned that the accuracy of these data, i.e., the agreement with experimental data, is not very good. It is advisable to restrain to energy (or enthalpy) differences of similar structures. Alternatively, approximate enthalpies of formation can be obtained more easily from semiempirical calculations (such as MNDO, AM1, or PM3).
Energies of some important free atoms (UB3LYP/6-311G(d,p) in Hartree): H, −0.502155930031; C, −37.8559889346; N, −54.5985431427; O, −75.0853856058; F, −99.7538096003; P, −341.280503655; S, −398.132082447; and Cl, −460.166160487.
Hence, the following energy of formation from the atoms is obtained for p-methylstyrene, C9H10, ΔEf atomic = −348.887196−9 × (−37.8559889346)−10 × (−0.502155930031) = −3.161736288290 Hartree = −8301.139 kJ/mol.
Enthalpies required to generate free atoms from the elements in the standard state (kJ/mol) : H, 218.00; C, 716.67; N, 472.68; O, 249.17; F, 78.4; P, 314.55; S, 276.98; and Cl, 121.29.
These values have to be added to the above-given atomic energy of formation (ignoring somewhat the difference between energy and enthalpy), yielding the following energy of formation for our example: −8301.14 + 9 × 716.67 + 10 × 218.00 = 328.89 kJ/mol.
The energy can be converted to the enthalpy by means of Eq. (1), assuming the validity of the ideal gas law; Δn is the change in the number of moles of gases:
In the example, 1 mole of product molecules (in the gas phase) is formed from 19 moles of atoms; therefore, Δn is −18 and the correction is −44.62 kJ/mol. Hence, the calculated enthalpy of formation of p-methylstyrene in the gas phase is ΔHf° = 284.27 kJ/mol.
The following experimental value for the enthalpy of formation of liquid p-methylstyrene is given in the literature: ΔHf°(l) = 114.6 kJ/mol ; adding the heat of vaporization of 47.6 kJ/mol , ΔHf°(g) = 162.2 kJ/mol for the gas phase. Thus, the calculated value deviates by about 120 kJ/mol.
By comparison, a semiempirical AM1 calculation yields an enthalpy of formation of ΔHf° = 140.3 kJ/mol, in better agreement with experiment.
3.1. Vibrational spectroscopy: infrared and Raman
Vibrational frequencies and hence infrared (IR) and Raman spectra can be calculated (Gaussian keyword “Freq”). In Gaussian 09, the infrared intensities are calculated by default, but the Raman intensities can also be obtained (keyword “Freq = Raman”). The calculated frequencies can be assigned to the respective molecular motions. The visualization of vibrations is easily achieved by tools such as molden.
As an example, p-cyanobenzaldehyde will be considered (Figure 1). The experimental IR data have been taken from the SDBS database . The two most prominent features are the C═O valence vibration at 1788 (exp. 1708) and the CN valence vibration at 2340 (2230) cm−1.
3.2. Nuclear magnetic resonance (NMR)
The calculation of isotropic chemical shielding values (σ) or shielding tensors requires some kind of “scaling” of the orbitals, for instance, the use of gauge-invariant atomic orbitals (GIAOs) . In Gaussian 09, the keyword “NMR” automatically invokes the use of GIAOs, and isotropic shielding values, anisotropies, and shielding tensors are calculated. In NMR experiments, however, not the shielding values are measured, but chemical shifts, which refer to some standard. For 1H, 13C, and 29Si NMR, tetramethylsilane (TMS) is used as a standard, and the respective chemical shifts of the three kinds of magnetic nuclei are set to zero (δH-1(TMS) = 0 ppm, δC-13(TMS) = 0 ppm, and δSi-29(TMS) = 0 ppm). Using the hybrid method B3LYP/6-311(d,p), the following average shielding values are obtained for protons and for 13C nuclei: σ(1H) = 31.3919 ppm and σ(13C) = 179.7024 ppm. The chemical shifts are then simply obtained by subtraction:
Using the abovementioned reference value for protons, the calculated chemical shifts are generally too small by about 0.5 ppm. In a survey of 21 natural products, a better fit for 13C nuclei, on the average, was obtained by using a reference value of 177.0 ppm instead .
For the example molecule p-methylstyrene, the following 1H and 13C chemical shifts (δ in ppm) were calculated (using the above-given calculated shielding values for TMS as reference); experimental 13C chemical shifts were taken from the NMRSHIFTDB database (Figure 2 and Table 1) .
|Pos.||C calc||C exp||H calc||H pred|
It is also possible to calculate NMR spin-spin coupling constants, i.e., J [Hz] (Gaussian keyword “NMR = SpinSpin”), but at the expense of computational time. The results obtained with the B3LYP hybrid functional are much better than those of HF ab initio calculations.
For p-methylstyrene, the following proton-proton spin-spin coupling constants J [Hz] were calculated: J13,14 = 6.86, J15,16 = 7.09 (o), J14,15 = 1.29, J13,16 = 1.26 (m), J13,15 = 0.41, J14,16 = 0.28 (p), J17,19 = 14.72 (trans), J17,18 = 10.11 (cis), and J18,19 = −0.53 (gem). These values are in accordance with those found in similar systems.
It should be mentioned that NMR data provide an excellent and sensitive test for the accuracy of quantum-chemical calculations.
3.3. Electron paramagnetic resonance (EPR)
In the case of free radicals, unrestricted calculations have to be performed in which different orbitals are assigned to α and β spins. Whereas unrestricted Hartree-Fock (UHF) calculations yield poor results for hyperfine couplings (HFC) because of serious problems due to spin contamination, calculations with the UB3LYP hybrid functional yield fairly acceptable results . The calculations yield Mulliken spin densities (better designated as spin populations), isotropic HFC (Fermi contact coupling constants), and anisotropic hyperfine tensors. g values and g tensors can also be calculated (in the Gaussian system, this requires the “NMR” keyword). The g value is a dimensionless proportionality factor which relates the magnetic moment to the angular momentum; the value for the free electron is ge = 2.00232, and only the electron spin is involved. In molecules, contributions from orbital momentum have to be taken into account, and the phenomenon becomes anisotropic. The calculated HFC and g values may be compared with experimental data from EPR spectroscopy (electron spin resonance, also called electron paramagnetic resonance) .
The method will be illustrated using the ubisemiquinone-Q1 radical anion as example, which serves as a model compound for coenzyme Q10.
Figure 3 shows the calculated Mulliken spin densities and the calculated proton HFC of this radical anion. The rotation of the long side chain is hindered; therefore, the two methylene protons are inequivalent. Comparison with experimental data (ethanol, 230 K, in parentheses) : 6.44 (5.84), methyl protons, and 3.56 (3.68) and 2.11 (2.17) MHz, methylene protons. The g tensor has been measured by high-field EPR experiments ; again, the experimental values are given in parentheses: gxx = 2.00826 (2.00646), gyy = 2.00601 (2.00542), gzz = 2.00207 (2.00222), and giso = 2.00545 (2.00470).
3.4. Electron spectroscopy (UV/vis)
In Hückel molecular orbital (HMO) theory, electronic excitation may be viewed as excitation of an electron from an occupied to an unoccupied orbital. The transition with the lowest energy, i.e., the longest wavelength, involves the excitation from the highest occupied molecular orbital (HOMO) to the lowest unoccupied molecular orbital (LUMO), although this transition might be forbidden.
In DFT, however, the Kohn-Sham orbitals are not suitable for this procedure, and a time-dependent Hamiltonian has to be used in the calculation (Gaussian keyword “TD”). The calculation gives the energies and the wavelengths of the excitations, the oscillator strengths f, and reports the orbitals which are involved. The vibrational fine structure and the consequences of the Franck-Condon principle are not taken into account.
In the case of the symmetrical crystal violet cation, the HOMO is represented by two degenerate orbitals, and two excitations have the same wavelength, calculated as 504.7 nm (f = 0.806); experimental data for the absorption maxima: 591.0 and 540.5 nm.
For further examples, see .
4. Substituent effects
4.1. Hammett σ parameters
The Hammett σ parameters refer to the acidities of substituted benzoic acids which will be considered in Section 4.3. Regarding electrophilic aromatic substitution (see Section 5.3), a modified set has to be used, at least for the para positions (σ+ parameters). The Hammett equation is
where ρ is the reaction parameter and k0 and ki are the rate constants for the unsubstituted and substituted compounds, respectively. The σ/σ+ parameters for electrophilic aromatic substitution have been determined from the relative stabilities of the σ complexes as averages for the following four reactions: protonation, bromination, nitration, and alkylation (by ethyl groups). A linear fit of σ/σ+ (literature data ) versus calculated σ (DFT) was determined for 17 substituents both in meta and in para positions, yielding a squared correlation coefficient of r2 = 0.932 (see Figure 4 and Table 2).
|Substituent||σ meta DFT||σ meta lit.||σ+ para DFT||σ+ para lit.|
4.2. Estimating inductive and mesomeric effects by virtual 19F NMR
The relative contributions of inductive (I) and mesomeric (M) effects might be inferred from a comparison of the Hammett σ/σ+ parameters (see Section 4.1) for the para (I + M) and meta (I + M/3) positions. Yet, a different approach is taken here. Using DFT, the obvious target to look for should be the charge density distribution. However, it turned out that the Mulliken charges, at least, did not yield satisfying results. Therefore, calculated isotropic 19F shielding values were used as a probe of local charge density. As a suitable system, 4-substituted (E,E)-1-fluoro-1,3-butadienes were chosen, in two conformations (Figure 5).
The geometries were optimized for the planar conformation, and the shielding values were calculated for this conformation (0°) and for the orthogonal conformation with a dihedral angle of 90° for the central single bond (see Figure 5). The relative shielding values σrel(19F) at 90° should be proportional to the inductive (I) effect, and the differences of the relative shielding values at 0° and at 90°, σrel(19F)0°−σrel(19F)90°, should be proportional to the mesomeric (M) effect. The reference compound is, of course, the unsubstituted compound (R = H). These data were calibrated against the Hammett σ/σ+ parameters, yielding a slope of −20.447, i.e., the data have to be divided by this factor. According to the sign convention of the Hammett σ/σ+ parameters, electron-withdrawing groups (EWG, −I, −M) have a positive sign (e.g., nitro and cyano), whereas electron-releasing groups (ERG, +I, +M) have a negative sign (e.g., alkyl groups) (see Table 3).
4.3. Acids and bases: pKa values
The calculation of gas-phase acidities is straightforward, but they do not correlate well with experimental pKa values . This situation is only partially improved by taking the bulk properties of the solvent (water) into account (Gaussian keyword “SCRF = (Solvent = Water)”). A much better approximation is obtained when additionally a few water molecules are taken into account explicitly, e.g., two water molecules in the case of carboxylic acids (see Figure 6).
Thus, Gibbs free energies for benzoic acid and a series of substituted benzoic acids (15 substituents both in meta and in para positions) as well as the respective anions were calculated. ΔG° = G°(anion) − G°(acid) was converted from Hartree to kJ/mol (see Appendix A), and log Ka was calculated according to
The relative pKa values were found to be quite reasonable but have to be scaled. Since the difference log Ka (substituted benzoic acid)-log Ka (benzoic acid) should be equal to the Hammett σ parameter, a linear fit of σ (literature data ) versus calculated Δ(log Ka) (DFT) was determined, yielding a slope of 0.3437 and a squared correlation coefficient of r2 = 0.967 (Figure 7). That means, the substituent effect on the calculated pKa value is overestimated by about a factor of 3.
For the calculation of absolute pKa values, the considerable Gibbs free solvation energy of the proton has to be taken into account (ΔG° = −1120.39 kJ/mol). The calculated Gibbs free energy of the dissociation of benzoic acid in the gas phase is ?G° = 1446.04 kJ/mol. For the hydration (dihydrate model, vide supra), ΔG° = −47.35 kJ/mol and ΔG° = −302.76 kJ/mol are obtained for benzoic acid and benzoate anion, respectively; in total, ΔG° = 1375.80 kJ/mol. Hence, the estimate for the Gibbs free energy of the dissociation of benzoic acid in aqueous solution is ΔG° = 1446.04 − 1375.80 = 70.24 kJ/mol corresponding to a pKa value of 12.31 (exp. 4.19). To put it differently, the model applied here accounts for about 96.7% of the true Gibbs free energy of solvation.
5. Reaction mechanisms
5.1. Pericyclic reactions
In a pericyclic reaction, σ or π bonds change concertedly (“simultaneously”) along a perimeter, i.e., a cycle. They have first been studied theoretically by Woodward and Hoffmann (“the conservation of orbital symmetry”) . Typical examples are sigmatropic reactions such as the Cope rearrangement, cycloadditions such as the Diels-Alder addition, or electrocyclic reactions (ring closures or openings).
Pericyclic reactions are particularly well amenable to DFT calculations; the transition structure can usually be obtained. The transition structure is a saddle point in the energy hyperspace, i.e., the energy has a maximum along the reaction path (the reaction coordinate), but is minimized with respect to all other coordinates. This can be checked by a frequency calculation. Exactly one frequency should be imaginary, namely, the one pertaining to the reaction coordinate. Thus, the reaction dynamics can be visualized by looking at that vibration.
5.1.1. Cope rearrangement
The Cope rearrangement is a [3,3] sigmatropic reaction. As an example, the degenerate Cope reaction of 1,5-hexadiene is shown (see Figures 8 and 9). The calculated activation energy is 129 kJ/mol (DFT).
5.1.2. Diels-Alder addition
The Diels-Alder addition is a [4 + 2] cycloaddition, a diene reacts with a dienophile to form a (substituted) cyclohexene. As an example, the Diels-Alder addition of acrylonitrile to cyclopentadiene leading to the endo-product is shown (see Figures 10 and 11). The calculated activation energy is 40 kJ/mol, and the calculated reaction energy is −132 kJ/mol (DFT).
5.1.3. Electrocyclic reactions
In an electrocyclic reaction, an unsaturated cycloalkane is formed from a conjugated polyene, or the reverse reaction occurs. Here, only thermally allowed electrocyclic reactions will be considered. For instance, cyclobutene is opened in a conrotatory manner to form 1,3-butadiene. (The calculated activation energy is 149 kJ/mol, and the calculated reaction energy is −39 kJ/mol, assuming that the most stable conformation of 1,3-butadiene is formed.) The example shown here is the disrotatory ring closure of 1,3,5-hexatriene to form 1,3-cyclohexadiene (see Figures 12 and 13). Starting with the most stable conformer of 1,3,5-hexatriene, the calculated activation energy is 252 kJ/mol, and the calculated reaction energy is −64 kJ/mol (DFT).
5.2. Nucleophilic aliphatic substitutions
The most important mechanisms for nucleophilic aliphatic substitutions are the single-step SN2 mechanism with backside attack of the nucleophile and a trigonal-bipyramidal transition state (for primary or secondary substrates) and the two-step SN1 mechanism with a carbenium ion intermediate (for secondary or tertiary substrates). Since these are ionic reactions, the progress in the gas phase may differ considerably from that in a polar solvent.
Considering first the degenerate SN2 reaction of fluoromethane with fluoride anion in the gas phase, the most stable species is a cluster of these two particles, which is formed in an exothermic reaction and calculated reaction energy −106 kJ/mol. The formation of the symmetric trigonal-bipyramidal transition state from this cluster requires an activation energy of 29 kJ/mol.
In the gas-phase reaction of chloromethane with fluoride anion (see Figures 14 and 15), the calculated reaction energy for the formation of fluoromethane and chloride anion at infinite distance is −198 kJ/mol. There is no activation energy for the forward reaction; the energy of the trigonal-bipyramidal transition state is lower by 19 kJ/mol than that of the cluster of chloromethane with fluoride. The most stable species is the cluster of fluoromethane with chloride anion, and the activation energy of the reverse reaction, starting with this cluster, would be 132 kJ/mol.
A typical example for an SN1 reaction is the reaction between tert-butanol and hydrogen chloride, yielding tert-butyl chloride and water (or the reverse reaction; see Figures 16 and 17). In the gas phase, this reaction is slightly exothermic with a calculated reaction energy of −8 kJ/mol (DFT). The concurrent elimination reaction (E1), yielding isobutene, water, and hydrogen chloride, is endothermic with a calculated reaction energy of 47 kJ/mol. For the reaction to proceed, it is necessary to form the protonated alcohol. Only the formation of an ion pair of chloride and protonated tert-butanol is conceivable, requiring an estimated activation energy of about 165 kJ/mol. The formation of a free tert-butyl carbenium ion is not feasible, because the energy required would be about 682 kJ/mol. The reaction should rather proceed by a backside attack similar to the SN2 reaction and stop at a cluster of tert-butyl chloride and water (reaction energy −19 kJ/mol).
A true SN1 mechanism requires an efficient solvation of the intermediate carbenium ion by polar solvent molecules.
5.3. Electrophilic aromatic substitutions
In gas-phase reactions of benzene with a reactive cationic electrophile such as H+, Br+, NO2+, or CH3CH2+ (cf. Section 4.1), the reaction proceeds without any energy barrier to the σ complex and stops there. Some kind of π complex is formed on the reaction path, but it is not a true intermediate because it is not characterized by a local energy minimum.
In a more realistic scenario, the electrophile is a less reactive complex, e.g., of a halogen, an alkyl, or an acyl chloride, with a Lewis acid such as aluminum chloride or iron(III) bromide. Now, the reaction will usually stop at the π complex stage. In order to force the reaction to proceed to the σ complex, a strongly activating substituent such as oxido (i.e., phenolate anion) was introduced. After removal or replacement of this substituent, the optimization procedure allowed the study of either the backward reaction or the forward reaction to the products, possibly after a modification of the arrangement of the reaction partners.
As an example, the chlorination of benzene catalyzed by aluminum chloride will be considered in detail (see Figures 18 and 19). The overall reaction in the gas phase, yielding chlorobenzene and hydrogen chloride, is exothermic with a calculated reaction energy of −131 kJ/mol and a Gibbs free reaction energy of ΔrG° = −143 kJ/mol (DFT). The reaction involves three intermediates, π complex 1, σ complex, and π complex 2, which were calculated as local minima and two transition structures, which could not yet be identified unambiguously. The crucial energy barrier is most likely the transition state leading to the σ complex. Taking the energies of the separated reactants as reference, the relative energies are as follows: −45 kJ/mol for π complex 1, approximately 208 kJ/mol for transition state 1, −80 kJ/mol for the σ complex, approximately 27 kJ/mol for transition state 2, −169 kJ/mol for π complex 2 (global minimum), and −131 kJ/mol for the separated products.
5.4. Additions and eliminations
5.4.1. Electrophilic addition
As an example, the addition of bromine to cyclohexene will be considered, yielding trans-1,2-dibromocyclohexane (see Figures 20 and 21). First, a π complex is formed (relative energy −22 kJ/mol with respect to the separated reactants) and, next, a bicyclic bromonium ion, which is more stable than the respective carbenium ion (by roughly 100 kJ/mol). Finally, the diaxial conformer of trans-1,2-dibromocyclohexane is formed, which is actually more stable than the diequatorial conformer by 7 kJ/mol. This finding is somewhat surprising and stands in contrast to previous assumptions. Apparently, electrostatic repulsion favors the diaxial form, whereas in monosubstituted cyclohexanes, the substituent prefers the equatorial position. The reaction energy is −108 kJ/mol, and the Gibbs free reaction energy is −50 kJ/mol; this is due to the unfavorable reaction entropy.
6. Conclusions and outlook
In the field of molecular chemistry, the use of DFT in combination with efficient software and modern computer equipment allows the development of “virtual chemistry,” i.e., the prediction of essentially all molecular properties and of reaction paths. To a certain extent, supramolecular chemistry is also accessible to this method; molecular clusters and microdroplets of solvents can be simulated. It stands to a reason, however, that computational time increases heavily with molecular size (or cluster size). In the case of ab initio calculations, the proportionality is to the fourth power of the size of the basis set; in DFT, the situation might be somewhat better; computational time is proportional to roughly the third power of the number of orbitals involved, judging from NMR calculations.
It should be pointed out that present-day DFT is only an approximate theory. Therefore, it is necessary to check the quality of the computational results against experimental data.
The following conversion factors have been used in this study: 1 Hartree (a.u.) = 2625.50 kJ/mol, 1 cal = 4.184 J, and pK = ΔG° [kJ/mol] /5.708008 (at T = 298.15 K).
Computational details. For all computations in Chapters 2 to 4, Gaussian 09 was used B3LYP/6–311(d,p) . For most computations in Chapter 5, deMon2k was used . For the preparation of input files, i.e., the generation of Z matrices, and the visualization of the results, molden  was of great help.
Sample Gaussian input file: methane
|%Chk = methane|
|h||1 hc2||2 hch3|
|h||1 hc2||3 hch3||2 dih4|
|h||1 hc2||4 hch3||3 dih5|
(The first line specifies the checkpoint file; the second the method and the task, in this case the geometry optimization; the fourth the title, the sixth the total charge, here 0; and the multiplicity, usually 1. Then, the Z matrix follows immediately.)