SIRIUS 7.5.0
Electronic structure library and applications
Public Member Functions | Private Attributes | List of all members
sirius::Atom Class Reference

Data and methods specific to the actual atom in the unit cell. More...

#include <atom.hpp>

Public Member Functions

 Atom (Atom_type const &type__, r3::vector< double > position__, r3::vector< double > vector_field__)
 Constructor. More...
 
void init ()
 Initialize atom. More...
 
void generate_radial_integrals (sddk::device_t pu__, mpi::Communicator const &comm__)
 Generate radial Hamiltonian and effective magnetic field integrals. More...
 
Atom_type const & type () const
 Return const reference to corresponding atom type object. More...
 
Atom_symmetry_classsymmetry_class ()
 Return reference to corresponding atom symmetry class. More...
 
Atom_symmetry_class const & symmetry_class () const
 Return const referenced to atom symmetry class. More...
 
int type_id () const
 Return atom type id. More...
 
r3::vector< double > const & position () const
 Return atom position in fractional coordinates. More...
 
void set_position (r3::vector< double > position__)
 Set atom position in fractional coordinates. More...
 
auto vector_field () const
 Return vector field. More...
 
int symmetry_class_id () const
 Return id of the symmetry class. More...
 
void set_symmetry_class (std::shared_ptr< Atom_symmetry_class > symmetry_class__)
 Set symmetry class of the atom. More...
 
void set_nonspherical_potential (double *veff__, double *beff__[3])
 Set muffin-tin potential and magnetic field. More...
 
void sync_radial_integrals (mpi::Communicator const &comm__, int const rank__)
 
void sync_occupation_matrix (mpi::Communicator const &comm__, int const rank__)
 
double const * h_radial_integrals (int idxrf1, int idxrf2) const
 
double * h_radial_integrals (int idxrf1, int idxrf2)
 
double const * b_radial_integrals (int idxrf1, int idxrf2, int x) const
 
template<spin_block_t sblock>
std::complex< double > radial_integrals_sum_L3 (int idxrf1__, int idxrf2__, std::vector< gaunt_L3< std::complex< double > > > const &gnt__) const
 
int num_mt_points () const
 
Radial_grid< double > const & radial_grid () const
 
double radial_grid (int idx) const
 
double mt_radius () const
 
int zn () const
 
int mt_basis_size () const
 
int mt_aw_basis_size () const
 
int mt_lo_basis_size () const
 
void set_occupation_matrix (const std::complex< double > *source)
 
void get_occupation_matrix (std::complex< double > *destination)
 
void set_uj_correction_matrix (const int l, const std::complex< double > *source)
 
bool apply_uj_correction ()
 
int uj_correction_l ()
 
auto uj_correction_matrix (int lm1, int lm2, int ispn1, int ispn2)
 
double & d_mtrx (int xi1, int xi2, int iv)
 
double const & d_mtrx (int xi1, int xi2, int iv) const
 
auto const & d_mtrx () const
 
auto & d_mtrx ()
 

Private Attributes

Atom_type const & type_
 Type of the given atom. More...
 
std::shared_ptr< Atom_symmetry_classsymmetry_class_
 Symmetry class of the given atom. More...
 
r3::vector< double > position_
 Position in fractional coordinates. More...
 
r3::vector< double > vector_field_
 Vector field associated with the current site. More...
 
sddk::mdarray< double, 2 > veff_
 Muffin-tin potential. More...
 
sddk::mdarray< double, 3 > h_radial_integrals_
 Radial integrals of the Hamiltonian. More...
 
sddk::mdarray< double, 2 > beff_ [3]
 Muffin-tin magnetic field. More...
 
sddk::mdarray< double, 4 > b_radial_integrals_
 Radial integrals of the effective magnetic field. More...
 
int lmax_pot_ {-1}
 Maximum l for potential and magnetic field. More...
 
sddk::mdarray< std::complex< double >, 4 > occupation_matrix_
 Unsymmetrized (sampled over IBZ) occupation matrix of the L(S)DA+U method. More...
 
sddk::mdarray< std::complex< double >, 4 > uj_correction_matrix_
 U,J correction matrix of the L(S)DA+U method. More...
 
bool apply_uj_correction_ {false}
 True if UJ correction is applied for the current atom. More...
 
int uj_correction_l_ {-1}
 Orbital quantum number for UJ correction. More...
 
sddk::mdarray< double, 3 > d_mtrx_
 Auxiliary form of the D_{ij} operator matrix of the pseudo-potential method. More...
 

Detailed Description

Data and methods specific to the actual atom in the unit cell.

Definition at line 36 of file atom.hpp.

Constructor & Destructor Documentation

◆ Atom()

sirius::Atom::Atom ( Atom_type const &  type__,
r3::vector< double >  position__,
r3::vector< double >  vector_field__ 
)
inline

Constructor.

Definition at line 90 of file atom.hpp.

Member Function Documentation

◆ init()

void sirius::Atom::init ( )
inline

Initialize atom.

Definition at line 98 of file atom.hpp.

◆ generate_radial_integrals()

void sirius::Atom::generate_radial_integrals ( sddk::device_t  pu__,
mpi::Communicator const &  comm__ 
)
inline

Generate radial Hamiltonian and effective magnetic field integrals.

Hamiltonian operator has the following representation inside muffin-tins:

\[ \hat H = -\frac{1}{2}\nabla^2 + \sum_{\ell m} V_{\ell m}(r) R_{\ell m}(\hat {\bf r}) = \underbrace{-\frac{1}{2} \nabla^2+V_{00}(r)R_{00}}_{H_{s}(r)} +\sum_{\ell=1} \sum_{m=-\ell}^{\ell} V_{\ell m}(r) R_{\ell m}(\hat {\bf r}) = \sum_{\ell m} \widetilde V_{\ell m}(r) R_{\ell m}(\hat {\bf r}) \]

where

\[ \widetilde V_{\ell m}(r) = \left\{ \begin{array}{ll} \frac{H_{s}(r)}{R_{00}} & \ell = 0 \\ V_{\ell m}(r) & \ell > 0 \end{array} \right. \]

Definition at line 141 of file atom.hpp.

◆ type()

Atom_type const & sirius::Atom::type ( ) const
inline

Return const reference to corresponding atom type object.

Definition at line 318 of file atom.hpp.

◆ symmetry_class() [1/2]

Atom_symmetry_class & sirius::Atom::symmetry_class ( )
inline

Return reference to corresponding atom symmetry class.

Definition at line 324 of file atom.hpp.

◆ symmetry_class() [2/2]

Atom_symmetry_class const & sirius::Atom::symmetry_class ( ) const
inline

Return const referenced to atom symmetry class.

Definition at line 330 of file atom.hpp.

◆ type_id()

int sirius::Atom::type_id ( ) const
inline

Return atom type id.

Definition at line 336 of file atom.hpp.

◆ position()

r3::vector< double > const & sirius::Atom::position ( ) const
inline

Return atom position in fractional coordinates.

Definition at line 342 of file atom.hpp.

◆ set_position()

void sirius::Atom::set_position ( r3::vector< double >  position__)
inline

Set atom position in fractional coordinates.

Definition at line 348 of file atom.hpp.

◆ vector_field()

auto sirius::Atom::vector_field ( ) const
inline

Return vector field.

Definition at line 354 of file atom.hpp.

◆ symmetry_class_id()

int sirius::Atom::symmetry_class_id ( ) const
inline

Return id of the symmetry class.

Definition at line 360 of file atom.hpp.

◆ set_symmetry_class()

void sirius::Atom::set_symmetry_class ( std::shared_ptr< Atom_symmetry_class symmetry_class__)
inline

Set symmetry class of the atom.

Definition at line 369 of file atom.hpp.

◆ set_nonspherical_potential()

void sirius::Atom::set_nonspherical_potential ( double *  veff__,
double *  beff__[3] 
)
inline

Set muffin-tin potential and magnetic field.

Definition at line 375 of file atom.hpp.

◆ sync_radial_integrals()

void sirius::Atom::sync_radial_integrals ( mpi::Communicator const &  comm__,
int const  rank__ 
)
inline

Definition at line 383 of file atom.hpp.

◆ sync_occupation_matrix()

void sirius::Atom::sync_occupation_matrix ( mpi::Communicator const &  comm__,
int const  rank__ 
)
inline

Definition at line 391 of file atom.hpp.

◆ h_radial_integrals() [1/2]

double const * sirius::Atom::h_radial_integrals ( int  idxrf1,
int  idxrf2 
) const
inline

Definition at line 396 of file atom.hpp.

◆ h_radial_integrals() [2/2]

double * sirius::Atom::h_radial_integrals ( int  idxrf1,
int  idxrf2 
)
inline

Definition at line 401 of file atom.hpp.

◆ b_radial_integrals()

double const * sirius::Atom::b_radial_integrals ( int  idxrf1,
int  idxrf2,
int  x 
) const
inline

Definition at line 406 of file atom.hpp.

◆ radial_integrals_sum_L3()

template<spin_block_t sblock>
std::complex< double > sirius::Atom::radial_integrals_sum_L3 ( int  idxrf1__,
int  idxrf2__,
std::vector< gaunt_L3< std::complex< double > > > const &  gnt__ 
) const
inline

Compute the following kinds of sums for different spin-blocks of the Hamiltonian:

\[ \sum_{L_3} \langle Y_{L_1} u_{\ell_1 \nu_1} | R_{L_3} h_{L_3} | Y_{L_2} u_{\ell_2 \nu_2} \rangle = \sum_{L_3} \langle u_{\ell_1 \nu_1} | h_{L_3} | u_{\ell_2 \nu_2} \rangle \langle Y_{L_1} | R_{L_3} | Y_{L_2} \rangle \]

Definition at line 420 of file atom.hpp.

◆ num_mt_points()

int sirius::Atom::num_mt_points ( ) const
inline

Definition at line 460 of file atom.hpp.

◆ radial_grid() [1/2]

Radial_grid< double > const & sirius::Atom::radial_grid ( ) const
inline

Definition at line 465 of file atom.hpp.

◆ radial_grid() [2/2]

double sirius::Atom::radial_grid ( int  idx) const
inline

Definition at line 470 of file atom.hpp.

◆ mt_radius()

double sirius::Atom::mt_radius ( ) const
inline

Definition at line 475 of file atom.hpp.

◆ zn()

int sirius::Atom::zn ( ) const
inline

Definition at line 480 of file atom.hpp.

◆ mt_basis_size()

int sirius::Atom::mt_basis_size ( ) const
inline

Definition at line 485 of file atom.hpp.

◆ mt_aw_basis_size()

int sirius::Atom::mt_aw_basis_size ( ) const
inline

Definition at line 490 of file atom.hpp.

◆ mt_lo_basis_size()

int sirius::Atom::mt_lo_basis_size ( ) const
inline

Definition at line 495 of file atom.hpp.

◆ set_occupation_matrix()

void sirius::Atom::set_occupation_matrix ( const std::complex< double > *  source)
inline

Definition at line 500 of file atom.hpp.

◆ get_occupation_matrix()

void sirius::Atom::get_occupation_matrix ( std::complex< double > *  destination)
inline

Definition at line 506 of file atom.hpp.

◆ set_uj_correction_matrix()

void sirius::Atom::set_uj_correction_matrix ( const int  l,
const std::complex< double > *  source 
)
inline

Definition at line 511 of file atom.hpp.

◆ apply_uj_correction()

bool sirius::Atom::apply_uj_correction ( )
inline

Definition at line 518 of file atom.hpp.

◆ uj_correction_l()

int sirius::Atom::uj_correction_l ( )
inline

Definition at line 523 of file atom.hpp.

◆ uj_correction_matrix()

auto sirius::Atom::uj_correction_matrix ( int  lm1,
int  lm2,
int  ispn1,
int  ispn2 
)
inline

Definition at line 528 of file atom.hpp.

◆ d_mtrx() [1/4]

double & sirius::Atom::d_mtrx ( int  xi1,
int  xi2,
int  iv 
)
inline

Definition at line 533 of file atom.hpp.

◆ d_mtrx() [2/4]

double const & sirius::Atom::d_mtrx ( int  xi1,
int  xi2,
int  iv 
) const
inline

Definition at line 538 of file atom.hpp.

◆ d_mtrx() [3/4]

auto const & sirius::Atom::d_mtrx ( ) const
inline

Definition at line 543 of file atom.hpp.

◆ d_mtrx() [4/4]

auto & sirius::Atom::d_mtrx ( )
inline

Definition at line 548 of file atom.hpp.

Member Data Documentation

◆ type_

Atom_type const& sirius::Atom::type_
private

Type of the given atom.

Definition at line 40 of file atom.hpp.

◆ symmetry_class_

std::shared_ptr<Atom_symmetry_class> sirius::Atom::symmetry_class_
private

Symmetry class of the given atom.

Definition at line 43 of file atom.hpp.

◆ position_

r3::vector<double> sirius::Atom::position_
private

Position in fractional coordinates.

Definition at line 46 of file atom.hpp.

◆ vector_field_

r3::vector<double> sirius::Atom::vector_field_
private

Vector field associated with the current site.

Definition at line 49 of file atom.hpp.

◆ veff_

sddk::mdarray<double, 2> sirius::Atom::veff_
private

Muffin-tin potential.

Definition at line 52 of file atom.hpp.

◆ h_radial_integrals_

sddk::mdarray<double, 3> sirius::Atom::h_radial_integrals_
private

Radial integrals of the Hamiltonian.

Definition at line 55 of file atom.hpp.

◆ beff_

sddk::mdarray<double, 2> sirius::Atom::beff_[3]
private

Muffin-tin magnetic field.

Definition at line 58 of file atom.hpp.

◆ b_radial_integrals_

sddk::mdarray<double, 4> sirius::Atom::b_radial_integrals_
private

Radial integrals of the effective magnetic field.

Definition at line 61 of file atom.hpp.

◆ lmax_pot_

int sirius::Atom::lmax_pot_ {-1}
private

Maximum l for potential and magnetic field.

Definition at line 64 of file atom.hpp.

◆ occupation_matrix_

sddk::mdarray<std::complex<double>, 4> sirius::Atom::occupation_matrix_
private

Unsymmetrized (sampled over IBZ) occupation matrix of the L(S)DA+U method.

Definition at line 67 of file atom.hpp.

◆ uj_correction_matrix_

sddk::mdarray<std::complex<double>, 4> sirius::Atom::uj_correction_matrix_
private

U,J correction matrix of the L(S)DA+U method.

Definition at line 70 of file atom.hpp.

◆ apply_uj_correction_

bool sirius::Atom::apply_uj_correction_ {false}
private

True if UJ correction is applied for the current atom.

Definition at line 73 of file atom.hpp.

◆ uj_correction_l_

int sirius::Atom::uj_correction_l_ {-1}
private

Orbital quantum number for UJ correction.

Definition at line 76 of file atom.hpp.

◆ d_mtrx_

sddk::mdarray<double, 3> sirius::Atom::d_mtrx_
private

Auxiliary form of the D_{ij} operator matrix of the pseudo-potential method.

The matrix is calculated for the scalar and vector effective fields (thus, it is real and symmetric).

\[ D_{\xi \xi'}^{\alpha} = \int V({\bf r}) Q_{\xi \xi'}^{\alpha}({\bf r}) d{\bf r} \]

The ionic part of the D-operator matrix is added in the D_operator class, when it is initialized.

Definition at line 86 of file atom.hpp.


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