2.19. MedeA MOPAC: Semiempirical Quantum Chemistry
Contents
| download: | pdf |
|---|
MOPAC only handles periodic systems without symmetry that are “sufficiently large” - typically on the order of 8 \({\mathring{\mathrm{A}}}\) in any direction. Please ensure your cell dimension is large enough. If you have a non-periodic system, you can right-click >> Create a periodic copy to create box, then right-click >> Edit cell to adjust the box size.
The user interface to MOPAC 2009/2012/2016 is based on flowcharts.
Once you have an active system, either a molecule or a periodic system, you bring up the interface by selecting Jobs >> New Job.
2.19.1. MOPAC: Semiempirical QM Stage
In the flowchart interface, clicking MOPAC: Semiempirical QM under the Methods section adds a MOPAC stage, you can double-click on the MOPAC box to open the flowchart. Once you have opened the MOPAC stage, you will find the following stages for its flowchart:
2.19.1.1. Initialization
Hamiltonian
Use this optional stage to select which Hamiltonian to use in MOPAC for the rest of the calculation. The choices are AM1, MNDO, MNDOD, PM3, PM6, PM7, and RM1, with the default of PM7. See the theory section for a discussion of the available Hamiltonians.
2.19.1.2. Energy
This section provides energy-based stages, such as those for calculating the energy or optimizing the structure.
Single Point Energy
Hartree-Fock wavefunction: automatic, RHF and UHF. This parameter defines the type of wavefunction to use in the SCF. automatic is the default, using an RHF wavefunction if there is an even number of electrons; UHF, for an odd number unless the type of calculation. RHF (restricted Hartree Fock) and UHF(unrestricted Hartree Fock) override the default, explicitly requesting such wavefunction. Note that a UHF wavefunction may not represent a pure spin state, e.g. a triplet or a singlet, whereas RHF, by construction, does.
SCF Convergence Scheme: automatic, Pulay or Camp-King. This option controls the method applied to converge the SCF equations. The automatic converger uses different approaches based on issues with the convergence. The Pulay or Camp-King convergers can be faster than the automatic approach, but may not be as robust. Also, they are not applicable in large systems when using the MOZYME functionality in MOPAC.
SCF Convergence: default 0.0001 kcal/mol The SCF is considered successfully converged when the energy difference between two successive iterations is smaller than this amount. If the requested value is too tight, i.e. too small, the SCF may never be considered to be converged. If too large, the results will not be accurate. This stage is also used in the Phonon module. See MedeA Phonon 2.0 section for more details.
Single Point Forces
This stage is used in the Phonon module. See MedeA Phonon 2.0 section for more details. The parameters in this stage are the same as in the Single Point Energy stage.
Structure Optimization
This stage handles optimizing the structure by minimizing the energy. The dialog for this stage has two main tabs: Optimization Control and SCF Control. The SCF Control tab has the same options described above for the Single Point Energy calculation. Optimization Control has the following options.
Minimizer: automatic, EF, BFGS, LBFGS and NLLSQ. automatic requests the EF (Eigenvector Following) UV/Vis Control minimizer for systems with up to 2000 degrees of freedom, and L-BFGS (limited memory Broyden-Fletcher-Goldfarb-Shanno) minimizer for larger systems. The other options explicitly request a given minimizer. EF explicitly requests the Eigenvector Following minimizer and can be useful to force its use on systems with more than 2000 degrees of freedom. BFGS requests the Broyden-Fletcher-Goldfarb-Shanno minimizer, which uses a quasi-Newton method, storing the inverse Hessian matrix in memory. Thus, the amount of memory required for larger systems may be prohibitive. LBFGS requests the limited memory Broyden-Fletcher-Goldfarb-Shanno minimizer, which recalculates the Hessian matrix as needed, so avoids the memory that the standard BFGS method needs. LBFGS is the default method for systems with more than 2000 degrees of freedom. NLLSQ invokes the Non-Linear Least Squares gradient minimization method, which will locate a minimum, transition state, or other points on the surface where the gradient norm is a minimum - which may not be a minimum or a transition state. The automatic and EF minimizers use two parameters to control how the Hessian is initialized and recalculated.
Starting Hessian: approximate and calculated. approximate, the default, uses an approximate Hessian matrix to start the minimization. For calculated the Hessian matrix is calculated using single-sided derivatives. This can be expensive for larger systems.
Recalculate hessian every #n steps in EF method: the default is 0. This option instructs MOPAC to recalculate the Hessian every so often during the minimization. If this parameter is set to 0 or left blank the Hessian will not be recalculated. As noted above, calculating the Hessian can be time-consuming, but may be useful for difficult minimizations.
For all the minimization methods, the following options are available:
Convergence Criterion: normal, precise or manual. normal, as its name implies, requests the normal convergence criterion for the minimization, which varies from minimizer to UV/Vis Control minimizer. precise tightens up the convergence criterion by a factor of 100, and also tightens up the associated SCF criteria appropriately, unless they are overridden by other parameters. manual allows you to specify the criterion, which defaults to 1.0 kcal/mol/\({\mathring{\mathrm{A}}}\). You may set the value to anything from 0.0 up; by default values less than 0.01 will be rounded up to 0.01 kcal/mol/\({\mathring{\mathrm{A}}}\).
Maximum iterations: default 10,000. This parameter limits the maximum number of steps that the minimization will allow, which is useful to prevent it from iterating forever if the convergence criteria are too tight.
2.19.1.3. Properties
IR/Raman Frequencies
This stage calculates the vibrational - infrared and Raman - frequencies and intensities for the molecule. The results can be graphically analyzed using the Analysis >> IR / Raman Spectra command.
This stage provides the following options:
IR/Raman Control: normal and precise. normal uses standard values for convergence criteria and employs single-sided differences to obtain the force constant matrix. precise tightens up the convergence criteria and uses two-sided differences to calculate the force constant matrix. This is more expensive and seldom needed.
You can check the Do not use symmetry box to prevent symmetry analysis of the molecule in the force constant calculations. This option is for expert users.
Note
Some artificial imaginary frequency vibrations may disappear when using this option. This leads however to more complex IR/Raman peak assignments, when trying to compare computed and experimental spectra.
The SCF Control options are identical to those for the Single Point Energy.
UV/Vis
This stage computes UV/Vis frequencies and intensities for molecules. The results can be graphically analyzed using the Analysis >> IR / Raman Spectra command. The dialog for this stage has two main tabs: UV/Vis Control and SCF Control
The UV/Vis Control stage provides the following options: CI level: CIS and CISD. CIS requests to include ground state and microstates from single-electron excitation in the UV/Vis spectrum calculation. CISD requests to include ground state and single and double electron excitations Number of orbitals in CI: allows you to indicate the number of configurations that will be used to compute the UV/Vis spectra. You can check the Ignore degenerate orbitals split by space box, while the default is to consider them (expert users). Precision and SCF Control have the same meanings and options than before.
Thermodynamics
This stage computes the force constants and derives the thermochemical properties from them: heats of formation, constant pressure heat capacities, and entropies at the requested temperature(s) are printed in Job.out. When computing these properties, the structure must be properly optimized before through a Structure Optimization stage: this is further checked when the mention System is a ground state is reported in Job.out of this stage or a prior IR/Raman Frequencies stage.
The temperatures for which the thermochemical properties are reported are declared in Thermo Control, by providing Initial temperature, Final temperature, and Temperature step. These temperatures can be given in different units. Other calculation parameters are the precision (precise or normal), whether or not to enforce the symmetry (Do not use symmetry box), and SCF Control. They have the same meaning as in the IR/Raman Frequencies stage.
2.19.1.4. Custom
This stage and the Extra input are for expert users. Custom Stage permits to edit the full MOPAC command line. All stages in Energy and Properties have also an Extra input banner. It permits to add extra keyword(s) to the existing command line of the corresponding stage. Be careful that it does not allow you to edit the full command line: extra input keywords are not always compatible with those of the stage. For instance, you cannot add keyword requesting structure optimization in a Single Point Energy stage, which explicitly requests a single point energy calculation. Use Custom Stage instead. Note that additional computed properties not accessible in the Energy and Properties stages are not reported in Job.out but in Stage_x/mopac.out file.
2.19.1.5. IR/Raman example
The flowcharts illustrate a simple calculation using Mopac, the IR/Raman spectrum of 1,4,5,8-Naphthalenetetracarboxylic dianhydride, called NTCDA for short.
Since this is a simple calculation, the main flowchart to the left is also very simple. The flowchart for the MOPAC, on the right, is slightly more complex. The Hamiltonian is set to PM6. The structure is then minimized carefully, using the precise criterion.
Finally, the IR/Raman spectrum is computed using default parameters, which are quite reasonable for most calculations.
After the calculation completes, the main information is contained in Job.out:
Stage 1: This Mopac calculation has 3 stages:
Stage 1.1: Using PM6 Hamiltonian for a subsequent stage
Stage 1.2: Minimize using PM6 Hamiltonian
Automatically choosing the minimizer (EF for small systems,
LBFGS for large) to precise convergence.
SCF convergence: 0.0001 kcal/mol
Empirical Formula Formula
C7O3H2 (C7O3H2)2
--------------- ---------------
Mopac Energy -1758.124 -3516.248 eV
Heat of formation -356.711 -713.422 kJ/mol
Dipole moment 0.001 Debye
Ionization potential 11.051 eV
HOMO energy -11.051 eV
LUMO energy -3.168 eV
COSMO area 240.380 Ang^2
COSMO volume 262.740 Ang^3
Point group D2h
This stage took 7.589 s for 410 SCFs
Stage 1.3: IR/Raman frequencies using PM6 Hamiltonian
and normal precision.
SCF convergence: 0.0001 kcal/mol
Mode Frequency (1/cm) Symmetry Intensity
------ ------------------ ---------- -----------
1 67.1 1Au 0.000
2 92.6 1B3u 0.859
3 96.8 1B1g 0.000
4 150.1 1B2g 0.000
5 150.5 2B3u 0.036
6 208.2 1B2u 1.564
7 219.6 3B3u 0.249
8 257.9 2B2g 0.000
In this example, this information includes not only thermodynamic data such as the heat of formation at 298 K, but also the vibrational modes along with their symmetry and intensity.
As was mentioned above, the results, the IR/Raman spectrum can also be analyzed using Analysis >> IR/Raman Spectra command. This brings up a list of previously calculated spectra. After you select a spectrum, it is brought up in a new window in MedeA:
The molecule is in the left pane. You can adjust the size of the panes using the small square button near the bottom of the “sash”. In the right pane is the calculated spectrum, in red. You can import experimental data using the Import… button. In this example, an imported experimental spectrum is shown with the black line. At the bottom of the window are three sliders, Scale factor, Line width and Intensity scaling. These allow you to adjust the calculated spectrum to match the experiment. Scale factor uniformly scales the frequency. Since Mopac frequencies are typically a little bit too high, scaling the frequencies by 0.9 to 1.0 will improve the fit. Likewise, the Line width slider is used to scale the width of the calculated lines. They are arbitrarily broadened with either a Gaussian or Lorentzian function, with a line width given by the slider. To change the broadening function used, click the More… button. The Intensity scaling slider adjusts the intensity of the calculated spectrum, typically in arbitrary units.
The legend at the middle left side of the graph is active. Clicking on the items will turn the display of the corresponding elements on and off. By default, a “stick” spectrum is not displayed, and the broadened line spectrum is. Finally, across the top of the graph are short blue lines that give all the fundamental frequencies, whether they are active or not. If you hover the mouse over one of the short lines, the frequency and intensity are display just below the graph. Right-clicking on one of the short lines provides a menu, which includes, among other things, the ability to animate the mode. The molecule in the left pane will vibrate appropriately for the mode; however, the amplitude is arbitrary and intentionally made quite large to help you visualize the mode.
| download: | pdf |
|---|