Exploring Potential Energy Surface with external forces - Journal of

Subscriber access provided by The University of Texas at El Paso (UTEP)

Quantum Electronic Structure

Exploring Potential Energy Surface with external forces Krzysztof Wolinski J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.8b00885 • Publication Date (Web): 31 Oct 2018 Downloaded from http://pubs.acs.org on November 4, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Exploring Potential Energy Surface with external forces

Krzysztof Wolinski

Department of Theoretical Chemistry

Maria Curie-Sklodowska University

pl. M.C. Sklodowska 3, 20-031 Lublin, Poland

e-mail [email protected]

1

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract We present a new method which allows automatic location of multiple local minima on the molecular Potential Energy Surface (PES). Standard geometry optimization procedures converge to a local minimum on the PES. They can be forced to penetrate other regions of the potential surface by applying external forces to nuclei in a molecule. One specific choice of such external forces is proposed in this paper, and preliminary results demonstrating its usefulness are presented.

Introduction A standard geometry optimizer locates one stationary point on the Potential Energy Surface (PES), usually the one which is in some sense closest to the starting geometry. To locate the absolute minimum is computationally hard but it is frequently less important than finding minima that are accessible with a relatively low energy barrier from a given initial structure. For smaller molecules this is usually done by repeating the geometry optimization starting with different initial structures, a procedure which is not suitable for the automatic location of possible reaction paths and products. Local minima on PES are separated by potential energy barriers. To overcome the barriers some work has to be done. This work can be done by external forces.

2


Page 2 of 47

Page 3 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


External forces have been widely used in mechanochemistry for the last forty years. Experimental techniques involving external forces have been exhaustively reviewed by Beyer and Clausen-Schaumann1. In the last twenty years a significant progress has been made also in the field of quantum mechanical description of mechanochemistry. There are theoretical methods which allow studies of molecular deformations caused by mechanical stress2–4. There are also methods for theoretical analysis of the stress-energy distribution in a molecule5–7. Artificial forces have been also used in order to induce chemical reactions8,9. External forces are also used in the Molecular Dynamics10 and similarly, different extra potentials are applied in Monte Carlo11 simulations in order to explore PES.

In 2009 we have proposed the Enforced Geometry Optimization (EGO)12,13 method for prediction of structural changes in molecules due to mechanical stress caused by external forces. At first this was done to simulate mechanochemical experiments14,15. The EGO method is closely related to three other models which have been proposed around the same time2–4. All these methods have been recently reviewed by Stauch and Dreuw16.

In the EGO method we perform standard geometry optimizations but in the presence of external forces acting on selected nuclei in a molecule. In principle, we can apply external forces to an arbitrary number of nuclei in arbitrary directions. In practice, we have applied, so far, external forces to one or more pairs of atoms pulling them apart or pushing them together. The EGO method turned out to be very succesful12,13,17–20 in describing conformational transitions in the pyranose ring as observed in AFM experiments14,15. Moreover, using the EGO method we have found new, previously unknown, high-energy but stable isomers of azobenzene21 and stilbene22. 3



Page 4 of 47

In those applications we have started the EGO calculations from a molecular geometry obtained with a standard geometry optimization (SGO). Thus, being at the local minimum on PES and using external forces in EGO we were able to leave that local minimum and obtain stressed structure which after relaxation (i.e. re-optimization without external forces) became a new local minimum on the original PES. However, the choice of applied external forces was not unique. We just put forces of a given magnitude on two atoms and observed the outcome. In this paper we propose an unique definition of external forces needed to move a molecule from one local minimum and cross a potential energy barrier separating it from another local minimum.

Theory The total energy of a molecule in the presence of constant external forces acting on nuclei can be written as 12,13

𝐸𝑡𝑜𝑡(𝒓) = 𝐸(𝒓) ― ∑𝐴𝒇𝑇𝐴 ⋅ 𝒓𝐴

(1)

where 𝐸(𝐫) is the energy of a system at the current geometry 𝐫 without external forces and 𝐟𝐴 denotes an external force applied to a nucleus 𝐴 at 𝐫𝐴 and 𝐫 = (𝐫1,𝐫2,…,𝐫𝑁) collects coordinates of all 𝑁 nuclei. The first derivative of the total energy (1) with respect to atomic coordinates (with a negative sign) gives the total force acting on atom 𝐴

𝒇𝑡𝑜𝑡 𝐴 =―

∂𝐸𝑡𝑜𝑡(𝒓) ∂𝒓𝐴

=―

∂𝐸(𝒓) ∂𝒓𝐴

(2)

+ 𝒇𝐴

4


Page 5 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


The total force 𝐟𝑡𝑜𝑡 𝐴 acting on nucleus 𝐴 is a sum of the internal force (negative energy gradient) and external force 𝐟𝐴. The potential energy surface modified by the external forces is called3 the force modified PES (FMPES). A molecular structure for which all total forces 𝐟𝑡𝑜𝑡 𝐴 are zero, corresponds to the local minimum (stationary point) on FMPES (1) but not on the original PES. Re-optimizing the geometry without external forces may give a new local minimum on the original PES.

The main issue here is the choice (definition) of external forces { 𝐟𝐴 } in (1) and (2) which will move us from the current minimum at 𝐫0 toward the potential energy barrier separating 𝐫0 from another local minimum on the original PES. To achieve this goal the chosen external forces should have the appropriate orientation (direction) and, in order to overcome the barrier, they should have the right magnitude. In contrast to (1) and (2) these external forces may not be, in general, constant.

The external forces that would move us from a current minimum toward an energy barrier should be adjusted somehow to the shape i.e. curvature of a current minimum. This kind of information is, of course, contained in the molecular Hessian 𝐇(𝐫0). What external forces could we derive from 𝐇(𝐫0)?

Such external forces should be very similar to the internal forces which will appear on atoms when we move just a little bit on PES from the current local minimum 𝐫0 to the closely located structure at 𝐫 = 𝐫0 + 𝐡. We have discussed this similarity or equivalence in our previous paper12. When a molecule is elongated from its equilibrium geometry (some bonds longer and fixed in 5



Page 6 of 47

Constrained Geometry Optimization (CGO)) then the internal forces are generated on nuclei. If these forces are applied at the equilibrium geometry as the external forces in EGO then the resulting structure is identical to that obtained in CGO. This is the main idea behind the COnstrained Geometries simulate External Force (COGEF) method2.

The internal forces at 𝐫 = 𝐫0 + 𝐡 can be estimated from the Taylor expansion around 𝐫0

1

𝐸(𝒓) = 𝐸(𝒓0) + 𝒈𝑇(𝒓0) ⋅ (𝒓 ― 𝒓0) + 2(𝒓 ― 𝒓0)𝑇 ⋅ 𝑯(𝒓0) ⋅ (𝒓 ― 𝒓0) + …

(3)

where 𝐠 and 𝐇 denote the energy gradient and Hessian, respectively, and the energy gradient at 𝐫0 is zero. In the harmonic approximation, i.e. neglecting cubic and higher order terms, the internal forces on atoms (negative energy gradient) at 𝐫 = 𝐫0 + 𝐡 can be in general written as

𝒇(𝒓) = ― 𝑯(𝒓0) ⋅ (𝒓 ― 𝒓0) = ― 𝑯(𝒓0) ⋅ 𝒉

(4)

Now we want to use the above formula to generate external forces. In principle this formula can be applied at every 𝐫 except 𝐫0where the displacement vector 𝐡 = 0. Therefore, we replace 𝐡 = 𝐫 ― 𝐫0 with an arbitrary displacement vector 𝐭 depending on 𝐫

𝒇(𝒓) = ― 𝐇(𝒓0) ⋅ 𝐭

(5)

It follows that applying external forces of +𝐇(𝐫0) ⋅ 𝐭 to nuclei at local minimum 𝐫0 shifts the minimum on a quadratic PES from 𝒓0 to the new minimum at 𝐫0 + 𝐭 on FMPES. Indeed, in the

6


Page 7 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


standard Newton-Raphson algorithm23 which converges in a single step on a quadratic surface, the new minimum is at

𝒓 = 𝒓0 + 𝑯 ―1(𝒓0) ⋅ [ + 𝑯(𝒓0) ⋅ 𝒕 ] = 𝒓0 + 𝒕

(6)

At the new minimum 𝐫0 + 𝐭 on FMPES the internal forces of ― 𝐇(𝐫0) ⋅ 𝐭 are canceled by the external forces of + 𝐇(𝐫0) ⋅ 𝐭. The situation is similar when we apply at 𝐫0 the external forces of ― 𝐇(𝐫0) ⋅ 𝐭. In this case we move from 𝐫0 to 𝐫0 ― 𝐭 which is also a minimum on another FMPES.

Thus, using at 𝐫0 external forces given by (5) we leave the local minimum as desired. However, the new issue now is the choice of the displacement vector 𝐭 which defines the external forces in (5). For a given, particular case one probably could figure out how to move atoms to achieve desired goal, however it would not be general. For more general case we may express the displacement vector 𝐭 in terms of the rescaled geometry vector 𝐫0 relative to the reference vector 𝐫𝑟𝑒𝑓

(7)

𝒕 = 𝑺0 ⋅ (𝒓0 ― 𝒓𝑟𝑒𝑓)

where 𝐒0 denotes the (3 ⋅ 𝑁 × 3 ⋅ 𝑁 )scaling (diagonal) matrix and 𝐫𝑟𝑒𝑓 is 3 ⋅ 𝑁 long the reference vector. The reference vector should be uniquely defined for a given structure and it should move with a molecule. The obvious choice is to refer coordinates of all atoms 𝐫 = (𝐫1,𝐫2 ,…,𝐫𝑁) to the center of mass 𝐫𝑐𝑜𝑚 so the relative geometry vector can be written as

7



Page 8 of 47

𝐫 ― 𝐫𝑟𝑒𝑓 = (𝐫1 ― 𝐫𝑐𝑜𝑚,𝒓2 ― 𝐫𝑐𝑜𝑚,…,𝐫𝑁– 𝐫𝑐𝑜𝑚)

In our case a molecule is always reoriented so its center of mass is placed at the origin. Hence, we have 𝐫𝑐𝑜𝑚 = (0, 0, 0) and 𝒓𝑟𝑒𝑓 = 0 and it will be omitted latter on.

Now the displacement vector is

(8)

𝒕 = 𝐒0 ⋅ 𝐫0

and external forces (5) at 𝒓0 can be written as

𝒇(𝒓0) = ― 𝑯(𝒓0) ⋅ 𝑺0 ⋅ 𝒓0

(9)

The diagonal scaling matrix 𝐒0 may be viewed as build out of 𝑁 atomic diagonal (3 × 3) matrices 𝐬𝑘 with s𝑘𝑥𝑥 , s𝑘𝑦𝑦 , s𝑘𝑧𝑧 diagonal elements for each atom 𝑘. In principle, one can imagine to use different 𝐬𝑘 matrices for each atom or different matrices for each type of atom or even the same matrix for all atoms. Thus there is some freedom here but the choice of the scaling matrix 𝐒0 may be crucial for the performance of the proposed method. It defines together with 𝐇(𝐫0), the orientation of the external forces (9). However, it is not so straight forward to define an unique scaling matrix for a general case. On the other hand, there is one unique vector which has an orientation specific for a given minimum 𝐫0 and it also represents forces. This is the 𝐇(𝐫0) ⋅ 𝐫0 vector. The external forces (9) which we want to apply at 𝐫0 may have the same specific orientation. May be it is not the best choice but we think that it is very reasonable. This can be achieved with the simplest possible scaling where 𝐒0 in (8) and (9) is just the unit matrix rescaled 8


Page 9 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


by a scalar parameter 𝑠: 𝐒0 = 𝑠 ⋅ 𝐈. With this setup the displacement vector 𝐭 is parallel to the geometry vector 𝐫0

(10)

𝒕 = 𝑠 ⋅ 𝒓0

and the external forces (9) given as

𝒇(𝒓0) = ―𝑠 ⋅ 𝑯(𝒓0) ⋅ 𝒓0

(11)

are parallel to 𝐇(𝐫0) ⋅ 𝐫0. This very simple scaling may limit the PES region accessible to external forces but it has an advantage of being quite general i.e. it can be apply in any case. Note that the orientation of the external forces (11) at a given minimum 𝐫0 depends only on its curvature as defined by the hessian matrix 𝐇(𝐫0).

These external forces will be applied at 𝐫0 in the first EGO optimization cycle. In every next EGO cycle at arbitrary 𝐫 we might use external forces generated the same way as

𝒇(𝐫) = ―𝑠 ⋅ 𝐇(𝐫0) ⋅ 𝐫

(12)

which is consistent with the total energy expression (FMPES)

1

𝐸𝑡𝑜𝑡(𝒓) = 𝐸(𝒓) + 2 ⋅ 𝑠 ⋅ 𝒓𝑇 ⋅ 𝑯(𝒓0) ⋅ 𝒓

9


(13)


Page 10 of 47

Note that in (12) and (13) for a structure at arbitrary 𝐫 the Hessian matrix 𝐇(𝐫0) is taken at the local minimum 𝐫0.

With the external forces (12) applied to atoms at the local minimum 𝐫0 the new geometry 𝐫1 obtained from 𝒓0in the first EGO cycle is

𝒓1 = 𝒓0 + 𝑯 ―1(𝒓0) ⋅ 𝒇(𝒓0) = 𝒓0 + 𝑯 ―1(𝒓0) ⋅ [ ―𝑠 ⋅ 𝑯(𝒓0) ⋅ 𝒓0] = (1 ― 𝑠) ⋅ 𝒓0

(14)

Note that two geometry vectors which differ only in sign (𝐫1 and ― 𝐫1), have all nuclear coordinates with opposite sign which corresponds to the inversion. Physically, they represent the same system and the energy at both geometries is the same. Thus, it is enough to consider in (14) only the case with (1 ― 𝑠) > 0 i.e. 𝑠 < 1. Now with s positive a molecule will shrink while with s negative will elongate.

Changing the sign of s will result in the change of the orientation of external forces (12) to opposite. It will also increase (with 𝑠 > 0) or decrease (𝑠 < 0) the total molecular energy on FMPES as compared to PES because the second term 𝐫𝑇 ⋅ 𝐇(𝐫0) ⋅ 𝐫 in (13) is always positive with the Hessian calculated at the minimum 𝐫0 (see the Supplementary Information for the proof).

It should be emphasized that in the EGO optimization, the Hessian matrix in (12) and (13) is the one calculated at the starting local minimum 𝒓0 on PES. This Hessian remains constant until the EGO optimization is done. The external forces (12) generated with this 𝐇(𝐫0) Hessian have a

10


Page 11 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


right orientation at the beginning of the EGO optimization when new molecular geometries are still close to the starting 𝐫0. However, when we move away from 𝐫0 this Hessian maybe less appropriate to generate correctly oriented external forces which may disturb the EGO optimization convergence.

The full Hessian H might be considered as a block build atomic pairs (3 × 3) matrices 𝐡𝐴𝐵 and then the external force acting on nucleus 𝐴 can be explicitly written as 𝐟𝐴 = ( 𝑓𝐴𝑥, 𝑓𝐴𝑦, 𝑓𝐴𝑧 ) with

𝑓𝐴𝑥

𝑓𝐴𝑦

= ―𝑠 ⋅

= ―𝑠 ⋅

𝑁

3

𝐵

𝑗

𝑁

3

𝐵

𝑗

∑ ∑𝐡 ∑ ∑𝐡 𝑁

𝐴𝐵 𝐵 𝑥𝑗 𝐫𝑗

𝐴𝐵 𝐵 𝑦𝑗 𝐫𝑗

(15)

3

𝐵 𝑓𝐴𝑧 = ―𝑠 ⋅ ∑𝐵 ∑𝑗 𝐡𝐴𝐵 𝑧𝑗 𝐫𝑗

where now the upper scripts denote nuclei and lower scripts denote Cartesian components 𝑥,𝑦,𝑧. It should be noted that the external force 𝐟𝐴 on a given nucleus 𝐴 depends on the location of ALL atoms so a force on each atom is correlated with positions of other atoms. In general, external forces on all atoms are different and they change in every optimization cycle.

Combined geometry optimization Let us now describe the entire procedure which will automatically locate a number of local minima on a molecular potential energy surface. Usually we construct preliminary molecular 11



structure using our GUI builder, the PQSmol program24. Using this structure in the standard geometry optimization (SGO) the first local minimum on PES is located. To conform that this is indeed a minimum we perform a vibrational analysis. The Hessian matrix from these calculations is used in the following EGO optimization according to (12) and (13). Once the EGO optimization is done and a local minimum on FMPES is found, external forces are switched off and SGO begins. The stressed EGO structure is relaxed and (hopefully) the next (second) local minimum on PES is found. For this structure the new Hessian matrix is calculated and then used in the following EGO optimization.

This procedure might be repeated as many times as one wants. Since the number of local minima on a given PES is not a priori known, a user may try to locate more local minima than actually exist on PES. In such a case we expect to obtain the same minimum more than once. If we begin our search from the first local minimum (found in ordinary SGO) then we need two consecutive geometry optimizations (EGO and SGO) to locate the next local minimum. Thus, we have to perform a sequence of SGO, EGO, SGO,…. geometry optimizations. Such combined optimization might be named the Standard and Enforced Geometry Optimization (SEGO) procedure.

As mentioned before, according to equation (14) with 0 < 𝑠 < 1, in the first EGO optimization cycle atomic coordinates get shorter and the molecule shrinks. The situation is different in all consecutive cycles:

12


Page 12 of 47

Page 13 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


∂𝐸(𝐫) [ { ∂𝐫 ] ― 𝑠 𝐇 ⋅ 𝐫 }}

𝐫𝑛 + 1 = 𝐫𝑛 + 𝐇𝑛―1 ⋅ ―

0

𝑛

𝑛

where the inverted Hessian 𝐇𝒏―1 (usually updated) corresponds to the current geometry 𝐫𝑛 while for external forces the starting 𝐇0(𝐫0) is still used. Moreover, the internal forces are also present.

Computational details All calculations carried out for these studies have been performed at the Density Functional Theory (DFT)25,26 level with the B3LYP27 exchange-correlation potential, using the 6-31g-dp basis set28. All calculations have been done in a parallel mode with the PQS program package24 and obtained molecular structures have been analyzed (visualized) using our Graphical User Interface, PQSMol.

In all EGO calculations the external forces were applied according to equation (12) and the maximum forces allowed on a given nucleus were limited to 0.15 - 0.20 au (1 au = 82400 pN). This last restriction defines the scaling factor 𝑠 in (12) and (13) as

𝑠 =

maximum force allowed maximum 𝐇(𝐫0) ⋅ 𝐫0 force on atom

where the maximum force allowed is set up by the user. All results reported in this paper have been obtained with the scaling factor 𝑠 being positive which means that molecules were shrunk in EGO (at least at the beginning). All our attempts to use negative 𝑠 (i.e. elongate molecules) failed: with small negative external forces elongated structures returned to the starting local 13



minimum after relaxation while with forces big enough some bonds were ruptured and molecules were destroyed .

Since we are interested in local minima on the original PES the extra energy term in equation (13) (arising from interactions of external forces with a molecule) was not included so we give the original PES energies in all tables, plots and figures.

To illustrate performance of the proposed method three molecules have been chosen:



butene (C4H8) ( five isomers),



dioxane(C4O2H8),



2,5-dihydroxyoxane (HO-[C5OH8]-OH).

Lewis structures of all these molecules are shown in the Supplementary Information.

Results and discussion For all three molecules considered here we will present the “history” of the combined geometry optimization (SEGO). The optimization history shows the energy changes during geometry optimization. Such a plot might be considered as a kind of “algorithmic reaction path” showing how the energy on ordinary PES changes when going from one to the next structure (optimization cycles). On the SEGO optimization history plot all local minima found can be seen as well as the energy barriers separating each two local minima. 14


Page 14 of 47

Page 15 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


As we mentioned before, one may try to locate a number of local minima, even more than exist for a given system. This information concerning the number of local minima to be found is given to the program through the input deck and the search for all requested local minima is done in one run. Thus, in some cases the same local minimum is found more than once. For all local minima found in one SEGO run, the energy and the lowest and the highest vibrational frequencies will be collected in separate tables.

We have considered here five isomers of butene : 1-butene in gauche conformation, cis-2-butene, isobutene, 1-butene in trans conformation and trans-2-butene. They are shown in Figures 1, 2, 3, 4 and 5 respectively. All of them are stable and correspond to local minima on PES except the trans conformer of 1-butene which is a saddle point on PES with one imaginary frequency. Energetically, isobutene is the lowest and 1-butene (trans) the highest isomer. In the case of butene we tried to locate 6 minima on the ordinary PES.

In Plot 1 the history of the SEGO optimization, starting from the structure close to gauche 1butene, is shown. The search for local minimum begins with SGO and ends at cycle 2 where the first local minimum (1-butene gauche) was found. Then there are consecutive EGO (external forces on) and SGO (external forces off) optimizations which begin and end as follows:

EGO cycles: 3-106, SGO cycles: 107-123,


EGO cycles: 249-346, SGO cycles: 347-377, 15




EGO cycles: 502-601, SGO cycles: 602-624

Six local minima corresponding to the stable structures (all vibrational frequencies real) have been found in this SEGO search no. 1 at

cycle

2

123

248

377

502

624

structure

1

2

3

4

5

6

These local minima are separated by five energy barriers at cycles: 13, 145, 271, 399 and 525.

A closer inspection shows that the structures 2,4 and 6 are the same as well as structures 3 and 5. Thus, we have found here 3 unique structures 1, 2(4,6) and 3(5). Moreover, the structures 2,4,6 are mirror images of the structures 3,5. All structures are shown in Figure 1.

We have performed the next SEGO search for multiple local minima starting from cis-2-butene (search no. 2). The optimization history is shown in Plot 2. Six local minima can be seen there at

cycle

2

137

250

377

492

619

structure

1

2

3

4

5

6

16


Page 16 of 47

Page 17 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


The structures corresponding to these local minima are shown in Figure 2. Among these, the structures number 1 and 4 are the same as well as number 3, 5 and 6. In this case we have found also 3 unique, stable structures corresponding to local minimum number 1(4), 2, and 3(5,6).

Plot 3 shows the SEGO geometry optimization history when the search for multiple minima was started from isobutene (search no. 3). The local minima have been found at

cycle

2

104

238

395

523

637

structure

1

2

3

4

5

6

and corresponding structures are shown in Figure 3. The structures number 4 and 6 are the same. Thus, in this case, 5 unique structures have been found: 1, 2, 3, 4(6) and 5. Note that the last structure number 5 is the cis-2-butene isomer.

Plot 4 presents the geometry optimization history for the SEGO search number 4 which started from the structure close to 1-butene in trans conformation. Six stationary points on the ordinary PES are there at

cycle

4

124

242

365

489

606

structure

1

2

3

4

5

6

The first structure is trans 1-butene (Cs symmetry) and it is not a local minimum. It is a saddle point with one imaginary frequency. All five remaining structures correspond to local minima. 17



The second structure is the gauche 1-butene isomer. Structures number 3 and 5 are the same as well as structures 4 and 6 and 3(5) are mirror images of 4(6). In this case the SEGO search has found one saddle point (1) and three unique local minima 2, 3(5) and 4(6). All structures are shown in Figure 4.

In the last SEGO search shown in Plot 5, the trans-2-butene isomer (C2h symmetry) has been found first. All six structures corresponding to the local minima are unique. They are in Plot 5 at

cycle

6

125

234

365

477

599

structure

1

2

3

4

5

6

The second structure is isobutene and the sixth structure is the gauche 1-butene conformer. All structures are shown in Figure 5.

To summarize: in all five SEGO searches for multiple local minima the following stable structures have been found:

searching from 1-butene gauche : 1, 2, 3 Figure 1.

searching from cis-2-butene : 1, 2, 3 Figure 2.

searching from isobutene : 1, 2, 3, 4, 5 Figure 3.

searching from 1-butene trans : 1, 2, 3 Figure 4.

18


Page 18 of 47

Page 19 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


searching from trans-2-butene : 1, 2, 3, 4, 5, 6 Figure 5.

The energies and the lowest and highest vibrational frequencies of all stable structures found in five SEGO searches are collected in Table 1. In this table we use two digit numbering: the first digit – search number, the second – structure number in this search. In Table 1 unique structures in a given SEGO search are separated by the dot lines. It can be seen that some of the stable structures found in these SEGO searches are the same (see the Supplementary Information for details). Overall in the case of butene, in all five SEGO searches, we have found ten unique stationary points on PES:

1.1, 1.2, 1.3 (Figure 1) – local minimum

2.1, 2.2, 2.3 (Figure 2) – local minimum

3.1, 3.2 (Figure 3) – local minimum

4.1 (Figure 4) – saddle point

5.1 (Figure 5) – local minimum

Among these, seven structures retain the butene character with one double C=C bond. Three structures number 2.2, 2.3 (Figure 2) and 3.2 (Figure 3) do not have any double bond anymore. The structure 2.3 is the 1-methyl-cyclopropane while 2.2 and 3.2 contain the divalent carbon (carbenes) and, to our best knowledge, they are not known. These two strange structures are the highest in energy (-157.108353 and -157.121683 au, respectively) and according to our results 19



they are predicted to be stable at the B3LYP level (Table 1). To verify stability of these isomers we have run the MP2 calculations as well. It turned out that the 3.2 structure is indeed stable but the 2.2 structure does not exist at the MP2 level. The geometry optimization at the MP2/6-31gdp level brings it to 1-methyl-cyclopropane (structure 2.3). In spite of the MP2 results for the 3.2 strange structure it is still possible that it is an artifact of the single-reference method.

In the case of dioxane with three known conformers we tried to locate six local minima. The SEGO optimization history is shown in Plot 6. Six local minima are present there at

cycle

4

108

198

291

384

477

structure

1

2

3

4

5

6

The structures corresponding to these minima are shown in Figure 6. The structures 2,4 and 6 are the same as well as structures 3 and 5. Thus, we have found here 3 unique structures 1, 2(4,6) and 3(5). Moreover, the structures 2,4,6 are mirror images of the structures 3,5. The energies and the lowest and highest vibrational frequencies are given in Table 2. Structures 1 and 3 are those which we expected to find : chair and 1,4 twisted boat. However, we did not find the 2,5-twisted boat.

The last molecular system considered in this study was the 2,5-dihydroxyoxane molecule, the pyranose ring in the chair conformation with two hydroxy groups in the axial orientation (HO[C5OH8]-OH). It is known from AFM experiments14,15 and theoretical studies12,13,19,20 that such a

20


Page 20 of 47

Page 21 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


system undergoes inversion into another chair with hydroxy groups in equatorial orientation when a molecule is elongated along a line connecting two linkage oxygens.

In this case we tried to locate ten local minima staring from the structure close to the chair conformation with two -OH groups in axial positions. The SEGO optimization history is shown in Plot 7. Ten local minima are present there at

cycle

9

115

257

389

526

596

661

769

835

891

structure

1

2

3

4

5

6

7

8

9

10

The first local minimum obtained in ordinary SGO corresponds to the chair conformation with two -OH groups in axial orientation. The second structure is the 1,4 twisted boat. The third one is an inverted chair with equatorial -oriented -OH groups. The next structure is an another twisted boat. Going from this local minimum (no 4) to the next one (no 5) the external forces cause the ring opening and all subsequent structures from 5 to 10 are different zig-zag type chains

Among these ten local minima only the structures corresponding to minimum no 8 and 10 are the same. All others are different. They are shown in Figure 7. The energies and the lowest and highest vibrational frequencies are given in Table 3.

21



Conclusions In this paper we proposed a new method for automatic location of multiple local minima on the Potential Energy Surface (PES). This method uses external forces to move over PES from one local minimum to another. The essence of this method is to perform a sequence of the standard (SGO) and enforced (EGO) geometry optimizations. The first one (SGO) locates a minimum on ordinary PES. Starting from that minimum the second (EGO) optimization with external forces finds a local minimum on force-modified PES (FMPES). The structure corresponding to this minimum on FMPES is then projected on ordinary PES (i.e. external forces are switched off) and the next SGO optimization hopefully finds a new local minimum on ordinary PES. The only condition which has to be satisfied is that external forces used in the EGO step are large enough to cross the potential energy barrier separating two local minima on PES. Since the proposed approach is based on a sequence of the standard (SGO) and enforced (EGO) geometry optimizations it has been named SEGO.

The crucial issue for performance of the SEGO method is the choice of the external forces used to penetrate PES. External forces are generated in SEGO using the energy Hessian calculated at the current minimum 𝐫0 on PES. At arbitrary molecular geometry 𝐫 these forces are given as a product of the Hessian matrix 𝐇(𝐫0) ⋅ and the geometry displacement vector 𝐡 = 𝐫 ― 𝐫0 or 𝐭 = 𝑠 ⋅ 𝐫0 . In principle, one could try to use both geometry displacements to generate external forces (except the very first EGO cycle starting from 𝐫0 where 𝐡 = 0). Of course, external forces

22


Page 22 of 47

Page 23 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


defined as 𝐟(𝐫) = ― 𝐇(𝐫0) ⋅ 𝐡 or as 𝐟(𝐫) = ― 𝐇(𝐫0) ⋅ 𝐭 are different in magnitude and orientation. In this study both approaches have been tested.

There was no problem with EGO optimizations when external forces of 𝐟(𝐫) = ― 𝐇(𝐫0) ⋅ 𝐭 were applied. However, there was a big problem in EGO when we tried to apply external forces of 𝐟 (𝐫) = ― 𝐇(𝐫0) ⋅ 𝐡. With the displacement vector 𝐡 = 𝐫 ― 𝐫0 we have started the EGO optimization from 𝐫0 using first 𝐭 = 𝑠 ⋅ 𝐫0 and then switched to 𝐡 after a few cycles. It turned out that external forces calculated with the displacement vector 𝐡 = 𝐫 ― 𝐫0 are NOT appropriate for location of multiple local minima. If they are big enough to do something they destroy (break apart) a molecule! It looks like the 𝐡 = 𝐫 ― 𝐫0 vector reorients the external forces 𝐟(𝐫) = ― 𝐇

(𝐫0) ⋅ 𝐡 in a wrong direction. A similar situation occurred when we tried to use external forces generated with a negative scaling factor 𝑠 in (12) i.e. a molecule was first elongated instead of shrunk. As mentioned before external forces applied in a “pulling” mode also destroy a molecule. This is perhaps not a surprise because it is easier to break bonds by pulling nuclei apart than to push them together.

Thus, all results reported in this paper have been obtained with external forces (12) generated using the displacement vector 𝐭 being the rescaled geometry vector (10) and with a positive scaling parameter 𝑠.

The preliminary results presented in this studies show that the SEGO method is quite effective. A number of unique local minima can be found automatically in one SEGO search (done in one

23



run). Some of the structures corresponding to these minima were unexpected and probably previously unknown (e.g.: 2.3 in Figure 2, 3.2 in Figure 3). However, it should be emphasized that not all expected local minima have been found which means that our external forces penetrate a limited part of PES.

Results of the SEGO search depend on the starting structure. Beginning the SEGO calculations with five different isomers of butene, multiple local minima were found in different order, some of them the same and some of them the different, which is not a surprise.

The SEGO outcome depends first of all on the magnitude and orientation of external forces used. External forces have to do work against the potential energy barrier (PEB) which separates two local minima on PES. If external forces are too small then PEB is not crossed in the EGO step and after switching them off, following SGO returns to the starting structure. Moreover, with two different in magnitude but big enough external forces the SEGO search may go along different paths.

With a given Hessian 𝐇(𝐫0) the magnitude of external forces is determined by the scaling parameter 𝑠. So far, the scaling in our calculations has been done in the simplest manner according to equations (10) and (12): all three coordinates of all atoms are rescaled by the same scalar factor 𝑠. In such a case the orientation of external forces (12) is defined only by the Hessian matrix. However, one can imagine more sophisticated approach where, for instance the 𝑥 ― and 𝑦 ― coordinates are scaled less than 𝑧 ― coordinate with 𝑠𝑥𝑥 = 𝑠𝑦𝑦 < 𝑠𝑧𝑧 . This could be used for a long chain molecule oriented along 𝑧 axis. Such scaling would change also the

24


Page 24 of 47

Page 25 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


orientation of external forces and might result in sampling a larger region of PES. We are going to try this approach in a near future.

In order to modify both, the magnitude and orientation of external forces, one may also play with the Hessian matrix. The results presented in this paper were obtained with external forces calculated according to equation (12) with the full Cartesian Hessian 𝐇(𝐫0). However, we have also implemented and tried to use atom-diagonal and full diagonal Hessian. The very preliminary results obtained this way are quite promising.

The efficiency of the SEGO optimizations as presented in this paper is not great. In most cases over 120 cycles were needed to go from one local minimum to another. The SEGO search consists of two consecutive EGO and SGO steps which needed about 100 and 20 cycles, respectively (out of 120). Thus, overall low efficiency of SEGO is due to quite slow convergence in the EGO step with external forces. Indeed, in the optimization history Plots 1-7 almost in all cases one can see a long energy plateau at the end of EGO (for instance EGO cycles 550-610 in Plot 1). This is more likely caused the fact that far away from an initial minimum at 𝐫0, the external forces generated with 𝐇(𝐫0) (12) do not have the right orientation anymore. However, this EGO efficiency can be significantly improved by ending the EGO optimization after the energy barrier is crossed i.e. at the beginning of a plateau, before a minimum on FMPES is found. It can be done by limiting the number of allowed optimization cycles. We have tested such an approach on the 1-butene isomer. We have limited the number of EGO cycles to 60 and we have obtained exactly the same local minima as presented before in Plot 1. However, before the whole SEGO search needed 624 optimization cycles while now it was reduced to 418. We 25



included the corresponding SEGO optimization history plot in the Supplementary Information section.

Finally, we would like to rise up one more issue concerning location of saddle points on PES. We have already pointed out in our previous papers12,13,22,22 that the EGO method is also very useful in location of transition states. A molecular structure which corresponds to the maximum energy on the EGO optimization history plot provides usually a very good starting guess for a transition state search. Within the SEGO method the situation is even better because the optimization history plots contain a number of local minima separated by energy barriers. Each structure corresponding to the energy maximum can be used as a guess in location of the transition state connecting two consecutive local minima. In principle, in one SEGO run it should be possible to find, let say, 𝑛 local minima and 𝑛 ― 1 transition states. Moreover, it could be done in full automatic manner. The work along this line is in progress.

Acknowledgments This research was supported in part by the U. S. National Science Foundation (NSF), Grant No. DMR-1609650 and the Arkansas Biosciences Institute, the major research component of the Arkansas Tobacco Settlement Proceeds Act of 2000.

References (1) Beyer, M. K.; Clausen-Schaumann, H. Mechanochemistry: The Mechanical Activation of Covalent Bonds. Chem. Rev. 2005, 105, 2921–2948.

26


Page 26 of 47

Page 27 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(2) Beyer, M. K. The Mechanical Strength of a Covalent Bond Calculated by Density Functional Theory. J. Chem. Phys. 2000, 112, 7307–7312. (3) Ong, M. T.; Leiding, J.; Tao, H.; Virshup, A. M.; Martínez, T. J. First Principles Dynamics and Minimum Energy Pathways for Mechanochemical Ring Opening of Cyclobutene. J. Am. Chem. Soc. 2009, 131, 6377–6379. (4) Ribas-Arino, J.; Shiga, M.; Marx, D. Understanding Covalent Mechanochemistry. Angew. Chem. Int. Ed. 2009, 48, 4190–4193. (5) Stauch, T.; Dreuw, A. A Quantitative Quantum-Chemical Analysis Tool for the Distribution of Mechanical Force in Molecules. J. Chem. Phys. 2014, 140, 134107. (6) Grandbois, M.; Beyer, M. Rief, M. How Strong Is a Covalent Bond? Science 1999, 283, 1727–1730. (7) Bofill, J. M.; Ribas-Ariño, J.; García, S. P.; Quapp, W. An Algorithm to Locate Optimal Bond Breaking Points on a Potential Energy Surface for Applications in Mechanochemistry and Catalysis. J. Chem. Phys. 2017, 147, 152710. (8) Maeda, S.; Morokuma, K. Communications: A Systematic Method for Locating Transition Structures of A+B→X Type Reactions. J. Chem. Phys. 2010, 132, 241102. (9) Maeda, S.; Morokuma, K. Finding Reaction Pathways of Type A + B → X: Toward Systematic Prediction of Reaction Mechanisms. J. Chem. Theory Comput. 2011, 7, 2335– 2345. (10) Suan Li, M.; Khanh Mai, B. Steered Molecular Dynamics-A Promising Tool for Drug Design. Curr. Bioinforma. 2012, 7, 342–351. (11) Shang, C.; Liu, Z.-P. Stochastic Surface Walking Method for Structure Prediction and Pathway Searching. J. Chem. Theory Comput. 2013, 9, 1838–1845. (12) Wolinski, K.; Baker, J. Theoretical Predictions of Enforced Structural Changes in Molecules. Mol. Phys. 2009, 107, 2403–2417. (13) Wolinski, K.; Baker, J. Geometry Optimization in the Presence of External Forces: A Theoretical Model for Enforced Structural Changes in Molecules. Mol. Phys. 2010, 108, 1845–1856.

27



(14) Marszalek, P. E.; Oberhauser, A. F.; Pang, Y.-P.; Fernandez, J. M. Polysaccharide Elasticity Governed by Chair-Boat Transitions of the Glucopyranose Ring. Nature 1998, 396, 661–664. (15) Marszalek, P. E.; Pang, Y.-P.; Li, H.; Yazal, J. E.; Oberhauser, A. F.; Fernandez, J. M. Atomic Levers Control Pyranose Ring Conformations. Proc. Natl. Acad. Sci. U. S. A. 1999, 96, 7894–7898. (16) Stauch, T.; Dreuw, A. Advances in Quantum Mechanochemistry: Electronic Structure Methods and Force Analysis. Chem. Rev. 2016, 116, 14137–14180. (17) Cybulska, J.; Brzyska, A.; Zdunek, A.; Woliński, K. Simulation of Force Spectroscopy Experiments on Galacturonic Acid Oligomers. PLoS ONE 2014, 9, e107896. (18) Brzyska, A.; Wolinski, K. Enforced Conformational Changes in the Structural Units of Glycosaminoglycan (Non-Sulfated Heparin-Based Oligosaccharides). RSC Adv. 2014, 4, 36640–36648. (19) Brzyska, A.; Płaziński, W.; Woliński, K. Force-Induced Structural Changes in NonSulfated Carrageenan Based Oligosaccharides – a Theoretical Study. Soft Matter 2018, 14, 6264–6277. (20) Wolinski, K.; Brzyska, A. Theoretical Studies of the Pyranose Ring under Mechanical Stress. Carbohydr. Res. 2018. (21) Baker, J.; Wolinski, K. Isomerization of Stilbene Using Enforced Geometry Optimization. J. Comput. Chem. 2011, 32, 43–53. (22) Baker, J.; Wolinski, K. Kinetically Stable High-Energy Isomers of C14H12 and C12H10N2 Derived from Cis-Stilbene and Cis-Azobenzene. J. Mol. Model. 2010, 17, 1335–1342. (23) Nocedal, J.; Wright, S. J. Numerical Optimization, 2nd ed.; Springer series in operations research; Springer: New York, 2006. (24) Baker, J.; Wolinski, K.; Janowski, T.; Saebo, S.; Pulay, P. PQS Version 4.0, Parallel Quantum Solutions,; Parallel Quantum Solutions P.O. Box 293 Fayetteville Arkansas 72702-0293 U.S.A. (25) Hohenberg, P.; Kohn, W. Inhomogeneous Electron Gas. Phys. Rev. 1964, 136, B864–B871.

28


Page 28 of 47

Page 29 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(26) Kohn, W.; Sham, L. J. Self-Consistent Equations Including Exchange and Correlation Effects. Phys. Rev. 1965, 140, A1133–A1138. (27) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648. (28) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. A. Self-Consistent Molecular Orbital Methods. XX. A Basis Set for Correlated Wave Functions. J. Chem. Phys. 1980, 72, 650.

1.1

1.2

1.3

1.4

29



1.5

Page 30 of 47

1.6

Figure 1. Structures from 1-butene (gauche conformation).

2.1

2.2

2.3

2.4

30


Page 31 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


2.5

2.6

Figure 2. Structures from cis-2-butene.

3.1

3.2

3.3

3.4

31



3.5

Page 32 of 47

3.6

Figure 3. Structures from isobutene.

4.1

4.2

4.3

4.4

32


Page 33 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


4.5 Figure 4. Structures from 1-butene (trans conformation).

4.6

5.1

5.2

5.3

5.4

33



Figure 5. Structures from 5.5 trans-2-butene.

34


Page 34 of 47

5.6

Page 35 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

2

3

4

5

6

35



Page 36 of 47

Figure 6. Structures from dioxane.

1

2

3

4

5

6

36


Page 37 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


7

8

9

10

Figure 7. Structures from pyranose ring.

37



Plot 1. Geometry optimization history for 1-butene (gauche conformation).

Plot 2. Geometry optimization history for cis-2-butene.

38


Page 38 of 47

Page 39 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Plot 3. Geometry optimization history for isobutene.

Plot 4. Geometry optimization history for 1-butene trans.

39



Plot 5. Geometry optimization history for trans-2-butene.

40


Page 40 of 47

Page 41 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Plot 6. Geometry optimization history for dioxane.

Plot 7. Geometry optimization history for pyranose ring.

41



42


Page 42 of 47

Page 43 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 1. Local minima on the Potential Energy Surface of butene found in the SEGO search with maximum external force allowed 0.175 au. Energy [au] and vibrational frequencies [cm-1] Structure

Energy

number

Vibrational

Frequencies

lowest

highest

SEGO search no 1 starting from 1-butene gauche 1.1

-157.232200728

163.55

3239.83

1.2

-157.232800953

108.38

3231.18

1.4

-157.232800447

109.34

3231.52

1.6

-157.232800214

108.61

3231.28

1.3

-157.232800983

109.08

3231.32

1.5

-157.232803495

109.40

3231.22

2, 4, 6 are mirror images of 3.,5

SEGO search no 2 starting from cis-2-butene 2.1

-157.236411768

134.50

3156.60

2.4

-157.236410613

132.99

3156.65

2.2

-157.108353217

226.45

3196.25

2.3

-157.223716831

225.21

3223.77

2.5

-157.223719566

224.93

3223.93

2.6

-157.223715260

223.67

3223.98

176.83

3231.50

SEGO search no 3 starting from isobutene 3.1

-157.238840556

43



Page 44 of 47

3.2

-157.121683302

101.65

3146.90

3.3

-157.108356341

226.12

3196.49

3.4

-157.223717334

224.37

3223.81

3.6

-157.223718048

225.85

3223.90

3.5a

-157.236408541

133.69

3156.64

SEGO search no 4 starting from 1-butene trans 4.1b

-157.229478186

-127.49

3230.21

4.2c

-157.232198738

162.92

3239.77

4.3

-157.232801948

109.11

3231.08

4.5

-157.232801205

109.23

3231.33

4.4

-157.232802605

109.04

3231.38

4.6

-157.232803806

109.13

3231.28

4, 6 are mirror images of 3, 5 ; structure no1 has Cs symmetry

SEGO search no 5 starting from trans-2-butene 5.1d

-157.238525479

179.95

3136.63

5.2e

-157.238840990

177.18

3231.48

5.3

-157.121684631

102.44

3146.73

5.4

-157.108355240

226.76

3196.63

5.5

-157.223718687

223.84

3223.77

a

b c d e cis-2-butane; 1 butene (trans); 1-butene(gauche); trans-2-butene; isobutene

44


Page 45 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 2. Local minima on the Potential Energy Surface of dioxane found in the SEGO search with maximum external force allowed 0.150 au. Energy [au] and vibrational frequencies [cm-1]

Structure

Energy

number

Vibrational

Frequencies

lowest

highest

1

-307.668706903

253.68

3099.33

2

-307.657705516

73.83

3081.38

4

-307.657709347

75.40

3081.18

6

-307.657709397

75.43

3081.19

3

-307.657704627

73.90

3081.41

5

-307.657702933

72.11

3081.16

2, 4, 6 are mirror images of 3, 5

45



Table 3. Local minima on the Potential Energy Surface of 2,5-dihydroxyoxane (pyranose ring) found in the SEGO search with maximum external force allowed 0.200 au. Energy [au] and vibrational frequencies [cm-1]

Structure

Energy

number

Vibrational

Frequencies

lowest

highest

1

-422.224333605

126.83

3799.64

2

-422.211394181

55.07

3807.16

3

-422.219331854

110.61

3816.27

4

-422.214089994

84.27

3818.73

5

-422.221492834

34.91

3799.58

6

-422.224503666

46.98

3814.12

7

-422.223044652

25.38

3802.21

8

-422.221441670

51.56

3797.69

10

-422.221439833

51.92

3797.83

9

-422.227958188

56.15

3666.75

46


Page 46 of 47

Page 47 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47



Exploring Potential Energy Surface with external forces - Journal of

Recommend Documents