ABACUS develop
Atomic-orbital Based Ab-initio Computation at UStc
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Private Types | Private Member Functions | Private Attributes | List of all members
hsolver::DiagoBPCG< T, Device > Class Template Reference

A class for diagonalization using the Blocked-PCG method. More...

#include <diago_bpcg.h>

Collaboration diagram for hsolver::DiagoBPCG< T, Device >:

Public Types

using HPsiFunc = std::function< void(T *, T *, const int, const int)>
 

Public Member Functions

 DiagoBPCG (const Real *precondition)
 Constructor for DiagoBPCG class.
 
 ~DiagoBPCG ()
 Destructor for DiagoBPCG class.
 
void init_iter (const int nband, const int nband_l, const int nbasis, const int ndim)
 Initialize the class before diagonalization.
 
void diag (const HPsiFunc &hpsi_func, T *psi_in, Real *eigenvalue_in, const std::vector< double > &ethr_band)
 Diagonalize the Hamiltonian using the BPCG method.
 

Private Types

using Real = typename GetTypeReal< T >::type
 
using ct_Device = typename ct::PsiToContainer< Device >::type
 
using setmem_var_op = ct::kernels::set_memory< Real, ct_Device >
 
using resmem_var_op = ct::kernels::resize_memory< Real, ct_Device >
 
using delmem_var_op = ct::kernels::delete_memory< Real, ct_Device >
 
using syncmem_var_h2d_op = ct::kernels::synchronize_memory< Real, ct_Device, ct::DEVICE_CPU >
 
using syncmem_var_d2h_op = ct::kernels::synchronize_memory< Real, ct::DEVICE_CPU, ct_Device >
 
using setmem_complex_op = ct::kernels::set_memory< T, ct_Device >
 
using delmem_complex_op = ct::kernels::delete_memory< T, ct_Device >
 
using resmem_complex_op = ct::kernels::resize_memory< T, ct_Device >
 
using syncmem_complex_op = ct::kernels::synchronize_memory< T, ct_Device, ct_Device >
 
using gemm_op = ModuleBase::gemm_op< T, Device >
 

Private Member Functions

void calc_prec ()
 Update the precondition array.
 
void calc_hpsi_with_block (const HPsiFunc &hpsi_func, T *psi_in, ct::Tensor &hpsi_out)
 Apply the H operator to psi and obtain the hpsi matrix.
 
void diag_hsub (const ct::Tensor &psi_in, const ct::Tensor &hpsi_in, ct::Tensor &hsub_out, ct::Tensor &eigenvalue_out)
 Diagonalization of the subspace matrix.
 
void rotate_wf (const ct::Tensor &hsub_in, ct::Tensor &psi_out, ct::Tensor &workspace_in)
 Inplace matrix multiplication to obtain the initial guessed wavefunction.
 
void calc_grad_with_block (const ct::Tensor &prec_in, ct::Tensor &err_out, ct::Tensor &beta_out, ct::Tensor &psi_in, ct::Tensor &hpsi_in, ct::Tensor &grad_out, ct::Tensor &grad_old_out)
 Calculate the gradient for all bands used in CG method.
 
void calc_hsub_with_block (const HPsiFunc &hpsi_func, T *psi_in, ct::Tensor &psi_out, ct::Tensor &hpsi_out, ct::Tensor &hsub_out, ct::Tensor &workspace_in, ct::Tensor &eigenvalue_out)
 Apply the Hamiltonian operator to psi and obtain the hpsi matrix.
 
void calc_hsub_with_block_exit (ct::Tensor &psi_out, ct::Tensor &hpsi_out, ct::Tensor &hsub_out, ct::Tensor &workspace_in, ct::Tensor &eigenvalue_out)
 Apply the Hamiltonian operator to psi and obtain the hpsi matrix.
 
void orth_projection (const ct::Tensor &psi_in, ct::Tensor &hsub_in, ct::Tensor &grad_out)
 Orthogonalize column vectors in grad to column vectors in psi.
 
void line_minimize (ct::Tensor &grad_in, ct::Tensor &hgrad_in, ct::Tensor &psi_out, ct::Tensor &hpsi_out)
 Optimize psi as well as the hpsi.
 
void orth_cholesky (ct::Tensor &workspace_in, ct::Tensor &psi_out, ct::Tensor &hpsi_out, ct::Tensor &hsub_out)
 Orthogonalize and normalize the column vectors in psi_out using Cholesky decomposition.
 
bool test_error (const ct::Tensor &err_in, const std::vector< double > &ethr_band)
 Checks if the error satisfies the given threshold.
 

Private Attributes

int n_band = 0
 the number of bands of all processes
 
int n_band_l = 0
 the number of bands of current process
 
int n_basis = 0
 the number of cols of the input psi
 
int n_dim = 0
 valid dimension of psi
 
int nline = 4
 max iter steps for all-band cg loop
 
ModuleBase::PGemmCN< T, Device > pmmcn
 parallel matrix multiplication
 
PLinearTransform< T, Device > plintrans
 
ct::DataType r_type = ct::DataType::DT_INVALID
 
ct::DataType t_type = ct::DataType::DT_INVALID
 
ct::DeviceType device_type = ct::DeviceType::UnKnown
 
ct::Tensor prec = {}
 
ct::Tensor h_prec = {}
 
ct::Tensor beta = {}
 The coefficient for mixing the current and previous step gradients, used in iterative methods.
 
ct::Tensor err_st = {}
 Error state value, if it is smaller than the given threshold, then exit the iteration.
 
ct::Tensor eigen = {}
 Calculated eigen.
 
ct::Tensor psi = {}
 
ct::Tensor hpsi = {}
 
ct::Tensor hsub = {}
 
ct::Tensor grad = {}
 
ct::Tensor hgrad = {}
 
ct::Tensor grad_old = {}
 
ct::Tensor work = {}
 work for some calculations within this class, including rotate_wf call
 
Device * ctx = {}
 ctx is nothing but the devices used in gemm_op (Device * ctx = nullptr;),
 
const Tone = nullptr
 
const Tzero = nullptr
 
const Tneg_one = nullptr
 
const T one_ = static_cast<T>(1.0)
 
const T zero_ = static_cast<T>(0.0)
 
const T neg_one_ = static_cast<T>(-1.0)
 

Detailed Description

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
class hsolver::DiagoBPCG< T, Device >

A class for diagonalization using the Blocked-PCG method.

Template Parameters
TThe floating-point type used for calculations.
DeviceThe device used for calculations (e.g., cpu or gpu).

Member Typedef Documentation

◆ ct_Device

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::ct_Device = typename ct::PsiToContainer<Device>::type
private

◆ delmem_complex_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::delmem_complex_op = ct::kernels::delete_memory<T, ct_Device>
private

◆ delmem_var_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::delmem_var_op = ct::kernels::delete_memory<Real, ct_Device>
private

◆ gemm_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::gemm_op = ModuleBase::gemm_op<T, Device>
private

◆ HPsiFunc

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::HPsiFunc = std::function<void(T*, T*, const int, const int)>

◆ Real

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::Real = typename GetTypeReal<T>::type
private

◆ resmem_complex_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::resmem_complex_op = ct::kernels::resize_memory<T, ct_Device>
private

◆ resmem_var_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::resmem_var_op = ct::kernels::resize_memory<Real, ct_Device>
private

◆ setmem_complex_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::setmem_complex_op = ct::kernels::set_memory<T, ct_Device>
private

◆ setmem_var_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::setmem_var_op = ct::kernels::set_memory<Real, ct_Device>
private

◆ syncmem_complex_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::syncmem_complex_op = ct::kernels::synchronize_memory<T, ct_Device, ct_Device>
private

◆ syncmem_var_d2h_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::syncmem_var_d2h_op = ct::kernels::synchronize_memory<Real, ct::DEVICE_CPU, ct_Device>
private

◆ syncmem_var_h2d_op

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
using hsolver::DiagoBPCG< T, Device >::syncmem_var_h2d_op = ct::kernels::synchronize_memory<Real, ct_Device, ct::DEVICE_CPU>
private

Constructor & Destructor Documentation

◆ DiagoBPCG()

template<typename T , typename Device >
hsolver::DiagoBPCG< T, Device >::DiagoBPCG ( const Real precondition)
explicit

Constructor for DiagoBPCG class.

Parameters
preconditionprecondition data passed by the "Hamilt_PW" class.

◆ ~DiagoBPCG()

template<typename T , typename Device >
hsolver::DiagoBPCG< T, Device >::~DiagoBPCG ( )

Destructor for DiagoBPCG class.

Member Function Documentation

◆ calc_grad_with_block()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::calc_grad_with_block ( const ct::Tensor prec_in,
ct::Tensor err_out,
ct::Tensor beta_out,
ct::Tensor psi_in,
ct::Tensor hpsi_in,
ct::Tensor grad_out,
ct::Tensor grad_old_out 
)
private

Calculate the gradient for all bands used in CG method.

prec_in[dim: n_basis_max, column major], err_out[dim: n_band, column major], beta_out[dim: n_band, column major], psi_in[dim: n_basis x n_band, column major, lda = n_basis_max], hpsi_in[dim: n_basis x n_band, column major, lda = n_basis_max], grad_out[dim: n_basis x n_band, column major, lda = n_basis_max], grad_old_out[dim: n_basis x n_band, column major, lda = n_basis_max],

Parameters
prec_inInput preconditioner.
err_outOutput error state value. If it is smaller than a given threshold, exit the iteration.
beta_outOutput beta coefficient.
psi_inInput wavefunction matrix.
hpsi_inProduct of psi_in and Hamiltonian.
grad_outOutput gradient matrix.
grad_old_outPrevious gradient matrix.
Note
The steps involved in optimization are:
  1. normalize psi
  2. calculate the epsilo
  3. calculate the gradient by hpsi - epsilo * psi
  4. gradient mix with the previous gradient
  5. Do precondition
Here is the call graph for this function:

◆ calc_hpsi_with_block()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::calc_hpsi_with_block ( const HPsiFunc hpsi_func,
T psi_in,
ct::Tensor hpsi_out 
)
private

Apply the H operator to psi and obtain the hpsi matrix.

This function calculates the matrix product of the Hamiltonian operator (H) and the input wavefunction (psi). The resulting matrix is stored in the output array hpsi_out.

Note
hpsi = H|psi>;

psi_in[dim: n_basis x n_band, column major, lda = n_basis_max], hpsi_out[dim: n_basis x n_band, column major, lda = n_basis_max].

Parameters
hpsi_funcA function computing the product of the Hamiltonian matrix H and a wavefunction blockvector X.
psi_inThe input wavefunction psi.
hpsi_outPointer to the array where the resulting hpsi matrix will be stored.
Here is the call graph for this function:

◆ calc_hsub_with_block()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::calc_hsub_with_block ( const HPsiFunc hpsi_func,
T psi_in,
ct::Tensor psi_out,
ct::Tensor hpsi_out,
ct::Tensor hsub_out,
ct::Tensor workspace_in,
ct::Tensor eigenvalue_out 
)
private

Apply the Hamiltonian operator to psi and obtain the hpsi matrix.

psi_out[dim: n_basis x n_band, column major, lda = n_basis_max], hpsi_out[dim: n_basis x n_band, column major, lda = n_basis_max], hsub_out[dim: n_band x n_band, column major, lda = n_band], eigenvalue_out[dim: n_basis_max, column major].

Parameters
hpsi_funcA function computing the product of matrix H and wavefunction blockvector X.
psi_inInput wavefunction pointer.
psi_outOutput wavefunction.
hpsi_outProduct of psi_out and Hamiltonian.
hsub_outSubspace matrix output.
eigenvalue_outComputed eigen.

◆ calc_hsub_with_block_exit()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::calc_hsub_with_block_exit ( ct::Tensor psi_out,
ct::Tensor hpsi_out,
ct::Tensor hsub_out,
ct::Tensor workspace_in,
ct::Tensor eigenvalue_out 
)
private

Apply the Hamiltonian operator to psi and obtain the hpsi matrix.

psi_out[dim: n_basis x n_band, column major, lda = n_basis_max], hpsi_out[dim: n_basis x n_band, column major, lda = n_basis_max], hsub_out[dim: n_band x n_band, column major, lda = n_band], eigenvalue_out[dim: n_basis_max, column major].

Parameters
hamilt_inPointer to the Hamiltonian object.
psi_inInput wavefunction.
psi_outOutput wavefunction.
hpsi_outProduct of psi_out and Hamiltonian.
hsub_outSubspace matrix output.
eigenvalue_outComputed eigen.

◆ calc_prec()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::calc_prec ( )
private

Update the precondition array.

This function updates the precondition array by copying the host precondition to the device in a 'gpu' runtime environment. The address of the precondition array is passed by the constructor function called by hsolver::HSolverPW::initDiagh. The precondition will be updated before the DiagoBPCG<T, Device>::diag call.

Note
prec[dim: n_band]
Parameters
devReference to the base_device::AbacusDevice_t object, speciy which device used in the calc_prec function.
precPointer to the host precondition array with [dim: n_band, column major]
h_precPointer to the host precondition array with [dim: n_band, column major].
d_precPointer to the device precondition array with [dim: n_band, column major].

◆ diag()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::diag ( const HPsiFunc hpsi_func,
T psi_in,
Real eigenvalue_in,
const std::vector< double > &  ethr_band 
)

Diagonalize the Hamiltonian using the BPCG method.

This function is called by the HsolverPW::solve() function.

Parameters
hpsi_funcA function computing the product of the Hamiltonian matrix H and a wavefunction blockvector X.
psi_inPointer to input wavefunction psi matrix with [dim: n_basis x n_band, column major].
eigenvalue_inPointer to the eigen array with [dim: n_band, column major].
Here is the caller graph for this function:

◆ diag_hsub()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::diag_hsub ( const ct::Tensor psi_in,
const ct::Tensor hpsi_in,
ct::Tensor hsub_out,
ct::Tensor eigenvalue_out 
)
private

Diagonalization of the subspace matrix.

All the matrix used in this function are stored and used as the column major. psi_in[dim: n_basis x n_band, column major, lda = n_basis_max], hpsi_in[dim: n_basis x n_band, column major, lda = n_basis_max], hpsi_out[dim: n_basis x n_band, column major, lda = n_basis_max], eigenvalue_out[dim: n_basis_max, column major].

Parameters
psi_inInput wavefunction matrix with [dim: n_basis x n_band, column major].
hpsi_inH|psi> matrix with [dim: n_basis x n_band, column major].
hsub_outOutput Hamiltonian subtracted matrix with [dim: n_band x n_band, column major]
eigenvalue_outComputed eigen array with [dim: n_band]
Here is the call graph for this function:

◆ init_iter()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::init_iter ( const int  nband,
const int  nband_l,
const int  nbasis,
const int  ndim 
)

Initialize the class before diagonalization.

This function allocates all the related variables, such as hpsi, hsub, before the diag call. It is called by the HsolverPW::initDiagh() function.

Parameters
nbandThe number of bands of all processes.
nband_lThe number of bands of current process.
nbasisThe number of basis functions. Leading dimension of psi.
ndimThe number of valid dimension of psi.
Here is the caller graph for this function:

◆ line_minimize()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::line_minimize ( ct::Tensor grad_in,
ct::Tensor hgrad_in,
ct::Tensor psi_out,
ct::Tensor hpsi_out 
)
private

Optimize psi as well as the hpsi.

Parameters
grad_inInput gradient array, [dim: n_basis x n_band, column major, lda = n_basis_max].
hgrad_inProduct of grad_in and Hamiltonian, [dim: n_basis x n_band, column major, lda = n_basis_max].
psi_outInput and output wavefunction array, [dim: n_basis x n_band, column major, lda = n_basis_max].
hpsi_outProduct of psi_out and Hamiltonian, [dim: n_basis x n_band, column major, lda = n_basis_max].
Note
The steps involved in optimization are:
  1. Normalize the gradient.
  2. Calculate theta.
  3. Update psi as well as hpsi.
Here is the call graph for this function:

◆ orth_cholesky()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::orth_cholesky ( ct::Tensor workspace_in,
ct::Tensor psi_out,
ct::Tensor hpsi_out,
ct::Tensor hsub_out 
)
private

Orthogonalize and normalize the column vectors in psi_out using Cholesky decomposition.

Parameters
workspace_inWorkspace memory, [dim: n_basis x n_band, column major, lda = n_basis_max]..
psi_outInput and output wavefunction array. [dim: n_basis x n_band, column major, lda = n_basis_max].
hpsi_outInput and output hpsi array. [dim: n_basis x n_band, column major, lda = n_basis_max].
hsub_outInput Hamiltonian product array. [dim: n_band x n_band, column major, lda = n_band].
Here is the call graph for this function:

◆ orth_projection()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::orth_projection ( const ct::Tensor psi_in,
ct::Tensor hsub_in,
ct::Tensor grad_out 
)
private

Orthogonalize column vectors in grad to column vectors in psi.

hsub_in and workspace_in are only used to store intermediate variables of the gemm operator.

Parameters
psi_inInput wavefunction array, [dim: n_basis x n_band, column major, lda = n_basis_max].
hsub_inSubspace matrix input, [dim: n_band x n_band, column major, lda = n_band].
grad_outInput and output gradient array, [dim: n_basis x n_band, column major, lda = n_basis_max]..
Note
This function is a member of the DiagoBPCG class.
Here is the call graph for this function:

◆ rotate_wf()

template<typename T , typename Device >
void hsolver::DiagoBPCG< T, Device >::rotate_wf ( const ct::Tensor hsub_in,
ct::Tensor psi_out,
ct::Tensor workspace_in 
)
private

Inplace matrix multiplication to obtain the initial guessed wavefunction.

hsub_in[dim: n_band x n_band, column major, lda = n_band], workspace_in[dim: n_basis x n_band, column major, lda = n_basis_max], psi_out[dim: n_basis x n_band, column major, lda = n_basis_max],

Parameters
hsub_inSubspace matrix input, dim [n_basis, n_band] with column major.
psi_outoutput wavefunction matrix with dim [n_basis, n_band], column major.
Here is the call graph for this function:

◆ test_error()

template<typename T , typename Device >
bool hsolver::DiagoBPCG< T, Device >::test_error ( const ct::Tensor err_in,
const std::vector< double > &  ethr_band 
)
private

Checks if the error satisfies the given threshold.

Parameters
err_inPointer to the error array.[dim: n_band, column major]
thr_inThe threshold.
Returns
Returns true if all error values are less than or equal to the threshold, false otherwise.
Here is the call graph for this function:

Member Data Documentation

◆ beta

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::beta = {}
private

The coefficient for mixing the current and previous step gradients, used in iterative methods.

◆ ctx

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
Device* hsolver::DiagoBPCG< T, Device >::ctx = {}
private

ctx is nothing but the devices used in gemm_op (Device * ctx = nullptr;),

◆ device_type

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::DeviceType hsolver::DiagoBPCG< T, Device >::device_type = ct::DeviceType::UnKnown
private

◆ eigen

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::eigen = {}
private

Calculated eigen.

◆ err_st

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::err_st = {}
private

Error state value, if it is smaller than the given threshold, then exit the iteration.

◆ grad

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::grad = {}
private

H|psi> - epsilo * psi, grad of the given problem. Dim: n_basis * n_band, column major, lda = n_basis_max.

◆ grad_old

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::grad_old = {}
private

◆ h_prec

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::h_prec = {}
private

◆ hgrad

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::hgrad = {}
private

◆ hpsi

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::hpsi = {}
private

◆ hsub

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::hsub = {}
private

◆ n_band

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
int hsolver::DiagoBPCG< T, Device >::n_band = 0
private

the number of bands of all processes

◆ n_band_l

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
int hsolver::DiagoBPCG< T, Device >::n_band_l = 0
private

the number of bands of current process

◆ n_basis

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
int hsolver::DiagoBPCG< T, Device >::n_basis = 0
private

the number of cols of the input psi

◆ n_dim

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
int hsolver::DiagoBPCG< T, Device >::n_dim = 0
private

valid dimension of psi

◆ neg_one

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
const T * hsolver::DiagoBPCG< T, Device >::neg_one = nullptr
private

◆ neg_one_

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
const T hsolver::DiagoBPCG< T, Device >::neg_one_ = static_cast<T>(-1.0)
private

◆ nline

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
int hsolver::DiagoBPCG< T, Device >::nline = 4
private

max iter steps for all-band cg loop

◆ one

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
const T* hsolver::DiagoBPCG< T, Device >::one = nullptr
private

◆ one_

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
const T hsolver::DiagoBPCG< T, Device >::one_ = static_cast<T>(1.0)
private

◆ plintrans

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
PLinearTransform<T, Device> hsolver::DiagoBPCG< T, Device >::plintrans
private

◆ pmmcn

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ModuleBase::PGemmCN<T, Device> hsolver::DiagoBPCG< T, Device >::pmmcn
private

parallel matrix multiplication

◆ prec

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::prec = {}
private

◆ psi

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::psi = {}
private

Pointer to the input wavefunction. Note: this pointer does not own memory, instead it ref the psi_in object. H|psi> matrix.

◆ r_type

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::DataType hsolver::DiagoBPCG< T, Device >::r_type = ct::DataType::DT_INVALID
private

◆ t_type

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::DataType hsolver::DiagoBPCG< T, Device >::t_type = ct::DataType::DT_INVALID
private

◆ work

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
ct::Tensor hsolver::DiagoBPCG< T, Device >::work = {}
private

work for some calculations within this class, including rotate_wf call

◆ zero

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
const T * hsolver::DiagoBPCG< T, Device >::zero = nullptr
private

◆ zero_

template<typename T = std::complex<double>, typename Device = base_device::DEVICE_CPU>
const T hsolver::DiagoBPCG< T, Device >::zero_ = static_cast<T>(0.0)
private

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