From 43ca34caf592018efab5c41b3103f70e6a84fabc Mon Sep 17 00:00:00 2001 From: dzzz2001 Date: Thu, 27 Mar 2025 16:51:08 +0800 Subject: [PATCH 1/2] reduce the memory consumption of cal_gint_vlocal --- source/module_hamilt_lcao/module_gint/gint.h | 1 + .../module_gint/gint_vl.cpp | 13 ++-- .../module_gint/gint_vl_cpu_interface.cpp | 63 +++---------------- .../module_gint/temp_gint/gint_vl.cpp | 8 +-- .../module_gint/temp_gint/gint_vl_metagga.cpp | 14 ++--- .../temp_gint/gint_vl_metagga_nspin4.cpp | 19 ++---- .../module_gint/temp_gint/gint_vl_nspin4.cpp | 13 +--- .../module_gint/temp_gint/phi_operator.cpp | 17 +++-- .../module_gint/temp_gint/phi_operator.h | 3 +- .../module_hcontainer/base_matrix.h | 12 ++++ 10 files changed, 55 insertions(+), 108 deletions(-) diff --git a/source/module_hamilt_lcao/module_gint/gint.h b/source/module_hamilt_lcao/module_gint/gint.h index 22d75e650b9..9f6dda6a96a 100644 --- a/source/module_hamilt_lcao/module_gint/gint.h +++ b/source/module_hamilt_lcao/module_gint/gint.h @@ -140,6 +140,7 @@ class Gint { //! psir_ylm: dim is [bxyz][LD_pool] //! psir_vlbr3: dim is [bxyz][LD_pool] //! hR: HContainer for storing the matrix elements + //! cal_meshball_vlocal is thread-safe! void cal_meshball_vlocal( const int na_grid, const int LD_pool, diff --git a/source/module_hamilt_lcao/module_gint/gint_vl.cpp b/source/module_hamilt_lcao/module_gint/gint_vl.cpp index da33737c3cb..4da4ff7fd28 100644 --- a/source/module_hamilt_lcao/module_gint/gint_vl.cpp +++ b/source/module_hamilt_lcao/module_gint/gint_vl.cpp @@ -19,7 +19,7 @@ #include #endif - +// this is a thread-safe function void Gint::cal_meshball_vlocal( const int na_grid, // how many atoms on this (i,j,k) grid const int LD_pool, @@ -36,6 +36,7 @@ void Gint::cal_meshball_vlocal( const int lgd_now = this->gridt->lgd; const int mcell_index = this->gridt->bcell_start[grid_index]; + std::vector hr_tmp; for(int ia1=0; ia1get_row_size(); const int n = tmp_matrix->get_col_size(); - + hr_tmp.resize(m * n); + ModuleBase::GlobalFunc::ZEROS(hr_tmp.data(), m*n); + int cal_pair_num=0; for(int ib=first_ib;ibget_pointer(), &n); + &beta, hr_tmp.data(), &n); + tmp_matrix->add_array_ts(hr_tmp.data()); } else { @@ -103,7 +107,8 @@ void Gint::cal_meshball_vlocal( dgemm_(&transa, &transb, &n, &m, &k, &alpha, &psir_vlbr3[ib][block_index[ia2]], &LD_pool, &psir_ylm[ib][block_index[ia1]], &LD_pool, - &beta, tmp_matrix->get_pointer(), &n); + &beta, hr_tmp.data(), &n); + tmp_matrix->add_array_ts(hr_tmp.data()); } } } diff --git a/source/module_hamilt_lcao/module_gint/gint_vl_cpu_interface.cpp b/source/module_hamilt_lcao/module_gint/gint_vl_cpu_interface.cpp index 5bbab01fc08..0eb684ad5fa 100644 --- a/source/module_hamilt_lcao/module_gint/gint_vl_cpu_interface.cpp +++ b/source/module_hamilt_lcao/module_gint/gint_vl_cpu_interface.cpp @@ -19,7 +19,6 @@ void Gint::gint_kernel_vlocal(Gint_inout* inout) { { /** * @brief When in OpenMP, it points to a newly allocated memory, */ - hamilt::HContainer hRGint_thread(*hRGint_kernel); std::vector block_iw(max_size,0); std::vector block_index(max_size+1,0); std::vector block_size(max_size,0); @@ -74,19 +73,8 @@ void Gint::gint_kernel_vlocal(Gint_inout* inout) { this->cal_meshball_vlocal( na_grid, LD_pool, block_size.data(), block_index.data(), grid_index, cal_flag.get_ptr_2D(),psir_ylm.get_ptr_2D(), psir_vlbr3.get_ptr_2D(), - &hRGint_thread); + hRGint_kernel); } - - #pragma omp critical - { - BlasConnector::axpy(hRGint_thread.get_nnr(), - 1.0, - hRGint_thread.get_wrapper(), - 1, - hRGint_kernel->get_wrapper(), - 1); - } - ModuleBase::TITLE("Gint_interface", "cal_gint_vlocal"); ModuleBase::timer::tick("Gint_interface", "cal_gint_vlocal"); } @@ -112,9 +100,6 @@ void Gint::gint_kernel_dvlocal(Gint_inout* inout) { #pragma omp parallel { - hamilt::HContainer pvdpRx_thread(pvdpRx_reduced[inout->ispin]); - hamilt::HContainer pvdpRy_thread(pvdpRy_reduced[inout->ispin]); - hamilt::HContainer pvdpRz_thread(pvdpRz_reduced[inout->ispin]); std::vector block_iw(max_size,0); std::vector block_index(max_size+1,0); std::vector block_size(max_size,0); @@ -159,34 +144,13 @@ void Gint::gint_kernel_dvlocal(Gint_inout* inout) { //and accumulates to the corresponding element in Hamiltonian this->cal_meshball_vlocal(na_grid, LD_pool, block_size.data(), block_index.data(), grid_index, cal_flag.get_ptr_2D(),psir_vlbr3.get_ptr_2D(), - dpsir_ylm_x.get_ptr_2D(), &pvdpRx_thread); + dpsir_ylm_x.get_ptr_2D(), &this->pvdpRx_reduced[inout->ispin]); this->cal_meshball_vlocal(na_grid, LD_pool, block_size.data(), block_index.data(), grid_index, cal_flag.get_ptr_2D(),psir_vlbr3.get_ptr_2D(), - dpsir_ylm_y.get_ptr_2D(), &pvdpRy_thread); + dpsir_ylm_y.get_ptr_2D(), &this->pvdpRy_reduced[inout->ispin]); this->cal_meshball_vlocal(na_grid, LD_pool, block_size.data(), block_index.data(), grid_index, cal_flag.get_ptr_2D(),psir_vlbr3.get_ptr_2D(), - dpsir_ylm_z.get_ptr_2D(), &pvdpRz_thread); - } - #pragma omp critical(gint_k) - { - BlasConnector::axpy(nnrg, - 1.0, - pvdpRx_thread.get_wrapper(), - 1, - this->pvdpRx_reduced[inout->ispin].get_wrapper(), - 1); - BlasConnector::axpy(nnrg, - 1.0, - pvdpRy_thread.get_wrapper(), - 1, - this->pvdpRy_reduced[inout->ispin].get_wrapper(), - 1); - BlasConnector::axpy(nnrg, - 1.0, - pvdpRz_thread.get_wrapper(), - 1, - this->pvdpRz_reduced[inout->ispin].get_wrapper(), - 1); + dpsir_ylm_z.get_ptr_2D(), &this->pvdpRz_reduced[inout->ispin]); } } ModuleBase::TITLE("Gint_interface", "cal_gint_dvlocal"); @@ -210,7 +174,6 @@ void Gint::gint_kernel_vlocal_meta(Gint_inout* inout) { { // define HContainer here to reference. //Under the condition of gamma_only, hRGint will be instantiated. - hamilt::HContainer hRGint_thread(*hRGint_kernel); std::vector block_iw(max_size,0); std::vector block_index(max_size+1,0); std::vector block_size(max_size,0); @@ -282,28 +245,18 @@ void Gint::gint_kernel_vlocal_meta(Gint_inout* inout) { //and accumulates to the corresponding element in Hamiltonian this->cal_meshball_vlocal( na_grid, LD_pool, block_size.data(), block_index.data(), grid_index, cal_flag.get_ptr_2D(), - psir_ylm.get_ptr_2D(), psir_vlbr3.get_ptr_2D(), &hRGint_thread); + psir_ylm.get_ptr_2D(), psir_vlbr3.get_ptr_2D(), hRGint_kernel); //integrate (d/dx_i psi_mu*vk(r)*dv) * (d/dx_i psi_nu) on grid (x_i=x,y,z) //and accumulates to the corresponding element in Hamiltonian this->cal_meshball_vlocal( na_grid, LD_pool, block_size.data(), block_index.data(), grid_index, cal_flag.get_ptr_2D(), - dpsir_ylm_x.get_ptr_2D(), dpsix_vlbr3.get_ptr_2D(), &hRGint_thread); + dpsir_ylm_x.get_ptr_2D(), dpsix_vlbr3.get_ptr_2D(), hRGint_kernel); this->cal_meshball_vlocal( na_grid, LD_pool, block_size.data(), block_index.data(), grid_index, cal_flag.get_ptr_2D(), - dpsir_ylm_y.get_ptr_2D(), dpsiy_vlbr3.get_ptr_2D(), &hRGint_thread); + dpsir_ylm_y.get_ptr_2D(), dpsiy_vlbr3.get_ptr_2D(), hRGint_kernel); this->cal_meshball_vlocal( na_grid, LD_pool, block_size.data(), block_index.data(), grid_index, cal_flag.get_ptr_2D(), - dpsir_ylm_z.get_ptr_2D(), dpsiz_vlbr3.get_ptr_2D(), &hRGint_thread); - } - -#pragma omp critical - { - BlasConnector::axpy(nnrg, - 1.0, - hRGint_thread.get_wrapper(), - 1, - hRGint_kernel->get_wrapper(), - 1); + dpsir_ylm_z.get_ptr_2D(), dpsiz_vlbr3.get_ptr_2D(), hRGint_kernel); } } diff --git a/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl.cpp b/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl.cpp index 7a28bc94fd4..0861549a3c6 100644 --- a/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl.cpp +++ b/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl.cpp @@ -33,7 +33,6 @@ void Gint_vl::cal_hr_gint_() PhiOperator phi_op; std::vector phi; std::vector phi_vldr3; - HContainer hr_gint_local(*hr_gint_); #pragma omp for schedule(dynamic) for(const auto& biggrid: gint_info_->get_biggrids()) { @@ -47,12 +46,7 @@ void Gint_vl::cal_hr_gint_() phi_vldr3.resize(phi_len); phi_op.set_phi(phi.data()); phi_op.phi_mul_vldr3(vr_eff_, dr3_, phi.data(), phi_vldr3.data()); - phi_op.phi_mul_phi_vldr3(phi.data(), phi_vldr3.data(), &hr_gint_local); - } -#pragma omp critical - { - BlasConnector::axpy(hr_gint_local.get_nnr(), 1.0, hr_gint_local.get_wrapper(), - 1, hr_gint_->get_wrapper(), 1); + phi_op.phi_mul_phi_vldr3(phi.data(), phi_vldr3.data(), *hr_gint_); } } } diff --git a/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_metagga.cpp b/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_metagga.cpp index 81b4bf14982..e59db503632 100644 --- a/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_metagga.cpp +++ b/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_metagga.cpp @@ -39,7 +39,6 @@ void Gint_vl_metagga::cal_hr_gint_() std::vector dphi_x_vldr3; std::vector dphi_y_vldr3; std::vector dphi_z_vldr3; - HContainer hr_gint_local(*hr_gint_); #pragma omp for schedule(dynamic) for(const auto& biggrid: gint_info_->get_biggrids()) { @@ -62,15 +61,10 @@ void Gint_vl_metagga::cal_hr_gint_() phi_op.phi_mul_vldr3(vofk_, dr3_, dphi_x.data(), dphi_x_vldr3.data()); phi_op.phi_mul_vldr3(vofk_, dr3_, dphi_y.data(), dphi_y_vldr3.data()); phi_op.phi_mul_vldr3(vofk_, dr3_, dphi_z.data(), dphi_z_vldr3.data()); - phi_op.phi_mul_phi_vldr3(phi.data(), phi_vldr3.data(), &hr_gint_local); - phi_op.phi_mul_phi_vldr3(dphi_x.data(), dphi_x_vldr3.data(), &hr_gint_local); - phi_op.phi_mul_phi_vldr3(dphi_y.data(), dphi_y_vldr3.data(), &hr_gint_local); - phi_op.phi_mul_phi_vldr3(dphi_z.data(), dphi_z_vldr3.data(), &hr_gint_local); - } -#pragma omp critical - { - BlasConnector::axpy(hr_gint_local.get_nnr(), 1.0, hr_gint_local.get_wrapper(), - 1, hr_gint_->get_wrapper(), 1); + phi_op.phi_mul_phi_vldr3(phi.data(), phi_vldr3.data(), *hr_gint_); + phi_op.phi_mul_phi_vldr3(dphi_x.data(), dphi_x_vldr3.data(), *hr_gint_); + phi_op.phi_mul_phi_vldr3(dphi_y.data(), dphi_y_vldr3.data(), *hr_gint_); + phi_op.phi_mul_phi_vldr3(dphi_z.data(), dphi_z_vldr3.data(), *hr_gint_); } } } diff --git a/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_metagga_nspin4.cpp b/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_metagga_nspin4.cpp index 9734605192f..1b727b97cc6 100644 --- a/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_metagga_nspin4.cpp +++ b/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_metagga_nspin4.cpp @@ -41,7 +41,6 @@ void Gint_vl_metagga_nspin4::cal_hr_gint_() std::vector dphi_x_vldr3; std::vector dphi_y_vldr3; std::vector dphi_z_vldr3; - std::vector> hr_gint_part_thread(nspin_, *hr_gint_part_[0]); #pragma omp for schedule(dynamic) for(const auto& biggrid: gint_info_->get_biggrids()) { @@ -66,20 +65,10 @@ void Gint_vl_metagga_nspin4::cal_hr_gint_() phi_op.phi_mul_vldr3(vofk_[is], dr3_, dphi_x.data(), dphi_x_vldr3.data()); phi_op.phi_mul_vldr3(vofk_[is], dr3_, dphi_y.data(), dphi_y_vldr3.data()); phi_op.phi_mul_vldr3(vofk_[is], dr3_, dphi_z.data(), dphi_z_vldr3.data()); - phi_op.phi_mul_phi_vldr3(phi.data(), phi_vldr3.data(), &hr_gint_part_thread[is]); - phi_op.phi_mul_phi_vldr3(dphi_x.data(), dphi_x_vldr3.data(), &hr_gint_part_thread[is]); - phi_op.phi_mul_phi_vldr3(dphi_y.data(), dphi_y_vldr3.data(), &hr_gint_part_thread[is]); - phi_op.phi_mul_phi_vldr3(dphi_z.data(), dphi_z_vldr3.data(), &hr_gint_part_thread[is]); - } - } -#pragma omp critical - { - for(int is = 0; is < nspin_; is++) - { - { - BlasConnector::axpy(hr_gint_part_thread[is].get_nnr(), 1.0, hr_gint_part_thread[is].get_wrapper(), - 1, hr_gint_part_[is]->get_wrapper(), 1); - } + phi_op.phi_mul_phi_vldr3(phi.data(), phi_vldr3.data(), *hr_gint_part_[is]); + phi_op.phi_mul_phi_vldr3(dphi_x.data(), dphi_x_vldr3.data(), *hr_gint_part_[is]); + phi_op.phi_mul_phi_vldr3(dphi_y.data(), dphi_y_vldr3.data(), *hr_gint_part_[is]); + phi_op.phi_mul_phi_vldr3(dphi_z.data(), dphi_z_vldr3.data(), *hr_gint_part_[is]); } } } diff --git a/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_nspin4.cpp b/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_nspin4.cpp index 25f582f8c64..694fcebfcca 100644 --- a/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_nspin4.cpp +++ b/source/module_hamilt_lcao/module_gint/temp_gint/gint_vl_nspin4.cpp @@ -34,7 +34,6 @@ void Gint_vl_nspin4::cal_hr_gint_() PhiOperator phi_op; std::vector phi; std::vector phi_vldr3; - std::vector> hr_gint_part_thread(nspin_, *hr_gint_part_[0]); #pragma omp for schedule(dynamic) for(const auto& biggrid: gint_info_->get_biggrids()) { @@ -50,17 +49,7 @@ void Gint_vl_nspin4::cal_hr_gint_() for(int is = 0; is < nspin_; is++) { phi_op.phi_mul_vldr3(vr_eff_[is], dr3_, phi.data(), phi_vldr3.data()); - phi_op.phi_mul_phi_vldr3(phi.data(), phi_vldr3.data(), &hr_gint_part_thread[is]); - } - } -#pragma omp critical - { - for(int is = 0; is < nspin_; is++) - { - { - BlasConnector::axpy(hr_gint_part_thread[is].get_nnr(), 1.0, hr_gint_part_thread[is].get_wrapper(), - 1, hr_gint_part_[is]->get_wrapper(), 1); - } + phi_op.phi_mul_phi_vldr3(phi.data(), phi_vldr3.data(), *hr_gint_part_[is]); } } } diff --git a/source/module_hamilt_lcao/module_gint/temp_gint/phi_operator.cpp b/source/module_hamilt_lcao/module_gint/temp_gint/phi_operator.cpp index 01b3fc69de3..d5a923c88a8 100644 --- a/source/module_hamilt_lcao/module_gint/temp_gint/phi_operator.cpp +++ b/source/module_hamilt_lcao/module_gint/temp_gint/phi_operator.cpp @@ -148,25 +148,29 @@ void PhiOperator::phi_mul_vldr3(const double* vl, const double dr3, const double } } +// this is a thread-safe function void PhiOperator::phi_mul_phi_vldr3( const double* phi, const double* phi_vldr3, - HContainer* hr) const + HContainer& hr) const { const char transa='N', transb='T'; const double alpha=1, beta=1; + std::vector tmp_hr; for(int i = 0; i < biggrid_->get_atoms_num(); ++i) { const auto atom_i = biggrid_->get_atom(i); const auto& r_i = atom_i->get_R(); const int iat_i = atom_i->get_iat(); + const int m = atoms_phi_len_[i]; for(int j = 0; j < biggrid_->get_atoms_num(); ++j) { const auto atom_j = biggrid_->get_atom(j); const auto& r_j = atom_j->get_R(); const int iat_j = atom_j->get_iat(); + const int n = atoms_phi_len_[j]; // only calculate the upper triangle matrix if(iat_i > iat_j) @@ -174,8 +178,11 @@ void PhiOperator::phi_mul_phi_vldr3( continue; } + tmp_hr.resize(m * n); + ModuleBase::GlobalFunc::ZEROS(tmp_hr.data(), m*n); + // FIXME may be r = r_j - r_i - const auto result = hr->find_matrix(iat_i, iat_j, r_i-r_j); + const auto result = hr.find_matrix(iat_i, iat_j, r_i-r_j); if(result == nullptr) { @@ -191,8 +198,10 @@ void PhiOperator::phi_mul_phi_vldr3( continue; } - dgemm_(&transa, &transb, &atoms_phi_len_[j], &atoms_phi_len_[i], &len, &alpha, &phi_vldr3[start_idx * cols_ + atoms_startidx_[j]], - &cols_,&phi[start_idx * cols_ + atoms_startidx_[i]], &cols_, &beta, result->get_pointer(), &atoms_phi_len_[j]); + dgemm_(&transa, &transb, &n, &m, &len, &alpha, &phi_vldr3[start_idx * cols_ + atoms_startidx_[j]], + &cols_,&phi[start_idx * cols_ + atoms_startidx_[i]], &cols_, &beta, tmp_hr.data(), &n); + + result->add_array_ts(tmp_hr.data()); } } } diff --git a/source/module_hamilt_lcao/module_gint/temp_gint/phi_operator.h b/source/module_hamilt_lcao/module_gint/temp_gint/phi_operator.h index ea15cb9804d..fafac776c39 100644 --- a/source/module_hamilt_lcao/module_gint/temp_gint/phi_operator.h +++ b/source/module_hamilt_lcao/module_gint/temp_gint/phi_operator.h @@ -55,10 +55,11 @@ class PhiOperator const double* phi, double* result) const; + // this is a thread-safe function void phi_mul_phi_vldr3( const double* phi, const double* phi_vldr3, - HContainer* hr) const; + HContainer& hr) const; void phi_dot_phi_dm( const double* phi, diff --git a/source/module_hamilt_lcao/module_hcontainer/base_matrix.h b/source/module_hamilt_lcao/module_hcontainer/base_matrix.h index 61ecb8c5994..db56ac46bfb 100644 --- a/source/module_hamilt_lcao/module_hcontainer/base_matrix.h +++ b/source/module_hamilt_lcao/module_hcontainer/base_matrix.h @@ -2,6 +2,7 @@ #define BASE_MATRIX_H #include +#include namespace hamilt { @@ -107,6 +108,15 @@ class BaseMatrix */ void set_size(const int& col_size_in, const int& row_size_in); + void add_array_ts(T* array) + { + std::lock_guard lock(mtx); + for (int i = 0; i < nrow_local * ncol_local; ++i) + { + value_begin[i] += array[i]; + } + } + private: bool allocated = false; @@ -118,6 +128,8 @@ class BaseMatrix // int current_multiple = 0; + // for thread safe + mutable std::mutex mtx; // number of rows and columns int nrow_local = 0; int ncol_local = 0; From a8b25040c4e89c350163c2af0cad0707f0aeabfd Mon Sep 17 00:00:00 2001 From: dzzz2001 Date: Fri, 28 Mar 2025 13:42:17 +0800 Subject: [PATCH 2/2] fix a bug --- .../module_gint/gint_vl.cpp | 33 +++---------------- 1 file changed, 5 insertions(+), 28 deletions(-) diff --git a/source/module_hamilt_lcao/module_gint/gint_vl.cpp b/source/module_hamilt_lcao/module_gint/gint_vl.cpp index 4da4ff7fd28..043f80225ad 100644 --- a/source/module_hamilt_lcao/module_gint/gint_vl.cpp +++ b/source/module_hamilt_lcao/module_gint/gint_vl.cpp @@ -84,34 +84,11 @@ void Gint::cal_meshball_vlocal( hr_tmp.resize(m * n); ModuleBase::GlobalFunc::ZEROS(hr_tmp.data(), m*n); - int cal_pair_num=0; - for(int ib=first_ib;ibib_length/4) - { - dgemm_(&transa, &transb, &n, &m, &ib_length, &alpha, - &psir_vlbr3[first_ib][block_index[ia2]], &LD_pool, - &psir_ylm[first_ib][block_index[ia1]], &LD_pool, - &beta, hr_tmp.data(), &n); - tmp_matrix->add_array_ts(hr_tmp.data()); - } - else - { - for(int ib=first_ib; ibadd_array_ts(hr_tmp.data()); - } - } - } + dgemm_(&transa, &transb, &n, &m, &ib_length, &alpha, + &psir_vlbr3[first_ib][block_index[ia2]], &LD_pool, + &psir_ylm[first_ib][block_index[ia1]], &LD_pool, + &beta, hr_tmp.data(), &n); + tmp_matrix->add_array_ts(hr_tmp.data()); } } }