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

#include <gint.h>

Inheritance diagram for Gint:
Collaboration diagram for Gint:

Public Types

using T_psir_func = std::function< const ModuleBase::Array_Pool< double > &(const ModuleBase::Array_Pool< double > &psir_ylm, const Grid_Technique &gt, const int grid_index, const int is, const std::vector< int > &block_iw, const std::vector< int > &block_size, const std::vector< int > &block_index, const ModuleBase::Array_Pool< bool > &cal_flag)>
 

Public Member Functions

 ~Gint ()
 
Gintoperator= (Gint &&rhs)
 move operator for the next ESolver to directly use its infomation
 
hamilt::HContainer< double > * get_hRGint () const
 
std::vector< hamilt::HContainer< double > * > get_DMRGint () const
 
int get_ncxyz () const
 
void cal_gint (Gint_inout *inout)
 the unified interface to grid integration
 
void prep_grid (const Grid_Technique &gt, const int &nbx_in, const int &nby_in, const int &nbz_in, const int &nbz_start_in, const int &ncxyz_in, const int &bx_in, const int &by_in, const int &bz_in, const int &bxyz_in, const int &nbxx_in, const int &ny_in, const int &nplane_in, const int &startz_current_in, const UnitCell *ucell_in, const LCAO_Orbitals *orb_in)
 preparing FFT grid
 
void initialize_pvpR (const UnitCell &unitcell, const Grid_Driver *gd, const int &nspin)
 calculate the neighbor atoms of each atom in this processor size of BaseMatrix with be the non-parallel version
 
void reset_DMRGint (const int &nspin)
 resize DMRGint to nspin and reallocate the memory
 
void transfer_DM2DtoGrid (std::vector< hamilt::HContainer< double > * > DM2D)
 transfer DMR (2D para) to DMR (Grid para) in elecstate_lcao.cpp
 

Public Attributes

const Grid_Techniquegridt = nullptr
 
const UnitCellucell
 
T_psir_func psir_func_1 = nullptr
 
T_psir_func psir_func_2 = nullptr
 

Protected Member Functions

void gpu_vlocal_interface (Gint_inout *inout)
 in cal_gint_gpu.cpp
 
void gpu_rho_interface (Gint_inout *inout)
 
void gpu_force_interface (Gint_inout *inout)
 
void gint_kernel_vlocal (Gint_inout *inout)
 in cal_gint_cpu.cpp
 
void gint_kernel_dvlocal (Gint_inout *inout)
 calculate H_mu_nu(local)=<phi_0|vlocal|dphi_R>
 
void gint_kernel_vlocal_meta (Gint_inout *inout)
 calculate vlocal in meta-GGA functionals
 
void gint_kernel_rho (Gint_inout *inout)
 calculate charge density rho(r)=\int D_munu \phi_mu \phi_nu
 
void gint_kernel_tau (Gint_inout *inout)
 used in meta-GGA functional
 
void gint_kernel_force (Gint_inout *inout)
 compute forces
 
void gint_kernel_force_meta (Gint_inout *inout)
 compute forces related to meta-GGA functionals
 
void cal_meshball_vlocal (const int na_grid, const int LD_pool, const int *const block_size, const int *const block_index, const int grid_index, const bool *const *const cal_flag, const double *const *const psir_ylm, const double *const *const psir_vlbr3, hamilt::HContainer< double > *hR)
 
void gint_kernel_force (const int na_grid, const int grid_index, const double delta_r, double *vldr3, const int is, const bool isforce, const bool isstress, ModuleBase::matrix *fvl_dphi, ModuleBase::matrix *svl_dphi, const UnitCell &ucell)
 
void gint_kernel_force_meta (const int na_grid, const int grid_index, const double delta_r, double *vldr3, double *vkdr3, const int is, const bool isforce, const bool isstress, ModuleBase::matrix *fvl_dphi, ModuleBase::matrix *svl_dphi, const UnitCell &ucell)
 
void cal_meshball_force (const int grid_index, const int na_grid, const int *const block_size, const int *const block_index, const double *const *const psir_vlbr3_DMR, const double *const *const dpsir_x, const double *const *const dpsir_y, const double *const *const dpsir_z, ModuleBase::matrix *force)
 
void cal_meshball_stress (const int na_grid, const int *const block_index, const double *const psir_vlbr3_DMR, const double *const dpsirr, ModuleBase::matrix *stress)
 
void gint_kernel_rho (const int na_grid, const int grid_index, const double delta_r, int *vindex, const int LD_pool, const UnitCell &ucell, Gint_inout *inout)
 
void cal_meshball_rho (const int na_grid, const int *const block_index, const int *const vindex, const double *const *const psir_ylm, const double *const *const psir_DMR, double *const rho)
 Use grid integrals to compute charge density in a meshball.
 
void gint_kernel_tau (const int na_grid, const int grid_index, const double delta_r, int *vindex, const int LD_pool, Gint_inout *inout, const UnitCell &ucell)
 Use grid integrals to compute kinetic energy density tau.
 
void cal_meshball_tau (const int na_grid, int *block_index, int *vindex, double **dpsix, double **dpsiy, double **dpsiz, double **dpsix_dm, double **dpsiy_dm, double **dpsiz_dm, double *rho)
 Use grid integrals to compute kinetic energy density tau.
 

Protected Attributes

int nbx
 variables related to FFT grid
 
int nby
 
int nbz
 
int ncxyz
 
int nbz_start
 
int bx
 
int by
 
int bz
 
int bxyz
 
int nbxx
 
int ny
 
int nplane
 
int startz_current
 
hamilt::HContainer< double > * hRGint = nullptr
 
std::vector< hamilt::HContainer< double > * > hRGint_tmp
 size of vec is 4, only used when nspin = 4
 
hamilt::HContainer< std::complex< double > > * hRGintCd = nullptr
 stores Hamiltonian in sparse format
 
std::vector< hamilt::HContainer< double > * > DMRGint
 stores DMR in sparse format
 
hamilt::HContainer< double > * DMRGint_full = nullptr
 tmp tools used in transfer_DM2DtoGrid
 
std::vector< hamilt::HContainer< double > > pvdpRx_reduced
 
std::vector< hamilt::HContainer< double > > pvdpRy_reduced
 
std::vector< hamilt::HContainer< double > > pvdpRz_reduced
 

Detailed Description

namely Gint_Gamma and Gint_k, which contain specific operations for gamma point/multi-k calculations

Member Typedef Documentation

◆ T_psir_func

using Gint::T_psir_func = std::function< const ModuleBase::Array_Pool<double>&( const ModuleBase::Array_Pool<double> &psir_ylm, const Grid_Technique &gt, const int grid_index, const int is, const std::vector<int> &block_iw, const std::vector<int> &block_size, const std::vector<int> &block_index, const ModuleBase::Array_Pool<bool> &cal_flag)>

Constructor & Destructor Documentation

◆ ~Gint()

Gint::~Gint ( )

Member Function Documentation

◆ cal_gint()

void Gint::cal_gint ( Gint_inout inout)

the unified interface to grid integration

Here is the call graph for this function:
Here is the caller graph for this function:

◆ cal_meshball_force()

void Gint::cal_meshball_force ( const int  grid_index,
const int  na_grid,
const int *const  block_size,
const int *const  block_index,
const double *const *const  psir_vlbr3_DMR,
const double *const *const  dpsir_x,
const double *const *const  dpsir_y,
const double *const *const  dpsir_z,
ModuleBase::matrix force 
)
protected

Use grid integrals to compute the atomic force contributions na_grid: how many atoms on this (i,j,k) grid block_size: dim is [na_grid], number of columns of a band block_index: dim is [na_grid+1], total number of atomis orbitals psir_vlbr3_DMR: dim is [bxyz][LD_pool] dpsir_x: dim is [bxyz][LD_pool] dpsir_y: dim is [bxyz][LD_pool] dpsir_z: dim is [bxyz][LD_pool]

Here is the caller graph for this function:

◆ cal_meshball_rho()

void Gint::cal_meshball_rho ( const int  na_grid,
const int *const  block_index,
const int *const  vindex,
const double *const *const  psir_ylm,
const double *const *const  psir_DMR,
double *const  rho 
)
protected

Use grid integrals to compute charge density in a meshball.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ cal_meshball_stress()

void Gint::cal_meshball_stress ( const int  na_grid,
const int *const  block_index,
const double *const  psir_vlbr3_DMR,
const double *const  dpsirr,
ModuleBase::matrix stress 
)
protected

Use grid integrals to compute the stress contributions na_grid: how many atoms on this (i,j,k) grid block_index: dim is [na_grid+1], total number of atomis orbitals

Here is the caller graph for this function:

◆ cal_meshball_tau()

void Gint::cal_meshball_tau ( const int  na_grid,
int *  block_index,
int *  vindex,
double **  dpsix,
double **  dpsiy,
double **  dpsiz,
double **  dpsix_dm,
double **  dpsiy_dm,
double **  dpsiz_dm,
double *  rho 
)
protected

Use grid integrals to compute kinetic energy density tau.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ cal_meshball_vlocal()

void Gint::cal_meshball_vlocal ( const int  na_grid,
const int  LD_pool,
const int *const  block_size,
const int *const  block_index,
const int  grid_index,
const bool *const *const  cal_flag,
const double *const *const  psir_ylm,
const double *const *const  psir_vlbr3,
hamilt::HContainer< double > *  hR 
)
protected

calculate local potential contribution to the Hamiltonian na_grid: how many atoms on this (i,j,k) grid block_size: dim is [block_size], number of columns of a band block_index: dim is [na_grid+1], total number of atomic orbitals grid_index: index of grid group, for tracing iat cal_flag: dim is [bxyz][na_grid], whether the atom-grid distance is larger than cutoff psir_ylm: dim is [bxyz][LD_pool] psir_vlbr3: dim is [bxyz][LD_pool] hR: HContainer for storing the <phi_0|V|phi_R> matrix elements cal_meshball_vlocal is thread-safe!

Here is the call graph for this function:
Here is the caller graph for this function:

◆ get_DMRGint()

std::vector< hamilt::HContainer< double > * > Gint::get_DMRGint ( ) const
inline

◆ get_hRGint()

hamilt::HContainer< double > * Gint::get_hRGint ( ) const
inline

◆ get_ncxyz()

int Gint::get_ncxyz ( ) const
inline

◆ gint_kernel_dvlocal()

void Gint::gint_kernel_dvlocal ( Gint_inout inout)
protected

calculate H_mu_nu(local)=<phi_0|vlocal|dphi_R>

Here is the call graph for this function:
Here is the caller graph for this function:

◆ gint_kernel_force() [1/2]

void Gint::gint_kernel_force ( const int  na_grid,
const int  grid_index,
const double  delta_r,
double *  vldr3,
const int  is,
const bool  isforce,
const bool  isstress,
ModuleBase::matrix fvl_dphi,
ModuleBase::matrix svl_dphi,
const UnitCell ucell 
)
protected

in gint_fvl.cpp calculate vl contributuion to force & stress via grid integrals

◆ gint_kernel_force() [2/2]

void Gint::gint_kernel_force ( Gint_inout inout)
protected

compute forces

Here is the call graph for this function:
Here is the caller graph for this function:

◆ gint_kernel_force_meta() [1/2]

void Gint::gint_kernel_force_meta ( const int  na_grid,
const int  grid_index,
const double  delta_r,
double *  vldr3,
double *  vkdr3,
const int  is,
const bool  isforce,
const bool  isstress,
ModuleBase::matrix fvl_dphi,
ModuleBase::matrix svl_dphi,
const UnitCell ucell 
)
protected

in gint_fvl.cpp calculate vl contributuion to force & stress via grid integrals used in meta-GGA calculations

◆ gint_kernel_force_meta() [2/2]

void Gint::gint_kernel_force_meta ( Gint_inout inout)
protected

compute forces related to meta-GGA functionals

Here is the call graph for this function:
Here is the caller graph for this function:

◆ gint_kernel_rho() [1/2]

void Gint::gint_kernel_rho ( const int  na_grid,
const int  grid_index,
const double  delta_r,
int *  vindex,
const int  LD_pool,
const UnitCell ucell,
Gint_inout inout 
)
protected

Use grid integrals to compute charge density in gint_k_rho.cpp calculate the charge density & kinetic energy density (tau) via grid integrals

◆ gint_kernel_rho() [2/2]

void Gint::gint_kernel_rho ( Gint_inout inout)
protected

calculate charge density rho(r)=\int D_munu \phi_mu \phi_nu

Here is the call graph for this function:
Here is the caller graph for this function:

◆ gint_kernel_tau() [1/2]

void Gint::gint_kernel_tau ( const int  na_grid,
const int  grid_index,
const double  delta_r,
int *  vindex,
const int  LD_pool,
Gint_inout inout,
const UnitCell ucell 
)
protected

Use grid integrals to compute kinetic energy density tau.

◆ gint_kernel_tau() [2/2]

void Gint::gint_kernel_tau ( Gint_inout inout)
protected

used in meta-GGA functional

Here is the call graph for this function:
Here is the caller graph for this function:

◆ gint_kernel_vlocal()

void Gint::gint_kernel_vlocal ( Gint_inout inout)
protected

in cal_gint_cpu.cpp

When in OpenMP, it points to a newly allocated memory,

Prepare block information

Evaluate psi and dpsi on grids

Here is the call graph for this function:
Here is the caller graph for this function:

◆ gint_kernel_vlocal_meta()

void Gint::gint_kernel_vlocal_meta ( Gint_inout inout)
protected

calculate vlocal in meta-GGA functionals

Here is the call graph for this function:
Here is the caller graph for this function:

◆ gpu_force_interface()

void Gint::gpu_force_interface ( Gint_inout inout)
protected
Here is the call graph for this function:
Here is the caller graph for this function:

◆ gpu_rho_interface()

void Gint::gpu_rho_interface ( Gint_inout inout)
protected
Here is the call graph for this function:
Here is the caller graph for this function:

◆ gpu_vlocal_interface()

void Gint::gpu_vlocal_interface ( Gint_inout inout)
protected

in cal_gint_gpu.cpp

Here is the call graph for this function:
Here is the caller graph for this function:

◆ initialize_pvpR()

void Gint::initialize_pvpR ( const UnitCell unitcell,
const Grid_Driver gd,
const int &  nspin 
)

calculate the neighbor atoms of each atom in this processor size of BaseMatrix with be the non-parallel version

Here is the call graph for this function:
Here is the caller graph for this function:

◆ operator=()

Gint & Gint::operator= ( Gint &&  rhs)

move operator for the next ESolver to directly use its infomation

Here is the caller graph for this function:

◆ prep_grid()

void Gint::prep_grid ( const Grid_Technique gt,
const int &  nbx_in,
const int &  nby_in,
const int &  nbz_in,
const int &  nbz_start_in,
const int &  ncxyz_in,
const int &  bx_in,
const int &  by_in,
const int &  bz_in,
const int &  bxyz_in,
const int &  nbxx_in,
const int &  ny_in,
const int &  nplane_in,
const int &  startz_current_in,
const UnitCell ucell_in,
const LCAO_Orbitals orb_in 
)

preparing FFT grid

Here is the call graph for this function:
Here is the caller graph for this function:

◆ reset_DMRGint()

void Gint::reset_DMRGint ( const int &  nspin)

resize DMRGint to nspin and reallocate the memory

Here is the call graph for this function:

◆ transfer_DM2DtoGrid()

void Gint::transfer_DM2DtoGrid ( std::vector< hamilt::HContainer< double > * >  DM2D)

transfer DMR (2D para) to DMR (Grid para) in elecstate_lcao.cpp

Here is the call graph for this function:
Here is the caller graph for this function:

Member Data Documentation

◆ bx

int Gint::bx
protected

◆ bxyz

int Gint::bxyz
protected

◆ by

int Gint::by
protected

◆ bz

int Gint::bz
protected

◆ DMRGint

std::vector<hamilt::HContainer<double>*> Gint::DMRGint
protected

stores DMR in sparse format

◆ DMRGint_full

hamilt::HContainer<double>* Gint::DMRGint_full = nullptr
protected

tmp tools used in transfer_DM2DtoGrid

◆ gridt

const Grid_Technique* Gint::gridt = nullptr

◆ hRGint

hamilt::HContainer<double>* Gint::hRGint = nullptr
protected

save the < phi_0i | V | phi_Rj > in sparse H matrix. stores Hamiltonian in sparse format

◆ hRGint_tmp

std::vector<hamilt::HContainer<double>*> Gint::hRGint_tmp
protected

size of vec is 4, only used when nspin = 4

◆ hRGintCd

hamilt::HContainer<std::complex<double> >* Gint::hRGintCd = nullptr
protected

stores Hamiltonian in sparse format

◆ nbx

int Gint::nbx
protected

variables related to FFT grid

◆ nbxx

int Gint::nbxx
protected

◆ nby

int Gint::nby
protected

◆ nbz

int Gint::nbz
protected

◆ nbz_start

int Gint::nbz_start
protected

◆ ncxyz

int Gint::ncxyz
protected

◆ nplane

int Gint::nplane
protected

◆ ny

int Gint::ny
protected

◆ psir_func_1

T_psir_func Gint::psir_func_1 = nullptr

◆ psir_func_2

T_psir_func Gint::psir_func_2 = nullptr

◆ pvdpRx_reduced

std::vector<hamilt::HContainer<double> > Gint::pvdpRx_reduced
protected

◆ pvdpRy_reduced

std::vector<hamilt::HContainer<double> > Gint::pvdpRy_reduced
protected

◆ pvdpRz_reduced

std::vector<hamilt::HContainer<double> > Gint::pvdpRz_reduced
protected

◆ startz_current

int Gint::startz_current
protected

◆ ucell

const UnitCell* Gint::ucell

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