Calculating electrostatic potential using "ESPGrid" file

I want to calculate electrostatic potential for my systems using a manually provided grid information in an external file. According to the manual this can be done using setting the ESP_GRID =n, where n is the number of grid points in the ESPGrid file.

I encountered a possible bug where the grid printed in the plot.esp file does not match the grid provided in the ESPGrid file. Here is the input and output grid:

ESPGrid:

-1.87308530    0.54366160    2.60058023
-1.34808530   -0.36566506    2.60058023
-1.08190665    2.34380368    1.83192688
-2.01404996    1.91810810    1.83192688
-2.56807041    1.05603546    1.83192688
-2.56807041    0.03128774    1.83192688
-1.56775556    2.50719571    0.78192688
-2.39495787    1.93621918    0.78192688
-2.86206311    1.04622450    0.78192688
-2.86206311    0.04109871    0.78192688

plot.esp

Grid positions and esp values from SCF density (all in a.u.)

       X            Y            Z           GS

  -0.17718237E+01   0.52919838E+01  -0.77511187E+00  -0.61382639E-01
  -0.23984743E+01   0.44931714E+01  -0.50732899E+00  -0.49403268E-01
   0.32499777E+00   0.51178931E+01  -0.94928234E+00  -0.97358514E-01
  -0.39227236E+00   0.56349953E+01  -0.14671982E+01  -0.65959505E-01
  -0.13648564E+01   0.56733032E+01  -0.17876988E+01  -0.48467693E-01
  -0.22839658E+01   0.52206545E+01  -0.18090276E+01  -0.38633580E-01
   0.56652072E+00   0.50511552E+01  -0.20905355E+01  -0.60713087E-01
  -0.24328775E+00   0.54246994E+01  -0.25541656E+01  -0.45553206E-01
  -0.12096365E+01   0.53849249E+01  -0.28277819E+01  -0.35759615E-01
  -0.21111467E+01   0.49409435E+01  -0.28487022E+01  -0.28941210E-01

Input file:

$molecule
0 1
    N       0.84399998       0.45600000       0.92299998
    H       1.37199998       0.49500000       1.79799998
    C       1.64300001       1.12500000      -0.11300000
    H       2.00000000       2.09699988       0.24300000
    H       0.98699999       1.30700004      -0.97100002
    C       2.77900004       0.15300000      -0.50400001
    H       3.68899989       0.38400000       0.05900000
    H       3.03099990       0.22200000      -1.56599998
    C       2.23099995      -1.24000001      -0.11200000
    H       2.86599994      -1.70700002       0.64700001
    H       2.18499994      -1.92799997      -0.96200001
    C       0.82999998      -0.94099998       0.46900001
    H       0.54200000      -1.60899997       1.28600001
    H       0.06800000      -1.02699995      -0.31400001
$end

$rem
   METHOD     hf
   BASIS      6-31g**
   IANLTY     200
   ESP_GRID   10        !643
$end

$plots
   plot the electrostatic potential on a line
   1   0.0   0.0
   1   0.0   0.0
  1  0.0   0.0
   0  0  0  0
   0
$end

Please let me know if there is any mistakes in the input file or in the formatting of the ESPGrid file.

Please add the following to your $rem:

$rem
no_reorient true
sym_ignore true   ! or point_group_symmetry = false
....other stuff...
$end

What NO_REORIENT=TRUE will do is to make sure that Q-Chem uses the user’s coordinate system (in $molecule), rather than rotating into the standard nuclear orientation (which is the default behavior). Since your ESPGrid is presumably made in the user’s frame, that’s the only choice that makes sense. Turning off symmetry avoids an occasional problem when you turn off reorientation for a molecule with nontrivial point group symmetry.

This is pretty subtle and I will add a cautionary note to the manual.

That solved the problem, thank you so much.

For future reference, I checked in a fix for Q-Chem 6.2.2 so that one can set NO_REORIENT = TRUE without messing around with the symmetry keywords. And the cautionary note in the manual should also appear with v. 6.2.2.

I am able to run the ESP calculations now, but I have a question about the output.

Electrostatic potential will be computed using the SCF density
   ... ESP will be evaluated on a read-in grid points from file ESPGrid.
   ... The data will be saved in the file plot.esp
  --- DONE ESP SUMMARY ---
 Total job time:  6.68s(wall), 29.41s(cpu)
 Mon Nov 25 15:35:12 2024

The --- DONE ESP SUMMARY --- line suggests there should be some kind of summary for the ESP calculation, say number of grid points, but the output just ends. Is there a keyword that triggers the ESP summary printing or this is how it’s supposed to be?

That is how it is supposed to be. The output is in the file plot.esp and also in some cube files, if you request them.