Functional (or its derivative) value inflates when rho (or sigma/gamma) very close to zero

Intro

When density is very close to zero, functional or its derivative values could be extremely large.

In this case, only one grid is considered. It's density $\rho$ is 2e-26, and $\gamma$ is 2e-50. Using functional B88 exchange only. $\gamma$ is defined as $$\gamma = \nabla \rho \cdot \nabla \rho$$

Calculate $e, E\rho, E\gamma$, where $e$ (exchange kernel) is defined as $$E = \int e \rho d r$$ and $$E_\gamma = \frac{ \partial (e \rho) }{ \partial \gamma }$$

Expected Behavior

Possibly in general case, such low density (and its derivative) should yield zero functional value (and its derivative).

Current Behavior

However when I'm using PySCF with xcfun, they could be as large as 1e+16 sometimes. In this specified case,

$e = -7955883.1346$ (marked as exc)
$E_\rho = -2.12157e-05$ (marked as vxc_rho, seems to be somehow correct)
$E_\gamma = -2.45617e+16$ (marked as vxc_gamma)

Possible Solution

Currently I have no idea (0.0) I'm actually new to xcfun, and not planning using xcfun directly in projects currently.

Steps to Reproduce (for bugs)

This problem could be reproduced by inserting the following c++ code to examples/CXX_host folder and compile:

#include <iostream>
#include <XCFun/xcfun.h>
using namespace std;

int main()
{
    double output[100];
    // Define a new functional object instance
    xcfun_t * fun = xcfun_new();
    // Set xc code as `Becke 88` only
    xcfun_set(fun, "BECKEX", 1);
    // Set only one grid, where `rho[0]` density `rho[1]` gamma
    double rho[2] = { 2e-26, 2e-50 };
    // Setup? (actually I does not know meaning of this line 0w0)
    xcfun_eval_setup(fun, XC_N_GNN, XC_PARTIAL_DERIVATIVES, 1);
    // Evaluate exc, vxc_rho, vxc_gamma
    xcfun_eval(fun, rho, output);
    // ... the following line is only used to compare result from PySCF, where PySCF convention is $E = \int e * rho$
    output[0] /= rho[0] + 1e-150;
    cout << "exc      : " << output[0] << endl;  # -7.95588e+06
    cout << "vxc_rho  : " << output[1] << endl;  # -2.12157e-05
    cout << "vxc_gamma: " << output[2] << endl;  # -2.45617e+16
}

This problem is originally found in PySCF. The code above is somehow mocking PySCF source code xcfun_itrf.c. Same values could be reproduced by the following python code:

import numpy as np
from pyscf import gto, dft
from pyscf.dft import libxc, xcfun

if __name__ == '__main__':
    ni = dft.numint.NumInt()
    ni.libxc = xcfun
    rho_01 = np.array([[2e-26], [8.16496580927726e-26], [8.16496580927726e-26], [8.16496580927726e-26]])
    exc, vxc, _, _ = ni.eval_xc("B88", rho_01, deriv=1)
    print("exc      : ", exc)  # [-7955883.13460294]
    print("vxc_rho  : ", vxc[0])  # [-2.12156884e-05]
    print("vxc_gamma: ", vxc[1])  # [-2.4561749e+16]

For the python code above, if replace ni.libxc = xcfun with ni.libxc = libxc, then it will give almost (or exactly) zero values, which seems more confident to me.

Context

So actually I didn't find real troubles, since in most cases, functional grid values should be multiplied with density grid values when calculating energy or its derivative. So when density grid values are extremely small, effect could be small. Actually I use xcfun (via PySCF interface) in my project Py_xDH, and all tests passed, so currently no obvious problem. But once seen these extreme functional values, my feeling is somehow subtle.

Your Environment

XCFun github master (b1dd9d74ef874cf577db503a64bdcd7d2cebb465) PySCF version 1.7.2.post2 with Python 3.8.1 g++ 7.5.0, Windows Subsystem of Linux (Build 18363)

PS

Line 19 file examples/CXX_host/CMakeLists.txt should be changed to v2.0.1.tar.gz, otherwise compilation could fail (T.T)

Dear Zhenyu Zhu, if I remember correctly the Becke88 functional is not differentiable at (0,0). Plotting it now I see that the gamma derivative blows up. To stabilize the calculation XCFun will modify the input density to be small but bigger than what you put in (something like (1e-14, 0)).

Regards, Ulf

Den tors 4 juni 2020 kl 14:31 skrev Zhenyu Zhu ajz34 < notifications@github.com>:

Intro

When density is very close to zero, functional or its derivative values could be extremely large.

In this case, only one grid is considered. It's density $\rho$ is 2e-26, and $\gamma$ is 2e-50. Using functional B88 exchange only. $\gamma$ is defined as $$\gamma = \nabla \rho \cdot \nabla \rho$$

Calculate $e, E\rho, E\gamma$, where $e$ (exchange kernel) is defined as $$E = \int e \rho d r$$ and $$E_\gamma = \frac{ \partial (e \rho) }{ \partial \gamma }$$ Expected Behavior

Possibly in general case, such low density (and its derivative) should yield zero functional value (and its derivative). Current Behavior

However when I'm using PySCF with xcfun, they could be as large as 1e+16 sometimes. In this specified case,

$e = -7955883.1346$ (marked as exc)

$E_\rho = -2.12157e-05$ (marked as vxc_rho, seems to be somehow correct)

$E_\gamma = -2.45617e+16$ (marked as vxc_gamma)

Possible Solution

Currently I have no idea (0.0) I'm actually new to xcfun, and not planning using xcfun directly in projects currently. Steps to Reproduce (for bugs)

This problem could be reproduced by inserting the following c++ code to examples/CXX_host https://github.com/dftlibs/xcfun/tree/master/examples/CXX_host folder and compile:

include

include <XCFun/xcfun.h>using namespace std;

int main() { double output[100]; // Define a new functional object instance xcfun_t fun = xcfun_new(); // Set xc code as Becke 88 only xcfun_set(fun, "BECKEX", 1); // Set only one grid, where rho[0] density rho[1] gamma double rho[2] = { 2e-26, 2e-50 }; // Setup? (actually I does not know meaning of this line 0w0) xcfun_eval_setup(fun, XC_N_GNN, XC_PARTIAL_DERIVATIVES, 1); // Evaluate exc, vxc_rho, vxc_gamma xcfun_eval(fun, rho, output); // ... the following line is only used to compare result from PySCF, where PySCF convention is $E = \int e rho$ output[0] /= rho[0] + 1e-150; cout << "exc : " << output[0] << endl; # -7.95588e+06 cout << "vxc_rho : " << output[1] << endl; # -2.12157e-05 cout << "vxc_gamma: " << output[2] << endl; # -2.45617e+16 }

This problem is originally found in PySCF. The code above is somehow mocking PySCF source code xcfun_itrf.c https://github.com/pyscf/pyscf/blob/master/pyscf/lib/dft/xcfun_itrf.c. Same values could be reproduced by the following python code:

import numpy as npfrom pyscf import gto, dftfrom pyscf.dft import libxc, xcfun if name == 'main': ni = dft.numint.NumInt() ni.libxc = xcfun rho01 = np.array([[2e-26], [8.16496580927726e-26], [8.16496580927726e-26], [8.16496580927726e-26]]) exc, vxc, , _ = ni.eval_xc("B88", rho_01, deriv=1) print("exc : ", exc) # [-7955883.13460294] print("vxc_rho : ", vxc[0]) # [-2.12156884e-05] print("vxc_gamma: ", vxc[1]) # [-2.4561749e+16]

For the python code above, if replace ni.libxc = xcfun with ni.libxc = libxc, then it will give almost (or exactly) zero values, which seems more confident to me. Context

So actually I didn't find real troubles, since in most cases, functional grid values should be multiplied with density grid values when calculating energy or its derivative. So when density grid values are extremely small, effect could be small. Actually I use xcfun (via PySCF interface) in my project Py_xDH https://github.com/ajz34/py_xdh, and all tests passed, so currently no obvious problem. But once seen these extreme functional values, my feeling is somehow subtle. Your Environment

XCFun github master (b1dd9d7 https://github.com/dftlibs/xcfun/commit/b1dd9d74ef874cf577db503a64bdcd7d2cebb465 ) PySCF version 1.7.2.post2 with Python 3.8.1 g++ 7.5.0, Windows Subsystem of Linux (Build 18363) PS

Line 19 file examples/CXX_host/CMakeLists.txt https://github.com/dftlibs/xcfun/blob/master/examples/CXX_host/CMakeLists.txt#L19 should be changed to v2.0.1.tar.gz, otherwise compilation could fail (T.T)

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/dftlibs/xcfun/issues/131, or unsubscribe https://github.com/notifications/unsubscribe-auth/ADD25GKOMZI4J34BJXI2KMDRU6H2VANCNFSM4NSSSA4Q .

dftlibs / xcfun