ABACUS develop
Atomic-orbital Based Ab-initio Computation at UStc
|
A class for diagonalization using the Blocked-PCG method. More...
#include <diago_bpcg.h>
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 > ðr_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 > ðr_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 T * | one = nullptr |
const T * | zero = nullptr |
const T * | neg_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) |
A class for diagonalization using the Blocked-PCG method.
T | The floating-point type used for calculations. |
Device | The device used for calculations (e.g., cpu or gpu). |
|
private |
|
private |
|
private |
|
private |
using hsolver::DiagoBPCG< T, Device >::HPsiFunc = std::function<void(T*, T*, const int, const int)> |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
explicit |
Constructor for DiagoBPCG class.
precondition | precondition data passed by the "Hamilt_PW" class. |
hsolver::DiagoBPCG< T, Device >::~DiagoBPCG | ( | ) |
Destructor for DiagoBPCG class.
|
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],
prec_in | Input preconditioner. |
err_out | Output error state value. If it is smaller than a given threshold, exit the iteration. |
beta_out | Output beta coefficient. |
psi_in | Input wavefunction matrix. |
hpsi_in | Product of psi_in and Hamiltonian. |
grad_out | Output gradient matrix. |
grad_old_out | Previous gradient matrix. |
|
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.
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].
hpsi_func | A function computing the product of the Hamiltonian matrix H and a wavefunction blockvector X. |
psi_in | The input wavefunction psi. |
hpsi_out | Pointer to the array where the resulting hpsi matrix will be stored. |
|
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].
hpsi_func | A function computing the product of matrix H and wavefunction blockvector X. |
psi_in | Input wavefunction pointer. |
psi_out | Output wavefunction. |
hpsi_out | Product of psi_out and Hamiltonian. |
hsub_out | Subspace matrix output. |
eigenvalue_out | Computed eigen. |
|
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].
hamilt_in | Pointer to the Hamiltonian object. |
psi_in | Input wavefunction. |
psi_out | Output wavefunction. |
hpsi_out | Product of psi_out and Hamiltonian. |
hsub_out | Subspace matrix output. |
eigenvalue_out | Computed eigen. |
|
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.
dev | Reference to the base_device::AbacusDevice_t object, speciy which device used in the calc_prec function. |
prec | Pointer to the host precondition array with [dim: n_band, column major] |
h_prec | Pointer to the host precondition array with [dim: n_band, column major]. |
d_prec | Pointer to the device precondition array with [dim: n_band, column major]. |
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.
hpsi_func | A function computing the product of the Hamiltonian matrix H and a wavefunction blockvector X. |
psi_in | Pointer to input wavefunction psi matrix with [dim: n_basis x n_band, column major]. |
eigenvalue_in | Pointer to the eigen array with [dim: n_band, column major]. |
|
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].
psi_in | Input wavefunction matrix with [dim: n_basis x n_band, column major]. |
hpsi_in | H|psi> matrix with [dim: n_basis x n_band, column major]. |
hsub_out | Output Hamiltonian subtracted matrix with [dim: n_band x n_band, column major] |
eigenvalue_out | Computed eigen array with [dim: n_band] |
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.
nband | The number of bands of all processes. |
nband_l | The number of bands of current process. |
nbasis | The number of basis functions. Leading dimension of psi. |
ndim | The number of valid dimension of psi. |
|
private |
Optimize psi as well as the hpsi.
grad_in | Input gradient array, [dim: n_basis x n_band, column major, lda = n_basis_max]. |
hgrad_in | Product of grad_in and Hamiltonian, [dim: n_basis x n_band, column major, lda = n_basis_max]. |
psi_out | Input and output wavefunction array, [dim: n_basis x n_band, column major, lda = n_basis_max]. |
hpsi_out | Product of psi_out and Hamiltonian, [dim: n_basis x n_band, column major, lda = n_basis_max]. |
|
private |
Orthogonalize and normalize the column vectors in psi_out using Cholesky decomposition.
workspace_in | Workspace memory, [dim: n_basis x n_band, column major, lda = n_basis_max].. |
psi_out | Input and output wavefunction array. [dim: n_basis x n_band, column major, lda = n_basis_max]. |
hpsi_out | Input and output hpsi array. [dim: n_basis x n_band, column major, lda = n_basis_max]. |
hsub_out | Input Hamiltonian product array. [dim: n_band x n_band, column major, lda = n_band]. |
|
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.
psi_in | Input wavefunction array, [dim: n_basis x n_band, column major, lda = n_basis_max]. |
hsub_in | Subspace matrix input, [dim: n_band x n_band, column major, lda = n_band]. |
grad_out | Input and output gradient array, [dim: n_basis x n_band, column major, lda = n_basis_max].. |
|
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],
hsub_in | Subspace matrix input, dim [n_basis, n_band] with column major. |
psi_out | output wavefunction matrix with dim [n_basis, n_band], column major. |
|
private |
Checks if the error satisfies the given threshold.
err_in | Pointer to the error array.[dim: n_band, column major] |
thr_in | The threshold. |
|
private |
The coefficient for mixing the current and previous step gradients, used in iterative methods.
|
private |
ctx is nothing but the devices used in gemm_op (Device * ctx = nullptr;),
|
private |
|
private |
Calculated eigen.
|
private |
Error state value, if it is smaller than the given threshold, then exit the iteration.
|
private |
H|psi> - epsilo * psi, grad of the given problem. Dim: n_basis * n_band, column major, lda = n_basis_max.
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
the number of bands of all processes
|
private |
the number of bands of current process
|
private |
the number of cols of the input psi
|
private |
valid dimension of psi
|
private |
|
private |
|
private |
max iter steps for all-band cg loop
|
private |
|
private |
|
private |
|
private |
parallel matrix multiplication
|
private |
|
private |
Pointer to the input wavefunction. Note: this pointer does not own memory, instead it ref the psi_in object. H|psi> matrix.
|
private |
|
private |
|
private |
work for some calculations within this class, including rotate_wf call
|
private |
|
private |