Geometry optimization by using LIBXC caused no SCF iteration with error NaN

[Jingyun Ye]

Hi, all

I am trying to use M06-2X functional in CP2K (PBE potential)  to compare the binding energy of CO2  with Gaussian results( also use M06-2X functional, same basis sets). 

But the calculation stops  at the SCF every time with the following information. But if I change the functional M06-2X to PBE, the job runs well. 

Anyone know what's the problem? Any suggestion will be help. 

Thanks very much.

 

 SCF WAVEFUNCTION OPTIMIZATION

  ----------------------------------- OT ---------------------------------------

  Allowing for rotations:  F

  Optimizing orbital energies:  F

  Minimizer      : CG                  : conjugate gradient

  Preconditioner : FULL_ALL            : diagonalization, state selective

  Precond_solver : DEFAULT

  Line search    : 2PNT                : 2 energies, one gradient

  stepsize       :    0.15000000

  energy_gap     :    0.00100000

  eps_taylor     :   0.10000E-15

  max_taylor     :             4

  mixed_precision    : F

  ----------------------------------- OT ---------------------------------------

  Step     Update method      Time    Convergence         Total energy    Change

  ------------------------------------------------------------------------------

  Trace(PS):                                  209.9999999999

  Electronic density on regular grids:                   NaN                 NaN

  Core density on regular grids:              210.0000000000        0.0000000000

  Total charge density on r-space grids:                 NaN

  Total charge density g-space grids:                    NaN

[Marco]

Hi,

Something that comes to mind is the inclusion of HF exchange in going from the PBE functional to the M06-2X functional. Since your input file parameters were not shown, I want to mention that you should probably include a &HF section in your input file (http://manual.cp2k.org/cp2k-2_5-branch/CP2K_INPUT/FORCE_EVAL/DFT/XC/HF.html) when using the M06-2X functional. If you need examples of the &HF section, I would look into the regtest files for sample inputs. I usually grep for the concepts of interest. For example, to get examples of the use of the &HF section, I would go into the QS directory:

cd /cp2k-2.5.1/tests/QS/

and then grep for &HF

grep -i '&HF' */*.inp

The grep will provide a list of files which make use of the &HF section.

Something like the following may be useful for specifying the HF exchange parameters in calculations that make use of HF exchange such as the M06-2X functional:
 &XC
    FUNCTIONAL_ROUTINE NEW
    DENSITY_CUTOFF 1.0E-10
    GRADIENT_CUTOFF 1.0E-10
    TAU_CUTOFF 1.0E-10
    &HF
      PW_HFX_BLOCKSIZE 50
      &INTERACTION_POTENTIAL
        CUTOFF_RADIUS 3.0
        POTENTIAL_TYPE TRUNCATED
        T_C_G_DATA ./t_c_g.dat
      &END INTERACTION_POTENTIAL
      &MEMORY
        EPS_STORAGE_SCALING 1.0E-01
        MAX_DISK_SPACE YYYY
        MAX_MEMORY YYYY
        STORAGE_LOCATION ./STORAGE/
      &END MEMORY
      &PERIODIC
        NUMBER_OF_SHELLS -1
      &END PERIODIC
      &SCREENING
        EPS_SCHWARZ 1.0E-08
        EPS_SCHWARZ_FORCES 1.0E-06
        SCREEN_P_FORCES TRUE
      &END SCREENING
    &END HF
    &XC_FUNCTIONAL B3LYP
    &END XC_FUNCTIONAL
    &XC_GRID
      USE_FINER_GRID
      XC_DERIV SPLINE2_SMOOTH
      XC_SMOOTH_RHO NN10
    &END XC_GRID
  &END XC

Best regards,
Marco

[Jingyun Ye]

Hi, Marco

I tried to combine HF with M06-2X to run the job. But the same error was obtained. Why we need to include HF to M06-2X? 

I attached input and output files, you may find something important information in it. 

Thanks very much. 

[Marco]

Hello,

I looked at your input file and everything seemed okay. The reason you need to include a &HF section in your input file for M06-2X is that the M06-2X functional is a global hybrid functional which incorporates 54% HF exchange. It is my understanding that the minimum required input specification in the &HF section for a calculation employing a hybrid functional is the &SCREENING keyword. Take a look at the input files provided in the regtest-libxc directory (for example the "H2O-hybrid-b3lyp_libxc_uks.inp" file).

~/input/cp2k/cp2k-2.4.0/tests/QS/regtest-libxc/

I looked at some of the source code files related to the calling of the XC functionalities. In particular, the xc_libxc.f file states the following as a comment at the beginning of the file.

!>      WARNING: In the subroutine libxc_lsd_calc, it could be that the
!>      ordering for the 1st index of v2lapltau, v2rholapl, v2rhotau,
!>      v2sigmalapl and v2sigmatau is not correct. For the moment it does not
!>      matter since the calculation of the 2nd derivatives for meta-GGA
!>      functionals is not implemented in CP2K.

If this is the case, then it is possible that a meta-GGA functional (which generally is a function of the laplacian of the density) such as M06-2X may not work properly in CP2K. Sorry I could not be of more help.

Best regards,
Marco

[Geoffrey Wood]

 

I was wondering if this problem is resolvable or will be addressed in version 2.6?  I can run Hybrid-GGA functionals (eg B3LYP and PBE0) using the libxc library but hybrid-meta-GGAs all fail at the SCF (i.e not a geometry optimization just a single-point). I've compiled CP2k version 2.5.1 svn:13632 and libxc version 2.01.  I've tried to run M06 and BR89 both of with stop at the SCF initiation. However the TPSS functional, which doesn't have an HF section as given in the CP2K tests does run (i.e  H2O-tpssx_libxc.inp). As a side question the libxc library contains both exchange and correlation parts of these functionals, i.e. XC_MGGA_X_M06 and XC_MGGA_C_M06 are these suppossd to be combined somehow?

Thanks

 

 

On Sunday, April 13, 2014 2:33:23 PM UTC-4, Jingyun Ye wrote:

[Marco]

Hi Geoffrey,

The libxc v2.1 and v2.2 were added to the SVN version of cp2k (https://groups.google.com/forum/#!searchin/cp2k/Marco/cp2k/kyPFdrTV4uc/-skBy_pEMYUJ). The xc_libxc.F was updated. Using the latest libxc versions may help. In my previous post regarding meta-GGA functionals, I may have made an incorrect statement. It may not be possible to perform cp2k calculations with meta-GGAs which require 2nd derivatives but single-point and geometry optimization calculations should run. Can you compile CP2K-2.5.1 with libxc version 2.1.0 and test your calculations again. The XC_MGGA_X_M06 and XC_MGGA_C_M06 components should be combined.

Best regards,
Marco

[Geoffrey Wood]

Thanks Marco,

I'll compile, test and give you my feedback.

Geoff.

[Geoffrey Wood]

Hi Marco,

 

I haven't had much success with this. I have compiled everything from scratch including gcc (v 4.8.3) lapack (v 3.1.4) atlas (v 3.10.1) scalapack (installer v 1.0.2)  libint (v 1.1.5) libxc (both v 2.1.0 and v 2.2.0).

 

CP2K (SVN 14141) compiles and links but I get a warning at the end saying:

 

/usr/bin/ld: Warning: size of symbol `build_deriv1_eri' changed from 10368 in /home/woodg07/work/software/cp2k_trunk/cp2k/lib/Linux-x86-64-gfortran/popt/libcp2k.a(hfx_libint_wrapper.o) to 2048 in /home/woodg07/work/software//libintBuild/lib/libderiv.a(init_libderiv.o)
/usr/bin/ld: Warning: size of symbol `build_eri' changed from 19208 in /home/woodg07/work/software/cp2k_trunk/cp2k/lib/Linux-x86-64-gfortran/popt/libcp2k.a(hfx_libint_wrapper.o) to 5000 in /home/woodg07/work/software//libintBuild/lib/libint.a(init_libint.o)
/usr/bin/ld: Warning: size of symbol `build_deriv1_eri' changed from 10368 in /home/woodg07/work/software/cp2k_trunk/cp2k/lib/Linux-x86-64-gfortran/popt/libcp2k.a(hfx_libint_wrapper.o) to 2048 in /home/woodg07/work/software//libintBuild/lib/libderiv.a(init_libderiv.o)
/usr/bin/ld: Warning: size of symbol `build_eri' changed from 19208 in /home/woodg07/work/software/cp2k_trunk/cp2k/lib/Linux-x86-64-gfortran/popt/libcp2k.a(hfx_libint_wrapper.o) to 5000 in /home/woodg07/work/software//libintBuild/lib/libint.a(init_libint.o)

 

 

If I run the tests for libxc (i.e. B3LYP and PBE0) these fail at the initial SCF iteration giving the error:

 

Program received signal SIGSEGV: Segmentation fault - invalid memory reference.

Backtrace for this error:

Program received signal SIGSEGV: Segmentation fault - invalid memory reference.

Backtrace for this error:

Program received signal SIGSEGV: Segmentation fault - invalid memory reference.

Backtrace for this error:

Program received signal SIGSEGV: Segmentation fault - invalid memory reference.

Backtrace for this error:
#0  0x2BD271D in _gfortran_backtrace at backtrace.c:258
#1  0x2BAF3E0 in _gfortrani_backtrace_handler at compile_options.c:129
#2  0x3C7BA3291F
#3  0x0
#0  0x2BD271D in _gfortran_backtrace at backtrace.c:258
#1  0x2BAF3E0 in _gfortrani_backtrace_handler at compile_options.c:129
#2  0x3C7BA3291F
#3  0x0

However if I just run a standard BLYP calculation then there are no problems and the job completes.

[Matthias Krack]

Dear Geoffrey,

you write that you compiled libint v1.1.5, but CP2K is only validated for libint v1.1.4.

Matthias

[Michael Banck]

On Thu, Jul 24, 2014 at 07:43:03AM -0700, Matthias Krack wrote:
> you write that you compiled libint v1.1.5, but CP2K is only validated for
> libint v1.1.4.

As a data point, Debian/Ubuntu ships libint-1.1.4 and we have not seen
any problems with that, including building cp2k-2.5.1 against it.

Michael

[Geoffrey Wood]

Hi Matthias and Michael..

 

Thanks for your comments. I have previously run cp2k 2.5.1 with libint 1.1.5 and tested it without any problems i.e. correct energies, gradients 2nd derivatives, BOMD using GGAs, hyrbrid GGAs, with and without libxc.  It is only with the current trunk that I'm running into problems.  I will re-link with libint 1.1.4 and provide feedback.

 

Thanks again.

 

Geoff.

[Geoffrey Wood]

Hello-

 

Thanks for everyone's comments. I discovered that part of the problem I was having was an error in my "arch" file. With that fixed I have compiled and linked CP2K (SVN 14141) libint (1.1.4) and libxc (2.1.0) and tested the XC libraries.  These are the results:

 

system:

CUTOFF=1000

 

    O   0.000000    0.000000   -0.065587
    H   0.000000   -0.757136    0.520545
    H   0.000000    0.757136    0.520545

    &KIND H
      BASIS_SET 6-31Gxx
      POTENTIAL ALL
    &END KIND
    &KIND O
      BASIS_SET 6-31Gxx
      POTENTIAL ALL
    &END KIND

 

XC_HYB_GGA_XC_PBEH (+25% HF x),  CP2k:  -76.33577878183063  another popular QM code:  -76.3357457

XC_HYB_GGA_XC_B3LYP (+20% HF x), CP2k: -76.41801937267151  another popular QM code: -76.4180549

XC_MGGA_C_M06 (+27% HF x), CP2K: -70.00722600780877

XC_MGGA_X_M06 (+27% HF x) CP2K: -76.02248099520743

XC_MGGA_X_M06 XC_MGGA_C_M06 (+27% HF x): CP2K:  -76.42254621238931  another popular QM code: -76.3842635

XC_HYB_MGGA_XC_M06 (+27% HF x): CP2K (crashes on SCF)

XC_HYB_MGGA_XC_M06 (no HF x included): CP2K (crashes on SCF)

 

So there still seems to be some problems with M06. I.e. the hybrid-meta XC doesn't seem to work and trying to combine the X and C parts together gives the incorrect energy.  Thoughts?

 

Geoff

[hut...@chem.uzh.ch]

Hi

Are you running GPW or GAPW? With that cutoff one never knows.

M06L might be another good test.

Juerg
--------------------------------------------------------------
Juerg Hutter                         Phone : ++41 44 635 4491
Institut für Chemie C                FAX   : ++41 44 635 6838
Universität Zürich                   E-mail: hut...@chem.uzh.ch
Winterthurerstrasse 190
CH-8057 Zürich, Switzerland
---------------------------------------------------------------

-----cp...@googlegroups.com wrote: -----
To: cp...@googlegroups.com
From: Geoffrey Wood
Sent by: cp...@googlegroups.com
Date: 07/28/2014 03:27PM
Subject: Re: [CP2K:5563] Re: Geometry optimization by using LIBXC caused no SCF iteration with error NaN

Hello-
 
Thanks for everyone's comments. I discovered that part of the problem I was having was an error in my "arch" file. With that fixed I have compiled and linked CP2K (SVN 14141) libint (1.1.4) and libxc (2.1.0) and tested the XC libraries.  These are the results:
 
system:

CUTOFF=10000

 
    O   0.000000    0.000000   -0.065587
    H   0.000000   -0.757136    0.520545
    H   0.000000    0.757136    0.520545
    &KIND H
      BASIS_SET 6-31Gxx
      POTENTIAL ALL
    &END KIND
    &KIND O
      BASIS_SET 6-31Gxx
      POTENTIAL ALL
    &END KIND

 
XC_HYB_GGA_XC_PBEH (+25% HF x),  CP2k:  -76.33577878183063  another popular QM code:  -76.3357457
XC_HYB_GGA_XC_B3LYP (+20% HF x), CP2k: -76.41801937267151  another popular QM code: -76.4180549
XC_MGGA_C_M06 (+27% HF x), CP2K: -70.00722600780877
XC_MGGA_X_M06 (+27% HF x) CP2K: -76.02248099520743
XC_MGGA_X_M06 XC_MGGA_C_M06 (+27% HF x): CP2K:  -76.42254621238931  another popular QM code: -76.3842635
XC_HYB_MGGA_XC_M06 (+27% HF x): CP2K (crashes on SCF)
XC_HYB_MGGA_XC_M06 (no HF x included): CP2K (crashes on SCF)
 
So there still seems to be some problems with M06. I.e. the hybrid-meta XC doesn't seem to work and trying to combine the X and C parts together gives the incorrect energy.  Thoughts?
 
Geoff
 
 

--
You received this message because you are subscribed to the Google Groups "cp2k" group.
To unsubscribe from this group and stop receiving emails from it, send an email to cp2k+uns...@googlegroups.com.
To post to this group, send email to cp...@googlegroups.com.
Visit this group at http://groups.google.com/group/cp2k.
For more options, visit https://groups.google.com/d/optout.

[Geoffrey Wood]

 

Hello

 

I am using GAPW and the cutoff is 1000 (not 10000, sorry for the error).  I used such a large cutoff because I was testing the convergence behavior of the energy wrt cutoff (I tried 200, 400, 800, 1000 and 4000) in the end I reported the 1000 numbers above as it seemed to give convergence to ~10^-5 Hartree, At this level the difference between the B3LYP and PBE0 energies between the two programs is very small.  Knowing that GAPW is a different approximation to using the Gaussian basis sets this difference is to be expected (as reported in the original paper). Using the same settings the difference between the M06 energies seems to be rather large, which is the concern I have. ll try the M06L functional and provide feedback.

 

Thanks, Geoff.

[Geoffrey Wood]

I ran some others these are the results:

 

Key Words		%HF	CP2K	G09	diff (kJ/mol)
XC_GGA_XC_HCTH_407		0	-76.40895	-76.40918	-0.6
XC_HYB_GGA_XC_B97_1		0.21	-76.39399	-76.39374	0.7
XC_MGGA_C_M06_L		0	-67.68914
XC_MGGA_X_M06_L		0	-76.04416
XC_MGGA_C_M06_L	XC_MGGA_X_M06_L	0	-76.43029	-76.41007	53.1
XC_HYB_MGGA_XC_MPWB1K		0.44	SCF crash	-76.38142

 

 

The MGGA seems to be causing the problems.

[hut...@chem.uzh.ch]

Hi

so this is the GAPW bug for meta-functionals. As I have mentioned
before, there is an inconsistency in the way the kinetic energy density
is handled in the GAPW one-center terms.

regards

Juerg
--------------------------------------------------------------
Juerg Hutter                         Phone : ++41 44 635 4491
Institut für Chemie C                FAX   : ++41 44 635 6838
Universität Zürich                   E-mail: hut...@chem.uzh.ch
Winterthurerstrasse 190
CH-8057 Zürich, Switzerland
---------------------------------------------------------------

-----cp...@googlegroups.com wrote: -----
To: cp...@googlegroups.com
From: Geoffrey Wood
Sent by: cp...@googlegroups.com

Date: 07/28/2014 10:43PM
Subject: Re: [CP2K:5567] Re: Geometry optimization by using LIBXC caused no SCF iteration with error NaN

[Geoffrey Wood]

Hi Juerg,

 

Thanks, I didn't realize this was a known bug.   This makes the comparisons I was trying to make difficult, as GPW does not handle all-electron calculations. Do you have any suggestions to get around this?

 

Cheers, Geoff.