Suggested procedure for locating bugs

Users are supplied with the source code for MOPAC, and, while the original code is fairly bug-free, after it has been modified there is a possibility that bugs may have been introduced. In these circumstances the author of the changes is obviously responsible for removing the offending bug, and the following ideas might prove useful in this context.

First of all, and most important, before any modifications are done a back-up copy of the standard MOPAC should be made. This will prove invaluable in pinpointing deviations from the standard working. This point cannot be over-emphasized--make a back-up before modifying MOPAC!

Clearly, a bug can occur almost anywhere, and a logical search sequence is necessary in order to minimize the time taken to locate it.

If possible, perform the debugging with a small molecule, in order to save time (debugging is, of necessity, time consuming) and to minimize output.

The two sets of subroutines in MOPAC, those involved with the electronics and those involved in the geometrics, are kept strictly separate, so the first question to be answered is which set contains the bug. If the heats of formation, derivatives, I.P.s, and charges, etc., are correct, the bug lies in the geometrics; if faulty, in the electronics.

 

Bugs discovered in the electronics subroutines

Use formaldehyde for this test. Use keywords 1SCF, DEBUG, and any others necessary.

The main steps are:

1.
Check the starting one-electron matrix and two-electron integral string, using the keyword HCORE. It is normally sufficient to verify that the two hydrogen atoms are equivalent, and that the π system involves only pz on oxygen and carbon. Note that numerical values are not checked, but only relative values.

If an error is found, use MOLDAT to verify the orbital character, etc.

If faulty the error lies in READMO, GETGEO or MOLDAT.

Otherwise the error lies in HCORE, H1ELEC or ROTATE.

If the starting matrices are correct, go on to step (2).

2.
Check the density or Fock matrix on every iteration, with the words FOCK or DENSITY. Check the equivalence of the two hydrogen atoms, and the π system, as in (1).

If an error is found, check the first Fock matrix. If faulty, the bug lies in ITER, probably in the Fock subroutines FOCK1 or FOCK2, or in the (guessed) density matrix (MOLDAT). An exception is in the UHF closed-shell calculation, where a small asymmetry is introduced to initiate the separation of the α and β UHF wavefunctions.

If no error is found, check the second Fock matrix. If faulty, the error lies in the density matrix DENSIT, or the diagonalization RSP.

If the Fock matrix is acceptable, check all the Fock matrices. If the error starts in iterations 2 to 4, the error probably lies in CNVG, if after that, in PULAY, if used.

If SCF is achieved, and the heat of formation is faulty, check HELECT. If C.I. was used check MECI.

If the derivatives are faulty, use DCART to verify the Cartesian derivatives. If these are faulty, check DCART and DHC. If they are correct, or not calculated, check the DERIV finite difference calculation. If the wavefunction is non-variationally optimized, check DERNVO.

 

Bugs discovered in the geometric subroutines

If the geometric calculation is faulty, use FLEPO or PRNT=5 to monitor the optimization, DERIV may also be useful here.

For the FORCE calculation, DCART or DERIV are useful for variationally optimized functions, COMPFG for non-variationally optimized functions.

For reaction paths, verify that FLEPO is working correctly; if so, then PATHS is faulty.

For saddle-point calculations, verify that FLEPO is working correctly; if so, then REACT1 is faulty.

Keep in mind the fact that MOPAC is a large calculation, and, while intended to be versatile, many combinations of options have not been tested. If a bug is found in the original code, please communicate details to Dr. James J. P. Stewart, Stewart Computational Chemistry, 15210 Paddington Circle, Colorado Springs, CO 80921-2512; E-mail: openmopac@gmail.com.