Libxc 4.0.x
Compatibility notice
This manual is for Libxc 4.0.x.
See the Change Log for a complete list of changes made to the library API over the different Libxc versions.
An example
Probably the best way to explain the usage of Libxc is through an example. The following small program calculates the xc energy for a given functional for several values of the density; the available C bindings can be found in header file xc.h. When compiling this example it is necessary to ink with libxc.a or libxc.so, depending if you compiled the static or shared versions of the library.
#include <stdio.h>
#include <xc.h>
int main()
{
xc_func_type func;
double rho[5] = {0.1, 0.2, 0.3, 0.4, 0.5};
double sigma[5] = {0.2, 0.3, 0.4, 0.5, 0.6};
double exc[5];
int i, vmajor, vminor, vmicro, func_id = 1;
xc_version(&vmajor, &vminor, &vmicro);
printf("Libxc version: %d.%d.%d\n", vmajor, vminor, vmicro);
if(xc_func_init(&func, func_id, XC_UNPOLARIZED) != 0){
fprintf(stderr, "Functional '%d' not found\n", func_id);
return 1;
}
switch(func.info->family)
{
case XC_FAMILY_LDA:
xc_lda_exc(&func, 5, rho, exc);
break;
case XC_FAMILY_GGA:
case XC_FAMILY_HYB_GGA:
xc_gga_exc(&func, 5, rho, sigma, exc);
break;
}
for(i=0; i<5; i+=1){
printf("%lf %lf\n", rho[i], exc[i]);
}
xc_func_end(&func);
}
The functionals are divided in families (LDA, GGA, etc.). Given a functional identifier xc.func_id, the functional is initialized by xc_func_init, and evaluated by xc_XXX_exc, which returns the energy per particle exc. Finally the function xc_func_end cleans up.
The output of the above example should look something like this:
Libxc version: 4.0.0
0.100000 -0.342809
0.200000 -0.431912
0.300000 -0.494416
0.400000 -0.544175
0.500000 -0.586194
Fortran 90 bindings are also included in Libxc. These can be found in the file [libxc_master.F90] (https://gitlab.com/libxc/libxc/blob/master/src/libxc_master.F90). In general, calling Libxc from Fortran is as simple as from C. Note that you also need to link with libxcf90.a or libxcf90.so when using the Fortran 90 interface. Here is the previous example in Fortran:
program lxctest
use xc_f90_types_m
use xc_f90_lib_m
implicit none
TYPE(xc_f90_pointer_t) :: xc_func
TYPE(xc_f90_pointer_t) :: xc_info
real(8) :: rho(5) = (/0.1, 0.2, 0.3, 0.4, 0.5/)
real(8) :: sigma(5) = (/0.2, 0.3, 0.4, 0.5, 0.6/)
real(8) :: exc(5)
integer :: i, vmajor, vminor, vmicro, func_id = 1
call xc_f90_version(vmajor, vminor, vmicro)
write(*,'("Libxc version: ",I1,".",I1,".",I1)') vmajor, vminor, vmicro
call xc_f90_func_init(xc_func, xc_info, func_id, XC_UNPOLARIZED)
select case (xc_f90_info_family(xc_info))
case(XC_FAMILY_LDA)
call xc_f90_lda_exc(xc_func, 5, rho(1), exc(1))
case(XC_FAMILY_GGA, XC_FAMILY_HYB_GGA)
call xc_f90_gga_exc(xc_func, 5, rho(1), sigma(1), exc(1))
end select
do i = 1, 5
write(*,"(F8.6,1X,F9.6)") rho(i), exc(i)
end do
call xc_f90_func_end(xc_func)
end program lxctest
Since Libxc 3.0 a Fortran 2003 interface is also provided. This interface is found in libxc_master.F03 and it follows more closely the C interface. Note that to use it you need to link with the libxcf03.a or libxcf03.so libraries. Here is the previous example using this interface:
program lxctest
use xc_f03_lib_m
implicit none
TYPE(xc_f03_func_t) :: xc_func
TYPE(xc_f03_func_info_t) :: xc_info
real(8) :: rho(5) = (/0.1, 0.2, 0.3, 0.4, 0.5/)
real(8) :: sigma(5) = (/0.2, 0.3, 0.4, 0.5, 0.6/)
real(8) :: exc(5)
integer :: i, vmajor, vminor, vmicro, func_id = 1
call xc_f03_version(vmajor, vminor, vmicro)
write(*,'("Libxc version: ",I1,".",I1,".",I1)') vmajor, vminor, vmicro
call xc_f03_func_init(xc_func, func_id, XC_UNPOLARIZED)
xc_info = xc_f03_func_get_info(xc_func)
select case (xc_f03_func_info_get_family(xc_info))
case(XC_FAMILY_LDA)
call xc_f03_lda_exc(xc_func, 5, rho(1), exc(1))
case(XC_FAMILY_GGA, XC_FAMILY_HYB_GGA)
call xc_f03_gga_exc(xc_func, 5, rho(1), sigma(1), exc(1))
end select
do i = 1, 5
write(*,"(F8.6,1X,F9.6)") rho(i), exc(i)
end do
call xc_f03_func_end(xc_func)
end program lxctest
The info structure
For each functional there is a structure that holds several pieces of information. The relevant part of this structure for the end user is defined as:
#define XC_FLAGS_HAVE_EXC (1 << 0) /* 1 */
#define XC_FLAGS_HAVE_VXC (1 << 1) /* 2 */
#define XC_FLAGS_HAVE_FXC (1 << 2) /* 4 */
#define XC_FLAGS_HAVE_KXC (1 << 3) /* 8 */
#define XC_FLAGS_HAVE_LXC (1 << 4) /* 16 */
#define XC_FLAGS_1D (1 << 5) /* 32 */
#define XC_FLAGS_2D (1 << 6) /* 64 */
#define XC_FLAGS_3D (1 << 7) /* 128 */
#define XC_FLAGS_HYB_CAM (1 << 8) /* 256 */
#define XC_FLAGS_HYB_CAMY (1 << 9) /* 512 */
#define XC_FLAGS_VV10 (1 << 10) /* 1024 */
#define XC_FLAGS_HYB_LC (1 << 11) /* 2048 */
#define XC_FLAGS_HYB_LCY (1 << 12) /* 4096 */
#define XC_FLAGS_STABLE (1 << 13) /* 8192 */
#define XC_FLAGS_DEVELOPMENT (1 << 14) /* 16384 */
#define XC_FLAGS_NEEDS_LAPLACIAN (1 << 15) /* 32768 */
typedef struct{
int number; /* identifier number */
int kind; /* XC_EXCHANGE, XC_CORRELATION, XC_EXCHANGE_CORRELATION, XC_KINETIC */
char *name; /* name of the functional, e.g. "PBE" */
int family; /* type of the functional, e.g. XC_FAMILY_GGA */
func_reference_type *refs[XC_MAX_REFERENCES]; /* index of the references */
int flags; /* see above for a list of possible flags */
double dens_threshold;
/* this allows to have external parameters in the functional */
int n_ext_params;
const func_params_type *ext_params;
void (*set_ext_params)(struct xc_func_type *p, const double *ext_params);
...
} xc_func_info_type;
For example, for the Slater exchange functional, this structure is defined as
const xc_func_info_type xc_func_info_lda_x = {
XC_LDA_X,
XC_EXCHANGE,
"Slater exchange",
XC_FAMILY_LDA,
{&xc_ref_Dirac1930_376, &xc_ref_Bloch1929_545, NULL, NULL, NULL},
XC_FLAGS_3D | XC_FLAGS_HAVE_EXC | XC_FLAGS_HAVE_VXC | XC_FLAGS_HAVE_FXC | XC_FLAGS_HAVE_KXC,
1e-24,
0, NULL, NULL,
...
};
The user of the library can access the information in the following way:
#include <stdio.h>
#include <xc.h>
int main()
{
int ii;
xc_func_type func;
xc_func_init(&func, XC_GGA_X_B88, XC_UNPOLARIZED);
printf("The functional '%s' is ", func.info->name);
switch (func.info->kind) {
case (XC_EXCHANGE):
printf("an exchange functional");
break;
case (XC_CORRELATION):
printf("a correlation functional");
break;
case (XC_EXCHANGE_CORRELATION):
printf("an exchange-correlation functional");
break;
case (XC_KINETIC):
printf("a kinetic energy functional");
break;
default:
printf("of unknown kind");
break;
}
printf(", it belongs to the '", func.info->name);
switch (func.info->family) {
case (XC_FAMILY_LDA):
printf("LDA");
break;
case (XC_FAMILY_GGA):
printf("GGA");
break;
case (XC_FAMILY_HYB_GGA):
printf("Hybrid GGA");
break;
case (XC_FAMILY_MGGA):
printf("MGGA");
break;
case (XC_FAMILY_HYB_MGGA):
printf("Hybrid MGGA");
break;
default:
printf("unknown");
break;
}
printf("' family and is defined in the reference(s):\n");
for(ii = 0; func.info->refs[ii] != NULL; ii++){
printf("[%d] %s\n", ii+1, func.info->refs[ii]->ref);
}
xc_func_end(&func);
}
The output of this example should look like this:
The functional 'Becke 88' belongs to the 'GGA' family and is defined in the reference(s):
[1] A. D. Becke, Phys. Rev. A 38, 3098 (1988)
Several routines are available to access the information about the functional from Fortran. Here is the previous example using the Fortran 90 interface:
program xcinfo
use xc_f90_types_m
use xc_f90_lib_m
implicit none
type(xc_f90_pointer_t) :: xc_func
type(xc_f90_pointer_t) :: xc_info
integer :: i
character(len=120) :: name, kind, family, ref
call xc_f90_func_init(xc_func, xc_info, XC_GGA_X_B88, XC_UNPOLARIZED)
call xc_f90_info_name(xc_info, name)
select case(xc_f90_info_kind(xc_info))
case (XC_EXCHANGE)
write(kind, '(a)') 'an exchange functional'
case (XC_CORRELATION)
write(kind, '(a)') 'a correlation functional'
case (XC_EXCHANGE_CORRELATION)
write(kind, '(a)') 'an exchange-correlation functional'
case (XC_KINETIC)
write(kind, '(a)') 'a kinetic energy functional'
case default
write(kind, '(a)') 'of unknown kind'
end select
select case (xc_f90_info_family(xc_info))
case (XC_FAMILY_LDA);
write(family,'(a)') "LDA"
case (XC_FAMILY_GGA);
write(family,'(a)') "GGA"
case (XC_FAMILY_HYB_GGA);
write(family,'(a)') "Hybrid GGA"
case (XC_FAMILY_MGGA);
write(family,'(a)') "MGGA"
case (XC_FAMILY_HYB_MGGA);
write(family,'(a)') "Hybrid MGGA"
case default;
write(family,'(a)') "unknown"
end select
write(*,'("The functional ''", a, "'' is ", a, ", it belongs to the ''", a, "'' family and is defined in the reference(s):")') &
trim(name), trim(kind), trim(family)
i = 0
call xc_f90_info_refs(xc_info, i, ref)
do while(i >= 0)
write(*, '(a,i1,2a)') '[', i, '] ', trim(ref)
call xc_f90_info_refs(xc_info, i, ref)
end do
call xc_f90_func_end(xc_func)
end program xcinfo
and using the Fortran 2003 interface:
program xcinfo
use xc_f03_lib_m
implicit none
type(xc_f03_func_t) :: xc_func
type(xc_f03_func_info_t) :: xc_info
integer :: i
character(len=120) :: kind, family
call xc_f03_func_init(xc_func, XC_GGA_X_B88, XC_UNPOLARIZED)
xc_info = xc_f03_func_get_info(xc_func)
select case(xc_f03_func_info_get_kind(xc_info))
case (XC_EXCHANGE)
write(kind, '(a)') 'an exchange functional'
case (XC_CORRELATION)
write(kind, '(a)') 'a correlation functional'
case (XC_EXCHANGE_CORRELATION)
write(kind, '(a)') 'an exchange-correlation functional'
case (XC_KINETIC)
write(kind, '(a)') 'a kinetic energy functional'
case default
write(kind, '(a)') 'of unknown kind'
end select
select case (xc_f03_func_info_get_family(xc_info))
case (XC_FAMILY_LDA);
write(family,'(a)') "LDA"
case (XC_FAMILY_GGA);
write(family,'(a)') "GGA"
case (XC_FAMILY_HYB_GGA);
write(family,'(a)') "Hybrid GGA"
case (XC_FAMILY_MGGA);
write(family,'(a)') "MGGA"
case (XC_FAMILY_HYB_MGGA);
write(family,'(a)') "Hybrid MGGA"
case default;
write(family,'(a)') "unknown"
end select
write(*,'("The functional ''", a, "'' is ", a, ", it belongs to the ''", a, "'' family and is defined in the reference(s):")') &
trim(xc_f03_func_info_get_name(xc_info)), trim(kind), trim(family)
i = 0
do while(i >= 0)
write(*, '(a,i1,2a)') '[', i+1, '] ', trim(xc_f03_func_reference_get_ref(xc_f03_func_info_get_references(xc_info, i)))
end do
call xc_f03_func_end(xc_func)
end program xcinfo
Initialization
A functional is initialized by the routine
int xc_func_init(xc_func_type *p, int functional, int nspin);
input:
functional: which functional do we want?
nspin: either XC_UNPOLARIZED or XC_POLARIZED
output:
p: structure that holds our functional
returns: 0 (OK) or -1 (ERROR)
External parameters
There are some functionals that require extra information, like the number of particles or the temperature. The number of parameters and a short description of their meaning is included in the info structure
typedef struct{
double value;
char *description;
} func_params_type;
typedef struct{
...
/* this allows to have external parameters in the functional */
int n_ext_params;
const func_params_type *ext_params;
void (*set_ext_params)(struct xc_func_type *p, const double *ext_params);
...
} xc_func_info_type;
Note that several functions are provided in the Fortran 2003 interface to access this information, but not in the Fortran 90 interface.
These parameters are set by default to some value and changing them is done using the xc_func_set_ext_params function
xc_func_set_ext_params(xc_func_type *p, double *ext_params);
input:
p: structure obtained from calling xc_func_init
ext_params: array of external parameters
Evaluation
The routines that should be used depend on the functional family. Note: The routines expect a nonnegative density.
LDA
This is done by the routines
void xc_lda_exc(const xc_lda_type *p, int np, double *rho, double *exc);
void xc_lda_vxc(const xc_lda_type *p, int np, double *rho, double *vxc);
void xc_lda_exc_vxc(const xc_lda_type *p, inp np, double *rho, double *exc, double *vxc);
void xc_lda_fxc(const xc_lda_type *p, int np, double *rho, double *fxc);
void xc_lda_kxc(const xc_lda_type *p, int np, double *rho, double *kxc);
void xc_lda(const xc_lda_type *p, int np, double *rho, double *exc, double *vxc, double *fxc, double *kxc);
input:
p: structure obtained from calling xc_func_init
np: number of points
rho[]: the density
output:
exc[]: the energy per unit particle
vxc[]: first derivative of the energy per unit volume
fxc[]: second derivative of the energy per unit volume
kxc[]: third derivative of the energy per unit volume
The derivatives are defined as
$$ \mathrm{vxc}_\alpha = \frac{d \epsilon}{d n_\alpha} \qquad \mathrm{fxc}_{\alpha\beta} = \frac{d^2 \epsilon}{d n_\alpha d n_\beta} \qquad \mathrm{kxc}_{\alpha\beta} = \frac{d^3 \epsilon}{d n_\alpha d n_\beta d n_\gamma} $$ where $\epsilon$ is the energy per unit volume.
If the functional was initialized with nspin=XC_UNPOLARIZED, the spin indices should be dropped from the previous expressions, resulting in arrays of size np. Otherwise, the parameters have dimensions rho[2*np], vxc[2*np], fxc[3*np] and kxc[4*np]. The energy per unit particle array always has dimensions exc[np].
The components of fxc are
$$ f_\mathrm{xc}[1] = f_{\mathrm{xc} \uparrow\uparrow} \qquad f_\mathrm{xc}[2] = f_{\mathrm{xc} \uparrow\downarrow} \qquad f_\mathrm{xc}[3] = f_{\mathrm{xc} \downarrow\downarrow} \qquad $$
and the ones of kxc are
$$ k_\mathrm{xc}[1] = k_{\mathrm{xc} \uparrow\uparrow\uparrow} \qquad k_\mathrm{xc}[2] = k_{\mathrm{xc} \uparrow\uparrow\downarrow} \qquad k_\mathrm{xc}[3] = k_{\mathrm{xc} \uparrow\downarrow\downarrow} \qquad k_\mathrm{xc}[4] = k_{\mathrm{xc} \downarrow\downarrow\downarrow} \qquad $$
GGA
This is done by the routines
void xc_gga_exc(const xc_func_type *p, int np, double *rho, double *sigma, double *exc)
void xc_gga_vxc(const xc_func_type *p, int np, double *rho, double *sigma, double double *vrho, double *vsigma)
void xc_gga_exc_vxc(const xc_func_type *p, int np, double *rho, double *sigma, double *exc, double *vrho, double *vsigma)
void xc_gga_fxc(const xc_func_type *p, int np, double *rho, double *sigma, double *v2rho2, double *v2rhosigma, double *v2sigma2)
void xc_gga(const xc_func_type *p, int np, double *rho, double *sigma, double *exc, double *vrho, double *vsigma, double *v2rho2, double *v2rhosigma, double *v2sigma2)
input:
p: structure obtained from calling xc_func_init
np: number of points
rho[]: the density
sigma[]: contracted gradients of the density
output:
exc[]: the energy per unit particle
vrho[]: first partial derivative of the energy per unit volume in terms of the density
vsigma[]: first partial derivative of the energy per unit volume in terms of sigma
v2rho2[]: second partial derivative of the energy per unit volume in terms of the density
v2rhosigma[]: second partial derivative of the energy per unit volume in terms of the density and sigma
v2sigma2[]: second partial derivative of the energy per unit volume in terms of sigma
The derivatives are defined as
$$ \mathrm{vrho}_{\alpha} = \frac{\partial \epsilon}{\partial n_\alpha} \qquad \mathrm{vsigma}_{\alpha} = \frac{\partial \epsilon}{\partial \sigma_\alpha} $$
$$ \mathrm{v2rho2}_{\alpha\beta} = \frac{d^2 \epsilon}{\partial n_\alpha \partial n_\beta} \qquad \mathrm{v2rhosigma}_{\alpha\beta} = \frac{\partial \epsilon}{\partial n_\alpha \partial \sigma_\beta} \qquad \mathrm{v2sigma2}_{\alpha\beta} = \frac{d^2 \epsilon}{\partial \sigma_\alpha \partial \sigma_\beta} $$
where $\epsilon$ is the energy per unit volume.
If the functional was initialized with nspin=XC_UNPOLARIZED, the spin indices should be dropped from the previous expressions, resulting in arrays of size np. Otherwise, the parameters have dimensions rho[2*np], vxc[2*np], sigma[3*np], vsigma[3*np], v2rho2[3*np], v2rhosigma[6*np], v2sigma2[6*np]. The energy per unit particle array always has dimensions exc[np]. The components of the contracted gradient sigma are
$$ \sigma[1] = \nabla n_\uparrow \cdot \nabla n_\uparrow \qquad \sigma[2] = \nabla n_\uparrow \cdot \nabla n_\downarrow \qquad \sigma[3] = \nabla n_\downarrow \cdot \nabla n_\downarrow \qquad $$
MetaGGA
This is done by the routines
void xc_mgga_exc(const xc_func_type *p, int np, double *rho, double *sigma, double *lapl, double *tau, double *exc)
void xc_mgga_vxc(const xc_func_type *p, int np, double *rho, double *sigma, double *lapl, double *tau, double *vrho, double *vsigma, double *vlapl, double *vtau)
void xc_mgga_exc_vxc(const xc_func_type *p, int np, double *rho, double *sigma, double *lapl, double *tau,
double * exc, double *vrho, double *vsigma, double *vlapl, double *vtau)
void xc_mgga_fxc(const xc_func_type *p, int np, double *rho, double *sigma, double *lapl, double *tau,
double *v2rho2, double *v2rhosigma, double *v2sigma2, double *v2rhotau, double *v2tausigma, double *v2tau2)
void xc_mgga(const xc_func_type *p, int np, double *rho, double *sigma, double *lapl, double *tau,
double *exc, double *vrho, double *vsigma, double *vlapl, double *vtau,
double *v2rho2, double *v2rhosigma, double *v2sigma2, double *v2rhotau, double *v2tausigma, double *v2tau2)
input:
p: structure obtained from calling xc_func_init
np: number of points
rho[]: the density
sigma[]: contracted gradients of the density
lapl[]: the laplacian of the density
tau[]: the kinetic energy density
output:
exc[]: energy density per unit particle
vrho[]: first partial derivative of the energy per unit volume in terms of the density
vsigma[]: first partial derivative of the energy per unit volume in terms of sigma
vlapl[]: first partial derivative of the energy per unit volume in terms of the laplacian of the density
vtau[]: first partial derivative of the energy per unit volume in terms of the kinetic energy density
The derivatives are defined as
$$ \mathrm{vrho}_{\alpha} = \frac{\partial \epsilon}{\partial n_\alpha} \qquad \mathrm{vsigma}_{\alpha} = \frac{\partial \epsilon}{\partial \sigma_\alpha} \qquad \mathrm{vlapl}_{\alpha} = \frac{\partial \epsilon}{\partial \nabla^2 n_\alpha} \qquad \mathrm{vtau}_{\alpha} = \frac{\partial \epsilon}{\partial \tau_\alpha} $$
Note that the kinetic energy density is defined with the $1/2$ factor (until Libxc version 1.2, tau was defined as twice the kinetic energy density):
$$ {\tau}_{\alpha} = \frac{1}{2}\sum_i | \nabla \psi_i |^2 $$
Accessor functions
Some functionals provide accessor functions that give parameters needed for calculating quantities outside of libxc:
The mixing parameter for global hybrids:
double xc_hyb_exx_coef(const xc_func_type *p);
Parameters for CAM-type hybrids, where omega is the range-separation parameter, alpha is the fraction of Hartree-Fock exchange, and beta is the fraction of short-range exchange:
void xc_hyb_cam_coef(const xc_func_type *p, double *omega, double *alpha, double *beta);
Non-local correlation parameters for functionals such as WB97X_V, VV10, LC_VV10, and WB97M_V:
void xc_nlc_coef(const xc_func_type *p, double *nlc_b, double *nlc_C);
Asymptotic parameter for AK13:
double xc_gga_ak13_get_asymptotic) (double homo);
Destruction
If you no longer need your xc functional, you should destroy it with the routine
void xc_func_end (xc_func_type *p);
input:
p: structure holding the functional to destroy
Other useful functions
xc_family_from_id
If you have the identifier of a functional (that was read, for example from an input file), and want to find out to which family it belongs to, you can use the routine
int xc_family_from_id(int functional);
input:
functional: identifier of the functional
returns:
XC_FAMILY_UNKNOWN: could not find the family
XC_FAMILY_LDA, XC_FAMILY_GGA, etc.
Analogously for Fortran:
integer function xc_family_from_id (functional)
xc_functional_get_name
If you know the identifier of a functional and want to know its name without having to call xc_func_init, you can use the function
char *xc_functional_get_name(int number);
input
number: identifier of the functional
returns:
The name of the functional
xc_functional_get_number
If you know the name of a functional and want to know its identifier, you can use the function
int xc_functional_get_number(const char *name);
input
name: then name of the functional in the form "xc_gga_x_pbe" or "gga_x_pbe". Note that the name is case insesitive.
returns:
The functional identifier
Density threshold
When the density is smaller than a certain threshold, the library will skip the evaluation of the functional and instead set all the output quantities to zero. This is to avoid underflows and/or overflows. Each functional is defined with a default threshold that should be adquate for most cases, but the range of densities in which a functional can be safely evaluated might depend on the computer architecture or the compiler and compiler flags used. Therefore, the library provides a function to change the threshold
void xc_func_set_dens_threshold(xc_func_type *p, double dens_threshold);
input
p: structure holding the functional
dens_threshold: the new density threshold
In order to help the user choose an appropriate value for the density threshold, the xc-threshold utility is provided to evaluate a functional for several different trial densities.
Version
The actual version of libxc being used can be obtained by calling
void xc_version (int *major, int *minor, int *micro);
output:
major: the major version number
minor: the minor version number
micro: the micro version number
This information is also available in the xc_version.h header file, where XC_VERSION, XC_MAJOR_VERSION, XC_MINOR_VERSION, and XC_MICRO_VERSION symbols are defined.
Linking via GNU autotools
An m4 macro for use in a configure script for a project linking to the Fortran 90 Libxc interface is available in the distribution, taken from Octopus.