I have a problem with geometry optimization in Gaussian software with M062X method. The structure was smoothly optimized using M062X/Gen. Nevertheless, I need very accurate results so afterwards the job was submitted once more, but this time with Opt=VeryTight and Int=UltraFine keywords. Nevertheless, I got the error:

Optimization stopped.
    -- Number of steps exceeded,  NStep= 126
    -- Flag reset to prevent archiving.

... and at the end:

Error termination request processed by link 9999.
Error termination via Lnk1e in /usr/local/gaussian/2016-C.01/g16/l9999.exe at Tue Aug 25 16:05:14 2020.

What can I do? The input is following and, as I told, geometry is taken directly from the previous step which was performed using the same method and basis but without Opt=VeryTight and Int=UltraFine keywords.

I'm stuck, don't know what to do. Any idea?

#p M062X/Gen
opt=(Cartesian,maxcycles=200,verytight,restart) freq
geom=check guess=read


0 1

H     0
S   3   1.00
     34.0613410              0.60251978D-02
      5.1235746              0.45021094D-01
      1.1646626              0.20189726
S   1   1.00
      0.32723041             1.0000000
S   1   1.00
      0.10307241             1.0000000
P   1   1.00
      0.8000000              1.0000000
P   1   1.00
      0.95774129632D-01      1.0000000
C     0
S   6   1.00
  13575.3496820              0.22245814352D-03
   2035.2333680              0.17232738252D-02
    463.22562359             0.89255715314D-02
    131.20019598             0.35727984502D-01
     42.853015891            0.11076259931
     15.584185766            0.24295627626
S   2   1.00
      6.2067138508           0.41440263448
      2.5764896527           0.23744968655
S   1   1.00
      0.57696339419          1.0000000
S   1   1.00
      0.22972831358          1.0000000
S   1   1.00
      0.95164440028D-01      1.0000000
S   1   1.00
      0.48475401370D-01      1.0000000
P   4   1.00
     34.697232244            0.53333657805D-02
      7.9582622826           0.35864109092D-01
      2.3780826883           0.14215873329
      0.81433208183          0.34270471845
P   1   1.00
      0.28887547253           .46445822433
P   1   1.00
      0.10056823671           .24955789874
D   1   1.00
      1.09700000             1.0000000
D   1   1.00
      0.31800000             1.0000000
D   1   1.00
      0.90985336424D-01      1.0000000
F   1   1.00
      0.76100000             1.0000000
O     0
S   6   1.00
  27032.3826310              0.21726302465D-03
   4052.3871392              0.16838662199D-02
    922.32722710             0.87395616265D-02
    261.24070989             0.35239968808D-01
     85.354641351            0.11153519115
     31.035035245            0.25588953961
S   2   1.00
     12.260860728            0.39768730901
      4.9987076005           0.24627849430
S   1   1.00
      1.1703108158           1.0000000
S   1   1.00
      0.46474740994          1.0000000
S   1   1.00
      0.18504536357          1.0000000
S   1   1.00
      0.70288026270D-01      1.0000000
P   4   1.00
     63.274954801            0.60685103418D-02
     14.627049379            0.41912575824D-01
      4.4501223456           0.16153841088
      1.5275799647           0.35706951311
P   1   1.00
      0.52935117943           .44794207502
P   1   1.00
      0.17478421270           .24446069663
P   1   1.00
      0.51112745706D-01      1.0000000
D   1   1.00
      2.31400000             1.0000000
D   1   1.00
      0.64500000             1.0000000
D   1   1.00
      0.14696477366          1.0000000
F   1   1.00
      1.42800000             1.0000000

Here is .xyz geometry:

O         -4.44721        1.17772       -0.00003
C         -3.68389        0.06486       -0.00000
O         -4.16201       -1.03886        0.00003
C         -2.24755        0.39366       -0.00001
C         -1.34975       -0.59256        0.00001
C          0.10468       -0.46469        0.00000
C          0.87630       -1.62398       -0.00001
C          2.26287       -1.55445       -0.00001
C          2.88570       -0.32236       -0.00001
O          4.23767       -0.13929       -0.00001
C          2.12513        0.85396        0.00001
O          2.73759        2.06167        0.00002
C          0.74994        0.77741        0.00001
H         -5.36979        0.88538       -0.00002
H         -1.98671        1.44280       -0.00003
H         -1.74093       -1.60620        0.00003
H          0.38777       -2.58968       -0.00001
H          2.86062       -2.45869       -0.00002
H          4.69415       -0.98608       -0.00002
H          3.69353        1.92694        0.00001
H          0.19248        1.70441        0.00002
  • $\begingroup$ A couple of things I'm noticing right off the bat. 1. You say the previous step was done with the same method and basis set, but your question mentions B3LYP/def2-TZVPD while the input file is for M06-2X/def2-TZVPD. 2. You have maxcycles=200, but your error message says it stopped after only 126. 3. You have opt=restart, but this would be to restart a failed job; you should only need geom=check to extract the successful geometry. Assuming those don't help to fix your issue, it would be useful to see how close the geometry came to converging (search for Item at the last step). $\endgroup$ – Tyberius Aug 25 '20 at 16:19
  • $\begingroup$ Oh, right... 1) Just copied .log of different file. Nevertheless, the main idea is the system converges at any DFT functional unless I specify keywords Opt=verytight & Int=ultrafine. Then, it fails. Changed to M06-2X to keep it consistent with .log 2) But 126th is step, not cycle. Steps are definied that way (from Gaussian doc): MaxStep=N ... Sets the maximum size for an optimization step (the initial trust radius) to 0.01N Bohr or radians. The default value for N is 30, isn't it? 3) Yes, I have been using geom=check and guess=check and it still fails. That's why I used opt=restart $\endgroup$ – farmaceut Aug 25 '20 at 17:49
  • $\begingroup$ For the example .log file, these are values of convergence: Maximum Force - 0.000006 - 0.000002 - NO.....RMS - Force - 0.000002 - 0.000001 - NO....Maximum Displacement - 0.000638 - 0.000006 - NO.......RMS - Displacement - 0.000129 - 0.000004 - NO.. As you can see they are rather close to converge.. $\endgroup$ – farmaceut Aug 25 '20 at 17:54
  • $\begingroup$ According to this researchgate.net/post/… "This means that the Gaussian job terminated abnormally in some fashion. Typically, it means that a geometry optimization has not converged." $\endgroup$ – Cody Aldaz Aug 25 '20 at 18:53
  • 1
    $\begingroup$ @Martin-マーチン Would it be wise to check unpruned grid, Int=(Grid=199974)? Well, that's some kind of result the functional doesn't work though it should. M06 are widely used when describing compounds I am studying. Isn't 500 cycles too much, wouldn't it pinpoint the geometry is definitely wrong? How exactly new guess influences results? No, it has not. $\endgroup$ – farmaceut Aug 25 '20 at 19:33

TL;DR: There really isn't a straightforward way to fix this error, or a way to guarantee success. This usually involves a very detailed hands-on troubleshooting session, because everything is dependent on everything and even tiny screws could lead to a significant change and often the underlying issues are worse than toggling a switch (or adding a keyword).

What does the error Error termination request processed by link 9999. actually mean?

Gaussian often isn't very straight forward with its errors. There are multiple explanation sites of common errors on the internet because of this. The description of it is actually given many lines above. If you look up link 9999 on the gaussian online manual, then you'll find:

  • L9999: Finalizes calculation and output

So that's the something-went-wrong-catch-all of Gaussian. That usually means Gaussian's internal procedures set a flag to prevent it from terminating normally.

In your specific case (and I experienced none other) it simply means the optimiser didn't find a stationary point within the allocated resources.

Anyone can reproduce this error quickly with the following minimal example:

#P PM6 opt(maxcycle=3) 

L9999 error producing input

O 0. 0. 0.
H 1. 0. 0.
H 0. 1. 0.

! Blank line at the end

Obviously, the above example forces the error for demonstration purposes only. BTW, the ! is the comment character; these lines will be ignored. (The last line is there because SE decides to remove crucial blank lines automagically. So the example is copy-pastable with the comment.)

Now, Gaussian isn't the best at communicating. You'll find (as you did) the following:

 Optimization stopped.
    -- Number of steps exceeded,  NStep=   3
    -- Flag reset to prevent archiving.
                           ! Non-Optimized Parameters !
                           ! (Angstroms and Degrees)  !
 --------------------------                            --------------------------

There it tells you that at NStep= 3 the number of steps is exceeded (more to that end later), and that'll reset some flag to prevent something. That's actually the internal procedure to cause the final message.

Steps and Cycles. Apparently Gaussian isn't very consistent in its terminology. Steps may be called cycles. So the keyword option maxcycle= sets the maximum NStep.

Unfortunately, the L9999 error isn't as much an error as it is a message to you: Basically the program is trying to tell you that it failed to reach what you have asked of it. The underlying issue is somewhere in the set-up of your calculation.

Fixing the 'error'

You need to perform a deep-dive into the nitty-gritty details of your calculation.

In the most common case, your starting guess simply was too far away from the stationary point. You can plot the convergence with each cycle to see what is going on. Most molecular viewers have that built in.
In the best case, you are a few cycles away from convergence. You can then simply restart the calculation. I recommend setting up a new calculation to prevent accidental overwriting. For the example above, this would be (in the same directory, requires > G09 D.01):

#P PM6 opt(maxcycle=300, restart) 

L9999 error restart

O 0. 0. 0.
H 1. 0. 0.
H 0. 1. 0.

! Blank line at the end

So increasing the maxcycle value might already fix the 'error'. Note that the default is 128, so you need to set something higher than that if you restart. Another common way is to simply extract the last coordinates and start a new calculation.
In most cases on tight optimisation criteria and a reasonably well-behaved molecule, the default cycles will suffice. You can always save time by running pre-optimisations first. I recommend starting with semi-empirics, then increase to a pure functional with a split valence basis set. The new Gaussian versions in conjunction with a queueing system can handle chaining those calculations easily (use %oldchk), so you don't waste time.

Everything beyond that is magic.

Often, when you get to know your molecules, you will also know whether the default number of cycles will suffice and you can adjust your default route section.

Your example is a bit more than the default as you are trying to converge to verytight. The first question you need to ask yourself is why you are doing this and whether this is really necessary. In my personal and very subjective experience, it takes about thrice the time (and cycles) to converge from tight to verytight as it takes from loose to tight.

In order to get to convergence in fewer cycles, tighten the integration grid. The default is int(ultrafine) so you do not need to specify this. If you are running on a smaller grid, then your calculations are not reliable (Get Ready to Recalculate by Derek Lowe). If you are doing benchmark tests, and you want them to be benchmarks for the method not the grid, then use int(superfine). It might also be a good idea to force tighter SCF convergence with scf(conver=9).

The estimation or calculation of the step and step direction is directly influenced by the quality of the wave function/ density matrix. If that's too coarse, you likely circle around the minimum.

You can adjust the step size to force it to go slower, opt(maxstep=10), the default is 30. Whenever you do that, you should increase the cycles appropriately. An okay rule of thumb is by halving the step size to double the cycles. This is usually a good idea when trying to converge transition states.

Sometimes, the estimation of the Hessian is not good enough. You might want to start your calculation with something better than a guess, then use opt(calcFC). Alternatively, you can run a frequency calculation of your guess structure and analyse that before doing anything else. This might help identifying lower lying issues in the molecular structure. If you do, opt(RCFC) will be your friend.
If you really, really want it, and you know exactly what you are doing (or you don't care anymore), then calculating force constants at every step is the brute-force method, i.e. opt(calcall). I do not recommend this, it is very, very wasteful. There are intermediate options to do these things, so study the manual for more.

I always recommend using redundant coordinates. Often enough the default is good enough. In rare cases they won't. The evil ones are near 0 or 180° angles, sometimes rings, cages, strongly delocalised parts of the molecule, long bonds, dispersive interactions, etc., the list is long. You can build your own redundant coordinates and that might help in those cases.
If you are not yet doing it, you should switch to generalised internal coordinates (requires G16): geom=GIC. If you need examples, browse on Chemistry.SE.

If you really, really want to add more brute to your force, you can try optimising in cartesian coordinates. Note that you will not take advantage of many features of the Berny algorithm, you'll introduce more constraints, and you will need more cycles. (For a well-behaved molecule I'd add about 25% to be safe, but that is again as subjective as it gets.)

If you are still convinced that what you are looking for exists, you can try a different algorithm. I have never used it, but I heard about people using it.

The final failure

In most cases, if you get the L9999 error and it cannot easily be fixed (and I would include increasing the grid here), there is a deep-lying issue. You should check whether other methods produce sensible results. If your optimisations with very tight criteria succeed at BP86, PBE, TPSS, B3LYP, PBE0, TPSSh, etc., but not on M06-2X, then the problem might be that the Minnesota functional doesn't describe your system well enough. (I've only picked that because of the example. B3LYP is much more likely to fail producing sensible results. Statistically speaking.)
The usual way to go about that is to look for benchmarks studies in that field. Many things have been done already and you don't have to reinvent the wheel. (If only I could remember some of the quality review articles. I'm sure [or at least hopeful] someone will comment them below.)

The post-final Hoorah!

There are alternative programs to Gaussian. It might not be a terrible idea to cross-validate.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.