From e6676228be5f40c1cadcacf72f398bc2c1700492 Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Thu, 28 Nov 2024 20:22:49 +0800 Subject: [PATCH 01/10] Refactor: Remove bpcg dependency on Psi and Hamilt --- source/module_hsolver/diago_bpcg.cpp | 53 +++++++++++++++++----------- source/module_hsolver/diago_bpcg.h | 21 +++++++---- source/module_hsolver/hsolver_pw.cpp | 24 +++++++++++-- 3 files changed, 69 insertions(+), 29 deletions(-) diff --git a/source/module_hsolver/diago_bpcg.cpp b/source/module_hsolver/diago_bpcg.cpp index 63d0956464..6a870b3461 100644 --- a/source/module_hsolver/diago_bpcg.cpp +++ b/source/module_hsolver/diago_bpcg.cpp @@ -27,14 +27,17 @@ DiagoBPCG::DiagoBPCG(const Real* precondition_in) template DiagoBPCG::~DiagoBPCG() { // Note, we do not need to free the h_prec and psi pointer as they are refs to the outside data - delete this->grad_wrapper; + // delete this->grad_wrapper; } template -void DiagoBPCG::init_iter(const psi::Psi &psi_in) { +void DiagoBPCG::init_iter(/*const T *psi_in,*/ const int nband, const int nbasis) { // Specify the problem size n_basis, n_band, while lda is n_basis - this->n_band = psi_in.get_nbands(); - this->n_basis = psi_in.get_nbasis(); + // this->n_band = psi_in.get_nbands(); + // this->n_basis = psi_in.get_nbasis(); + this->n_band = nband; + this->n_basis = nbasis; + // All column major tensors @@ -52,8 +55,9 @@ void DiagoBPCG::init_iter(const psi::Psi &psi_in) { this->prec = std::move(ct::Tensor(r_type, device_type, {this->n_basis})); //TODO: Remove class Psi, using ct::Tensor instead! - this->grad_wrapper = new psi::Psi(1, this->n_band, this->n_basis, psi_in.get_ngk_pointer()); - this->grad = std::move(ct::TensorMap(grad_wrapper->get_pointer(), t_type, device_type, {this->n_band, this->n_basis})); + // this->grad_wrapper = new psi::Psi(1, this->n_band, this->n_basis, psi_in.get_ngk_pointer()); + // this->grad = std::move(ct::TensorMap(grad_wrapper->get_pointer(), t_type, device_type, {this->n_band, this->n_basis})); + this->grad = std::move(ct::Tensor(t_type, device_type, {this->n_band, this->n_basis})); } template @@ -174,16 +178,19 @@ void DiagoBPCG::rotate_wf( template void DiagoBPCG::calc_hpsi_with_block( - hamilt::Hamilt* hamilt_in, - const psi::Psi& psi_in, + // hamilt::Hamilt* hamilt_in, + const HPsiFunc& hpsi_func, + // const psi::Psi& psi_in, + T *psi_in, ct::Tensor& hpsi_out) { // calculate all-band hpsi - psi::Range all_bands_range(1, psi_in.get_current_k(), 0, psi_in.get_nbands() - 1); - hpsi_info info(&psi_in, all_bands_range, hpsi_out.data()); - hamilt_in->ops->hPsi(info); + // psi::Range all_bands_range(1, psi_in.get_current_k(), 0, psi_in.get_nbands() - 1); + // hpsi_info info(&psi_in, all_bands_range, hpsi_out.data()); + // hamilt_in->ops->hPsi(info); + hpsi_func(psi_in, hpsi_out.data(), this->n_basis, this->n_band); - return; + // return; } template @@ -207,8 +214,10 @@ void DiagoBPCG::diag_hsub( template void DiagoBPCG::calc_hsub_with_block( - hamilt::Hamilt *hamilt_in, - const psi::Psi &psi_in, + // hamilt::Hamilt* hamilt_in, + const HPsiFunc& hpsi_func, + // const psi::Psi& psi_in, + T *psi_in, ct::Tensor& psi_out, ct::Tensor& hpsi_out, ct::Tensor& hsub_out, @@ -216,7 +225,7 @@ void DiagoBPCG::calc_hsub_with_block( ct::Tensor& eigenvalue_out) { // Apply the H operator to psi and obtain the hpsi matrix. - this->calc_hpsi_with_block(hamilt_in, psi_in, hpsi_out); + this->calc_hpsi_with_block(hpsi_func, psi_in, hpsi_out); // Diagonalization of the subspace matrix. this->diag_hsub(psi_out,hpsi_out, hsub_out, eigenvalue_out); @@ -250,19 +259,21 @@ void DiagoBPCG::calc_hsub_with_block_exit( template void DiagoBPCG::diag( - hamilt::Hamilt* hamilt_in, - psi::Psi& psi_in, + // hamilt::Hamilt* hamilt_in, + const HPsiFunc& hpsi_func, + // psi::Psi& psi_in, + T *psi_in, Real* eigenvalue_in) { const int current_scf_iter = hsolver::DiagoIterAssist::SCF_ITER; // Get the pointer of the input psi - this->psi = std::move(ct::TensorMap(psi_in.get_pointer(), t_type, device_type, {this->n_band, this->n_basis})); + this->psi = std::move(ct::TensorMap(psi_in /*psi_in.get_pointer()*/, t_type, device_type, {this->n_band, this->n_basis})); // Update the precondition array this->calc_prec(); // Improving the initial guess of the wave function psi through a subspace diagonalization. - this->calc_hsub_with_block(hamilt_in, psi_in, this->psi, this->hpsi, this->hsub, this->work, this->eigen); + this->calc_hsub_with_block(hpsi_func, psi_in, this->psi, this->hpsi, this->hsub, this->work, this->eigen); setmem_complex_op()(this->grad_old.template data(), 0, this->n_basis * this->n_band); @@ -293,7 +304,7 @@ void DiagoBPCG::diag( syncmem_complex_op()(this->grad_old.template data(), this->grad.template data(), n_basis * n_band); // Calculate H|grad> matrix - this->calc_hpsi_with_block(hamilt_in, this->grad_wrapper[0], this->hgrad); + this->calc_hpsi_with_block(hpsi_func, this->grad.data(), /*this->grad_wrapper[0],*/ this->hgrad); // optimize psi as well as the hpsi // 1. normalize grad @@ -305,7 +316,7 @@ void DiagoBPCG::diag( this->orth_cholesky(this->work, this->psi, this->hpsi, this->hsub); if (current_scf_iter == 1 && ntry % this->nline == 0) { - this->calc_hsub_with_block(hamilt_in, psi_in, this->psi, this->hpsi, this->hsub, this->work, this->eigen); + this->calc_hsub_with_block(hpsi_func, psi_in, this->psi, this->hpsi, this->hsub, this->work, this->eigen); } } while (ntry < max_iter && this->test_error(this->err_st, this->all_band_cg_thr)); diff --git a/source/module_hsolver/diago_bpcg.h b/source/module_hsolver/diago_bpcg.h index c14472de2e..37d7e9eddc 100644 --- a/source/module_hsolver/diago_bpcg.h +++ b/source/module_hsolver/diago_bpcg.h @@ -52,7 +52,9 @@ class DiagoBPCG * * @param psi_in The input wavefunction psi. */ - void init_iter(const psi::Psi &psi_in); + void init_iter(/*const T *psi_in,*/ const int nband, const int nbasis); + + using HPsiFunc = std::function; /** * @brief Diagonalize the Hamiltonian using the BPCG method. @@ -63,7 +65,10 @@ class DiagoBPCG * @param psi The input wavefunction psi matrix with [dim: n_basis x n_band, column major]. * @param eigenvalue_in Pointer to the eigen array with [dim: n_band, column major]. */ - void diag(hamilt::Hamilt *phm_in, psi::Psi &psi, Real *eigenvalue_in); + void diag(// hamilt::Hamilt *phm_in, + const HPsiFunc& hpsi_func, + // psi::Psi& psi_in, + T *psi_in, Real *eigenvalue_in); private: @@ -139,8 +144,10 @@ class DiagoBPCG * @param hpsi_out Pointer to the array where the resulting hpsi matrix will be stored. */ void calc_hpsi_with_block( - hamilt::Hamilt* hamilt_in, - const psi::Psi& psi_in, + // hamilt::Hamilt* hamilt_in, + const HPsiFunc& hpsi_func, + // const psi::Psi& psi_in, + T *psi_in, ct::Tensor& hpsi_out); /** @@ -228,8 +235,10 @@ class DiagoBPCG * @param eigenvalue_out Computed eigen. */ void calc_hsub_with_block( - hamilt::Hamilt* hamilt_in, - const psi::Psi& psi_in, + // hamilt::Hamilt* hamilt_in, + const HPsiFunc& hpsi_func, + // const psi::Psi& psi_in, + T *psi_in, ct::Tensor& psi_out, ct::Tensor& hpsi_out, ct::Tensor& hsub_out, ct::Tensor& workspace_in, ct::Tensor& eigenvalue_out); diff --git a/source/module_hsolver/hsolver_pw.cpp b/source/module_hsolver/hsolver_pw.cpp index 945a65f87d..b2d40db522 100644 --- a/source/module_hsolver/hsolver_pw.cpp +++ b/source/module_hsolver/hsolver_pw.cpp @@ -467,9 +467,29 @@ void HSolverPW::hamiltSolvePsiK(hamilt::Hamilt* hm, } else if (this->method == "bpcg") { + const int nband = psi.get_nbands(); + const int nbasis = psi.get_nbasis(); + auto ngk_pointer = psi.get_ngk_pointer(); + // hpsi_func (X, HX, ld, nvec) -> HX = H(X), X and HX blockvectors of size ld x nvec + auto hpsi_func = [hm, ngk_pointer](T* psi_in, T* hpsi_out, const int ld_psi, const int nvec) { + ModuleBase::timer::tick("DavSubspace", "hpsi_func"); + + // Convert "pointer data stucture" to a psi::Psi object + auto psi_iter_wrapper = psi::Psi(psi_in, 1, nvec, ld_psi, ngk_pointer); + + psi::Range bands_range(true, 0, 0, nvec - 1); + + using hpsi_info = typename hamilt::Operator::hpsi_info; + hpsi_info info(&psi_iter_wrapper, bands_range, hpsi_out); + hm->ops->hPsi(info); + + ModuleBase::timer::tick("DavSubspace", "hpsi_func"); + }; DiagoBPCG bpcg(pre_condition.data()); - bpcg.init_iter(psi); - bpcg.diag(hm, psi, eigenvalue); + // bpcg.init_iter(psi); + bpcg.init_iter(nband, nbasis); + // bpcg.diag(hm, psi, eigenvalue); + bpcg.diag(hpsi_func, psi.get_pointer(), eigenvalue); } else if (this->method == "dav_subspace") { From 8d549fdb86eb6fed861312a1fb0e12762d441a2b Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Thu, 28 Nov 2024 20:48:53 +0800 Subject: [PATCH 02/10] Test: change bpcg tests to fit new interface --- .../module_hsolver/test/diago_bpcg_test.cpp | 30 ++++++++++++++++--- 1 file changed, 26 insertions(+), 4 deletions(-) diff --git a/source/module_hsolver/test/diago_bpcg_test.cpp b/source/module_hsolver/test/diago_bpcg_test.cpp index 39eacf5db6..e723486b4e 100644 --- a/source/module_hsolver/test/diago_bpcg_test.cpp +++ b/source/module_hsolver/test/diago_bpcg_test.cpp @@ -130,10 +130,32 @@ class DiagoBPCGPrepare psi_local.fix_k(0); double start, end; start = MPI_Wtime(); - bpcg.init_iter(psi_local); - bpcg.diag(ha,psi_local,en); - bpcg.diag(ha,psi_local,en); - bpcg.diag(ha,psi_local,en); + using T = std::complex; + const int dim = DIAGOTEST::npw; + const std::vector &h_mat = DIAGOTEST::hmatrix; + auto hpsi_func = [h_mat, dim](T *psi_in, T *hpsi_out, + const int ld_psi, const int nvec) { + auto one = std::make_unique(1.0); + auto zero = std::make_unique(0.0); + const T *one_ = one.get(); + const T *zero_ = zero.get(); + + base_device::DEVICE_CPU *ctx = {}; + // hpsi_out(dim * nvec) = h_mat(dim * dim) * psi_in(dim * nvec) + hsolver::gemm_op()( + ctx, 'N', 'N', + dim, nvec, dim, + one_, + h_mat.data(), dim, + psi_in, ld_psi, + zero_, + hpsi_out, ld_psi); + }; + bpcg.init_iter(nband, npw); + bpcg.diag(hpsi_func, psi_local.get_pointer(), en); + // bpcg.diag(ha,psi_local,en); + // bpcg.diag(ha,psi_local,en); + // bpcg.diag(ha,psi_local,en); end = MPI_Wtime(); //if(mypnum == 0) printf("diago time:%7.3f\n",end-start); delete [] DIAGOTEST::npw_local; From 2d21ef90f071a337d2c839f472aac31ab754ed86 Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Thu, 28 Nov 2024 20:54:35 +0800 Subject: [PATCH 03/10] Test: make bpcg restart the same time as before in test --- source/module_hsolver/test/diago_bpcg_test.cpp | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/source/module_hsolver/test/diago_bpcg_test.cpp b/source/module_hsolver/test/diago_bpcg_test.cpp index e723486b4e..fa5f8b3c07 100644 --- a/source/module_hsolver/test/diago_bpcg_test.cpp +++ b/source/module_hsolver/test/diago_bpcg_test.cpp @@ -152,7 +152,9 @@ class DiagoBPCGPrepare hpsi_out, ld_psi); }; bpcg.init_iter(nband, npw); - bpcg.diag(hpsi_func, psi_local.get_pointer(), en); + bpcg.diag(hpsi_func, psi_local.get_pointer(), en); + bpcg.diag(hpsi_func, psi_local.get_pointer(), en); + bpcg.diag(hpsi_func, psi_local.get_pointer(), en); // bpcg.diag(ha,psi_local,en); // bpcg.diag(ha,psi_local,en); // bpcg.diag(ha,psi_local,en); From 2088b8122879e05083c63391b7cb9104dcc61b01 Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Thu, 28 Nov 2024 21:22:44 +0800 Subject: [PATCH 04/10] Tests: add the template disambiguator for dependent names for bpcg tests --- source/module_hsolver/diago_bpcg.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/source/module_hsolver/diago_bpcg.cpp b/source/module_hsolver/diago_bpcg.cpp index 6a870b3461..31a4adb829 100644 --- a/source/module_hsolver/diago_bpcg.cpp +++ b/source/module_hsolver/diago_bpcg.cpp @@ -304,7 +304,7 @@ void DiagoBPCG::diag( syncmem_complex_op()(this->grad_old.template data(), this->grad.template data(), n_basis * n_band); // Calculate H|grad> matrix - this->calc_hpsi_with_block(hpsi_func, this->grad.data(), /*this->grad_wrapper[0],*/ this->hgrad); + this->calc_hpsi_with_block(hpsi_func, this->grad.template data(), /*this->grad_wrapper[0],*/ this->hgrad); // optimize psi as well as the hpsi // 1. normalize grad From 585857db1c0d9d2716e258e2a3f8145df54b0b85 Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Fri, 29 Nov 2024 15:06:43 +0800 Subject: [PATCH 05/10] Refactor: clean up bpcg --- source/module_hsolver/diago_bpcg.h | 11 +--- .../module_hsolver/test/diago_bpcg_test.cpp | 51 +++++++++---------- 2 files changed, 27 insertions(+), 35 deletions(-) diff --git a/source/module_hsolver/diago_bpcg.h b/source/module_hsolver/diago_bpcg.h index 37d7e9eddc..b60b54cb39 100644 --- a/source/module_hsolver/diago_bpcg.h +++ b/source/module_hsolver/diago_bpcg.h @@ -52,7 +52,7 @@ class DiagoBPCG * * @param psi_in The input wavefunction psi. */ - void init_iter(/*const T *psi_in,*/ const int nband, const int nbasis); + void init_iter(const int nband, const int nbasis); using HPsiFunc = std::function; @@ -65,10 +65,7 @@ class DiagoBPCG * @param psi The input wavefunction psi matrix with [dim: n_basis x n_band, column major]. * @param eigenvalue_in Pointer to the eigen array with [dim: n_band, column major]. */ - void diag(// hamilt::Hamilt *phm_in, - const HPsiFunc& hpsi_func, - // psi::Psi& psi_in, - T *psi_in, Real *eigenvalue_in); + void diag(const HPsiFunc& hpsi_func, T *psi_in, Real *eigenvalue_in); private: @@ -144,9 +141,7 @@ class DiagoBPCG * @param hpsi_out Pointer to the array where the resulting hpsi matrix will be stored. */ void calc_hpsi_with_block( - // hamilt::Hamilt* hamilt_in, const HPsiFunc& hpsi_func, - // const psi::Psi& psi_in, T *psi_in, ct::Tensor& hpsi_out); @@ -235,9 +230,7 @@ class DiagoBPCG * @param eigenvalue_out Computed eigen. */ void calc_hsub_with_block( - // hamilt::Hamilt* hamilt_in, const HPsiFunc& hpsi_func, - // const psi::Psi& psi_in, T *psi_in, ct::Tensor& psi_out, ct::Tensor& hpsi_out, ct::Tensor& hsub_out, ct::Tensor& workspace_in, diff --git a/source/module_hsolver/test/diago_bpcg_test.cpp b/source/module_hsolver/test/diago_bpcg_test.cpp index fa5f8b3c07..4d3abe52ae 100644 --- a/source/module_hsolver/test/diago_bpcg_test.cpp +++ b/source/module_hsolver/test/diago_bpcg_test.cpp @@ -132,7 +132,7 @@ class DiagoBPCGPrepare start = MPI_Wtime(); using T = std::complex; const int dim = DIAGOTEST::npw; - const std::vector &h_mat = DIAGOTEST::hmatrix; + const std::vector &h_mat = DIAGOTEST::hmatrix_local; auto hpsi_func = [h_mat, dim](T *psi_in, T *hpsi_out, const int ld_psi, const int nvec) { auto one = std::make_unique(1.0); @@ -155,9 +155,7 @@ class DiagoBPCGPrepare bpcg.diag(hpsi_func, psi_local.get_pointer(), en); bpcg.diag(hpsi_func, psi_local.get_pointer(), en); bpcg.diag(hpsi_func, psi_local.get_pointer(), en); - // bpcg.diag(ha,psi_local,en); - // bpcg.diag(ha,psi_local,en); - // bpcg.diag(ha,psi_local,en); + bpcg.diag(hpsi_func, psi_local.get_pointer(), en); end = MPI_Wtime(); //if(mypnum == 0) printf("diago time:%7.3f\n",end-start); delete [] DIAGOTEST::npw_local; @@ -243,29 +241,30 @@ TEST(DiagoBPCGTest, Hamilt) } }*/ +// This test will not work for now! // bpcg for a 2x2 matrix -#ifdef __MPI -#else -TEST(DiagoBPCGTest, TwoByTwo) -{ - int dim = 2; - int nband = 2; - ModuleBase::ComplexMatrix hm(2, 2); - hm(0, 0) = std::complex{4.0, 0.0}; - hm(0, 1) = std::complex{1.0, 0.0}; - hm(1, 0) = std::complex{1.0, 0.0}; - hm(1, 1) = std::complex{3.0, 0.0}; - // nband, npw, sub, sparsity, reorder, eps, maxiter, threshold - DiagoBPCGPrepare dcp(nband, dim, 0, true, 1e-4, 50, 1e-10); - hsolver::DiagoIterAssist>::PW_DIAG_NMAX = dcp.maxiter; - hsolver::DiagoIterAssist>::PW_DIAG_THR = dcp.eps; - HPsi> hpsi; - hpsi.create(nband, dim); - DIAGOTEST::hmatrix = hm; - DIAGOTEST::npw = dim; - dcp.CompareEigen(hpsi.precond()); -} -#endif +// #ifdef __MPI +// #else +// TEST(DiagoBPCGTest, TwoByTwo) +// { +// int dim = 2; +// int nband = 2; +// ModuleBase::ComplexMatrix hm(2, 2); +// hm(0, 0) = std::complex{4.0, 0.0}; +// hm(0, 1) = std::complex{1.0, 0.0}; +// hm(1, 0) = std::complex{1.0, 0.0}; +// hm(1, 1) = std::complex{3.0, 0.0}; +// // nband, npw, sub, sparsity, reorder, eps, maxiter, threshold +// DiagoBPCGPrepare dcp(nband, dim, 0, true, 1e-4, 50, 1e-10); +// hsolver::DiagoIterAssist>::PW_DIAG_NMAX = dcp.maxiter; +// hsolver::DiagoIterAssist>::PW_DIAG_THR = dcp.eps; +// HPsi> hpsi; +// hpsi.create(nband, dim); +// DIAGOTEST::hmatrix = hm; +// DIAGOTEST::npw = dim; +// dcp.CompareEigen(hpsi.precond()); +// } +// #endif TEST(DiagoBPCGTest, readH) { From effef630a6c8ffc2591dac3a4984b63d138c5f5f Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Fri, 29 Nov 2024 15:20:20 +0800 Subject: [PATCH 06/10] Docs: new BPCG interface --- source/module_hsolver/diago_bpcg.h | 19 +++++++++++-------- source/module_hsolver/hsolver_pw.cpp | 4 ++-- 2 files changed, 13 insertions(+), 10 deletions(-) diff --git a/source/module_hsolver/diago_bpcg.h b/source/module_hsolver/diago_bpcg.h index b60b54cb39..353ec225a8 100644 --- a/source/module_hsolver/diago_bpcg.h +++ b/source/module_hsolver/diago_bpcg.h @@ -50,7 +50,8 @@ class DiagoBPCG * This function allocates all the related variables, such as hpsi, hsub, before the diag call. * It is called by the HsolverPW::initDiagh() function. * - * @param psi_in The input wavefunction psi. + * @param nband The number of bands. + * @param nbasis The number of basis functions. Leading dimension of psi. */ void init_iter(const int nband, const int nbasis); @@ -61,8 +62,9 @@ class DiagoBPCG * * This function is called by the HsolverPW::solve() function. * - * @param phm_in A pointer to the hamilt::Hamilt object representing the Hamiltonian operator. - * @param psi The input wavefunction psi matrix with [dim: n_basis x n_band, column major]. + * @param hpsi_func A function computing the product of the Hamiltonian matrix H + * and a wavefunction blockvector X. + * @param psi_in Pointer to input wavefunction psi matrix with [dim: n_basis x n_band, column major]. * @param eigenvalue_in Pointer to the eigen array with [dim: n_band, column major]. */ void diag(const HPsiFunc& hpsi_func, T *psi_in, Real *eigenvalue_in); @@ -105,7 +107,7 @@ class DiagoBPCG /// work for some calculations within this class, including rotate_wf call ct::Tensor work = {}; - psi::Psi* grad_wrapper; + // psi::Psi* grad_wrapper; /** * @brief Update the precondition array. * @@ -136,7 +138,8 @@ class DiagoBPCG * 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]. * - * @param hamilt_in A pointer to the hamilt::Hamilt object representing the Hamiltonian operator. + * @param hpsi_func A function computing the product of the Hamiltonian matrix H + * and a wavefunction blockvector X. * @param psi_in The input wavefunction psi. * @param hpsi_out Pointer to the array where the resulting hpsi matrix will be stored. */ @@ -222,8 +225,8 @@ class DiagoBPCG * hsub_out[dim: n_band x n_band, column major, lda = n_band], * eigenvalue_out[dim: n_basis_max, column major]. * - * @param hamilt_in Pointer to the Hamiltonian object. - * @param psi_in Input wavefunction. + * @param hpsi_func A function computing the product of matrix H and wavefunction blockvector X. + * @param psi_in Input wavefunction pointer. * @param psi_out Output wavefunction. * @param hpsi_out Product of psi_out and Hamiltonian. * @param hsub_out Subspace matrix output. @@ -316,7 +319,7 @@ class DiagoBPCG */ bool test_error(const ct::Tensor& err_in, Real thr_in); - using hpsi_info = typename hamilt::Operator::hpsi_info; + // using hpsi_info = typename hamilt::Operator::hpsi_info; using ct_Device = typename ct::PsiToContainer::type; using setmem_var_op = ct::kernels::set_memory; diff --git a/source/module_hsolver/hsolver_pw.cpp b/source/module_hsolver/hsolver_pw.cpp index b2d40db522..0ecc048631 100644 --- a/source/module_hsolver/hsolver_pw.cpp +++ b/source/module_hsolver/hsolver_pw.cpp @@ -467,8 +467,8 @@ void HSolverPW::hamiltSolvePsiK(hamilt::Hamilt* hm, } else if (this->method == "bpcg") { - const int nband = psi.get_nbands(); - const int nbasis = psi.get_nbasis(); + const int nband = psi.get_nbands(); + const int nbasis = psi.get_nbasis(); auto ngk_pointer = psi.get_ngk_pointer(); // hpsi_func (X, HX, ld, nvec) -> HX = H(X), X and HX blockvectors of size ld x nvec auto hpsi_func = [hm, ngk_pointer](T* psi_in, T* hpsi_out, const int ld_psi, const int nvec) { From eb61c5ceeac34631b095d643d7264bf6c1ce8e07 Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Fri, 29 Nov 2024 15:22:20 +0800 Subject: [PATCH 07/10] clean up bpcg code --- source/module_hsolver/diago_bpcg.cpp | 8 +------- 1 file changed, 1 insertion(+), 7 deletions(-) diff --git a/source/module_hsolver/diago_bpcg.cpp b/source/module_hsolver/diago_bpcg.cpp index 31a4adb829..edc7359ec9 100644 --- a/source/module_hsolver/diago_bpcg.cpp +++ b/source/module_hsolver/diago_bpcg.cpp @@ -31,7 +31,7 @@ DiagoBPCG::~DiagoBPCG() { } template -void DiagoBPCG::init_iter(/*const T *psi_in,*/ const int nband, const int nbasis) { +void DiagoBPCG::init_iter(const int nband, const int nbasis) { // Specify the problem size n_basis, n_band, while lda is n_basis // this->n_band = psi_in.get_nbands(); // this->n_basis = psi_in.get_nbasis(); @@ -178,9 +178,7 @@ void DiagoBPCG::rotate_wf( template void DiagoBPCG::calc_hpsi_with_block( - // hamilt::Hamilt* hamilt_in, const HPsiFunc& hpsi_func, - // const psi::Psi& psi_in, T *psi_in, ct::Tensor& hpsi_out) { @@ -214,9 +212,7 @@ void DiagoBPCG::diag_hsub( template void DiagoBPCG::calc_hsub_with_block( - // hamilt::Hamilt* hamilt_in, const HPsiFunc& hpsi_func, - // const psi::Psi& psi_in, T *psi_in, ct::Tensor& psi_out, ct::Tensor& hpsi_out, @@ -259,9 +255,7 @@ void DiagoBPCG::calc_hsub_with_block_exit( template void DiagoBPCG::diag( - // hamilt::Hamilt* hamilt_in, const HPsiFunc& hpsi_func, - // psi::Psi& psi_in, T *psi_in, Real* eigenvalue_in) { From c1e8046f98dd30bdd5b52fda57bfba6715c18d34 Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Fri, 29 Nov 2024 19:01:05 +0800 Subject: [PATCH 08/10] clean useless code --- source/module_hsolver/diago_bpcg.cpp | 10 -------- source/module_hsolver/hsolver_pw.cpp | 2 -- .../module_hsolver/test/diago_bpcg_test.cpp | 24 ------------------- 3 files changed, 36 deletions(-) diff --git a/source/module_hsolver/diago_bpcg.cpp b/source/module_hsolver/diago_bpcg.cpp index edc7359ec9..82b1cdf0a8 100644 --- a/source/module_hsolver/diago_bpcg.cpp +++ b/source/module_hsolver/diago_bpcg.cpp @@ -33,8 +33,6 @@ DiagoBPCG::~DiagoBPCG() { template void DiagoBPCG::init_iter(const int nband, const int nbasis) { // Specify the problem size n_basis, n_band, while lda is n_basis - // this->n_band = psi_in.get_nbands(); - // this->n_basis = psi_in.get_nbasis(); this->n_band = nband; this->n_basis = nbasis; @@ -54,9 +52,6 @@ void DiagoBPCG::init_iter(const int nband, const int nbasis) { this->prec = std::move(ct::Tensor(r_type, device_type, {this->n_basis})); - //TODO: Remove class Psi, using ct::Tensor instead! - // this->grad_wrapper = new psi::Psi(1, this->n_band, this->n_basis, psi_in.get_ngk_pointer()); - // this->grad = std::move(ct::TensorMap(grad_wrapper->get_pointer(), t_type, device_type, {this->n_band, this->n_basis})); this->grad = std::move(ct::Tensor(t_type, device_type, {this->n_band, this->n_basis})); } @@ -183,12 +178,7 @@ void DiagoBPCG::calc_hpsi_with_block( ct::Tensor& hpsi_out) { // calculate all-band hpsi - // psi::Range all_bands_range(1, psi_in.get_current_k(), 0, psi_in.get_nbands() - 1); - // hpsi_info info(&psi_in, all_bands_range, hpsi_out.data()); - // hamilt_in->ops->hPsi(info); hpsi_func(psi_in, hpsi_out.data(), this->n_basis, this->n_band); - - // return; } template diff --git a/source/module_hsolver/hsolver_pw.cpp b/source/module_hsolver/hsolver_pw.cpp index 0ecc048631..4572495ccd 100644 --- a/source/module_hsolver/hsolver_pw.cpp +++ b/source/module_hsolver/hsolver_pw.cpp @@ -486,9 +486,7 @@ void HSolverPW::hamiltSolvePsiK(hamilt::Hamilt* hm, ModuleBase::timer::tick("DavSubspace", "hpsi_func"); }; DiagoBPCG bpcg(pre_condition.data()); - // bpcg.init_iter(psi); bpcg.init_iter(nband, nbasis); - // bpcg.diag(hm, psi, eigenvalue); bpcg.diag(hpsi_func, psi.get_pointer(), eigenvalue); } else if (this->method == "dav_subspace") diff --git a/source/module_hsolver/test/diago_bpcg_test.cpp b/source/module_hsolver/test/diago_bpcg_test.cpp index 4d3abe52ae..6274f75d74 100644 --- a/source/module_hsolver/test/diago_bpcg_test.cpp +++ b/source/module_hsolver/test/diago_bpcg_test.cpp @@ -241,30 +241,6 @@ TEST(DiagoBPCGTest, Hamilt) } }*/ -// This test will not work for now! -// bpcg for a 2x2 matrix -// #ifdef __MPI -// #else -// TEST(DiagoBPCGTest, TwoByTwo) -// { -// int dim = 2; -// int nband = 2; -// ModuleBase::ComplexMatrix hm(2, 2); -// hm(0, 0) = std::complex{4.0, 0.0}; -// hm(0, 1) = std::complex{1.0, 0.0}; -// hm(1, 0) = std::complex{1.0, 0.0}; -// hm(1, 1) = std::complex{3.0, 0.0}; -// // nband, npw, sub, sparsity, reorder, eps, maxiter, threshold -// DiagoBPCGPrepare dcp(nband, dim, 0, true, 1e-4, 50, 1e-10); -// hsolver::DiagoIterAssist>::PW_DIAG_NMAX = dcp.maxiter; -// hsolver::DiagoIterAssist>::PW_DIAG_THR = dcp.eps; -// HPsi> hpsi; -// hpsi.create(nband, dim); -// DIAGOTEST::hmatrix = hm; -// DIAGOTEST::npw = dim; -// dcp.CompareEigen(hpsi.precond()); -// } -// #endif TEST(DiagoBPCGTest, readH) { From 4a74ea03ed98640e4e8a214046a4d34f00479656 Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Fri, 29 Nov 2024 19:10:36 +0800 Subject: [PATCH 09/10] clean useless code --- source/module_hsolver/diago_bpcg.cpp | 1 - source/module_hsolver/diago_bpcg.h | 2 -- 2 files changed, 3 deletions(-) diff --git a/source/module_hsolver/diago_bpcg.cpp b/source/module_hsolver/diago_bpcg.cpp index 82b1cdf0a8..aaa21877bc 100644 --- a/source/module_hsolver/diago_bpcg.cpp +++ b/source/module_hsolver/diago_bpcg.cpp @@ -27,7 +27,6 @@ DiagoBPCG::DiagoBPCG(const Real* precondition_in) template DiagoBPCG::~DiagoBPCG() { // Note, we do not need to free the h_prec and psi pointer as they are refs to the outside data - // delete this->grad_wrapper; } template diff --git a/source/module_hsolver/diago_bpcg.h b/source/module_hsolver/diago_bpcg.h index 353ec225a8..269fb696b5 100644 --- a/source/module_hsolver/diago_bpcg.h +++ b/source/module_hsolver/diago_bpcg.h @@ -319,8 +319,6 @@ class DiagoBPCG */ bool test_error(const ct::Tensor& err_in, Real thr_in); - // using hpsi_info = typename hamilt::Operator::hpsi_info; - using ct_Device = typename ct::PsiToContainer::type; using setmem_var_op = ct::kernels::set_memory; using resmem_var_op = ct::kernels::resize_memory; From e3d3a5936b16ca16938755fb3580d93f2b6f2644 Mon Sep 17 00:00:00 2001 From: Chen Nuo <49788094+Cstandardlib@users.noreply.github.com> Date: Fri, 29 Nov 2024 19:11:54 +0800 Subject: [PATCH 10/10] clean useless code --- source/module_hsolver/diago_bpcg.h | 1 - 1 file changed, 1 deletion(-) diff --git a/source/module_hsolver/diago_bpcg.h b/source/module_hsolver/diago_bpcg.h index 269fb696b5..2ca8167f9e 100644 --- a/source/module_hsolver/diago_bpcg.h +++ b/source/module_hsolver/diago_bpcg.h @@ -107,7 +107,6 @@ class DiagoBPCG /// work for some calculations within this class, including rotate_wf call ct::Tensor work = {}; - // psi::Psi* grad_wrapper; /** * @brief Update the precondition array. *