View source for MRD:bs517-201905

== Units and indexing ==
Throughout this document all energies are expressed in kcal mol<nowiki><sup>-1</sup></nowiki>, distances in Å, time is an arbitrary parametrisation of the trajectory. For each system with atoms A, B and C, r<sub>1</sub>=r<sub>AB</sub>, r<sub>2</sub>=r<sub>BC</sub>, p<sub>1</sub>=p<sub>AB</sub> and p<sub>2</sub>=p<sub>BC</sub>.

== Exercise 1: H + H<sub>2</sub> system ==

=== Transition states ===
'''On a potential energy surface diagram, how is the transition state mathematically defined? How can the transition state be identified, and how can it be distinguished from a local minimum of the potential energy surface?'''
* Maximum of the minimum energy path between the reactants and products
* Saddle point on the potential energy surface; <math>\text{2D: }\nabla i: \frac{\partial V_i}{\partial p_i}=0, D=[V_{xx}V_{yy}-\left(V_{xy}\right)^2]<0</math>.
*Can be located by running simulations starting from initial conditions close to, but not equal to it; these simulations either (for a sufficiently small change in initial conditions) oscillate around the transition state,  or more likely, end up on the reactant/product side. A tiny change in initial conditions might change the outcome of the "reaction" completely. Simulations with the exact transition state as initial conditions will stay at the transition state indefinitely (though this cannot be achieved in practice).
*Distinguishing from local minima: simulations starting from sufficiently close to a local minimum will always stay there, oscillating; <math>\text{2D: }\nabla i: \frac{\partial V_i}{\partial p_i}=0, D=[V_{xx}V_{yy}-\left(V_{xy}\right)^2]>0, V_{xx}>0</math>.

{{fontcolor1|red|This equation looks just copied from somewhere and is not explained, please define the variables. Also it takes the gradient of the potential energy surface with respect to momentum? This doesn't really make sence. A discussion of orthogonal directions having either positive or negative second derivatives of potential energy with respect to interatomic distance would be more appropriate. [[User:Pu12|Pu12]] ([[User talk:Pu12|talk]]) 12:59, 24 May 2019 (BST)}}

=== Locating the transition state ===
'''Report your best estimate of the transition state position (r<sub>ts</sub>) and explain your reasoning illustrating it with a “Internuclear Distances vs Time” plot for a relevant trajectory.'''

The transition state was approximated by simulations with initial conditions r<sub>1</sub>=r<sub>2</sub>=cst., p<sub>1</sub>=p<sub>2</sub>=0. For initial conditions with internuclear distances sufficiently close to that of the transition state, the atoms stay in the simulation area and an oscillatory behaviour of the atoms was observed. Contrary to the oscillations observed in simulations leading to reactions, these happen in one dimension, i.e. they go thru the transition state as their equilibrium point in every period.

The final approximation of the transition state was reached by changing the initial conditions  r<sub>1</sub>=r<sub>2</sub>=cst., so the amplitude of the vibrations observed in the "Internuclear Distances vs Time" plot approached zero. The smallest amplitude was observed at  r<sub>1</sub>=r<sub>2</sub>=0.9077425=r<sub>ts</sub> (Figure 1-3).
[[File:Bs517-201905.1.2.3 0.95.png|thumb|Figure 1: Internuclear Distances vs Time graph of H+H<sub>2</sub> from r<sub>1</sub>=r<sub>2</sub>=0.95, p<sub>1</sub>=p<sub>2</sub>=0. The stoms oscillate at a low amplitude, so we are close to the transition state, but not exactly on it.]]
[[File:Bs517-201905.1.2.1 0.9077.png|thumb|Figure 2: Internuclear Distances vs Time graph of H+H<sub>2</sub> from r<sub>1</sub>=r<sub>2</sub>=0.9077≈r<sub>ts</sub>, p<sub>1</sub>=p<sub>2</sub>=0. Almost no oscillation can be observed on the simulation, we almost exactly on the transition state.]]
[[File:bs517-201905.1.2.4_0.9077425.png|thumb|Figure 3: Internuclear Distances vs Time graph of H+H<sub>2</sub> from r<sub>1</sub>=r<sub>2</sub>=0.9077=r<sub>ts</sub> p<sub>1</sub>=p<sub>2</sub>=0. Vibrations can still be seen, even this close to the transition state, however, only at very high zoom.]]
{{fontcolor1|red|Good use of graphs. [[User:Pu12|Pu12]] ([[User talk:Pu12|talk]]) 16:03, 23 May 2019 (BST)}}
=== Minimum energy path ===

[[File:bs517-201905.1.3.1_0.9177425_mep.png|thumb|Figure 4: Contour plot of the ''mep'' simulation with initial conditions r<sub alt="r_1=t_{ts}+\delta, r_2=r_{ts}, p_1=p_2=0, \delta=0.01">1</sub>=r<sub>ts</sub>+0.01, r<sub>2</sub>=r<sub>ts</sub>, p<sub>1</sub>=p<sub>2</sub>=0.]]'''Comment on how the ''mep'' and the trajectory you just calculated differ.'''

As we can see, the trajectories in Figures 4 and 5 differ in their vibrational states. As ''mep''  simulations reset the momentum of each particle to zero after every step, the particles will always move in the direction of the fastest potential energy descent around them, i.e. <math>\nabla s \in \text{steps}: p_s\parallel -\nabla V</math>. However, when running a ''Dynamics'' simulation the momenta are not reset, so the particles will not travel in the direction of the fastest descent, they will only be accelerated in its direction.{{fontcolor1|red|Good. [[User:Pu12|Pu12]] ([[User talk:Pu12|talk]]) 16:03, 23 May 2019 (BST)}} This retention of their velocities causes them to overshoot the bottom of the well every time, forcing them to vibrate while moving along the surface.
[[File:bs517-201905.1.3.2_0.9177425_dyn.png|thumb|Figure 5: Contour plot of the ''Dynamics'' simulation with initial conditions r<sub alt="r_1=t_{ts}+\delta, r_2=r_{ts}, p_1=p_2=0, \delta=0.01">1</sub>=r<sub>ts</sub>+0.01, r<sub>2</sub>=r<sub>ts</sub>, p<sub>1</sub>=p<sub>2</sub>=0.]]


As we can see from Figures 6 and 7, B and C get closer together and farther apart from A as the reaction progresses. Running the simulation from the final positions of the above trajectory and the final momenta reversed gives us back our original initial conditions at the end of the simulation (if the simulation is done for the same number of frames with the same step size). This supports the theory that when only conservative forces are present all processes could run forwards or backwards in time as well.

Offsetting r<sub>2</sub> from the transition state instead of r<sub>1</sub> gives similar results, the only difference being that the reaction goes in the other direction towards the reactants.

[[File:bs517-201905.1.3.3_0.9177425_dyn_idvt.png|thumb|Figure 6: Internuclear Distances vs Time plot of the ''Dynamics'' simulation with initial conditions r<sub alt="r_1=t_{ts}+\delta, r_2=r_{ts}, p_1=p_2=0, \delta=0.01">1</sub>=r<sub>ts</sub>+0.01, r<sub>2</sub>=r<sub>ts</sub>, p<sub>1</sub>=p<sub>2</sub>=0.]]

[[File:bs517-201905.1.3.4_0.9177425_dyn_mvt.png|thumb|Figure 7: Momenta vs Time plot of the ''Dynamics'' simulation with initial conditions r<sub alt="r_1=t_{ts}+\delta, r_2=r_{ts}, p_1=p_2=0, \delta=0.01">1</sub>=r<sub>ts</sub>+0.01, r<sub>2</sub>=r<sub>ts</sub>, p<sub>1</sub>=p<sub>2</sub>=0.]]

=== Reactive and unreactive trajectories ===
'''Complete the table below by adding the total energy, whether the trajectory is reactive or unreactive, and provide a plot of the trajectory and a small description for what happens along the trajectory. What can you conclude from the table?'''

{| class="wikitable" border=1
! p<sub>1</sub> !! p<sub>2</sub> !! E<sub>tot</sub> !! Reactive? !! Description of the dynamics !! Trajectory
|-
| -1.25 || -2.5  ||-99.018||Yes||Weak B-C oscillations after bond has formed || [[File:Bs517-201905.1.3.5.1.png]]
|-
| -1.5  || -2.0  ||-100.456||No||Reaction does not reach TS || [[File:Bs517-201905.1.3.5.2.png]]
|-
| -1.5  || -2.5  ||-98.956||Yes||Weak B-C oscillations after bond has formed || [[File:Bs517-201905.1.3.5.3.png]]
|-
| -2.5  || -5.0  ||-84.956||No||B "bounces" voilently between A and C until finally settling with A again || [[File:Bs517-201905.1.3.5.4.png]]
|-
| -2.5  || -5.2  ||-83.416||Yes||B "bounces" between A and C, but settles on C || [[File:Bs517-201905.1.3.5.5.png]]
|}
{{fontcolor1|red|You should discuss the relative motions of the atoms more in the above descriptions and refer to the crossing of the transition state more. [[User:Pu12|Pu12]] ([[User talk:Pu12|talk]]) 16:03, 23 May 2019 (BST)}}
As we can see from the above table, the distribution of energies is extremely important when determining whether a reaction is going to be successful or not. In the second to last case for example the system enters the product region at a point, but it "bounces back", ultimately ending up as the starting materials again.

'''State what are the main assumptions of Transition State Theory. Given the results you have obtained, how will Transition State Theory predictions for reaction rate values compare with experimental values?'''

The assumptions of Transition State Theory are as follows:<ref><nowiki>T. Bligaard, J.K. Nørskov, in Chemical Bonding at Surfaces and Interfaces, 2008 (https://www.sciencedirect.com/topics/chemistry/transition-state-theory)</nowiki></ref>
* Born-Oppenheimer approximation
* quantum tunnelling effects are ignored
* Boltzmann distribution of energies in the reactants
* If transition state is reached with a velocity towards product state, system will not come back towards the starting materials
Experimental values are very close to those predicted by TST if all assumptions hold true. TST overestimates reaction rates if, for example, barrier recrossing is significant, and underestimates them if quantum tunnelling pathways are involved in the forwards direction.
{{fontcolor1|red|Good. [[User:Pu12|Pu12]] ([[User talk:Pu12|talk]]) 16:03, 23 May 2019 (BST)}}
== Exercise 2: F-H-H system ==

=== PES inspection ===
'''By inspecting the potential energy surfaces, classify the F + H<sub>2</sub> and H + HF reactions according to their energetics (endothermic or exothermic). How does this relate to the bond strength of the chemical species involved?'''

'''Locate the approximate position of the transition state.'''

The reaction between F + H<sub>2</sub> is exothermic and the reverse reaction is endothermic. We can see this from the potential energy surface, as the equipotential lines go lower on the side of the HF species. This is due to the strong H-F bond, that has a high ionic contribution, especially compared to the non-polar H-H bond.

The approximate transition state is r<sub>1</sub>=1.81069389, r<sub>2</sub>=0.7449985 (Figure 8, 9). This was found by doing ''Dynamics'' simulations with p<sub>1</sub>=p<sub>2</sub>=0, and changing r<sub>1</sub> and r<sub>2</sub> towards the F + H<sub>2</sub> side bit by bit if the reaction went to the HF + H state, and ''vice versa''.

[[File:Bs517-201905.2.1.1 1.81069389.png|thumb|Figure 8: Internuclear Distances vs Time plot of the F-H-H system from the initial conditions r<sub>1</sub>=1.81069389, r<sub>2</sub>=0.7449985]]
[[File:Bs517-201905.2.1.2 1.81069389.png|thumb|Figure 9: Contour plot of the F-H-H system from the initial conditions r<sub>1</sub>=1.81069389, r<sub>2</sub>=0.7449985]]'''Report the activation energy for both reactions.'''

Transition state energy: -103.752

HF + H energy: -133.910

H<sub>2</sub> + F energy: -103.973

E<sub>A</sub>(F+H<sub>2</sub>->HF+H)=0.221

E<sub>A</sub>(HF+H->H<sub>2</sub>+F)=30.158

{{fontcolor1|red|Correct. [[User:Pu12|Pu12]] ([[User talk:Pu12|talk]]) 16:03, 23 May 2019 (BST)}}
=== Reaction dynamics ===
'''In light of the fact that energy is conserved, discuss the mechanism of release of the reaction energy. Explain how this could be confirmed experimentally.'''

As we can see from Figure 10, the H-F bond at the end of the reaction vibrates much more, than the H-H bond did at the start of the reaction. This means that throughout the reaction potential energy is converted into heat energy, e.i. the vibration of the atoms. This could be confirmed by for example a simple calorimetry experiment in a bomb calorimeter.

{{fontcolor1|red|Also translational energy is converted into vibrational. [[User:Pu12|Pu12]] ([[User talk:Pu12|talk]]) 16:03, 23 May 2019 (BST)}}

[[File:Bs517-201905.2.1.3.png|thumb|Figure 10: Controur plot of the F+ H<sub>2</sub> system from the initial conditions r<sub>1</sub>=2.424212428049714, r<sub>2</sub>=0.7399019047604279, p<sub>1</sub>=-1, p<sub>2</sub>=-0.4]]

'''Setup a calculation starting on the side of the reactants of F + H<sub>2</sub>, at the bottom of the well r<sub>HH</sub> = 0.74, r<sub>FH</sub> = 2.4, with a momentum p<sub>FH</sub> = -0.5, and explore several values of p<sub>HH</sub> in the range -3 to 3 (explore values also close to these limits). What do you observe? Note that we are putting a significant amount of energy (much more than the activation energy) into the system on the H - H vibration.'''

{| class="wikitable" border=1
! p<sub>1</sub> !! p<sub>2</sub> !! E<sub>tot</sub> !! Reactive? !! Description of the dynamics || Trajectory
|-
| -0.5 || -3 || -96.279 || No || Barrier recrossing || [[File:Bs517-201905.2.2.1 -3.png]]
|-
| -0.5 || -2.8 || -96.816|| No|| Barrier recrossing and re-recrossing || [[File:Bs517-201905.2.2.1 -2.9.png]]
|-
| -0.5 || -2.8 || -97.339 || Yes  || Barrier recrossing and re-recrossing || [[File:Bs517-201905.2.2.1 -2.8.png]]
|-
| -0.5 || -2.6 || -98.319|| No|| Barrier recrossing || [[File:Bs517-201905.2.2.1 -2.6.png]]
|-
| -0.5 || -2.4 || 99.219|| Yes|| Barrier recrossing and re-recrossing || [[File:Bs517-201905.2.2.1 -2.4.png]]
|-
| -0.5 || -2.2 || -100.039|| No|| Does not reach TS || [[File:Bs517-201905.2.2.1 -2.2.png]]
|-
| -0.5 || -2.0 || -100.779|| No|| Does not reach TS || [[File:Bs517-201905.2.2.1 -2.0.png]]
|-
| -0.5 || -1.5 || -102.279|| No|| Does not reach TS || [[File:Bs517-201905.2.2.1 -1.5.png]]
|-
| -0.5 || -1 || -103.379|| No|| Does not reach TS|| [[File:Bs517-201905.2.2.1 -1.png]]
|-
| -0.5 || 0 || -103.779|| No|| Does not reach TS|| [[File:Bs517-201905.2.2.1 0.png]]
|-
| -0.5 || 1 || -102.279|| No|| Does not reach TS|| [[File:Bs517-201905.2.2.1 1.png]]
|-
| -0.5 || 1.5 || -100.779|| No|| Does not reach TS || [[File:Bs517-201905.2.2.1 1.5.png]]
|-
| -0.5 || 2 || -98.779|| No|| Barrier recrossing || [[File:Bs517-201905.2.2.1 -2.6.png]]
|-
| -0.5 || 2.2 || -97.839|| No|| Barrier recrossing || [[File:Bs517-201905.2.2.1 2.2.png]]
|-
| -0.5 || 2.4 || -98.819|| No|| Barrier recrossing || [[File:Bs517-201905.2.2.1 2.4.png]]
|-
| -0.5 || 2.6 || -95.719|| No||Barrier recrossing || [[File:Bs517-201905.2.2.1 2.6.png]]
|-
| -0.5 || 2.8 || -94.539|| No|| Multiple barrier recrossings and re-recrossings || [[File:Bs517-201905.2.2.1 2.8.png]]
|-
| -0.5 || 2.9 || -93.919|| No|| Barrier recrossing || [[File:Bs517-201905.2.2.1 2.9.png]]
|-
| -0.5 || 3 || -93.279|| Yes|| Multiple barrier recrossings and re-recrossings || [[File:Bs517-201905.2.2.1 3.png]]
|-
| -0.8 || 0.1 || -103.484|| Yes|| Only one barrier crossing, reaction is successful even though the overall energy is much lower, than in any of the other cases with a successful reaction --> better energy distribution || [[File:Bs517-201905.2.2.2.1.png]]
|}
We can see from the table that the initial energy distribution is even more important for this reaction, than for the reactoin between H and H<sub>2</sub>. Also, just putting more energy into the system may or may not yield a successful reaction.

Sometimes changing the initial conditions bit a tiny amount changes the outcome of the reaction completely, so though everything here is deterministic, there is an element of chaos theory involved.

'''The cases studied are an illustration of Polanyi's empirical rules. Discuss how the distribution of energy between different modes (translation and vibration) affect the efficiency of the reaction, and how this is influenced by the position of the transition state.'''

Polanyi's empirical rules are based on Hammond's postulate. They state that for endothermic reactions vibrational energy is preferred to translational energy on the reactant side, and the other way around for exothermic reactions. This is because to get up to the top of an endothermic barrier, a system with low vibrational excitation and high translational energy will not be able to pass the potential wall as it just "bounces off", usually a high vibational excitation is required to get up while going back and forth between the two walls. In the exothermic case the system is required to have high translational energy to, instead of just bouncing between the walls and going back to the reactant side, actually pass the barrier and descend to the potential well.<ref name=":0">Steinfeld, J.I., Francisco, Joseph S & Hase, William L, 1989. Chemical kinetics and dynamics, Englewood Cliffs, N.J.: Prentice Hall.</ref>

As we can see from Figures 11-15, for a reaction of HF + F with high initial viabrational energy and low initial translational energy (most of the energy of the reaction is stored in the vibrations of the F-H bond initially), the reaction succeeds even though the total energy of the system is very low (-102.463), meaning that this must be a very efficient trajectory. At first, the system bounces between the walls of the potential well (F-H bond vibrates). Next, the HF hydrogen gets passed back and forth between the F and other H atoms, finally settling with the other H atom. We can see on Figure 13, the H molecule also vibrates at the end of the reaction, but it does so significantly less violently, than the HF molecule did.

{{fontcolor1|red|Good report overall 5/5. [[User:Pu12|Pu12]] ([[User talk:Pu12|talk]]) 12:59, 24 May 2019 (BST)}}

[[File:Bs517-201905.2.3.2 0.925 2.25 -9 -1.5.png|thumb|Figure 11: Contour plot of the reaction of F + H<sub>2</sub> initial conditions r<sub>1</sub>=0.925, r<sub>2</sub>=2.25, p<sub>1</sub>=-9, p<sub>2</sub>=-1.5.]]
[[File:Bs517-201905.2.3.2 0.925 2.25 -9 -1.5_surf.png|thumb|Figure 12: Surface plot of the reaction of F + H<sub>2</sub> initial conditions r<sub>1</sub>=0.925, r<sub>2</sub>=2.25, p<sub>1</sub>=-9, p<sub>2</sub>=-1.5.]]
[[File:Bs517-201905.2.3.3 0.925 2.25 -9 -1.5 idvt.png|thumb|Figure 13: Internuclear distances vs Time plot of the reaction of F + H<sub>2</sub> initial conditions r<sub>1</sub>=0.925, r<sub>2</sub>=2.25, p<sub>1</sub>=-9, p<sub>2</sub>=-1.5.]]
[[File:Bs517-201905.2.3.3 0.925 2.25 -9 -1.5 evt.png|thumb|Figure 14: Energies vs Time plot of the reaction of F + H<sub>2</sub> initial conditions r<sub>1</sub>=0.925, r<sub>2</sub>=2.25, p<sub>1</sub>=-9, p<sub>2</sub>=-1.5.]]
[[File:Bs517-201905.2.3.3 0.925 2.25 -9 -1.5 mvt.png|thumb|Figure 15: Momenta vs Time plot of the reaction of F + H<sub>2</sub> initial conditions r<sub>1</sub>=0.925, r<sub>2</sub>=2.25, p<sub>1</sub>=-9, p<sub>2</sub>=-1.5.]]