ABACUS develop
Atomic-orbital Based Ab-initio Computation at UStc
Loading...
Searching...
No Matches
Public Member Functions | Private Member Functions | Private Attributes | List of all members
Bessel_Basis Class Reference

#include <bessel_basis.h>

Collaboration diagram for Bessel_Basis:

Public Member Functions

 Bessel_Basis ()
 
 ~Bessel_Basis ()
 
void init (const bool start_from_file, const double &ecutwfc, const int &ntype, const int &lmax_in, const bool &smooth, const double &sigma, const double &rcut_in, const double &tol_in, const UnitCell &ucell, const double &dk=0.01, const double &dr=0.01)
 Initialization of Bessel function related matrices.
 
const int & get_ecut_number () const
 return number of SBFs used for one chi (see details for more information)
 
double Polynomial_Interpolation (const int &it, const int &l, const int &ic, const double &gnorm) const
 Cubic spline interpolation for matrix Faln.
 
double Polynomial_Interpolation2 (const int &l, const int &ie, const double &gnorm) const
 Cubic spline interpolation for matrix TableOne.
 
const double & get_ecut () const
 get energy cutoff, which is used to truncate SBF Jlq.
 
const double & get_rcut () const
 cutoff radius of radial SBF Jlq.
 
const double & get_tolerence () const
 
const bool & get_smooth () const
 check if SBFs are smoothed (mohan add 2009-08-28)
 
const double & get_sigma () const
 get sigma the stddev (standard deviation) used in smooth function (Gaussian function)
 

Private Member Functions

void allocate_C4 (const int &ntype, const int &lmax, const int &nmax, const int &ecut_number, const UnitCell &ucell)
 Allocate memory for C4 matrix and initialize all elements to one.
 
void readin_C4 (const std::string &name, const int &ntype, const int &ecut, const int &rcut, const int &ecut_number, const double &tolerence, const UnitCell &ucell)
 Read C4 from external file. Presently an O(N^2) search algorithm is used. A HTML parser is needed in the future to improve performance.
 
void init_TableOne ()
 
void init_Faln (const int &ntype, const int &lmax, const int &nmax, const int &ecut_number, const UnitCell &ucell)
 calculate F_{aln}(it, il, in, ik) = sum_{ie}{C4(it, il, in, ie)*TableOne(il, ie, ik)}, where TableOne is overlap integral between two spherical bessel functions (jle(r) and jlk(r))
 
void init_TableOne (const bool smooth_in, const double &sigma_in, const double &ecut, const double &rcut, const double &dr, const double &dk, const int &lmax, const int &ecut_number, const double &tolerence)
 calculate element value of TableOne matrix
 

Private Attributes

ModuleBase::realArray Faln
 the most important array to calculate spillage, has dimension (ntype, lmax+1, max_n, nk)
 
ModuleBase::realArray C4
 Coefficients to be optimized!
 
ModuleBase::realArray TableOne
 matrix whose elements are int{dr r^2 j_l(qr)*j_l(kr)}, has dimension (lmax+1, nq, nk)
 
int kmesh =0
 mesh of k vector, k is in j_l(k*r)
 
double Dk
 grid of k
 
int Ecut_number
 number of q vector, q is in j_l(q*r)
 
double rcut =0.0
 Cutoff radius (in a.u.) of SBFs, for any SBF j_l(qr), r>=rcut, j_l(q*r) = 0 (if not smoothed)
 
double ecut =0.0
 energy cutoff for determining kmesh and number of SBFs
 
double tolerence =0.0
 
bool smooth =false
 whether smooth SBFs around cutoff radius, resulting in non-zero values. For importance of smooth of SBFs, see J. Phys.: Condens. Matter 22 (2010) 445501, eqn 6. (mohan add 2009-01-18)
 
double sigma =0.0
 stddev of smooth function (Gaussian function, centered at rcut)
 
int nwfc =0
 number of localized wave functions
 

Constructor & Destructor Documentation

◆ Bessel_Basis()

Bessel_Basis::Bessel_Basis ( )

◆ ~Bessel_Basis()

Bessel_Basis::~Bessel_Basis ( )

Member Function Documentation

◆ allocate_C4()

void Bessel_Basis::allocate_C4 ( const int &  ntype,
const int &  lmax,
const int &  nmax,
const int &  ecut_number,
const UnitCell ucell 
)
private

Allocate memory for C4 matrix and initialize all elements to one.

Parameters
ntypenumber of atom types
lmaxmaximal angular momentum of localized orbitals
nmaxmaximal principal quantum number of localized orbitals
ecut_numbernumber of SBFs
Here is the call graph for this function:
Here is the caller graph for this function:

◆ get_ecut()

const double & Bessel_Basis::get_ecut ( ) const
inline

get energy cutoff, which is used to truncate SBF Jlq.

Parameters
<br>
Returns
energy cutoff in Ry
Here is the caller graph for this function:

◆ get_ecut_number()

const int & Bessel_Basis::get_ecut_number ( ) const
inline

return number of SBFs used for one chi (see details for more information)

atomic orbital is constructed always with not only one set of SBFs. For different sets, they are marked with different chi(s), similar with concept of contracted GTOs. For one chi, it is 'q' the summation index, and q is in SBFs like: j_l(q*r), where l is the order of SBF.

Returns
number of SBFs
Here is the caller graph for this function:

◆ get_rcut()

const double & Bessel_Basis::get_rcut ( ) const
inline

cutoff radius of radial SBF Jlq.

Parameters
<br>
Returns
cutoff radius in a.u.
Here is the caller graph for this function:

◆ get_sigma()

const double & Bessel_Basis::get_sigma ( ) const
inline

get sigma the stddev (standard deviation) used in smooth function (Gaussian function)

Parameters
<br>
Returns
stddev of smooth function
Here is the caller graph for this function:

◆ get_smooth()

const bool & Bessel_Basis::get_smooth ( ) const
inline

check if SBFs are smoothed (mohan add 2009-08-28)

Attention
in this case, the Jlq are not the true Jlq.
Parameters
<br>
Returns
boolean whether SBFs are smoothed
Here is the caller graph for this function:

◆ get_tolerence()

const double & Bessel_Basis::get_tolerence ( ) const
inline
Here is the caller graph for this function:

◆ init()

void Bessel_Basis::init ( const bool  start_from_file,
const double &  ecutwfc,
const int &  ntype,
const int &  lmax_in,
const bool &  smooth,
const double &  sigma,
const double &  rcut_in,
const double &  tol_in,
const UnitCell ucell,
const double &  dk = 0.01,
const double &  dr = 0.01 
)

Initialization of Bessel function related matrices.

Used for a specific group of C4 coefficients. 2021-01-04, mohan added a new input parameter lmax_in, if we only generate numerical atomic orbitals based on spherical Bessel functions, lmax_in = ucell.lmax. However, if we want to generate Spherical Bessel functions (SBF) for descriptor, then the lmax_in is controlled by user.

Note
This function is called in source_io/numerical_basis.cpp and source_io/numerical_descriptor.cpp
Parameters
start_from_filewhether read C4 coefficients stored in external files
ecutwfccutoff for numerical atomic orbitals
ntypeatom types
lmax_inmaximal angular momentum for numerical orbitals
smoothwhether smooth SBFs when perform integration to calculate value of matrix element of TableOne. For details, see J. Phys.: Condens. Matter 22 (2010) 445501
sigmastddev of Gaussian function for smoothing SBFs
rcut_incutoff radius for SBFs
tol_inaccurancy control for SBFs
dkkspace grid
drrealspace grid
ucellUnitCell class object, ucell.nmax will be used in this function
Here is the call graph for this function:
Here is the caller graph for this function:

◆ init_Faln()

void Bessel_Basis::init_Faln ( const int &  ntype,
const int &  lmax,
const int &  nmax,
const int &  ecut_number,
const UnitCell ucell 
)
private

calculate F_{aln}(it, il, in, ik) = sum_{ie}{C4(it, il, in, ie)*TableOne(il, ie, ik)}, where TableOne is overlap integral between two spherical bessel functions (jle(r) and jlk(r))

Parameters
ntypenumber of atomtype
lmaxmaximal angular momentum
nmaxmaximal chi
ecut_numbernumber of SBFs
Here is the call graph for this function:
Here is the caller graph for this function:

◆ init_TableOne() [1/2]

void Bessel_Basis::init_TableOne ( )
private
Here is the caller graph for this function:

◆ init_TableOne() [2/2]

void Bessel_Basis::init_TableOne ( const bool  smooth_in,
const double &  sigma_in,
const double &  ecut,
const double &  rcut,
const double &  dr,
const double &  dk,
const int &  lmax,
const int &  ecut_number,
const double &  tolerence 
)
private

calculate element value of TableOne matrix

(be called in Bessel_Basis::init(), used for outputing overlap Q matrix) initialize the table whose matrix element is the result of integral int{dr r^2 jle(r)*jlk(r)}, TableOne has three subscript (l, ie, ik), the first runs over orbitals' angular momentum and ie, ik run over ecut_number and kmesh SBFs

Parameters
smooth_inwhether jle(r) SBF is smoothed by a Gaussian function
sigma_instddev for controlling smearing of Gaussian function for smoothing jle(r)
ecutwfcplanewave kinetic energy cutoff for controlling kspace sampling
rcutcutoff radius of SBFs
drrealspace grid
dkkspace grid
lmaxmaximal angular momentum for SBFs
ecut_numbernumber of SBFs
tolerenceaccurancy of SBFs
Here is the call graph for this function:

◆ Polynomial_Interpolation()

double Bessel_Basis::Polynomial_Interpolation ( const int &  it,
const int &  l,
const int &  ic,
const double &  gnorm 
) const

Cubic spline interpolation for matrix Faln.

Parameters
itatom type index
langular momentum
icchi index
gnormnorm of G+k vector
Returns
interpolated value
Here is the caller graph for this function:

◆ Polynomial_Interpolation2()

double Bessel_Basis::Polynomial_Interpolation2 ( const int &  l,
const int &  ie,
const double &  gnorm 
) const

Cubic spline interpolation for matrix TableOne.

Parameters
langular momentum
ieq index (see explanation in note of function BesselBasis::get_ecut_number())
gnormnorm of G+k vector
Returns
interpolated value
Here is the caller graph for this function:

◆ readin_C4()

void Bessel_Basis::readin_C4 ( const std::string &  name,
const int &  ntype,
const int &  ecut,
const int &  rcut,
const int &  ecut_number,
const double &  tolerence,
const UnitCell ucell 
)
private

Read C4 from external file. Presently an O(N^2) search algorithm is used. A HTML parser is needed in the future to improve performance.

Parameters
namename of external file where C4-stored file information is contained
ntypenumber of atom types
ecutenergy cutoff
rcutcutoff radius
ecut_numbernumber of SBFs
tolerenceaccurancy of SBFs, here only used for consistency check
Here is the call graph for this function:
Here is the caller graph for this function:

Member Data Documentation

◆ C4

ModuleBase::realArray Bessel_Basis::C4
private

Coefficients to be optimized!

◆ Dk

double Bessel_Basis::Dk
private

grid of k

◆ ecut

double Bessel_Basis::ecut =0.0
private

energy cutoff for determining kmesh and number of SBFs

◆ Ecut_number

int Bessel_Basis::Ecut_number
private

number of q vector, q is in j_l(q*r)

◆ Faln

ModuleBase::realArray Bessel_Basis::Faln
private

the most important array to calculate spillage, has dimension (ntype, lmax+1, max_n, nk)

◆ kmesh

int Bessel_Basis::kmesh =0
private

mesh of k vector, k is in j_l(k*r)

◆ nwfc

int Bessel_Basis::nwfc =0
private

number of localized wave functions

◆ rcut

double Bessel_Basis::rcut =0.0
private

Cutoff radius (in a.u.) of SBFs, for any SBF j_l(qr), r>=rcut, j_l(q*r) = 0 (if not smoothed)

◆ sigma

double Bessel_Basis::sigma =0.0
private

stddev of smooth function (Gaussian function, centered at rcut)

◆ smooth

bool Bessel_Basis::smooth =false
private

whether smooth SBFs around cutoff radius, resulting in non-zero values. For importance of smooth of SBFs, see J. Phys.: Condens. Matter 22 (2010) 445501, eqn 6. (mohan add 2009-01-18)

◆ TableOne

ModuleBase::realArray Bessel_Basis::TableOne
private

matrix whose elements are int{dr r^2 j_l(qr)*j_l(kr)}, has dimension (lmax+1, nq, nk)

◆ tolerence

double Bessel_Basis::tolerence =0.0
private

The documentation for this class was generated from the following files: