Chapter 25. One-Dimensional Reaction Pathway

Table of Contents

Input File (rxnpath/rxnpath1.dat):
Archive File (rxnpath/rxnpath1.arc):
Output File (rxnpath/rxnpath1.out) :

This example illustrates the use of the reaction coordinate method for locating an approximate transition state for further refinement. It involves computation of the rotational barrier of ethane using the MINDO3 Hamiltonian. The energies and reaction coordinates could be plotted to produce a potential map such as that shown in Figure 4.1, “Potential energy along a reaction coordinate”. It should be noted in this connection that a reaction path may fail to pass through a saddle point (be discontinuous). Also, the intermediate points generated along the reaction path are not defined in a strict mathematical sense because one of the degrees of freedom (the reaction coordinate) is constrained and so not allowed to optimize.

Instructions for providing coordinate values in the input file can be found in the section called “Extra Input Data”.

Input File (rxnpath/rxnpath1.dat):

  mindo3 rhf singlet t=auto grad truste 1
Ethane rotational barrier
Reaction Path calculation
 C              0.000000  0    0.000000  0    0.000000  0    0    0    0
 C              1.544846  1    0.000000  0    0.000000  0    1    0    0
 H              1.116935  1  110.665270  1    0.000000  0    2    1    0
 H              1.116935  1  110.665213  1  120.000000  1    2    1    3
 H              1.116935  1  110.665213  1 -120.000000  1    2    1    3
 H              1.116935  1  110.665213  1   60.000000 -1    1    2    3  2
 H              1.116935  1  110.665213  1  120.000000  1    1    2    6
 H              1.116935  1  110.665213  1 -120.000000  1    1    2    6
 0              0.000000  0    0.000000  0    0.000000  0    0    0    0
$$ rxnpath - values along reaction path  3
  50 40 30 20 10 0 -10 -20 -30 -40 -50 -60   4
$$ end of extra data
	

1

A reaction path calculation will be performed if one and only one of optimization flags are marked with a -1, so there is no keyword need for a reaction path. TRUSTE has been specified to perform the constrained optimization at each step on the relaxation path. TRUSTE is used by default but BFGS, DFP, or EF may also be used instead.

2

The dihedral angle of this atom is marked with a -1 for the optimization flag. This is the reaction coordinate. The hydrogens after the one defining the reaction coordinate have their dihedrals defined with respect to it. A change in the dihedral of this atom will rotate the entire methyl group because of the way the geometry is defined. This is a common and useful way to define the geometries of methyl groups.

3

This is the extra input section marker for reaction path data. Note, that this marker can be shortened to $$ rxnpath. Details of these markers are found in the section called “Extra Input Data”

4

The line of data at the bottom of the file gives the values which the reaction coordinate is to take on in the course of the study. The first calculation will use 60° (from the input geometry) and subsequent calculations will be done at 50°, 40°, and so forth (from the extra data section). All geometric parameters marked for optimization (except the reaction coordinate itself) will be optimized at each of these points. The archive file and output file will contain each of these results. There is considerable savings in computational effort in using the reaction coordinate method instead of a series of individual calculations, in that the density matrix from the last step in the reaction coordinate is used as the starting point for the next step. This approximation can also cause difficulties if the wavefunction becomes undefined or does not evolve rapidly enough in relation to the reaction coordinate. (If this occurs, try NEWDEN which forces the density to be recomputed at each step.)

Archive File (rxnpath/rxnpath1.arc):

This file is somewhat abbreviated to conserve space. The marker ========= is used to show where some of the results were omitted. A complete copy of the file is found in the test suite results directory.

 Timestamp:  2011-08-31-12-51-45-00000014C8-win64
 User Info:  John Millam, Nahum, 

                     SUMMARY OF MINDO/3 CALCULATION
                                                       Aug-31-2011
                          AMPAC Version 10.0.1
                             Presented by:
  
                        Semichem, Inc.
                        www.semichem.com
  
 FORMULA: C2H6
 Ethane rotational barrier
 Reaction Path calculation

     GEOMETRY OPTIMIZED : ENERGY MINIMIZED
     SCF FIELD WAS ACHIEVED
 
          FINAL HEAT OF FORMATION   =       -19.856372 kcal
                                    =       -83.098917 kJ
          ELECTRONIC ENERGY         =      -918.458568 eV
          CORE-CORE REPULSION       =       575.250965 eV
          TOTAL ENERGY              =      -343.207603 eV
          GRADIENT NORM             =         0.097962 
          RMS GRADIENT NORM         =         0.023759 
          UNSTABLE MODE(S)          =         0 ( ESTIMATE  )
          FOR REACTION COORDINATE   =        60.000000 degrees    1
          REACTION GRADIENT         =         0.000000 kcal/radian    2
          IONIZATION POTENTIAL      =        11.758127 eV
          HOMO-LUMO GAP             =        14.407093 eV
          DIPOLE                    =         0.000000 debyes
          MOLECULAR WEIGHT          =        30.069400 
          MOLECULAR POINT GROUP     = D3d     0.100000
          NO. OF FILLED LEVELS      =         7 (OCC = 2)
          TOTAL NUMBER OF ORBITALS  =        14
          COMPUTATION TIME          =         0.00     SECONDS

          FINAL GEOMETRY OBTAINED                                       CHARGE
 MINDO3 RHF SINGLET T=AUTO GRAD TRUSTE
 Ethane rotational barrier
 Reaction Path calculation
  C     0.000000  0    0.000000  0    0.000000  0      0     0     0    0.0761
  C     1.477145  1    0.000000  0    0.000000  0      1     0     0    0.0761
  H     1.111407  1  113.134821  1    0.000000  0      2     1     0   -0.0254
  H     1.111407  1  113.134821  1  120.000009  1      2     1     3   -0.0254
  H     1.111407  1  113.134821  1 -120.000010  1      2     1     3   -0.0254
  H     1.111407  1  113.134822  1   60.000000 -1      1     2     3   -0.0254    3
  H     1.111407  1  113.134818  1  119.999997  1      1     2     6   -0.0254
  H     1.111407  1  113.134821  1 -119.999998  1      1     2     6   -0.0254
  0     0.000000  0    0.000000  0    0.000000  0      0     0     0
 Timestamp:  2011-08-31-12-51-45-00000014C8-win64
 User Info:  John Millam, Nahum, 

                     SUMMARY OF MINDO/3 CALCULATION
                                                       Aug-31-2011
                          AMPAC Version 10.0.1
                             Presented by:
  
                        Semichem, Inc.
                        www.semichem.com
  
 FORMULA: C2H6
 Ethane rotational barrier
 Reaction Path calculation

     GEOMETRY OPTIMIZED : ENERGY MINIMIZED
     SCF FIELD WAS ACHIEVED
 
          FINAL HEAT OF FORMATION   =       -19.795135 kcal
                                    =       -82.842640 kJ
          ELECTRONIC ENERGY         =      -918.428736 eV
          CORE-CORE REPULSION       =       575.223788 eV
          TOTAL ENERGY              =      -343.204948 eV
          GRADIENT NORM             =         0.136921 
          RMS GRADIENT NORM         =         0.033208 
          UNSTABLE MODE(S)          =         0 ( ESTIMATE  )
          FOR REACTION COORDINATE   =        50.000000 degrees    4
          REACTION GRADIENT         =        -0.688205 kcal/radian
          IONIZATION POTENTIAL      =        11.748697 eV
          HOMO-LUMO GAP             =        14.397865 eV
          DIPOLE                    =         0.001340 debyes
          MOLECULAR WEIGHT          =        30.069400 
          MOLECULAR POINT GROUP     = D3      0.100000
          NO. OF FILLED LEVELS      =         7 (OCC = 2)
          TOTAL NUMBER OF ORBITALS  =        14
          COMPUTATION TIME          =         0.02     SECONDS

          FINAL GEOMETRY OBTAINED                                       CHARGE
 MINDO3 RHF SINGLET T=AUTO GRAD TRUSTE
 Ethane rotational barrier
 Reaction Path calculation
  C     0.000000  0    0.000000  0    0.000000  0      0     0     0    0.0761
  C     1.477286  1    0.000000  0    0.000000  0      1     0     0    0.0761
  H     1.111390  1  113.172356  1    0.000000  0      2     1     0   -0.0253
  H     1.111257  1  113.219586  1  120.364741  1      2     1     3   -0.0250
  H     1.111540  1  113.092673  1 -119.643604  1      2     1     3   -0.0258
  H     1.111390  1  113.172362  1   50.000000 -1      1     2     3   -0.0253    5
  H     1.111257  1  113.219582  1  120.364729  1      1     2     6   -0.0250
  H     1.111540  1  113.092672  1 -119.643592  1      1     2     6   -0.0258
  0     0.000000  0    0.000000  0    0.000000  0      0     0     0

==========

 Timestamp:  2011-08-31-12-51-45-00000014C8-win64
 User Info:  John Millam, Nahum, 

                     SUMMARY OF MINDO/3 CALCULATION
                                                       Aug-31-2011
                          AMPAC Version 10.0.1
                             Presented by:
  
                        Semichem, Inc.
                        www.semichem.com
  
 FORMULA: C2H6
 Ethane rotational barrier
 Reaction Path calculation

     GEOMETRY OPTIMIZED : ENERGY MINIMIZED
     SCF FIELD WAS ACHIEVED
 
          FINAL HEAT OF FORMATION   =       -19.150443 kcal
                                    =       -80.144602 kJ
          ELECTRONIC ENERGY         =      -918.152857 eV
          CORE-CORE REPULSION       =       574.975865 eV
          TOTAL ENERGY              =      -343.176992 eV
          GRADIENT NORM             =         0.145858 
          RMS GRADIENT NORM         =         0.035376 
          UNSTABLE MODE(S)          =         0 ( ESTIMATE  )
          FOR REACTION COORDINATE   =        20.000000 degrees
          REACTION GRADIENT         =        -1.268164 kcal/radian    6
          IONIZATION POTENTIAL      =        11.734594 eV
          HOMO-LUMO GAP             =        14.386343 eV
          DIPOLE                    =         0.000933 debyes
          MOLECULAR WEIGHT          =        30.069400 
          MOLECULAR POINT GROUP     = D3      0.100000
          NO. OF FILLED LEVELS      =         7 (OCC = 2)
          TOTAL NUMBER OF ORBITALS  =        14
          COMPUTATION TIME          =         0.00     SECONDS

          FINAL GEOMETRY OBTAINED                                       CHARGE
 MINDO3 RHF SINGLET T=AUTO GRAD TRUSTE
 Ethane rotational barrier
 Reaction Path calculation
  C     0.000000  0    0.000000  0    0.000000  0      0     0     0    0.0760
  C     1.478837  1    0.000000  0    0.000000  0      1     0     0    0.0760
  H     1.111310  1  113.346051  1    0.000000  0      2     1     0   -0.0252
  H     1.110935  1  113.593357  1  120.736162  1      2     1     3   -0.0243
  H     1.111760  1  113.152213  1 -119.254457  1      2     1     3   -0.0264
  H     1.111310  1  113.346059  1   20.000000 -1      1     2     3   -0.0252
  H     1.110935  1  113.593353  1  120.736135  1      1     2     6   -0.0243
  H     1.111760  1  113.152210  1 -119.254430  1      1     2     6   -0.0264
  0     0.000000  0    0.000000  0    0.000000  0      0     0     0
 Timestamp:  2011-08-31-12-51-45-00000014C8-win64
 User Info:  John Millam, Nahum, 

                     SUMMARY OF MINDO/3 CALCULATION
                                                       Aug-31-2011
                          AMPAC Version 10.0.1
                             Presented by:
  
                        Semichem, Inc.
                        www.semichem.com
  
 FORMULA: C2H6
 Ethane rotational barrier
 Reaction Path calculation

     GEOMETRY OPTIMIZED : ENERGY MINIMIZED
     SCF FIELD WAS ACHIEVED
 
          FINAL HEAT OF FORMATION   =       -18.970678 kcal
                                    =       -79.392287 kJ
          ELECTRONIC ENERGY         =      -918.060403 eV
          CORE-CORE REPULSION       =       574.891207 eV
          TOTAL ENERGY              =      -343.169197 eV
          GRADIENT NORM             =         0.072150 
          RMS GRADIENT NORM         =         0.017499 
          UNSTABLE MODE(S)          =         0 ( ESTIMATE  )
          FOR REACTION COORDINATE   =        10.000000 degrees
          REACTION GRADIENT         =        -0.743838 kcal/radian
          IONIZATION POTENTIAL      =        11.735142 eV
          HOMO-LUMO GAP             =        14.387602 eV
          DIPOLE                    =         0.000343 debyes
          MOLECULAR WEIGHT          =        30.069400 
          MOLECULAR POINT GROUP     = D3      0.100000
          NO. OF FILLED LEVELS      =         7 (OCC = 2)
          TOTAL NUMBER OF ORBITALS  =        14
          COMPUTATION TIME          =         0.02     SECONDS

          FINAL GEOMETRY OBTAINED                                       CHARGE
 MINDO3 RHF SINGLET T=AUTO GRAD TRUSTE
 Ethane rotational barrier
 Reaction Path calculation
  C     0.000000  0    0.000000  0    0.000000  0      0     0     0    0.0759
  C     1.479327  1    0.000000  0    0.000000  0      1     0     0    0.0759
  H     1.111337  1  113.412828  1    0.000000  0      2     1     0   -0.0253
  H     1.111106  1  113.543014  1  120.413103  1      2     1     3   -0.0248
  H     1.111515  1  113.340628  1 -119.588183  1      2     1     3   -0.0259
  H     1.111337  1  113.412775  1   10.000000 -1      1     2     3   -0.0253
  H     1.111106  1  113.543050  1  120.413195  1      1     2     6   -0.0248
  H     1.111515  1  113.340649  1 -119.588278  1      1     2     6   -0.0259
  0     0.000000  0    0.000000  0    0.000000  0      0     0     0
 Timestamp:  2011-08-31-12-51-45-00000014C8-win64
 User Info:  John Millam, Nahum, 

                     SUMMARY OF MINDO/3 CALCULATION
                                                       Aug-31-2011
                          AMPAC Version 10.0.1
                             Presented by:
  
                        Semichem, Inc.
                        www.semichem.com
  
 FORMULA: C2H6
 Ethane rotational barrier
 Reaction Path calculation

     GEOMETRY OPTIMIZED : ENERGY MINIMIZED
     SCF FIELD WAS ACHIEVED
 
          FINAL HEAT OF FORMATION   =       -18.903827 kcal
                                    =       -79.112517 kJ
          ELECTRONIC ENERGY         =      -918.017449 eV
          CORE-CORE REPULSION       =       574.851151 eV
          TOTAL ENERGY              =      -343.166298 eV
          GRADIENT NORM             =         0.215401 
          RMS GRADIENT NORM         =         0.052243 
          UNSTABLE MODE(S)          =         0 ( ESTIMATE  )
          FOR REACTION COORDINATE   =         0.000000 degrees
          REACTION GRADIENT         =         0.005405 kcal/radian    7
          IONIZATION POTENTIAL      =        11.734696 eV
          HOMO-LUMO GAP             =        14.387297 eV
          DIPOLE                    =         0.000208 debyes
          MOLECULAR WEIGHT          =        30.069400 
          MOLECULAR POINT GROUP     = D3h     0.100000
          NO. OF FILLED LEVELS      =         7 (OCC = 2)
          TOTAL NUMBER OF ORBITALS  =        14
          COMPUTATION TIME          =         0.05     SECONDS

          FINAL GEOMETRY OBTAINED                                       CHARGE
 MINDO3 RHF SINGLET T=AUTO GRAD TRUSTE
 Ethane rotational barrier
 Reaction Path calculation
  C     0.000000  0    0.000000  0    0.000000  0      0     0     0    0.0759
  C     1.479559  1    0.000000  0    0.000000  0      1     0     0    0.0759
  H     1.111302  1  113.477840  1    0.000000  0      2     1     0   -0.0253
  H     1.111385  1  113.403539  1  119.926319  1      2     1     3   -0.0254
  H     1.111250  1  113.504329  1 -120.036232  1      2     1     3   -0.0252
  H     1.111301  1  113.477992  1    0.000000 -1      1     2     3   -0.0253
  H     1.111386  1  113.403447  1  119.926022  1      1     2     6   -0.0254
  H     1.111251  1  113.504273  1 -120.035935  1      1     2     6   -0.0252
  0     0.000000  0    0.000000  0    0.000000  0      0     0     0
 Timestamp:  2011-08-31-12-51-45-00000014C8-win64
 User Info:  John Millam, Nahum, 

                     SUMMARY OF MINDO/3 CALCULATION
                                                       Aug-31-2011
                          AMPAC Version 10.0.1
                             Presented by:
  
                        Semichem, Inc.
                        www.semichem.com
  
 FORMULA: C2H6
 Ethane rotational barrier
 Reaction Path calculation

     GEOMETRY OPTIMIZED : ENERGY MINIMIZED
     SCF FIELD WAS ACHIEVED
 
          FINAL HEAT OF FORMATION   =       -18.970662 kcal
                                    =       -79.392219 kJ
          ELECTRONIC ENERGY         =      -918.058220 eV
          CORE-CORE REPULSION       =       574.889023 eV
          TOTAL ENERGY              =      -343.169196 eV
          GRADIENT NORM             =         0.080598 
          RMS GRADIENT NORM         =         0.019548 
          UNSTABLE MODE(S)          =         0 ( ESTIMATE  )
          FOR REACTION COORDINATE   =       -10.000000 degrees
          REACTION GRADIENT         =         0.744132 kcal/radian    8
          IONIZATION POTENTIAL      =        11.734676 eV
          HOMO-LUMO GAP             =        14.387100 eV
          DIPOLE                    =         0.000029 debyes
          MOLECULAR WEIGHT          =        30.069400 
          MOLECULAR POINT GROUP     = D3      0.100000
          NO. OF FILLED LEVELS      =         7 (OCC = 2)
          TOTAL NUMBER OF ORBITALS  =        14
          COMPUTATION TIME          =         0.02     SECONDS

          FINAL GEOMETRY OBTAINED                                       CHARGE
 MINDO3 RHF SINGLET T=AUTO GRAD TRUSTE
 Ethane rotational barrier
 Reaction Path calculation
  C     0.000000  0    0.000000  0    0.000000  0      0     0     0    0.0760
  C     1.479337  1    0.000000  0    0.000000  0      1     0     0    0.0760
  H     1.111297  1  113.455775  1    0.000000  0      2     1     0   -0.0253
  H     1.111564  1  113.298462  1  119.596941  1      2     1     3   -0.0259
  H     1.111093  1  113.547580  1 -120.433148  1      2     1     3   -0.0248
  H     1.111297  1  113.455757  1  -10.000000 -1      1     2     3   -0.0253
  H     1.111564  1  113.298474  1  119.596983  1      1     2     6   -0.0259
  H     1.111093  1  113.547589  1 -120.433197  1      1     2     6   -0.0248
  0     0.000000  0    0.000000  0    0.000000  0      0     0     0

==========

 Timestamp:  2011-08-31-12-51-45-00000014C8-win64
 User Info:  John Millam, Nahum, 

                     SUMMARY OF MINDO/3 CALCULATION
                                                       Aug-31-2011
                          AMPAC Version 10.0.1
                             Presented by:
  
                        Semichem, Inc.
                        www.semichem.com
  
 FORMULA: C2H6
 Ethane rotational barrier
 Reaction Path calculation

     GEOMETRY OPTIMIZED : ENERGY MINIMIZED
     SCF FIELD WAS ACHIEVED
 
          FINAL HEAT OF FORMATION   =       -19.856323 kcal
                                    =       -83.098713 kJ
          ELECTRONIC ENERGY         =      -918.455084 eV
          CORE-CORE REPULSION       =       575.247482 eV
          TOTAL ENERGY              =      -343.207601 eV
          GRADIENT NORM             =         0.080077 
          RMS GRADIENT NORM         =         0.019421 
          UNSTABLE MODE(S)          =         0 ( ESTIMATE  )
          FOR REACTION COORDINATE   =       -60.000000 degrees       9
          REACTION GRADIENT         =        -0.001448 kcal/radian
          IONIZATION POTENTIAL      =        11.758166 eV
          HOMO-LUMO GAP             =        14.407233 eV
          DIPOLE                    =         0.000233 debyes
          MOLECULAR WEIGHT          =        30.069400 
          MOLECULAR POINT GROUP     = D3d     0.100000
          NO. OF FILLED LEVELS      =         7 (OCC = 2)
          TOTAL NUMBER OF ORBITALS  =        14
          COMPUTATION TIME          =         0.02     SECONDS

          FINAL GEOMETRY OBTAINED                                       CHARGE
 MINDO3 RHF SINGLET T=AUTO GRAD TRUSTE
 Ethane rotational barrier
 Reaction Path calculation
  C     0.000000  0    0.000000  0    0.000000  0      0     0     0    0.0761
  C     1.477194  1    0.000000  0    0.000000  0      1     0     0    0.0761
  H     1.111417  1  113.099625  1    0.000000  0      2     1     0   -0.0254
  H     1.111372  1  113.152439  1  119.968537  1      2     1     3   -0.0253
  H     1.111365  1  113.164965  1 -119.997299  1      2     1     3   -0.0253
  H     1.111415  1  113.100849  1  -60.000000 -1      1     2     3   -0.0254
  H     1.111373  1  113.151721  1  119.966091  1      1     2     6   -0.0253
  H     1.111366  1  113.164460  1 -119.994850  1      1     2     6   -0.0253
  0     0.000000  0    0.000000  0    0.000000  0      0     0     0
	

1

The value used in the input file is the first position to which the reaction coordinate is set, in this case 60°.

2

The reaction gradient is computed for each value of the reaction coordinate, and represents the slope of the energy profile along the reaction coordinate. The reaction gradient should approach zero and change sign when a transition state is crossed in the reaction coordinate. If a bond length is chosen as the reaction coordinate the units are in terms of kcal/Angstroms, otherwise they are in kcal/radian.

3

While the other values in the geometry have been optimized, the reaction coordinate is not altered. Due to this constraint, this is not a valid stationary point. The gnorm value, though very low for this geometry, does not reflect any contribution from unoptimized parameters such as the reaction coordinate.

4 5

The reaction coordinate for the next step corresponds to the first item in the list given in the extra data section.

6

Several points have been omitted, but the calculation is now at a 20° torsional angle. The reaction gradient is now quite significant and is negative.

7

The value of 0° for the reaction coordinate should theoretically be the very top in energy and correspond to the transition state for this process. The reaction gradient value is now almost zero, indicating that this point is very near the transition state. This is probably a good approximate geometry for the transition state, and the geometry here could now be gradient minimized using LTRD or another appropriate procedure.

8

The reaction coordinate has passed over the transition state and we are now proceeding downhill. Note that the reaction gradient has changed sign. In the present case, one of the reaction coordinate values was close to the transition state and a very small value was obtained. In some cases, a sudden change from a large negative to a large positive value can be observed in a reaction gradient between two steps. An approximate transition state geometry can be interpolated from this information or a reaction coordinate with smaller intervals could be attempted.

9

This is the final point of the reaction coordinate and the value corresponds to the last value in the list given in the extra data section.

Output File (rxnpath/rxnpath1.out) :

Due to its size, this file was not placed in the manual. It simply contains the full output as called for by the keywords in an .out file for each step in the reaction coordinate. A full copy may be found in the results directory of the test suite.