Intrinsic Reaction Coordinate/Minimum Energy Pathway (MEPGS) Calculation

Clicked A Few Times

12:24:38 PM PDT - Mon, Jun 3rd 2013

Hello Folks,

I am running into problems getting an IRC calculation to run and was hoping someone else might have some insight.

I am receiving an error that "The energy has increased" and the calculation dies before doing the first point of the MEP/IRC. I have tried decreasing the step size (stride) to no avail. I have tried tightening up the energy convergence criteria in the TS calculation. At this point I am considering tightening up the DFT grid size, but that has more ramifications since grids need to be consistent between calculations to compare energies.

Thanks for your help!

Input deck:

TITLE molecule

geometry noautosym
 C                    -1.56583532     0.89980740    -0.46106375
 C                    -0.63078461     1.82698720    -0.86157945
 C                    -0.50773429     2.24918622    -2.30919751
 H                    -1.44949302     2.07714793    -2.84402532
 H                     0.28275637     1.67258199    -2.80899191
 H                    -0.24479137     3.30808189    -2.40178582
 C                    -1.55906112     0.29368306     0.83927529
 H                    -2.21346146     0.46986517    -1.22690280
 S                     0.57452173     2.28434060     0.27251061
 C                    -0.11188030    -0.95552821     1.12771631
 H                    -2.39624750    -0.36907737     1.06390518
 H                    -1.28782623     0.97753389     1.66409184
 S                     1.38219886    -0.04057327     0.92765129
 C                    -0.32367076    -2.07715447     0.09224193
 C                    -0.39215468    -1.42587090     2.56874464
 H                    -0.25842731    -0.61046728     3.28498429
 H                     0.29935385    -2.23912992     2.84821336
 H                    -1.41677948    -1.81613178     2.66570304
 H                    -0.15719796    -1.71031957    -0.92379821
 H                    -1.34409850    -2.48456759     0.15829454
 H                     0.37792543    -2.90958221     0.27310725
end

basis spherical
 C library 6-31G*
 H library 6-31G*
 S library 6-31G*
END

driver
 maxiter 200
end

dft
 xc b3lyp
 iterations 2000
end

task dft freq

freq
 reuse molecule.hess
end

mepgs
 maxmep 300
 maxiter 200
 inhess 2
 noxyz
 evib 0.0005
 stride 0.01
 backward
 mswg
end

task dft mepgs

Last few lines of output file:

                 ----------------------------------------
                 |  Time  |  1-e(secs)   |  2-e(secs)   |
                 ----------------------------------------
                 |  CPU   |       0.01   |       2.91   |
                 ----------------------------------------
                 |  WALL  |       0.01   |       3.07   |
                 ----------------------------------------

          Expected decrease in energy              0.000500

          Obtained decrease in energy             -0.005967

 -------------------------
 The energy has increased
 -------------------------
 ------------------------------------------------------------------------
 dft mepgs failed                                                                        0
 ------------------------------------------------------------------------
 This type of error is most commonly associated with calculations not reaching convergence criteria
 ------------------------------------------------------------------------

Gets Around

7:18:33 PM PDT - Sun, Aug 4th 2013
Your molecule.hess has imaginary frequency mode? Try to set FORWARD direction.

Clicked A Few Times

1:08:55 PM PST - Tue, Feb 18th 2014

Hello Vladimir,

Here is the weird thing. The vibrational analysis says there is an imaginary mode (as you can see we run the frequency and pipe that into the IRC):

----------------------------------------------------------------------------
 Normal Eigenvalue ||           Projected Infra Red Intensities
  Mode   [cm**-1]  || [atomic units] [(debye/angs)**2] [(KM/mol)] [arbitrary]
 ------ ---------- || -------------- ----------------- ---------- -----------
    1     -295.810 ||    0.036177           0.835        35.267      20.733
    2       -0.000 ||    0.000044           0.001         0.043       0.025
    3       -0.000 ||    0.000015           0.000         0.014       0.009
    4       -0.000 ||    0.000013           0.000         0.012       0.007
    5        0.000 ||    0.000039           0.001         0.038       0.022
    6        0.000 ||    0.000037           0.001         0.036       0.021
    7        0.000 ||    0.000022           0.001         0.022       0.013
    8      124.486 ||    0.002547           0.059         2.483       1.460
    9      146.712 ||    0.000637           0.015         0.621       0.365
   10      185.619 ||    0.001857           0.043         1.810       1.064
   11      207.290 ||    0.000024           0.001         0.024       0.014
   12      236.593 ||    0.001386           0.032         1.351       0.795
   13      273.500 ||    0.003481           0.080         3.393       1.995
   14      308.056 ||    0.005062           0.117         4.935       2.901

As for setting it to forward direction, sometimes it gives me the same "energy increased" error, other times it works fine.

Gets Around

12:11:19 AM PDT - Sat, May 24th 2014
Can anybody point me to the documentation about mepgs. I can't find any. So mswg - mass weight coordinstes flip - ??? forward - forward irc direction backward - backward irc direction ircboth - ??? (both directions) evib - Desired decrease in energy when following the imaginary normal mode away from a saddle point. (default=0.0001 Hartree) ??? stride - Determines how far apart points on the reaction path will be. STRIDE is used to calculate the step taken. (default = 0.1 sqrt(amu)Bohr) ??? nircopt - number of optimization steps ??? maxmep - max optimization steps ??? nptopt - ??? eprec - ???? Sometimes I've got error: norm not preserved What's meaning of this error? What's meaning of 2-nd geometry in .bxyz and .fxyz files? There is no corresponding energies in `grep \& .nwo`. 1-st geometry is initial. 3-rd is a "point 0" of IRC The answer can be given in the terms of http://www.cup.uni-muenchen.de/ch/compchem/geom/irc1.html

Gets Around

3:47:24 AM PDT - Sat, May 24th 2014
evib = 0.0 is possible and 1st and 2nd geometry is equal in that case

Forum >> NWChem's corner >> NWChem functionality