Computational Labs

Module 1

Introduction to Force Field Methods

Force field (FF) methods, also referred to as molecular mechanics, expresses the electronic energy as a parametric function of the nuclear coordinates. The parameters are fitted to experimental or high level computational data. FF methods focus purely on the atoms rather than individual electrons. In addition, nuclear motion present in the system is ignored.

FF methods describe molecules in terms of balls and springs. Different atoms have different sizes and softness and different bonds have different lengths and stiffness. Also, FF calculations enable the relative energies and barriers for the inter-conversion of different conformations to be made. This information is useful when trying to assess which product is the major/minor product in a given reaction. There are many different types of force fields that can be used to calculate properties of a molecule however in this experiment the Allinger’s MM2 method will be used.

The Force Field Energy is the sum of the following energy terms:

The stretch energy: the energy function for stretching a bond between two atoms
The bend energy: the energy required to bend an angle
the torsion energy: the energy required to rotate around a bond
Van de Waals energy and electrostatic charge (dipole-dipole): the non bonded atom-atom interactions
Cross terms: describes coupling between the first three ¹

Stable molecules typically relate to the minimum point on the potential energy surface diagram. Stable molecules correspond to minima on the potential energy surface. This point can be obtained by reducing the FF energy as a function of the nuclear co-ordinates.

Practical considerations

Force field methods are only models of the actual quantum mechanical systems. The success of such calculations is largely dependent on the user since the user is responsible for deciding the initial geometry of the bonds, a start guess is made. The main advantage of FF methods of modelling is that calculations can be generated reasonably quickly, provided the molecule is not too large. The quality of the predicted geometries and total energies depends on the quality of the parameters used to calculate these values. A limitation associated with the FF method is that they are Zero dimensional, thus one cannot evaluate the probable error of a specific results within the method. To assess the accuracy of a given result one must compare the result to calculations on similar classes of molecules.

Key literature

1. F. Jensen, Second Edition 207, 'Introduction to Computational Chemistry', pp 20-40

1.2 Modelling Using Molecular Mechanics

The Hydrogenation of Cyclopentadiene Dimer

Background

The dimerisation of cyclopentadiene is an example of a Diels Alder reaction where one cyclopentadiene molecule acts as the 4 π-electron diene (s-cis conformation) and the other as a 2 π-electron dieneophile. The Diels – Alder reaction forms a new six membered ring which adopts a cage like structure containing two double bonds. The type of product dimer can be classified as either ‘exo’ or ‘endo’. To evaluate which type of dimer is dominant, one must consider whether the reaction is controlled by thermodynamics or kinetics.

Thermodynamic control favours the reaction route which produces the product with the lowest total energy hence greatest product stability. Conversely, kinetic control favours the reaction route with the lowest transition state energy hence the reaction route with the lowest activation energy barrier. Image 2 displays the energy level profile for the dimerization reaction involving two cyclopentadiene molecules, see discussion later for more detail.

This aim of this part of the investigation was to use ChemBio 3D to draw and model structures 1-4, see image 3. The geometry of each individual structure was optimised using the MM2 energy minimisation calculation, as the lowest energy form of any molecule is always preferred.

Results

**Table 1: Force Field Energy Components: Structures 1-4, all energies have units kcal/mol**
Structure	Stretch	Bend	Stretch-Ben	Torsion	Non- 1,4 VDW	1,4 VDW	Dipole/Dipole	Total Energy
1	1.2834	20.6175	-0.8384	7.6239	-1.4023	4.2394	0.3769	31.8994
2	1.2435	20.8455	-0.8321	9.5094	-1.5143	4.3092	0.4456	34.0068
3	1.2693	19.8617	-0.8324	10.7938	-1.2004	5.6432	0.1621	35.6973
4	1.0985	14.5364	-0.5456	12.5039	-1.0741	4.5015	0.1407	31.1613

Discussion

Structures 1 and 2

Structures 1 and 2 correspond to the exo and endo products respectively.

The total energy of 2 is 2.1074 kcal mol-1 greater than that of 1. This shows that the exo product is thermodynamically more stable than the endo product. To understand why structure 1 is lower in total energy than 2 the individual energy components must be analysed.

The energy component with the largest numerical difference between structures 1 and 2 is likely to be the main contributor to this stability indifference.

The torsional energy component of the endo dimer is 1.8855 kcal mol-1 higher than the exo dimer. Torsional energy describes the energy change associated with the rotation around a B-C bond in a four atom sequence A-B-C-D, where A-B, B-C and C-D are bonded. The torsion angle is defined as the angle formed by the A-B and the C-D bonds, looking down the B-C bond. For example, the energy cost for rotation around a double bond is greater than the cost for rotation around a single bond since there is more rigidity associated with double bonds. Although both structures contain the same number of double and single bonds, the the endo dimer may experience slightly more steric hindrance between the ‘upward’ facing hydrogen’s adjacent to the triangular shaped C-C-C atom bridge. Comparing 'like for like' atoms between structures 1 and 2 the torsion angle in the exo dimer [ H(13)-C(3)-H(12)-C(2) ] is -112.9 degrees whereas the torsion angle in the endo dimer [H(12)- C(2) -C(3)- H(13)] is -41.9 degrees.

Since the torsion angle in the endo dimer is more acute there is likely to be greater steric interactions hence greater instability. The larger the substituent groups present the more significant the endo/exo torsion energy difference is likely to be since larger sized groups experience greater steric interactions.

According to the MM2 theoretical calculations the exo dimer is expected to be the major product. However, empirical evidence shows that the thermal dimerisation of cyclopentadiene yields the endo product exclusively. This means that the reaction must be kinetically controlled. ²

According to literature MP3 calculations, the electronic activation energy of the Diels-Alder reaction with cyclopentadiene as the diene and dienophille is lower for the endo dimer (20.71 kcal mol-1 ) than the exo dimer (23.42 kcal mol-1 ). The lower the activation energy for a given reaction, the lower the energy of the corresponding transition state hence the faster the product formation rate. This provides evidence to support why the endo dimer is formed exclusively. The reaction energy level profile, image 2, illustrates the relative energy levels of the exo and endo transition states and the final dimer products.³

Frontier orbital interaction theory can be used to explain why the endo transition state has a relatively lower energy than the exo transition state. When two cyclopentadiene molecules come together the symmetry of the reacting alkene orbitals is correct for bond formation. In addition, at the back of the diene the symmetry of the orbitals is appropriate for bonding interactions to exist across the ‘space between orbitals’. Although these interactions do not lead to the formal bonds, the stereochemistry of the dimer is adjusted so that the favourable interactions can exist. This specific type of favourable orbital interaction is referred to as secondary orbital overlap (SOO). SSO is defined as the positive overlap of a nonactive frame in the frontier molecular orbitals of a pericyclic reaction.⁴

The carbon atoms interact strongly in the endo transition state but not in the exo transition state. It is difficult to quantify the magnitude of the stabilising interaction and it is also difficult to understand the importance of these secondary orbital interactions relative to other intermolecular interactions in the transition state.

Structures 3 and 4

Structure 3 and 4 result from the hydrogenation reaction of the cyclopentadiene endo dimer. The total energy of structure 3 is 4.536 kcal mol-1 greater than 4. This suggests that 4 is thermodynamically more stable than 3. The two largest energy components that contribute most to the total energy differences are: the bending and 1,4 Van Der Waal (VDW) energies.

The bending energy of structure 4 is 5.3253 kcal mol-1 lower than that of 3. The bending energy refers to the energy required for bending an angle formed by three atoms. The lower bending energy in 4 can be explained by the presence of less strain in the molecule overall. Generally, the larger the bond angles are the less strained the system is due to reduced steric hindrance. For example, (see image 6) the angle between the labelled carbon atoms 5-6-7 in both structures 3 and 4 are 97.0 degrees and 101.0 degrees respectively. These results illustrate that smaller bond angles contribute to higher energy systems. It is interesting to compare the bond angles in the alkene carbon regions as this provides further evidence of relative bond strain. In structure 3 the carbons 1-5-6 and 5-1-2 have bond angles of 107.8 and 107.7 degrees respectively. Whereas, the carbons 3-8-9 and 8-9-10 have bond angles of 112.4 and 113.0 degrees respectively. The relatively smaller bond angles compared to the idealised 120 degree bond angle associated with sp2 carbons, in structure 3 shows that this hydrogenated dimer is likely to be the minor product.The 1,4 VDW energy describes the repulsion and attraction between atoms that are not directly bonded. Structure 4 is 1.1417 kcal mol-1 lower than 3.

It is possible to assess the relative ease of hydrogenation of each double bond in the endo dimer 2 (see image 4 carbon pairs 5-1 and 8-9). From a thermodynamic perspective the hydrogenation of the double bond corresponding to carbons 5-1, dimer 2, produces a more stable product. This is confirmed by the lower stretch, bend and van der Waals energy components of 4 relative to 3. Analysis of the C(5)=C(1) and C(8)=C(9) bond lengths of dimer 2 show that the the C(8)=C(9) is slightly shorter than C(5)=C(1)( 0.001 Angstrom difference)hence the C(8)=C(9)is likely to be slightly stronger and less prone to hydrogenation. Finally, comparing the alkene bond regions in dimer 2:

C(5)-C(1)-C(2) bond angle is 107.7 degrees

whereas

C(3)-C(8)-C(9) is 112.3 degrees

Ideally, the sp2 carbon centres (8 and 1) should have bond angles of 120 degrees. Since alkene C(8)-(C)9 is closer to this ideal value it is less likely to be hydrogenated. Alkene C(5)-C(1) is further away from the ideal angle, in fact it is more similar to the sp3 hybridised ideal bond angle, 109.5 degrees. As a result, there is a lower energy cost associated with the C(5)-C(1) double bond.

Key literature

2. DOI:http://pubs.acs.org/doi/pdf/10.1021/ja01494a060

3. DOI:http://pubs.acs.org/doi/pdf/10.1021/ja00060a048

4. DOI:http://pubs.acs.org/doi/pdf/10.1021/jo00384a016

____________________________________________________________________________________________________________________________________

1.2.2 Modelling Using Molecular Mechanics

Stereochemistry of Nucleophilic additions to a pyridinium ring (NAD+ analogue)

Image 9:MM2 Calculation Reactant 5 geometry optimistation test 2

Reaction 1: reactant 5 to form product 6

Image 6 shows the reaction scheme whereby the optically active derivative of prolinol (5) which reacts with methyl magnesium iodide. This results in the alkylation of the pyridine ring at the 4-position, with the absolute stereochemistry shown in reactant 6, image 6.

First, the structure of reactant 5 was drawn out using ChemBio 3D. A MM2 energy minimisation calculation was run to achieve the optimum geometry. The initial total energy of the molecule was 43.1309 kcal/mol, see image 7.

To identify the geometry corresponding to the lowest total energy, the structure of reactant 5 was manually altered several times by changing individual atom positions . Given the time constraints of the investigation, it was only possible to conduct a few different structural variation tests(note: not all tests have been illustrated).

After several tests, the lowest total energy achieved was 43.113 kcal/mol, see image 8.

Literature studies suggest that N-methyl and N-benzyl derivatives undergo highly regioselective and steroselective addition with Grignard reagents perhaps due to attractions between the amide oxygen atom (oxygen contains two lone pairs, electron rich density) and the Grignard reagent (Mg atom is electron deficient since it has electropositive properties). It has also been suggested that the organic ligand (CH3) is delivered from the magnesium to the carbon atom (4) in a conjugative manner, producing an intermediate enolate like structure, 13 see image 10. According to Dreiding stereomodels the seven membered ring in reactant 5 is not very flexible thus only two conformational minima geometries are possible:

Conformer 1. the amide carbonyl oxygen to be located above the plane of the pyridine ring, positioned anti to the hydrogen atom

Conformer 2. the carbonyl group and the pyridine ring are both located in the same plane

To test Dreiding steromodel theory the oxygen carbonyl was moved so that the dihedral angle between the carbonyl C-O and the aromatic C-C , was as close as possible to 180 degrees, see image 11. However, when the MM2 energy minimisation calculation was run the dihedral angle changed to -169.7 degrees (image 11) and the total energy value was 43.1416 kcal/mol. After several attempts It was not possible to manually reduce the total energy any further below 43.113 kcal/mol.⁵

One of the difficulties associated with ChemBio3D manual modelling is that there is a vast number of structure variation that are possible as each of the atoms that constitute the whole structure can be moved independently. As a result, there are a vast number of different total energies. In order to get a better understanding of the true total energy, a algorithm could be created to automatically model all the possible structural outcomes. A graph showing all the possible total energy values could be computed, the extreme outliers could be removed and finally an average energy value could be calculated. However, for large molecules such an algorithm may require a lot of computer memory process time.

Studies have shown that the mechanism for the methyl group addition to the pyridine ring, carbon (4), goes via a six centred transition state, see image 10 for diagrammatic explanation. Base on this empirical evidence I tested distorting the position of the carbonyl oxygen significantly below that of the pyridine ring and then running the energy minimisation calculation. Interestingly, after ‘minimisation’ the carbonyl group and the pyridine ring automatically reverted to the almost coplanar position. After optimisation, the carbonyl oxygen atoms remained slightly up relative to the aromatic ring (see image 8) thus the favourable interactions between the Grignard reagent (specifically the Mg atom) and the carbonyl oxygen also direct the methyl group to the upward face of the aromatic ring, see product 6 image 6. To extent this investigation further one could observe whether the regio-steroselectively is significantly increase when bulkier R groups are attached to the Grignard reagents.

When the MM2 energy minimisation calculation was completed for reactant 5 and the Grignard reactant MeMgI an error message appeared. This may be because the MM2 minimisation calculated does not include parameters for Mg metal atoms. It is difficult to identify a suitable force field for compounds involving metal coordination as there are few good parameters available. This is mainly due to the large variation of bonding possible around metals compared to organic molecules. To improve the model it would be necessary to include data to produce the necessary parameters to model the geometry of the MeMgI reactant also.

The same energy minimisation process was repeated on reactant 5 using the MMFF94 minimisation method. The total energy value obtained was 57.4335 kcal/mol. The reason for this significant total energy difference could be due to the different force fields parameters used. It is not possible to compare the total energies of a molecule using different force fields.

Reaction 2: reactant 7 to product 8

Similarly, reactant 7 was draw and modelled using ChemBio3D. The first MM2 total energy minimisation test produced a structure with a total energy of 36.5385 kcal/mol, see image 12. The geometry was altered to achieve the lowest total energy structure possible. Given the time limitations, the lowest total energy structure produced was 36.1213 kcal/mol, see image 15. Image 15 shows that the carbonyl oxygen atom is downward facing relative to the aromatic ring, the opposite position compared to reactant 5.

The mechanism of reaction 2 is a SN2 type nucleophillic addition reaction, occurring at carbon 4 of the aromatic ring. The carbonyl oxygen is relatively bulky thus the PHNH2 group is more likely to attack from the opposite site to the carbonyl oxygen. The carbonyl oxygen faces downwards relative to the aromatic ring, see image 15, hence the PHNH2 group will attack from the top face. Also the PHNH2 group cannot form any favourable interactions with the electronegative oxygen atom since the nitrogen and the phenyl rings are electronegative and electron rich respectively thus only unfavourable repulsive forces may be present.

Click on the buttons below to to view reactants 5 and 7 in 3D

Key Literature 5. A. G. Shultz, L. Flood and J. P. Springer, J. Org. Chemistry, 1986, 51, 838 DOI:10.1021/jo00356a016

____________________________________________________________________________________________________________________________________

1.2.3 Modelling Using Molecular Mechanics

Stereochemistry and Reactivity of an Intermediate in the Synthesis of Taxol

According to literature studies the key intermediates in the total Taxol synthesis may have originated via a oxy-Cope rearrangement mechanism. Intermediates 9 and 10 are structural isomers determined by the conformation of the cyclohexane ring :

Intermediate 9 forms a twist boat cyclohexane structure, carbonyl group ‘up’

Intermediate 10 forms a chair cyclohexane structure, carbonyl group ‘down’

Both intermediate structures were initially modelled using ChemBio3D (no atom position changes were made at this stage) and the total energies were minimised using the MM2 force field method. The results showed that intermediate 10 has a lower total energy than 9 thus it is possible to conclude that 10 is more stable than 9. From a visual perspective, intermediate 9 contains the carbonyl group pointing ‘upwards’ thus it may experience unfavourable steric interaction between the bridging carbon group with two methyl groups attached leading to a higher energy.

Click on the buttons below to view the intermediates in the total synthesis of Taxol:

The cyclohexane ring can adopt many different conformations including: boat, twist boat and chair. Several variations of the ring conformation were tested for both intermediates and the lowest energies resulted in a twist boat, intermediate 9, and a chair, intermediate 10.

The alkene located adjacent to the carbon bridgehead reacts very slowly which is unusual since most alkene are very reactive towards electron deficient species. This unusual property may be caused by the highly substituted nature of this particular alkene bond preventing oncoming reactants from easily accessing this region. Although both intermediates are quite strained due to the two fused rings, the bridgehead olefins are less reactive that expected. This unusual property referred to as alkene hyper-stability. Olefin strain energy is based on the concept the difference between the strain energy of an olefin and its parent hydrocarbon indicates the reactivity and stability of bridgehead olefins. This concept is the opposite of Bredt’s Rule since the olefin energies are negative reflecting greater stability.

__________________________________________________________________________________________________________________________________

Modelling Using Semi-Empirical Molecular Orbital Theory

In this section the focus is quantum mechanical electronic properties of reactivity, how electrons influence bonds and the derived spectroscopic properties of molecules. Using ChemBio3D, the MM2 force field was initially used to optimise the geometry of compound 12, the result showed that the minimum total energy was 17.9022 kcal/mol .

The component energies include:

Stretch: 0.6182
Bend: 4.8034
Stretch-Bend: 0.0407
Torsion: 7.6289
Non-1,4 VDW: -1.0968
1,4 VDW: 5.7962
Dipole/Dipole: 0.1116
Total Energy: 17.9022 kcal/mol

Once the MM2 energy minimisation calculation was complete the MOPAC/PM6 force field was used to generate a valence-electron molecular wave-function description of the electron density distribution. Using these quantum mechanical results various molecular orbital diagrams were created including: HOMO, HOMO-1, the LUMO, LUMO+1 and LUMO+2 (See Diagrams to the right).

Compound 12 reacts readily with electrophilic reagents such as dichlorocarbene or peracid due to the alkene regions of high electron density. The alkene HOMO interacts with the electrophile LUMO thus to evaluate which of the two alkenes is more susceptible to attack one must consider:

the size of the two molecular orbitals form each alkene region
the relative difference in energy between the HOMO (alkene) and LUMO (electrophile)- as the smaller the energy gap (HOMO-LUMO Gap)the more energetically favourable the interactions are leading to an addition reaction.

It is also important to consider the direction of attack of the incoming electrophile relative to the cyclopropyl bridge.

According to literature studies it is thought that incoming electrophiles react regiospecifically with the alkene double bond located endo to the chlorine atom on the cyclopropyl bridge. Since the chlorine and hydrogen atoms are relative similar in size, steric hindrance does not play a significant role when determining which face of the molecule the electrophile enters from. Instead, the regioselectivity is mainly caused by orbital and electrostatic factors. Unfortunately, the modelling methods used in this experiment do not allow us to view specific factors such as the transition state geometry, which may provide useful information about the reaction mechanism.

To evaluate which of the two alkenes is subject to electrophilic attack we must focus on the HOMO diagram of compound 12. This shows that the electron density cloud over the alkene syn to the chlorine atom on the bridging cyclopropyl ring is greater in size compared to the anti alkene. Thus it is possible to conclude that the syn alkene is more susceptible to electrohpilic attack. Literature studies also confirm this conclusion , see image: Literature Calculated PM3 HOMO compound 12.

**Table 2: MM2 Force Field Energy Components: Structures: Compound 12 (diene) and Hydrogenated Compound 12 (monoalkene) all energies have units kcal/mol**
Structure	Stretch	Bend	Stretch-Ben	Torsion	Non- 1,4 VDW	1,4 VDW	Dipole/Dipole	Total Energy	Image
12 (diene)	0.6119	4.7619	0.0393	7.6409	-1.0820	5.8159	0.1118	17.8996
12 hydrogenated (monoalkene)	0.8988	4.6837	0.0163	10.7394	-1.0305	6.9624	0.0708	22.3410

Section 2

A density function method was used to measure influence of the Cl-C bond on the vibrational frequencies of compound molecule 12 and its' dihydro derivative (monoalkene).

__________________________________________________________________________________________________________________________

Structure based Mini project using DFT-based Molecular orbital methods

The Birch reduction of unsaturated cyclic ketone 4 using t-BuOH leads to the alcohol 6 in 87 % yield and the ketone 5 in a ratio 10:1 and further conversion of 5 into 6 was resulted in 92% yield using a Birch reduction under the same conditions: Li, NH3, t-BuOH, THF, -78 °C.

According to literature ⁶ the 1H and 13C NMR spectra were recorded on a Bruker Advance 500 spectrometer at 500 and 125 MHz, respectively. The solvent used was CDCl3 as solvent and all reactions were fun under an atomosphere of argon.

How to...

a) test whether or not the conversion of 5 into 6 was achieved in the Birch reduction reaction

Infra Red Spectra Analysis could be used to assess which functional groups are present. The key definitive stretch frequency to be confident that the reaction has worked would be evidence of the alcohol -OH strech at around 3550-3200 cm ^-1 . Presence of a stretch around 1870-1540 cm ^-1 would indicate that there is still some of the ketone 5 present ⁷

b)which stereoisomer was produced

The literature results from the optical rotation show Compound 6: [a]²⁰_D -11.6. Thus the computational optical rotation result can be compared to the literature result and the product 6 could also be experimentally tested in the laboratory.

Why was this stereoisomer obtained?

Optical rotation of Compound 6 is reported: Compound 6: [a]²⁰_D -11.6 Optical rotation of Compound 5: [a]²⁰_D -19.5

Key Literature

6. L. Castellanos et al. / Stereoselective synthesis of (–)-4-epiaxinyssamine Tetrahedron 63 (2007) 1544–1552(DOI:10.1016/j.tetasy.2005.02.012 7. R.M Silverstein, C.G Bassler, T.C Morrill, Spectrometric Identification of Organic Compounds, pp73-83

__________________________________________________________________________________________________________________________________

References