QM/MM optimization

We are trying to run large QM/MM calculation (8,000 atoms) and are only able to were only able to do that using Cartesian coordinates. Unless GEOM_OPT_COORDS=0, Q-Chem crashes with “integer overflow”. Why is that?

Also, when it does run, I see that the total time per optimization cycle is much longer that the wall time for energy and gradient calculation. Does anyone know where the time is spent? Is there something wrong with our setup?

This sounds like a problem that should be reported to Q-Chem support (support@q-chem.com). Please include the full input and output files in your report.

John Herbert explained to me what the issue is:

  1. To work with delocalized internal coordinates, one has to generate a giant matrix (Wilson B-matrix) and diagonalize the metric Bt*B. The dimension of the matrix is several times more than NAtoms. So this becomes impractical for large systems. Hence, one can only work with the Cartesian coordinates.
  2. Q-Chem’s default algorithm for optimization is a guasi-Newton algorithm that requires diagonalizing a Hessian-size matrix at each step. Obviously, this becomes impractical for very large systems. Hence, when NAtoms is very large, one should switch to the steepest descent L-BFGS algorithm by using GEOM_OPT_ALGORITHM = -1, courtesy of then Bhaskar Rana and John Herbert.

Bottom line, for our system with ~8,000 atoms, using
GEOM_OPT_COORDS=0
GEOM_OPT_ALGORITHM = -1
works. However, optimization does take quite a lot of steps, which is expected for these algorithms.

The options in Anna’s reply,
GEOM_OPT_COORDS=0
GEOM_OPT_ALGORITHM = -1
are steepest descent in Cartesian coordinates. It avoids the diagonalization step that is needed to generate delocalized internal coordinates, with is intractable when number of atoms is large (as in typical QM/MM). Alternative is limited-memory (L-) BFGS, still in Cartesians but quasi-Newton instead of steepest descent. (This is Bhaskar Rana’s contribution.)

A very simple example of the setup for L-BFGS QM/MM-optimization (water dimer with one QM water) is:

$molecule
0 1
O 1.3584299158 0.1073418692 -0.2758823010 101 2 3 0 0
H 0.3884464033 -0.0182409613 -0.1252830261 88 1 0 0 0
H 1.7255084907 -0.4596377323 0.4452211001 88 1 0 0 0
O -1.3286731559 -0.2344591932 0.1344130752 101 5 6 0 0
H -1.7784222051 0.6341523014 0.2664866218 88 4 0 0 0
H -1.7690342253 -0.5495694315 -0.6907120483 88 4 0 0 0
$end

$rem
basis sto-3g
exchange hf
qm_mm_interface janus
force_field charmm27
user_connect true
jobtype opt
no_reorient true
symmetry false
$end

$forceman
QMMM-LBFGS
LBFGS_M 10
$end

$qm_atoms
1 2 3
$end

1 Like