Terriable unphysical SCF iteration Energy with pure DFT functionals in MRCC

Dear all,

Recently, I have run some calculations with pure DFT functionals, I have found that there may be  seriously bug in MRC when run pure DFT functionals :

below is the input simple file:

basis=def2-SVP
scftype=RHF
calc=PBE
mem=25GB
scfalg=direct
charge=0
mult=1

unit=angs
geom=xyz
780

   O     -9.132801    0.787796    0.000000
   H     -5.172563   -8.452087    0.000000
   H     -8.393218   -1.749255    0.000000
   H     -6.251572   -6.388824    0.000000
   C     -9.892516   -0.186569    0.000000
   N     -6.213506   -8.456263    0.000000
   H     -8.092538   -4.336187    0.000000
   C     -9.457744   -1.569674    0.000000
   C     -6.878186   -7.275679    0.000000
   H    -11.591334    0.975663    0.000000
   O     -6.135848  -10.721773   -0.000000
........

The Output:

 ***KOHN-SHAM ENERGY IN STEP   1 IS        -22418.7290256274282001 [AU]
 ======================================================================
 ITERATION STEP   2
 CPU time [min]:    79.965                   Wall time [min]:    10.354
 
 Constructing Fock-matrix...
 
 First step of Coulomb...
   1% done.
  12% done.
  23% done.
  34% done.
  45% done.
  56% done.
  67% done.
  78% done.
  89% done.
 100% done.
 CPU time [min]:   106.827                   Wall time [min]:    13.964
 
 Second step of Coulomb...
   1% done.
  12% done.
  23% done.
  34% done.
  45% done.
  56% done.
  67% done.
  78% done.
  89% done.
 100% done.
 CPU time [min]:   127.778                   Wall time [min]:    16.583
 
 Calculating the exchange-correlation matrix...
 Exchange-correlation energy [au]:  -3151.086993903663
 CPU time [min]:   137.399                   Wall time [min]:    17.902
 
 RMS of [F,P]:                                     15425.43231431328968
 RMS of difference density:                            0.10065297386253
 Gap [au]:  0.00005612
 ALPHA OCC:1480 310
 BETA  OCC:1480 310
 
 ***KOHN-SHAM ENERGY IN STEP   2 IS       2150664.2358872508630157 [AU]
 ======================================================================
 ITERATION STEP   3
 CPU time [min]:   147.411                   Wall time [min]:    19.157

The SCF iteration  step 1 seems reasonable, but when the step 2, there is a  Terriable unphysical SCF iteration Energy,

change the DFT functional the any other pure DFT functionals(LDA GGA mete-GGA) will give similar error！

This example was run by the 2023 release with all newest patch.

Another example is below, this example can be reproduced by the original 2023 release without any patch, but it seems can not be reproduced by 2023 patch version:

input file:

basis=def2-TZVPP
calc=PBE
mem=15000MB
charge=0
scfiguess=sad

unit=angs
geom=xyz
168

O     -1.375483   -2.595682    3.766493
O     -0.372842   -0.530460    3.884040
O      2.330189    0.523153    1.799683
O      0.284227    2.513642   -0.248388
O      2.387040    3.300481    0.286055
O      3.927242    1.967703   -2.726166
O      2.219188   -1.067386   -2.033834
O      4.209258   -1.245902   -3.161601
O      1.762226   -1.219421   -5.851607
O      1.036363    3.109592   -3.708459
O      0.248679    0.966591   -3.882879
O     -2.558629   -0.090636   -1.936082
O     -0.328113   -2.389073    0.234426
O     -2.547790   -2.722917   -0.303139
O     -4.075727   -1.344773    2.683503
O     -2.386873    1.594877    2.030106
O     -4.417155    1.786823    3.078386
O     -1.879430    1.877978    5.676578
N     -3.058883    0.124765    4.751587
H     -3.721790   -0.211253    4.040066
......

and the output:

.....
 Calculating the exchange-correlation matrix...
 Exchange-correlation energy [au]:   -532.206056870897
 CPU time [min]:    42.509                   Wall time [min]:     5.415
 
 RMS of [F,P]:                                     47350.04503486274189
 RMS of difference density:                            0.10203989774679
 Gap [au]:  0.00707702
 ALPHA OCC: 300
 BETA  OCC: 300
 
 ***KOHN-SHAM ENERGY IN STEP   3 IS      11219036.1995280999690294 [AU]
 ======================================================================
.......

By the  way, only pure DFT functionals show this abnormality, no other  type DFT.

the attachment is the detail output files for further investigations and debug(outputs have all detail  inputs).


Qinqing


 

Dear Qinqing,

the smaller example looks like something that was solved in the patched version.
The 780 atom system is quite large and seems to have some delocalized electronic structure. You can see that the gap is almost zero in the first iteration. Vanishing gap is a known issue of pure DFT methods for such systems.

Could you, please, double check is pure DFT is meaningful for the large system? For example, you should try it with other programs or other SCF settings in MRCC.
Alternatively, I can recommend using hybrid DFT method with local density fitting via the excrad keyword.

Best wishes,
Peter

Dear professor Peter,

Yes, as I have mentioned, the second example can not be reproduced by the patched version.

However, For the bigger example, It was converged smoothly with ORCA and Gaussian16, and also the real gap calculated by ORCA and Gaussian16 are very similar, it is about 1.42 eV, however MRCC give 0.0015 eV? This is terribly too small and wrong gap.

We only need WFN for further AIM related analysis, the Pure DFT WFN is indeed good for the purpose, aother reason is that pure DFT will be the fastest one. all for this is because such big molecule we can only do a geometry optimization by a semiempirical method -- gfn-XTB, however the semiempirical  xtb WFN is not very good sometimes very poor for AIM purpose！So the Pure DFT is the good choice .

Greetings，
Qinqing

 

Dear Qinqing,
Thank you very much for reporting this issue. Please download the patches for integ and integ.f, which will fix this issue.

Dear Prof. Mihaly Kallay,

Thanks so much for your kindly help. Now It is OK, MRCC gives almost the same gap & energy as ORCA/Gaussian of my system.

By the way, it seems that the integral efficiency is slightly improved this time?

Greetings,
Qinqing

Dear Qinqing,
The integral code has not been modified. Probably compilers are better.