QMol_DFT_Vc_LDA_soft

Local-density approximation (LDA) correlation potential and functional for soft-Coulomb electron-electron interaction potentials.

Description

Use QMol_DFT_Vc_LDA_soft to describe the LDA correlation potential and functional for the soft-Coulomb electron-electron interaction potential [Helbig 2011]

$$ {\mathcal{V}}_{{\mathrm{e}\mathrm{e}}} (r)=\frac{1}{\sqrt{r^2 +1}},~~~~~~(1) $$

in DFT models. Note that the softening parameter is set to unity and cannot be adjusted. For an up- and down-spin one-body densities $\rho^{\uparrow }$ and $\rho^{\downarrow }$ , the LDA correlation energy functional is defined as

$$ H_{\textrm{c}} \left\lbrack \rho^{\uparrow } ,\rho^{\downarrow } \right\rbrack =\int \left(\rho^{\uparrow } \left(x\right)+\rho^{\downarrow } \left(x\right)\right)\varepsilon_{\textrm{c}} \left(\rho^{\uparrow } \left(x\right),\rho^{\downarrow } \left(x\right)\right)~dx=\int \frac{1}{2r_s \left(x\right)}\varepsilon_{\textrm{c}} \left(r_s \left(x\right),\xi \left(x\right)\right)~dx,~~~~~~(2) $$

with the exchange energy per particle $\varepsilon_{\textrm{c}}$ [Helbig 2011] and parametrization in terms of the Wigner-Seitz radius $r_s$ and spin polarization $\xi$ as

$$ \epsilon_c \left(\rho_{\uparrow } ,\rho_{\downarrow } \right) = \epsilon_c \left(r_s ,\xi \right) = \epsilon_c^0 \left(r_s \right) + \xi^2 \left\lbrack \epsilon_c^1 \left(r_s \right)-\epsilon_c^0 \left(r_s \right)\right\rbrack , ~~~~ r_s =\frac{1}{2\left(\rho_{\uparrow } +\rho_{\downarrow } \right)} ~~~~ \textrm{and} ~~~~ \xi =\frac{\rho_{\uparrow } -\rho_{\downarrow } }{\rho_{\uparrow } +\rho_{\downarrow } }. ~~~~~~(3) $$

This leads to the LDA correlation potential

$$ {\mathcal{V}} _{\textrm{c}}^{\uparrow } = \frac{\delta E _{\textrm{c}} }{\delta \rho^{\uparrow } } = \varepsilon _{\textrm{c}}^0 \left(r _s \right) + \xi \left(2-\xi \right)\left\lbrack \varepsilon _{\textrm{c}}^1 \left(r _s \right) - \varepsilon _{\textrm{c}}^0 \left(r _s \right) \right\rbrack - r _s \left\lbrace \partial _{r _s } \varepsilon _{\textrm{c}}^0 \left(r _s \right) + \xi^2 \left\lbrack \partial _{r _s } \varepsilon _{\textrm{c}}^1 \left(r _s \right) - \partial _{r _s } \varepsilon _{\textrm{c}}^0 \left(r _s \right) \right\rbrack \right\rbrace ,~~~~~~(4.1) $$

$$ {\mathcal{V}} _{\textrm{c}}^{\downarrow } = \frac{\delta E _{\textrm{c}} }{\delta \rho^{\downarrow } } = \varepsilon _{\textrm{c}}^0 \left(r _s \right) - \xi \left(2 + \xi \right)\left\lbrack \varepsilon _{\textrm{c}}^1 \left(r _s \right) - \varepsilon _{\textrm{c}}^0 \left(r _s \right) \right\rbrack - r _s \left\lbrace \partial _{r _s } \varepsilon _{\textrm{c}}^0 \left(r _s \right) + \xi^2 \left\lbrack \partial _{r _s } \varepsilon _{\textrm{c}}^1 \left(r _s \right) - \partial _{r _s } \varepsilon _{\textrm{c}}^0 \left(r _s \right) \right\rbrack \right\rbrace ,~~~~~~(4.2) $$

given that

$$ \frac{\partial r_s }{\partial \rho^{\uparrow } } = \frac{\partial r_s }{\partial \rho^{\downarrow } } = -2r_s^2 , ~~~~ \frac{\partial \xi }{\partial \rho^{\uparrow } } = \frac{2\rho^{\downarrow } }{{\left(\rho^{\uparrow } + \rho^{\downarrow } \right)}^2 } = \left(1-\xi \right)2r_s ~~~~ \textrm{and} ~~~~ \frac{\partial \xi }{\partial \rho^{\downarrow } } =\frac{2\rho^{\uparrow } }{{\left(\rho^{\uparrow } +\rho^{\downarrow } \right)}^2 } = -\left(1+\xi \right)2r_s . $$

As defined in [Helbig 2011], both $\epsilon_c^0$ and $\epsilon_c^1$ components share the same form of parametrization

$$ \epsilon_c^{0,1} \left(r_s \right)=-\frac{1}{2}\frac{r_s +Er_s^2 }{A+Br_s +Cr_s^2 +Dr_s^3 }\ln \left(1+\alpha r_s +\beta r_s^m \right),~~~~~~(5) $$

and its derivative reads

$$ \partial_{r_s } \epsilon_c^{0,1} \left(r_s \right)=-\frac{1}{2}\frac{A+2AEr_s +\left(BE-C\right)r_s^2 -2Dr_s^3 -DEr_s^4 }{{\left(A+Br_s +Cr_s^2 +Dr_s^3 \right)}^2 }\ln \left(1+\alpha r_s +\beta r_s^m \right)-\frac{1}{2}\frac{r_s +Er_s^2 }{A+Br_s +Cr_s^2 +Dr_s^3 }\frac{\alpha +\beta mr_s^{m-1} }{1+\alpha r_s +\beta r_s^m }. $$

QMol_DFT_Vc_LDA_soft supports average-density self-interaction correction (ADSIC) [Legrand 2002].

Class properties

isInitialized (isInit)

Whether the potential object is properly initialized. This is used throughout the QMol-grid package to check that the potential object holds meaningful information and is ready for use. Changing its isSpinPol may cause simulations to fail or produce erroneous results.

type

Flavor of DFT functional [ 'LDA_C' ]

• This is a constant property, which can be used by other components of the QMol-grid package to determine the flavor/type of functional a given object belongs to.

Class methods

Creation

constructor

Create a LDA correlation object.

obj = QMol_DFT_Vc_LDA_soft;

Changing class properties

reset

Reset the object by deleting/re-initializing all run-time properties of the class and updating the isInitialized flag to false.

obj.reset;

• This is the common reset method available to all classes throughout the QMol-grid package.

clear

Clear all class properties

obj.clear;

Initializing the object

initialize

Initialize a QMol_DFT_Vc_LDA_soft object and set the isInitialized flag to true.

obj.initialize(DFT);

• DFT is the DFT-model handle object, i.e., QMol_DFT_spinPol or QMol_DFT_spinRes, to which the LDA Slater exchange functional is attached.
• To avoid any mismatch in internal properties, initialize first reset the object before performing the initialization.

Initialize a QMol_DFT_Vc_LDA_soft object for a specific flavor of self-interaction correction (SIC) and set the isInitialized flag to true.

obj.initialize(DFT,SIC);

• Omitted, [] or 'none' input SIC disable SIC.
• 'ADSIC' input SIC uses ADSIC [Legrand 2002].

Run-time documentation

getMemoryProfile

Get an estimate of the memory held by a QMol_DFT_Vc_LDA_soft object with either

mem = obj.getMemoryProfile;
mem = obj.getMemoryProfile(false);

• The object must be properly initialized with a domain discretization.
• The estimate only includes the discretization of member components on the domain grid and ignores other (small) properties.
• The output mem is the estimated size in bytes.

Additionally display the detail of the memory footprint with

mem = obj.getMemoryProfile(true);

showDocumentation

Display the run-time documentation for the specific configuration of a QMol_DFT_Vc_LDA_soft object

ref = obj.showDocumentation;

• The output ref is a cell vector containing the list of references to be included in the bibliography.

LDA correlation functional

Before using any of its LDA correlation functional methods, the QMol_DFT_Vc_LDA_soft object must be properly initialized.

getEnergy

Get the LDA correlation energy for the parent DFT object

E = obj.getEnergy;

• This computes the correlation energy associated with the Kohn-Sham orbitals in the parent DFT model. To do so, it first compute the one-body density associated with the Kohn-Sham orbitals and their occupation parameters.
• The output scalar E contains the numerical evaluation of the exchange energy of Eq. (2).
• This is equivalent to, but more efficient than, obj.getEnergy(DFT.getDensity) with DFT being the same DFT-model handle object used to initialize the external-potential object.

Get the LDA correlation energy for a specific one-body density

E = obj.getEnergy(rho);

• Same as above but with using the input one-body-density object rho in the numerical computation of the energy.
• rho should be a QMol_DFT_density object defined over the same domain discretization as the parent DFT model.

getPotential

Get the LDA correlation potential for the parent DFT object with either

V = obj.getPotential;
V = obj.getPotential([]);

• This computes the correlation potential associated with the Kohn-Sham orbitals in the parent DFT model. To do so, it first compute the one-body density $\rho$ associated with the Kohn-Sham orbitals and their occupation parameters.
• The output Kohn-Sham potential object contains the numerical evaluation of the exchange potential of Eq. (4).
• This creates a new Kohn-Sham potential object V in which the exchange potential is stored.
• For spin-restricted models, access the discretization of the exchange potential with V.potential.
• For spin-restricted models, access the discretization of the up- and down-spin exchange potentials with V.potentialUp and V.potentialDown, respectively.
• Note that getPotential does not initialize the output potential object V.

Get the LDA correlation potential for a specific one-body density

V = obj.getPotential(rho);

• Same as above but with using the input one-body-density object rho in the numerical computation of the potential.
• rho should be a QMol_DFT_density object defined over the same domain discretization as the parent DFT model.

Overwrite the LDA correlation potential in an existing Kohn-Sham potential object with any of

obj.getPotential([],V);         % use parent DFT density
obj.getPotential([],V,false);
obj.getPotential(rho,V);        % supply the density
obj.getPotential(rho,V,false);

• This is similar to above without creating a new Kohn-Sham potential object.
• Any content in the input object V is erased before assigning the exchange potential to it.

Add the LDA correlation potential to an existing Kohn-Sham potential object

obj.getPotential([],V,true);    % use parent DFT density
obj.getPotential(rho,V,true);   % supply the density

• This is formally equivalent to the in-place addition $\mathcal{V}\gets \mathcal{V}+{\mathcal{V}}_{{\mathrm{X}}}^{{\mathrm{L}\mathrm{D}\mathrm{A}}}$ .

getPotentialDerivative

Get the LDA correlation potential gradient for the parent DFT object with either

DV = obj.getPotentialDerivative(1);
DV = obj.getPotentialDerivative(1,[]);

• This computes the correlation potential gradient associated with the Kohn-Sham orbitals in the parent DFT model. To do so, it first compute the one-body density $\rho$ associated with the Kohn-Sham orbitals and their occupation parameters.
• The output Kohn-Sham potential gradient object contains the numerical evaluation of the exchange potential gradient.
• In practice, getPotentialDerivative computes the correlation potential -- like in getPotential -- from which it determines the gradient via fast-Fourier transforms. This is made possible by the local nature of the LDA correlation potential which sould therefore have periodic or vanishing boundary conditions over the domain, like the one-body density and Kohn-Sham orbitals do.
• This creates a new Kohn-Sham potential gradient object DV in which the exchange potential gradient is stored.
• For spin-restricted models, access the discretization of the exchange potential gradient with DV.potentialGradient.
• For spin-restricted models, access the discretization of the up- and down-spin exchange potential gradients with DV.potentialGradientUp and DV.potentialGradientDown, respectively.
• Note that getPotentialDerivative does not initialize the output potential gradient object DV.
• Note that the first input 1 is required. This is to provide a uniform signature with higher dimension where the dimension along which the gradient component is applied must be specified.

Get the LDA correlation potential gradient for a specific one-body density

DV = obj.getPotentialDerivative(1,rho);

• Same as above but with using the input one-body-density object rho in the numerical computation of the potential gradient.
• rho should be a QMol_DFT_density object defined over the same domain discretization as the parent DFT model.

Overwrite the LDA correlation potential gradient in an existing Kohn-Sham potential gradient object with any of

obj.getPotentialDerivative(1,[],DV);        % use parent DFT density
obj.getPotentialDerivative(1,[],DV,false);
obj.getPotentialDerivative(1,rho,DV);       % supply the density
obj.getPotentialDerivative(1,rho,DV,false);

• This is similar to above without creating a new Kohn-Sham potential gradient object.
• Any content in the input object DV is erased before assigning the exchange potential gradient to it.

Add the LDA correlation potential gradient to an existing Kohn-Sham potential gradient object with either

obj.getPotentialDerivative(1,[],DV,true);   % use parent DFT density
obj.getPotentialDerivative(1,rho,DV,true);  % supply the density

• This is formally equivalent to the in-place addition $\nabla \mathcal{V}\gets \nabla \mathcal{V}+\nabla {\mathcal{V}}_{{\mathrm{C}}}^{{\mathrm{L}\mathrm{D}\mathrm{A}}}$ .

Examples

Create a discretization domain

disc = QMol_disc('xspan',-20:.1:25);

Create a LDA correlation functional object

V_C = QMol_DFT_Vc_LDA_soft;

Create a minimal DFT-model object required to initialize the exchange functional class and display the run-time documentation

DFT = QMol_DFT_spinRes('discretization',disc);
disc.initialize(DFT);
V_C.initialize(DFT);
V_C.showDocumentation;

yielding

  * Correlation functional               local-density approximation (LDA)
    for a one-dimensional (1D) soft-Coulomb electron-electron interaction
    potential of the form Vee(x) = 1 / sqrt( x^2 + 1 ) [Helbig 2011].
    V-01.21.001 (07/10/2024)                                     F. Mauger

Display the estimated memory footprint for the object

V_C.getMemoryProfile(true);

  * Correlation functional (LDA soft Coulomb)

QMol_DFT_Vc_LDA_soft objects do not store any large data and their memory footprint is approximated to 0 (which is why the previous command does not display any number).

Since we do not define a full DFT system, we need to generate the one-body density

% Create one-body density
rho = DFT.discretization.DFT_allocateDensity;
rho.set('density',3*exp(-(disc.x(:)+5).^2/4)+2*exp(-(disc.x(:)-5).^2/6));
rho.initialize(disc);

Plot the exchange potential and its gradient

% Get potential and gradient
V  = V_C.getPotential(rho);
DV = V_C.getPotentialDerivative(1,rho);

% Plot the results
figure; hold on
plot(disc.xspan,V.potential,'-','LineWidth',2,'DisplayName','V_{C}^{LDA}')
plot(disc.xspan,DV.potentialGradient,'-','LineWidth',2','DisplayName','{\nabla}V_{C}^{LDA}')
xlabel('position (a.u.)'); xlim(disc.xspan([1 end]));
ylabel('potential/gradient')
legend show

Test suite

Run the test suite for the class in normal or summary mode respectively with

QMol_test.test('DFT_Vc_LDA_soft');
QMol_test.test('-summary','DFT_Vc_LDA_soft');

For developers

Other hidden class properties

QMol_DFT_Vc_LDA_soft defines a handful of additional transient/constant and hidden properties to facilitate and speed up computations. These properties cannot be edited by any function outside of the object (SetAccess=private attribute).

DFT

DFT-model object [ [] (default) | QMol_DFT_spinPol handle object | QMol_DFT_spinRes handle object ]

• This is a copy of the DFT-model handle object passed to initialize.
• Un-initialized QMol_DFT_Vc_LDA_soft objects, i.e., isInitialized == false , have empty DFT.
• For practical reasons, DFT is editable by QMol_DFT classes.

SIC

Flavor of self-interaction correction (SIC) [ [] (default) | 0 | 1 ]

• SIC == 0 corresponds to no SIC.
• SIC == 1 corresponds to ADSIC [Legrand 2002].

tol

Density threshold [ nonnegative scalar (default 1e-10) ]

• Wherever the density falls below tol, the exchange-energy per particle and potential are forced to zero values. This is to avoid numerical artefacts associated with dividing with very small numbers in Eqs. (5).

References

[Helbig 2011] N. Helbig, J.I. Fuks, M. Casula, M.J. Verstraete, M.A.L. Marques, I.V. Tokatly, and A. Rubio, "Density functional theory beyond the linear regime: Validating an adiabatic local density approximation," Physical Review A 83, 032503 (2011).

[Legrand 2002] C. Legrand, E. Suraud, and P.-G. Reinhard, "Comparison of self-interaction-corrections for metal clusters," Journal of Physics B: Atomic, Molecular and Optical Physics 35, 1115 (2002).

Notes

The results displayed in this documentation page were generated using version 01.21 of the QMol-grid package.

• QMol_DFT_Vc_LDA_soft was introduced in version 01.00.
• getMemoryProfile was introduced in version 01.10.