QTAIM population analysis from CCSD(T) densities

[Nisha Mehta]

(This question is on behalf of my PI, Prof. Gershom (Jan M.L.) Martin.)

Dear MOLPRO gurus:

As part of a study on the sensitivity of various population analysis schemes to the basis set and electron correlation method, we were trying to compute QTAIM (a.k.a., “Bader”) partial charges. We are able to obtain those for various DFT functionals (including double hybrids), HF, and MP2 from certain other codes (such as the one named after Carl Friedrich Gauss that bans lots of people) by feeding .wfn, .wfx, or formatted checkpoint files to either Multiwfn or AIMALL. But obviously we wanted the “gold standard” CCSD(T) method as well, and MOLPRO exports none of these file formats. I tried .molden, but it seems I always get HF orbitals and hence HF QTAIM charges, no matter what I do.

 

Using the following simple example for carbon monoxide, we were able to generate fine-meshed density cubes (1.6GB apiece) that we then fed into Multiwfn.

 

gthresh,energy=1d-10

geom={c;o,1,r}

r=1.1314 ang

basis,avtz,h=vtz

{hf;orbital,2130.2}

{mp2;dm,2140.2}

{ccsd;dm,2150.2}

{ccsd(t);dm,2160.2}

{cube,hf_1,-1,400,400,800;density,2130.2}

{cube,mp2_1,-1,400,400,800;density,2140.2}

{cube,ccsd_1,-1,400,400,800;density,2150.2}

{cube,ccsdt_1,-1,400,400,800;density,2160.2}

 

The resulting HF and MP2 charges agree to about 3 decimal places with the values obtained from G****ian 16’s .wfx file using AIMALL, so we are reasonably confident that things work mechanically. For completeness, these are the annotated input menu options we fed into Multiwfn

17 # basin analysis

1 # generate basins and local attractors

2 # Generate the basins by using the grid data stored in memory

2 # Integrate real space functions in the basins

0 # The values of the grid data stored in memory

-10 # Exit menu

-10 # Exit program

 

For example, for the above example (basin 1=C, basin 2=O; partial charges are obviously 6-basin1 and 8-basin2)

 

HF:

   #Basin        Integral(a.u.)      Volume(a.u.^3)

       1          4.6269932148        882.12974564

       2          9.3728780326        976.13734046

compare with AIMALL from .wfx: 9.372524

 

MP2:

   #Basin        Integral(a.u.)      Volume(a.u.^3)

       1          4.8267488785        910.88545978

       2          9.1729876642        947.38162632

compare with AIMALL from .wfx: 9.171682

 

and then

CCSD:   #Basin        Integral(a.u.)      Volume(a.u.^3)       1          4.7692576607        902.12596702       2          9.2305276526        956.14111909CCSD(T):   #Basin        Integral(a.u.)      Volume(a.u.^3)       1          4.7964555831        907.08959726       2          9.2033035565        951.17748885

 

So this works, after a fashion, but is there a simpler/cleverer way that is more amenable to automation on a set of a few hundreds or thousands of molecules? Many thanks in advance for any suggestions.

 

Best wishes

Jan/Gershom

[Nisha Mehta]

(This question is on behalf of my PI, Prof Gershom (Jan M.L.) Martin: 

*** "Computational quantum chemistry and more" ********* MBP M1 ******Outlook****************
Dr. Gershom (Jan M.L.) Martin FRSC   |   Baroness Thatcher Professor of Chemistry
Department of Molecular Chemistry and Materials Science
Weizmann Institute of Science | Kimmelman Building, Room 361 | 7610001 Reḥovot, ISRAEL