EFP-EFP polarization damping result in polarization blow-up

Questions
#1

I am writing to update an discrepancy I came across when calculating energies for an ionic liquid-based dianion system. The QM region was modeled using wb97x-D/def2-SVPD, while 25 pairs of ionic liquid ions were modeled using the EFP approach. EFP parameters were computed as recommended in the Q-Chem and GAMESS manuals.

I observed unphysical polarization energy values in my system, which resulted in improper SCF convergence and incorrect molecular orbital eigenvalues. Despite trying several different approaches, including changing basis sets, altering initial guesses, and modifying the creation of EFP parameters, the underlying problem persisted.

The issue was traced back to the use of Gaussian-type polarization damping. According to literature, it is advisable to reduce the damping parameters (α and β) from the standard 0.6 for systems with ionic nature. However, in Q-Chem, these parameters are hardcoded.

To address this, I turned off the polarization damping. This led to quicker SCF convergence and correct molecular orbital calculations. The effect of damping was compared across different numbers of EFP pairs (1, 5, 10, and 25), showing that the impact of damping increases drastically with the number of EFP pairs.

Below is a comparison of the polarization energy with and without damping for 25 pairs of ionic liquid ions.

with EFP_POL_DAMP true

 Compute EFP polarization integrals 
   44  -22872.5430533622      1.38E-06 Convergence criterion met
 ---------------------------------------
 SCF time:  CPU 640396.65 s  wall 76330.02 s


    EFP ENERGY COMPONENTS (ATOMIC UNITS)

               ELECTROSTATIC ENERGY    -4.1780488440
                POLARIZATION ENERGY 52804.0805454944
                  DISPERSION ENERGY    -0.7397106719
          EXCHANGE-REPULSION ENERGY     9.0921589672
   OVERLAP-BASED CHARGE PENETRATION    -1.8540906348

    QM-NUC/EFP ELECTROSTATIC ENERGY   -57.2059254920
           QM/EFP DISPERSION ENERGY     0.0000000000
   QM/EFP EXCHANGE-REPULSION ENERGY     0.0000000000

                   TOTAL EFP ENERGY 52749.1949288189
       EFP CORRECTION TO SCF ENERGY   -54.8856166756


 SCF   energy = -22872.54305336
 Total energy = -22872.54305336
 
 --------------------------------------------------------------
 
                    Orbital Energies (a.u.)
 --------------------------------------------------------------
 
 Alpha MOs
 -- Occupied --
 -673.83778  -545.17285  -469.79332  -397.40559  -354.75212  -333.93041
 -318.93038  -306.57508  -298.05127  -279.46761  -270.53063  -263.54706
 -259.24247  -255.45852  -249.18515  -242.20049  -238.18967  -229.70989
 -227.23305  -221.39125  -219.97045  -213.47327  -212.23343  -209.91795
 -206.95313  -205.58516  -203.44355  -202.67708  -200.49162  -199.10233
 -197.14041  -195.90925  -194.65548  -192.30848  -192.16222  -188.02439
 -187.32685  -186.99037  -185.58730  -182.60919  -180.24186  -179.55960
 -177.04247  -175.89665  -174.94050  -172.76140  -170.94657  -169.45529
 -167.74927  -166.23647  -165.96357  -164.01991  -162.15817  -161.58569
 -160.17025  -159.67263  -158.91274  -157.52205  -157.00040  -155.76835
 -155.61514  -154.21065  -154.04698  -152.57195  -151.48568  -150.14844
 -149.95352  -148.70560  -147.78328  -147.49485  -146.99516  -146.61645
 -145.61395  -144.78558  -143.50388  -143.02028  -142.37614  -141.67867
 -140.48427  -140.15604  -139.29513  -138.83377  -138.57593  -137.38393
 -136.55325  -136.27050  -135.12164  -134.83167  -133.53131  -133.12865
 -132.67962  -132.54281  -132.01357  -131.40612  -130.85856  -129.88070
 -129.39048  -128.41683  -128.14975  -127.42676  -127.24884  -126.30083
 -125.39247  -125.00921  -123.95623  -123.33874  -123.17791  -122.77375
 -122.45697  -121.95809  -121.95030  -120.83411  -120.70580  -120.42082
 -119.82009  -119.35116  -119.17663  -118.63768  -118.44757  -117.60057
 -117.29418  -117.16854  -116.63224  -116.17088  -116.07994  -115.89663
 -115.45975  -114.82612  -114.59704  -114.19701  -113.69750  -113.38761
 -113.28695  -112.83132  -112.65047  -112.34064  -111.86567  -111.55104
 -111.05577  -110.96100  -110.48779  -109.96099  -109.63875  -109.11249
 -108.94568  -108.50007  -108.43995  -107.76076  -107.53973  -107.19465
 -106.58253  -105.94816  -105.64261  -105.53019  -105.19811  -104.85686
 -104.58387  -104.32349  -103.82586  -103.61029  -103.48507  -103.27297
 -103.12601  -102.78488  -102.40562  -101.85456  -101.69389  -101.19445
 -101.10477  -101.02269  -100.70027  -100.63426  -100.42966  -100.01755
  -99.75056   -99.40956   -99.12198   -98.95411   -98.70935   -98.57370
  -98.19326   -97.98472   -97.77851   -97.66823   -97.40885   -97.07757
  -96.90553   -96.57003   -96.32634   -96.21865   -96.11135   -95.94020
  -95.66347   -95.36258   -94.97658   -94.91968   -94.14062   -93.98190
  -93.69079   -93.56516   -93.43392   -93.26344   -93.16269   -93.07014
  -92.98614   -92.78358   -92.52152   -92.47047   -92.39809
 -- Virtual --
  -92.17736   -92.07575   -91.88818   -91.35308   -91.17899   -91.17838
  -91.00501   -90.92010   -90.64829   -90.33121   -90.31128   -90.20623
  -90.09498   -89.75840   -89.72369   -89.70264   -89.50781   -89.46771
  -89.21257   -89.10360   -88.90889   -88.88470   -88.72915   -88.57438
  -88.47252   -88.29749   -88.28996   -88.18746   -88.11980   -87.91107
  -87.75398   -87.60417   -87.41595   -87.35709   -87.01720   -86.94570
  -86.86775   -86.79036   -86.60535   -86.46818   -86.39568   -86.14400
  -86.06432   -85.97692   -85.87750   -85.58652   -85.34557   -85.16861
  -85.13157   -84.89086   -84.85803   -84.80168   -84.64142   -84.41147
  -84.21896   -84.19195   -84.09043   -84.01808   -83.68264   -83.53038

with EFP_POL_DAMP false

 Compute EFP polarization integrals 
   16   -3282.4883518808      9.15E-06 Convergence criterion met
 ---------------------------------------
 SCF time:  CPU 161699.89 s  wall 19962.06 s


    EFP ENERGY COMPONENTS (ATOMIC UNITS)

               ELECTROSTATIC ENERGY    -4.1780488440
                POLARIZATION ENERGY    13.0624874439
                  DISPERSION ENERGY    -0.7397106719
          EXCHANGE-REPULSION ENERGY     9.0795985206
   OVERLAP-BASED CHARGE PENETRATION    -1.8506834187

    QM-NUC/EFP ELECTROSTATIC ENERGY   -57.2059254920
           QM/EFP DISPERSION ENERGY     0.0000000000
   QM/EFP EXCHANGE-REPULSION ENERGY     0.0000000000

                   TOTAL EFP ENERGY   -41.8322824621
       EFP CORRECTION TO SCF ENERGY   -54.8947699061


 SCF   energy = -3282.48835188
 Total energy = -3282.48835188
 
 --------------------------------------------------------------
 
                    Orbital Energies (a.u.)
 --------------------------------------------------------------
 
 Alpha MOs
 -- Occupied --
-89.1462 -89.0424 -19.2204 -19.1816 -19.1470 -19.1432 -19.1396 -19.1381
-19.1292 -19.0692 -19.0421 -19.0321 -14.4608 -14.4515 -10.3052 -10.2916
-10.2892 -10.2839 -10.2819 -10.2787 -10.2762 -10.2726 -10.2703 -10.2617
-10.2584 -10.2533 -10.2493 -10.2415 -10.2397 -10.2370 -10.2316 -10.2308
-10.2302 -10.2265 -10.2265 -10.2257 -10.2242 -10.2240 -10.2219 -10.2131
-10.2108 -10.1980 -10.1948 -10.1851 -10.1849 -10.1716 -10.1697 -10.1610
-10.1556 -10.1556 -10.1526 -10.1494 -10.1439 -10.1308 -10.1284 -10.1187
 -8.1617  -8.0565  -6.1172  -6.1165  -6.1143  -6.0122  -6.0108  -6.0100
 -1.8436  -1.1399  -1.0537  -1.0524  -1.0432  -1.0339  -1.0291  -1.0051
 -1.0033  -0.9871  -0.9859  -0.8964  -0.8839  -0.8749  -0.8689  -0.8649
 -0.8562  -0.8424  -0.8360  -0.8064  -0.7928  -0.7902  -0.7750  -0.7682
 -0.7580  -0.7463  -0.7425  -0.7397  -0.7172  -0.7154  -0.7066  -0.6899
 -0.6846  -0.6804  -0.6770  -0.6625  -0.6499  -0.6401  -0.6277  -0.6159
 -0.6112  -0.6066  -0.5934  -0.5876  -0.5767  -0.5681  -0.5654  -0.5524
 -0.5484  -0.5431  -0.5354  -0.5312  -0.5194  -0.5180  -0.5069  -0.4991
 -0.4929  -0.4891  -0.4851  -0.4828  -0.4809  -0.4735  -0.4728  -0.4692
 -0.4678  -0.4671  -0.4625  -0.4612  -0.4566  -0.4532  -0.4471  -0.4427
 -0.4364  -0.4285  -0.4275  -0.4248  -0.4226  -0.4210  -0.4174  -0.4083
 -0.4072  -0.4015  -0.3993  -0.3960  -0.3904  -0.3902  -0.3855  -0.3847
 -0.3835  -0.3813  -0.3752  -0.3725  -0.3696  -0.3687  -0.3656  -0.3633
 -0.3619  -0.3601  -0.3597  -0.3569  -0.3528  -0.3506  -0.3491  -0.3457
 -0.3399  -0.3341  -0.3315  -0.3301  -0.3270  -0.3233  -0.3210  -0.3169
 -0.3111  -0.3098  -0.3082  -0.3027  -0.2983  -0.2918  -0.2906  -0.2854
 -0.2738  -0.2715  -0.2676  -0.2654  -0.2629  -0.2572  -0.2541  -0.2508
 -0.2491  -0.2474  -0.2416  -0.2408  -0.2344  -0.2319  -0.2223  -0.2130
 -0.2085  -0.2064  -0.1986  -0.1811  -0.1778  -0.1669  -0.1634  -0.1368
 -0.1182
 -- Virtual --
 -0.0500  -0.0058   0.0368   0.0585   0.0761   0.0885   0.1027   0.1110
  0.1222   0.1356   0.1479   0.1547   0.1594   0.1650   0.1696   0.1767
  0.1776   0.1794   0.1858   0.1923   0.2065   0.2102   0.2113   0.2219
  0.2235   0.2266   0.2273   0.2361   0.2366   0.2397   0.2440   0.2457
  0.2476   0.2482   0.2553   0.2563   0.2581   0.2609   0.2681   0.2714
  0.2725   0.2765   0.2768   0.2798   0.2836   0.2851   0.2899   0.2913
  0.2923   0.2935   0.2948   0.2975   0.2999   0.3004   0.3044   0.3072
  0.3093   0.3122   0.3144   0.3162   0.3185   0.3204   0.3216   0.3229
  0.3268   0.3276   0.3284   0.3295   0.3301   0.3315   0.3357   0.3375
  0.3382   0.3403   0.3440   0.3464   0.3500   0.3511   0.3521   0.3546
  0.3555   0.3567   0.3589   0.3593   0.3623   0.3635   0.3646   0.3674
#2

Turning off polarization damping seems like a bad idea. Might work okay at one geometry but can easily lead to “polarization catastrophe” if two molecules get close. (In an IL system with anions, I would expect this to be especially pervasive.) I suggest you submit a feature request to Q-Chem to make the polarization damping user-modifiable.

#3

Thank you for your suggestion.

In the meantime, I would like to make this modification myself until there’s a new release. Given that I have access to the source code through the developers, I tried to update the parse.c file in libefp to achieve this effect. However, I am unsure of the procedure to compile and install this specific addition without affecting the rest of the build. Could you walk me through the process?

#4

I think it’s great that you want to do it yourself. Can you please place a question in the Q-Chem developers’ forum on the Trac site? That way other developers (and Q-Chem staff scientists) may see it and respond.

#5

If you already have signed the required Q-Chem NDA for developers, you should have login access to Q-Chem Trac, where you can find all the required resources for compiling Q-Chem and developing Q-Chem code and external libraries. If you have not yet signed the Q-Chem NDA or if you have questions related to Q-Chem code development, please contact Q-Chem Developer Support (devsupport@q-chem.com).