#include <bessel_basis.h>
|
| 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)
|
|
|
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
|
|
|
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
|
|
◆ Bessel_Basis()
Bessel_Basis::Bessel_Basis |
( |
| ) |
|
◆ ~Bessel_Basis()
Bessel_Basis::~Bessel_Basis |
( |
| ) |
|
◆ 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
-
ntype | number of atom types |
lmax | maximal angular momentum of localized orbitals |
nmax | maximal principal quantum number of localized orbitals |
ecut_number | number of SBFs |
◆ get_ecut()
const double & Bessel_Basis::get_ecut |
( |
| ) |
const |
|
inline |
get energy cutoff, which is used to truncate SBF Jlq.
- Parameters
-
- Returns
- energy cutoff in Ry
◆ 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
◆ get_rcut()
const double & Bessel_Basis::get_rcut |
( |
| ) |
const |
|
inline |
cutoff radius of radial SBF Jlq.
- Parameters
-
- Returns
- cutoff radius in a.u.
◆ get_sigma()
const double & Bessel_Basis::get_sigma |
( |
| ) |
const |
|
inline |
get sigma the stddev (standard deviation) used in smooth function (Gaussian function)
- Parameters
-
- Returns
- stddev of smooth 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
-
- Returns
- boolean whether SBFs are smoothed
◆ get_tolerence()
const double & Bessel_Basis::get_tolerence |
( |
| ) |
const |
|
inline |
◆ 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_file | whether read C4 coefficients stored in external files |
ecutwfc | cutoff for numerical atomic orbitals |
ntype | atom types |
lmax_in | maximal angular momentum for numerical orbitals |
smooth | whether smooth SBFs when perform integration to calculate value of matrix element of TableOne. For details, see J. Phys.: Condens. Matter 22 (2010) 445501 |
sigma | stddev of Gaussian function for smoothing SBFs |
rcut_in | cutoff radius for SBFs |
tol_in | accurancy control for SBFs |
dk | kspace grid |
dr | realspace grid |
ucell | UnitCell class object, ucell.nmax will be used in 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
-
ntype | number of atomtype |
lmax | maximal angular momentum |
nmax | maximal chi |
ecut_number | number of SBFs |
◆ init_TableOne() [1/2]
void Bessel_Basis::init_TableOne |
( |
| ) |
|
|
private |
◆ 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_in | whether jle(r) SBF is smoothed by a Gaussian function |
sigma_in | stddev for controlling smearing of Gaussian function for smoothing jle(r) |
ecutwfc | planewave kinetic energy cutoff for controlling kspace sampling |
rcut | cutoff radius of SBFs |
dr | realspace grid |
dk | kspace grid |
lmax | maximal angular momentum for SBFs |
ecut_number | number of SBFs |
tolerence | accurancy of SBFs |
◆ 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
-
it | atom type index |
l | angular momentum |
ic | chi index |
gnorm | norm of G+k vector |
- Returns
- interpolated value
◆ Polynomial_Interpolation2()
double Bessel_Basis::Polynomial_Interpolation2 |
( |
const int & |
l, |
|
|
const int & |
ie, |
|
|
const double & |
gnorm |
|
) |
| const |
Cubic spline interpolation for matrix TableOne.
- Parameters
-
l | angular momentum |
ie | q index (see explanation in note of function BesselBasis::get_ecut_number()) |
gnorm | norm of G+k vector |
- Returns
- interpolated value
◆ 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
-
name | name of external file where C4-stored file information is contained |
ntype | number of atom types |
ecut | energy cutoff |
rcut | cutoff radius |
ecut_number | number of SBFs |
tolerence | accurancy of SBFs, here only used for consistency check |
◆ C4
Coefficients to be optimized!
◆ Dk
◆ 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
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
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:
- /home/runner/work/abacus-develop/abacus-develop/source/source_io/bessel_basis.h
- /home/runner/work/abacus-develop/abacus-develop/source/source_io/bessel_basis.cpp