diff --git a/source/Makefile.Objects b/source/Makefile.Objects index 04cbff62277..d8b6c391aec 100644 --- a/source/Makefile.Objects +++ b/source/Makefile.Objects @@ -677,7 +677,6 @@ OBJS_MODULE_RI=conv_coulomb_pot_k.o\ Matrix_Orbs22.o\ RI_2D_Comm.o\ Mix_DMk_2D.o\ - Mix_Matrix.o\ symmetry_rotation.o\ symmetry_rotation_output.o\ symmetry_irreducible_sector.o\ diff --git a/source/source_lcao/module_ri/CMakeLists.txt b/source/source_lcao/module_ri/CMakeLists.txt index fd3b167ac68..ec1344b1d58 100644 --- a/source/source_lcao/module_ri/CMakeLists.txt +++ b/source/source_lcao/module_ri/CMakeLists.txt @@ -8,7 +8,6 @@ if (ENABLE_LIBRI) Matrix_Orbs22.cpp RI_2D_Comm.cpp Mix_DMk_2D.cpp - Mix_Matrix.cpp ) if(ENABLE_LCAO) diff --git a/source/source_lcao/module_ri/Exx_LRI_interface.h b/source/source_lcao/module_ri/Exx_LRI_interface.h index 55ca9312daf..f4d4282b1c5 100644 --- a/source/source_lcao/module_ri/Exx_LRI_interface.h +++ b/source/source_lcao/module_ri/Exx_LRI_interface.h @@ -124,7 +124,8 @@ class Exx_LRI_Interface std::shared_ptr> exx_ptr; private: - Mix_DMk_2D mix_DMk_2D; + + Mix_DMk_2D mix_DMk_2D; bool exx_spacegroup_symmetry = false; ModuleSymmetry::Symmetry_rotation symrot_; diff --git a/source/source_lcao/module_ri/Exx_LRI_interface.hpp b/source/source_lcao/module_ri/Exx_LRI_interface.hpp index 5272cb7b10b..6a0c74aa091 100644 --- a/source/source_lcao/module_ri/Exx_LRI_interface.hpp +++ b/source/source_lcao/module_ri/Exx_LRI_interface.hpp @@ -35,9 +35,7 @@ void Exx_LRI_Interface::cal_exx_ions(const UnitCell& ucell, const bool ModuleBase::TITLE("Exx_LRI_Interface","cal_exx_ions"); if(!this->flag_finish.init) { throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } - this->exx_ptr->cal_exx_ions(ucell, write_cv); - this->flag_finish.ions = true; } @@ -49,13 +47,9 @@ void Exx_LRI_Interface::cal_exx_elec(const std::vectorflag_finish.init || !this->flag_finish.ions) - { - throw std::runtime_error("Exx init unfinished when " - +std::string(__FILE__)+" line "+std::to_string(__LINE__)); - } + { throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } this->exx_ptr->cal_exx_elec(Ds, ucell, pv, p_symrot); - this->flag_finish.elec = true; } @@ -64,17 +58,11 @@ void Exx_LRI_Interface::cal_exx_force(const int& nat) { ModuleBase::TITLE("Exx_LRI_Interface","cal_exx_force"); if(!this->flag_finish.init || !this->flag_finish.ions) - { - throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); - } + { throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } if(!this->flag_finish.elec) - { - throw std::runtime_error("Exx Hamiltonian unfinished when "+std::string(__FILE__) - +" line "+std::to_string(__LINE__)); - } + { throw std::runtime_error("Exx Hamiltonian unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } this->exx_ptr->cal_exx_force(nat); - this->flag_finish.force = true; } @@ -83,24 +71,17 @@ void Exx_LRI_Interface::cal_exx_stress(const double& omega, const doub { ModuleBase::TITLE("Exx_LRI_Interface","cal_exx_stress"); if(!this->flag_finish.init || !this->flag_finish.ions) - { - throw std::runtime_error("Exx init unfinished when " - +std::string(__FILE__)+" line "+std::to_string(__LINE__)); - } + { throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } if(!this->flag_finish.elec) - { - throw std::runtime_error("Exx Hamiltonian unfinished when " - +std::string(__FILE__)+" line "+std::to_string(__LINE__)); - } + { throw std::runtime_error("Exx Hamiltonian unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } this->exx_ptr->cal_exx_stress(omega, lat0); - this->flag_finish.stress = true; } template void Exx_LRI_Interface::exx_before_all_runners( - const K_Vectors& kv, + const K_Vectors& kv, const UnitCell& ucell, const Parallel_2D& pv) { @@ -147,14 +128,15 @@ void Exx_LRI_Interface::exx_beforescf(const int istep, if(GlobalC::exx_info.info_global.cal_exx) { if (this->exx_spacegroup_symmetry) - {this->mix_DMk_2D.set_nks(kv.get_nkstot_full() * (PARAM.inp.nspin == 2 ? 2 : 1), PARAM.globalv.gamma_only_local);} + { this->mix_DMk_2D.set_nks(kv.get_nkstot_full() * (PARAM.inp.nspin == 2 ? 2 : 1)); } else - {this->mix_DMk_2D.set_nks(kv.get_nks(), PARAM.globalv.gamma_only_local);} + { this->mix_DMk_2D.set_nks(kv.get_nks()); } - if(GlobalC::exx_info.info_global.separate_loop) - { this->mix_DMk_2D.set_mixing(nullptr); } + if (GlobalC::exx_info.info_global.separate_loop) + { this->mix_DMk_2D.set_mixing_plain(GlobalC::exx_info.info_global.mixing_beta_for_loop1); } else { this->mix_DMk_2D.set_mixing(chgmix.get_mixing()); } + // for exx two_level scf this->two_level_step = 0; } @@ -171,13 +153,13 @@ void Exx_LRI_Interface::exx_eachiterinit(const int istep, ModuleBase::TITLE("Exx_LRI_Interface","exx_eachiterinit"); if (GlobalC::exx_info.info_global.cal_exx) { - if (!GlobalC::exx_info.info_global.separate_loop - && (this->two_level_step - || istep > 0 + if (!GlobalC::exx_info.info_global.separate_loop + && (this->two_level_step + || istep > 0 || PARAM.inp.init_wfc == "file") // non separate loop case - || (GlobalC::exx_info.info_global.separate_loop - && PARAM.inp.init_wfc == "file" - && this->two_level_step == 0 + || (GlobalC::exx_info.info_global.separate_loop + && PARAM.inp.init_wfc == "file" + && this->two_level_step == 0 && iter == 1) ) // the first iter in separate loop case { @@ -186,33 +168,25 @@ void Exx_LRI_Interface::exx_eachiterinit(const int istep, auto cal = [this, &ucell,&kv, &flag_restart](const elecstate::DensityMatrix& dm_in) { if (this->exx_spacegroup_symmetry) - { this->mix_DMk_2D.mix(symrot_.restore_dm(kv,dm_in.get_DMK_vector(), *dm_in.get_paraV_pointer()), flag_restart); } + { this->mix_DMk_2D.mix(symrot_.restore_dm(kv, dm_in.get_DMK_vector(), *dm_in.get_paraV_pointer()), flag_restart); } else { this->mix_DMk_2D.mix(dm_in.get_DMK_vector(), flag_restart); } - const std::vector>>> - Ds = PARAM.globalv.gamma_only_local - ? RI_2D_Comm::split_m2D_ktoR( - ucell, - *this->exx_ptr->p_kv, - this->mix_DMk_2D.get_DMk_gamma_out(), - *dm_in.get_paraV_pointer(), - PARAM.inp.nspin) - : RI_2D_Comm::split_m2D_ktoR( - ucell, - *this->exx_ptr->p_kv, - this->mix_DMk_2D.get_DMk_k_out(), - *dm_in.get_paraV_pointer(), - PARAM.inp.nspin, - this->exx_spacegroup_symmetry); - - if (this->exx_spacegroup_symmetry && GlobalC::exx_info.info_ri.exx_symmetry_realspace) + const std::vector>>> Ds = + RI_2D_Comm::split_m2D_ktoR( + ucell, + *this->exx_ptr->p_kv, + this->mix_DMk_2D.get_DMk_out(), + *dm_in.get_paraV_pointer(), + PARAM.inp.nspin, + this->exx_spacegroup_symmetry); + if(this->exx_spacegroup_symmetry && GlobalC::exx_info.info_ri.exx_symmetry_realspace) { this->cal_exx_elec(Ds, ucell,*dm_in.get_paraV_pointer(), &this->symrot_); } else { this->cal_exx_elec(Ds, ucell,*dm_in.get_paraV_pointer()); } }; if(istep > 0 && flag_restart) - { cal(*dm_last_step); } + { cal(*this->dm_last_step); } else { cal(dm); } } @@ -387,21 +361,23 @@ bool Exx_LRI_Interface::exx_after_converge( // if init_wfc == "file", DM is calculated in the 1st iter of the 1st two-level step, so we mix it here const bool flag_restart = (this->two_level_step == 0 && PARAM.inp.init_wfc != "file") ? true : false; - if (this->exx_spacegroup_symmetry) - {this->mix_DMk_2D.mix(symrot_.restore_dm(kv, dm.get_DMK_vector(), *dm.get_paraV_pointer()), flag_restart);} + if(this->exx_spacegroup_symmetry) + { this->mix_DMk_2D.mix(symrot_.restore_dm(kv, dm.get_DMK_vector(), *dm.get_paraV_pointer()), flag_restart); } else - {this->mix_DMk_2D.mix(dm.get_DMK_vector(), flag_restart);} - - // GlobalC::exx_lcao.cal_exx_elec(p_esolver->LOC, p_esolver->LOWF.wfc_k_grid); - const std::vector>, RI::Tensor>>> - Ds = std::is_same::value //gamma_only_local - ? RI_2D_Comm::split_m2D_ktoR(ucell,*this->exx_ptr->p_kv, this->mix_DMk_2D.get_DMk_gamma_out(), *dm.get_paraV_pointer(), nspin) - : RI_2D_Comm::split_m2D_ktoR(ucell,*this->exx_ptr->p_kv, this->mix_DMk_2D.get_DMk_k_out(), *dm.get_paraV_pointer(), nspin, this->exx_spacegroup_symmetry); - - if (this->exx_spacegroup_symmetry && GlobalC::exx_info.info_ri.exx_symmetry_realspace) + { this->mix_DMk_2D.mix(dm.get_DMK_vector(), flag_restart); } + const std::vector>>> Ds = + RI_2D_Comm::split_m2D_ktoR( + ucell, + *this->exx_ptr->p_kv, + this->mix_DMk_2D.get_DMk_out(), + *dm.get_paraV_pointer(), + nspin, + this->exx_spacegroup_symmetry); + if(this->exx_spacegroup_symmetry && GlobalC::exx_info.info_ri.exx_symmetry_realspace) { this->cal_exx_elec(Ds, ucell, *dm.get_paraV_pointer(), &this->symrot_); } else { this->cal_exx_elec(Ds, ucell, *dm.get_paraV_pointer()); } // restore DM but not Hexx + iter = 0; this->two_level_step++; diff --git a/source/source_lcao/module_ri/Mix_DMk_2D.cpp b/source/source_lcao/module_ri/Mix_DMk_2D.cpp index 286075a5d66..15b3e8f7dd7 100644 --- a/source/source_lcao/module_ri/Mix_DMk_2D.cpp +++ b/source/source_lcao/module_ri/Mix_DMk_2D.cpp @@ -4,69 +4,85 @@ //======================= #include "Mix_DMk_2D.h" +#include "source_base/module_mixing/plain_mixing.h" #include "source_base/tool_title.h" -Mix_DMk_2D &Mix_DMk_2D::set_nks(const int nks, const bool gamma_only_in) +#include + +template +Mix_DMk_2D::~Mix_DMk_2D() +{ + if(this->flag_del_mixing) + delete this->mixing; +} + +template +void Mix_DMk_2D::set_nks(const int nks) { - ModuleBase::TITLE("Mix_DMk_2D", "set_nks"); - this->gamma_only = gamma_only_in; - if (this->gamma_only) - this->mix_DMk_gamma.resize(nks); - else - this->mix_DMk_k.resize(nks); - return *this; + this->mix_DMk.clear(); + this->mix_DMk.resize(nks); } -Mix_DMk_2D &Mix_DMk_2D::set_mixing(Base_Mixing::Mixing* mixing_in) +template +void Mix_DMk_2D::set_mixing(Base_Mixing::Mixing* mixing_in) { - ModuleBase::TITLE("Mix_DMk_2D","set_mixing"); - if(this->gamma_only) - for (Mix_Matrix>& mix_one : this->mix_DMk_gamma) - mix_one.init(mixing_in); - else - for (Mix_Matrix>>& mix_one : this->mix_DMk_k) - mix_one.init(mixing_in); - return *this; + if(this->flag_del_mixing) + delete this->mixing; + this->mixing = mixing_in; + this->flag_del_mixing = false; } -Mix_DMk_2D &Mix_DMk_2D::set_mixing_beta(const double mixing_beta) +template +void Mix_DMk_2D::set_mixing_plain(const double& mixing_beta) { - ModuleBase::TITLE("Mix_DMk_2D","set_mixing_beta"); - if(this->gamma_only) - for (Mix_Matrix>& mix_one : this->mix_DMk_gamma) - mix_one.mixing_beta = mixing_beta; - else - for (Mix_Matrix>>& mix_one : this->mix_DMk_k) - mix_one.mixing_beta = mixing_beta; - return *this; + if(this->flag_del_mixing) + delete this->mixing; + this->mixing = new Base_Mixing::Plain_Mixing(mixing_beta); + this->flag_del_mixing = true; } -void Mix_DMk_2D::mix(const std::vector>& dm, const bool flag_restart) +template +void Mix_DMk_2D::mix(const std::vector>& dm, const bool flag_restart) { - ModuleBase::TITLE("Mix_DMk_2D","mix"); - assert(this->mix_DMk_gamma.size() == dm.size()); - for(int ik=0; ikmix_DMk_gamma[ik].mix(dm[ik], flag_restart); + ModuleBase::TITLE("Mix_DMk_2D", "mix"); + if (flag_restart) + { this->restart_all(dm); } + else + { this->mix_all(dm); } } -void Mix_DMk_2D::mix(const std::vector>>& dm, const bool flag_restart) + +template +std::vector*> Mix_DMk_2D::get_DMk_out() const { - ModuleBase::TITLE("Mix_DMk_2D","mix"); - assert(this->mix_DMk_k.size() == dm.size()); - for(int ik=0; ikmix_DMk_k[ik].mix(dm[ik], flag_restart); + std::vector*> DMk_out(this->mix_DMk.size()); + for (int ik = 0; ik < this->mix_DMk.size(); ++ik) + { DMk_out[ik] = &this->mix_DMk[ik].data_out; } + return DMk_out; } -std::vector*> Mix_DMk_2D::get_DMk_gamma_out() const +template +void Mix_DMk_2D::restart_all(const std::vector>& data_in) { - std::vector*> DMk_out(this->mix_DMk_gamma.size()); - for(int ik=0; ikmix_DMk_gamma.size(); ++ik) - DMk_out[ik] = &this->mix_DMk_gamma[ik].get_data_out(); - return DMk_out; + assert(this->mix_DMk.size() == data_in.size()); + assert(this->mixing != nullptr); + for (int ik = 0; ik < data_in.size(); ++ik) + { + this->mix_DMk[ik].data_out = data_in[ik]; + this->mixing->init_mixing_data(this->mix_DMk[ik].mixing_data, data_in[ik].size(), sizeof(Tdata)); + } } -std::vector>*> Mix_DMk_2D::get_DMk_k_out() const + +template +void Mix_DMk_2D::mix_all(const std::vector>& data_in) { - std::vector>*> DMk_out(this->mix_DMk_k.size()); - for(int ik=0; ikmix_DMk_k.size(); ++ik) - DMk_out[ik] = &this->mix_DMk_k[ik].get_data_out(); - return DMk_out; -} \ No newline at end of file + assert(this->mix_DMk.size() == data_in.size()); + assert(this->mixing != nullptr); + for (int ik = 0; ik < data_in.size(); ++ik) + { + this->mixing->push_data(this->mix_DMk[ik].mixing_data, this->mix_DMk[ik].data_out.data(), data_in[ik].data(), nullptr, false); + this->mixing->mix_data(this->mix_DMk[ik].mixing_data, this->mix_DMk[ik].data_out.data()); + } +} + +template class Mix_DMk_2D; +template class Mix_DMk_2D>; diff --git a/source/source_lcao/module_ri/Mix_DMk_2D.h b/source/source_lcao/module_ri/Mix_DMk_2D.h index aefe1b297bd..a37348415e2 100644 --- a/source/source_lcao/module_ri/Mix_DMk_2D.h +++ b/source/source_lcao/module_ri/Mix_DMk_2D.h @@ -6,66 +6,62 @@ #ifndef MIX_DMK_2D_H #define MIX_DMK_2D_H -#include "Mix_Matrix.h" #include "source_base/module_mixing/mixing.h" -#include "source_base/matrix.h" +#include #include +template class Mix_DMk_2D { public: + ~Mix_DMk_2D(); + /** - * @brief Sets the number of k-points and gamma_only flag. + * @brief Sets the number of k-points. * @param nks Number of k-points. - * @param gamma_only_in Flag indicating if only gamma point is mixed. - * @return Reference to the current object. */ - Mix_DMk_2D &set_nks(const int nks, const bool gamma_only_in); + void set_nks(const int nks); /** * @brief Sets the mixing mode. * @param Mixing Mixing pointer. - * @return Reference to the current object. */ - Mix_DMk_2D &set_mixing(Base_Mixing::Mixing* mixing_in); + void set_mixing(Base_Mixing::Mixing* mixing_in); /** - * @brief Sets the mixing beta. - * @param mixing_beta Mixing beta. - * @return Reference to the current object. - */ - Mix_DMk_2D &set_mixing_beta(const double mixing_beta); - - /** - * @brief Mixes the double density matrix. - * @param dm Double Density matrix. - * @param flag_restart Flag indicating whether restart mixing. + * @brief Sets Base_Mixing::Plain_Mixing. + * @param mixing_beta mixing beta for plain mixing. */ - void mix(const std::vector>& dm, const bool flag_restart); + void set_mixing_plain(const double& mixing_beta); /** - * @brief Mixes the complex density matrix. - * @param dm Complex density matrix. + * @brief Mixes the density matrix. + * @param dm Density matrix. * @param flag_restart Flag indicating whether restart mixing. */ - void mix(const std::vector>>& dm, const bool flag_restart); + void mix(const std::vector>& dm, const bool flag_restart); /** - * @brief Returns the double density matrix. - * @return Double density matrices for each k-points. + * @brief Returns the density matrix. + * @return Density matrices for each k-points. */ - std::vector*> get_DMk_gamma_out() const; - /** - * @brief Returns the complex density matrix. - * @return Complex density matrices for each k-points. - */ - std::vector>*> get_DMk_k_out() const; + std::vector*> get_DMk_out() const; private: - std::vector>> mix_DMk_gamma; - std::vector>>> mix_DMk_k; - bool gamma_only; + struct DMk_Mix_Data + { + std::vector data_out; + Base_Mixing::Mixing_Data mixing_data; + }; + + void restart_all(const std::vector>& data_in); + + void mix_all(const std::vector>& data_in); + + std::vector mix_DMk; + Base_Mixing::Mixing* mixing = nullptr; + bool flag_del_mixing = false; }; -#endif \ No newline at end of file +#endif diff --git a/source/source_lcao/module_ri/Mix_Matrix.cpp b/source/source_lcao/module_ri/Mix_Matrix.cpp deleted file mode 100644 index bb7ad3c8404..00000000000 --- a/source/source_lcao/module_ri/Mix_Matrix.cpp +++ /dev/null @@ -1,164 +0,0 @@ -//======================= -// AUTHOR : Peize Lin -// DATE : 2023-05-09 -//======================= - -#ifndef MIX_DATA_HPP -#define MIX_DATA_HPP - -#include "Mix_Matrix.h" -#include "source_base/matrix.h" -#include "source_base/complexmatrix.h" -#include "source_base/tool_title.h" - -template<> -void Mix_Matrix::mix(const ModuleBase::matrix& data_in, const bool flag_restart) -{ - if(separate_loop) - { - this->mixing = new Base_Mixing::Plain_Mixing(this->mixing_beta); - } - - if(flag_restart) - { - this->data_out = data_in; - this->mixing->init_mixing_data(this->matrix_data, data_in.nc*data_in.nr, sizeof(*data_in.c)); - } - else - { - this->mixing->push_data(this->matrix_data, data_out.c, data_in.c, nullptr, false); - this->mixing->mix_data(this->matrix_data, data_out.c); - } - - if(separate_loop) - { - delete this->mixing; - this->mixing = nullptr; - } -} -template<> -void Mix_Matrix::mix(const ModuleBase::ComplexMatrix& data_in, const bool flag_restart) -{ - if (separate_loop) - { - this->mixing = new Base_Mixing::Plain_Mixing(this->mixing_beta); - } - - if (flag_restart) - { - this->data_out = data_in; - this->mixing->init_mixing_data(this->matrix_data, data_in.nc * data_in.nr, sizeof(*data_in.c)); - } - else - { - this->mixing->push_data(this->matrix_data, data_out.c, data_in.c, nullptr, false); - this->mixing->mix_data(this->matrix_data, data_out.c); - } - - if (separate_loop) - { - delete this->mixing; - this->mixing = nullptr; - } -} - -template<> -void Mix_Matrix>::mix(const std::vector& data_in, const bool flag_restart) -{ - if (separate_loop) - { - this->mixing = new Base_Mixing::Plain_Mixing(this->mixing_beta); - } - - if (flag_restart) - { - this->data_out = data_in; - this->mixing->init_mixing_data(this->matrix_data, data_in.size(), sizeof(*data_in.data())); - } - else - { - this->mixing->push_data(this->matrix_data, data_out.data(), data_in.data(), nullptr, false); - this->mixing->mix_data(this->matrix_data, data_out.data()); - } - - if (separate_loop) - { - delete this->mixing; - this->mixing = nullptr; - } -} -template<> -void Mix_Matrix>>::mix(const std::vector>& data_in, const bool flag_restart) -{ - if (separate_loop) - { - this->mixing = new Base_Mixing::Plain_Mixing(this->mixing_beta); - } - - if (flag_restart) - { - this->data_out = data_in; - this->mixing->init_mixing_data(this->matrix_data, data_in.size(), sizeof(*data_in.data())); - } - else - { - this->mixing->push_data(this->matrix_data, data_out.data(), data_in.data(), nullptr, false); - this->mixing->mix_data(this->matrix_data, data_out.data()); - } - - if (separate_loop) - { - delete this->mixing; - this->mixing = nullptr; - } -} - -/* -// mix for ct::Tensor if useful in the future. comment out for UT-coverage -template<> -void Mix_Matrix::mix(const ct::Tensor& data_in, const bool flag_restart) -{ - ModuleBase::TITLE("Mix_Matrix", "mix"); - if (separate_loop) - { - this->mixing = new Base_Mixing::Plain_Mixing(this->mixing_beta); - } - - if (flag_restart) - { - this->data_out = data_in; - this->mixing->init_mixing_data(this->matrix_data, data_in.NumElements(), sizeof(data_in.SizeOfType(data_in.data_type()))); - } - else - { - switch (data_in.data_type()) - { - case ct::DataType::DT_DOUBLE: - this->mixing->push_data(this->matrix_data, data_out.data(), data_in.data(), nullptr, false); - this->mixing->mix_data(this->matrix_data, data_out.data()); - break; - case ct::DataType::DT_COMPLEX_DOUBLE: - this->mixing->push_data(this->matrix_data, data_out.data>(), data_in.data>(), nullptr, false); - this->mixing->mix_data(this->matrix_data, data_out.data>()); - break; - default: - throw std::invalid_argument("Mix_Matrix: data type not supported"); - break; - } - } - - if (separate_loop) - { - delete this->mixing; - this->mixing = nullptr; - } -} -*/ - -template class Mix_Matrix; -template class Mix_Matrix; -template class Mix_Matrix>; -template class Mix_Matrix>>; -// template class Mix_Matrix; - -#endif \ No newline at end of file diff --git a/source/source_lcao/module_ri/Mix_Matrix.h b/source/source_lcao/module_ri/Mix_Matrix.h deleted file mode 100644 index fa11bbdb007..00000000000 --- a/source/source_lcao/module_ri/Mix_Matrix.h +++ /dev/null @@ -1,53 +0,0 @@ -//======================= -// AUTHOR : Peize Lin -// DATE : 2023-05-09 -//======================= - -#ifndef MIX_MATRIX_H -#define MIX_MATRIX_H - -#include - -#include "source_base/module_mixing/mixing.h" -#include "source_base/module_mixing/plain_mixing.h" -template -class Mix_Matrix -{ - public: - /** - * @brief init Mix_Matrix with Mixing pointer - * - * @param mixing_in Mixing pointer - */ - void init(Base_Mixing::Mixing* mixing_in) - { - this->mixing = mixing_in; - if (this->mixing == nullptr) - this->separate_loop = true; - } - /** - * @brief Get the data out object - * - * @return const Tdata& - */ - const Tdata& get_data_out() const - { - return this->data_out; - } - - /** - * @brief Mixes the input data according to the set mixing mode. - * @param data_in Input data to be mixed. - * @param flag_restart Flag indicating whether restart mixing. - */ - void mix(const Tdata& data_in, const bool flag_restart); - - double mixing_beta = 1.0; - - private: - Tdata data_out; - Base_Mixing::Mixing* mixing = nullptr; - Base_Mixing::Mixing_Data matrix_data; - bool separate_loop = false; -}; -#endif \ No newline at end of file diff --git a/source/source_lcao/module_ri/RPA_LRI.hpp b/source/source_lcao/module_ri/RPA_LRI.hpp index 075d9c2fccf..0323df970d0 100644 --- a/source/source_lcao/module_ri/RPA_LRI.hpp +++ b/source/source_lcao/module_ri/RPA_LRI.hpp @@ -113,14 +113,14 @@ void RPA_LRI::cal_postSCF_exx(const elecstate::DensityMatrix this->p_kv = &kv; this->orb_cutoff_ = orb.cutoffs(); - Mix_DMk_2D mix_DMk_2D; + Mix_DMk_2D mix_DMk_2D; bool exx_spacegroup_symmetry = (PARAM.inp.nspin < 4 && ModuleSymmetry::Symmetry::symm_flag == 1); if (exx_spacegroup_symmetry) - {mix_DMk_2D.set_nks(kv.get_nkstot_full() * (PARAM.inp.nspin == 2 ? 2 : 1), PARAM.globalv.gamma_only_local);} + {mix_DMk_2D.set_nks(kv.get_nkstot_full() * (PARAM.inp.nspin == 2 ? 2 : 1));} else - {mix_DMk_2D.set_nks(kv.get_nks(), PARAM.globalv.gamma_only_local);} + {mix_DMk_2D.set_nks(kv.get_nks());} - mix_DMk_2D.set_mixing(nullptr); + mix_DMk_2D.set_mixing_plain(1.0); ModuleSymmetry::Symmetry_rotation symrot; if (exx_spacegroup_symmetry) { @@ -138,9 +138,13 @@ void RPA_LRI::cal_postSCF_exx(const elecstate::DensityMatrix else { mix_DMk_2D.mix(dm.get_DMK_vector(), true); } const std::vector>>> - Ds = PARAM.globalv.gamma_only_local - ? RI_2D_Comm::split_m2D_ktoR(ucell,kv, mix_DMk_2D.get_DMk_gamma_out(), *dm.get_paraV_pointer(), PARAM.inp.nspin) - : RI_2D_Comm::split_m2D_ktoR(ucell,kv, mix_DMk_2D.get_DMk_k_out(), *dm.get_paraV_pointer(), PARAM.inp.nspin, exx_spacegroup_symmetry); + Ds = RI_2D_Comm::split_m2D_ktoR( + ucell, + kv, + mix_DMk_2D.get_DMk_out(), + *dm.get_paraV_pointer(), + PARAM.inp.nspin, + exx_spacegroup_symmetry); // set parameters for bare Coulomb potential GlobalC::exx_info.info_global.ccp_type = Conv_Coulomb_Pot_K::Ccp_Type::Hf; // not used now, Hf/Ccp -> singularity_correction, see conv_coulomb_pot_k.cpp diff --git a/source/source_lcao/module_ri/test/CMakeLists.txt b/source/source_lcao/module_ri/test/CMakeLists.txt index b1a0f3f90e6..0565ed6a73c 100644 --- a/source/source_lcao/module_ri/test/CMakeLists.txt +++ b/source/source_lcao/module_ri/test/CMakeLists.txt @@ -4,7 +4,7 @@ remove_definitions(-D__ROCM) AddTest( TARGET MODULE_RI_dm_mixing_test LIBS parameter base ${math_libs} device - SOURCES dm_mixing_test.cpp ../Mix_DMk_2D.cpp ../Mix_Matrix.cpp + SOURCES dm_mixing_test.cpp ../Mix_DMk_2D.cpp ) AddTest( TARGET MODULE_RI_ri_cv_io_test @@ -16,4 +16,4 @@ AddTest( LIBS parameter SOURCES abfs-vector3_order_test.cpp ) -install(DIRECTORY support DESTINATION ${CMAKE_CURRENT_BINARY_DIR}) \ No newline at end of file +install(DIRECTORY support DESTINATION ${CMAKE_CURRENT_BINARY_DIR}) diff --git a/source/source_lcao/module_ri/test/dm_mixing_test.cpp b/source/source_lcao/module_ri/test/dm_mixing_test.cpp index 3d8f3308a37..e63dc95afa4 100644 --- a/source/source_lcao/module_ri/test/dm_mixing_test.cpp +++ b/source/source_lcao/module_ri/test/dm_mixing_test.cpp @@ -1,9 +1,7 @@ #include "gmock/gmock.h" #include "gtest/gtest.h" -#include "source_base/complexmatrix.h" #include "source_base/module_mixing/broyden_mixing.h" #include "source_lcao/module_ri/Mix_DMk_2D.h" -#include "source_lcao/module_ri/Mix_Matrix.h" /************************************************ * unit test of charge_mixing.cpp & Mix_DMk_2D.cpp @@ -11,8 +9,6 @@ /** * - Tested Functions: - * - Mix_Matrix::mix: - * mix the matrix data according to the set mixing mode. * - Mix_DMk_2D::mix: * mix the density matrix data according to the set mixing mode. * @@ -24,30 +20,6 @@ class DM_Mixing_Test : public ::testing::Test DM_Mixing_Test() { mixing = new Base_Mixing::Broyden_Mixing(ndim, mixing_beta); - mix_data = std::vector(2); - mix_complexdata = std::vector(3); - mix_data[0].create(nr, nc); - mix_data[1].create(nr, nc); - for (int i = 0; i < nr; ++i) - { - for (int j = 0; j < nc; ++j) - { - mix_data[0](i, j) = i * nc + j; - mix_data[1](i, j) = i * nc + j + 0.2; - } - } - mix_complexdata[0].create(nr, nc); - mix_complexdata[1].create(nr, nc); - mix_complexdata[2].create(nr, nc); - for (int i = 0; i < nr; ++i) - { - for (int j = 0; j < nc; ++j) - { - mix_complexdata[0](i, j) = std::complex{double(i), double(j)}; - mix_complexdata[1](i, j) = std::complex{double(i), double(j) + 0.2}; - mix_complexdata[2](i, j) = std::complex{double(i) + 0.8, double(j)}; - } - } mix_data_vector = std::vector>(2); mix_complexdata_vector = std::vector>>(3); mix_data_vector[0].resize(nr * nc); @@ -84,99 +56,16 @@ class DM_Mixing_Test : public ::testing::Test const double mixing_beta = 0.3; protected: - std::vector mix_data; - std::vector mix_complexdata; std::vector> mix_data_vector; std::vector>> mix_complexdata_vector; }; -TEST_F(DM_Mixing_Test, Mix_Matrix) -{ - // Separete loop - // matrix - Mix_Matrix mix_matrix; - mix_matrix.init(nullptr); - mix_matrix.mixing_beta = mixing_beta; - for (int istep = 0; istep < 2; ++istep) - { - mix_matrix.mix(mix_data[istep], (istep == 0)); - } - - ModuleBase::matrix data_out = mix_matrix.get_data_out(); - for (int i = 0; i < nr; ++i) - { - for (int j = 0; j < nc; ++j) - { - EXPECT_DOUBLE_EQ(data_out(i, j), (1 - mixing_beta) * mix_data[0](i, j) + mixing_beta * mix_data[1](i, j)); - } - } - // vector - Mix_Matrix> mix_vector; - mix_vector.init(nullptr); - mix_vector.mixing_beta = mixing_beta; - for (int istep = 0; istep < 2; ++istep) - { - mix_vector.mix(mix_data_vector[istep], (istep == 0)); - } - - std::vector data_out_vector = mix_vector.get_data_out(); - for (int i = 0; i < nr; ++i) - { - for (int j = 0; j < nc; ++j) - { - EXPECT_DOUBLE_EQ(data_out_vector[i * nc + j], (1 - mixing_beta) * mix_data_vector[0][i * nc + j] + mixing_beta * mix_data_vector[1][i * nc + j]); - } - } - //vector > - Mix_Matrix>> mix_vector_complex; - mix_vector_complex.init(nullptr); - mix_vector_complex.mixing_beta = mixing_beta; - for (int istep = 0; istep < 2; ++istep) - { - mix_vector_complex.mix(mix_complexdata_vector[istep], (istep == 0)); - } - - std::vector> data_out_vector_complex = mix_vector_complex.get_data_out(); - for (int i = 0; i < nr; ++i) - { - for (int j = 0; j < nc; ++j) - { - EXPECT_DOUBLE_EQ(data_out_vector_complex[i * nc + j].real(), ((1 - mixing_beta) * mix_complexdata_vector[0][i * nc + j] + mixing_beta * mix_complexdata_vector[1][i * nc + j]).real()); - EXPECT_DOUBLE_EQ(data_out_vector_complex[i * nc + j].imag(), ((1 - mixing_beta) * mix_complexdata_vector[0][i * nc + j] + mixing_beta * mix_complexdata_vector[1][i * nc + j]).imag()); - } - } - - // Broyden mix - Mix_Matrix mix_complexmatrix; - mix_complexmatrix.init(mixing); - mixing->coef = { 1.1, -0.1 }; - for (int istep = 0; istep < 3; ++istep) - { - mix_complexmatrix.mix(mix_complexdata[istep], (istep == 0)); - } - ModuleBase::ComplexMatrix com_data_out = mix_complexmatrix.get_data_out(); - for (int i = 0; i < nr; ++i) - { - for (int j = 0; j < nc; ++j) - { - std::complex first_step_result - = (1 - mixing_beta) * mix_complexdata[0](i, j) + mixing_beta * mix_complexdata[1](i, j); - std::complex second_step_result - = (1 - mixing_beta) * first_step_result + mixing_beta * mix_complexdata[2](i, j); - std::complex ref = second_step_result * mixing->coef[1] + first_step_result * mixing->coef[0]; - EXPECT_DOUBLE_EQ(com_data_out(i, j).real(), ref.real()); - EXPECT_DOUBLE_EQ(com_data_out(i, j).imag(), ref.imag()); - } - } -} - TEST_F(DM_Mixing_Test, Mix_DMk_2D) { //Gamma only - Mix_DMk_2D mix_dmk_gamma; - mix_dmk_gamma.set_nks(1, true); - mix_dmk_gamma.set_mixing(nullptr); - mix_dmk_gamma.set_mixing_beta(mixing_beta); + Mix_DMk_2D mix_dmk_gamma; + mix_dmk_gamma.set_nks(1); + mix_dmk_gamma.set_mixing_plain(1.0); std::vector>> dm_gamma(2); dm_gamma[0] = std::vector>(1); dm_gamma[0][0] = mix_data_vector[0]; @@ -186,20 +75,19 @@ TEST_F(DM_Mixing_Test, Mix_DMk_2D) { mix_dmk_gamma.mix(dm_gamma[istep], (istep == 0)); } - std::vector*> dm_gamma_out = mix_dmk_gamma.get_DMk_gamma_out(); + std::vector*> dm_gamma_out = mix_dmk_gamma.get_DMk_out(); for (int i = 0; i < nr; ++i) { for (int j = 0; j < nc; ++j) { - EXPECT_DOUBLE_EQ(dm_gamma_out[0][0][i * nc + j], (1 - mixing_beta) * mix_data_vector[0][i * nc + j] + mixing_beta * mix_data_vector[1][i * nc + j]); + EXPECT_DOUBLE_EQ(dm_gamma_out[0][0][i * nc + j], mix_data_vector[1][i * nc + j]); } } // not Gamma only - Mix_DMk_2D mix_dmk; - mix_dmk.set_nks(1, false); - mix_dmk.set_mixing(nullptr); - mix_dmk.set_mixing_beta(mixing_beta); + Mix_DMk_2D> mix_dmk; + mix_dmk.set_nks(1); + mix_dmk.set_mixing_plain(1.0); std::vector>>> dm(2); dm[0] = std::vector>>(1); dm[0][0] = mix_complexdata_vector[0]; @@ -209,15 +97,41 @@ TEST_F(DM_Mixing_Test, Mix_DMk_2D) { mix_dmk.mix(dm[istep], (istep == 0)); } - std::vector>*> dm_out = mix_dmk.get_DMk_k_out(); + std::vector>*> dm_out = mix_dmk.get_DMk_out(); + for (int i = 0; i < nr; ++i) + { + for (int j = 0; j < nc; ++j) + { + EXPECT_DOUBLE_EQ(dm_out[0][0][i * nc + j].real(), mix_complexdata_vector[1][i * nc + j].real()); + EXPECT_DOUBLE_EQ(dm_out[0][0][i * nc + j].imag(), mix_complexdata_vector[1][i * nc + j].imag()); + } + } + + // Shared Broyden mix + Mix_DMk_2D> mix_dmk_broyden; + mix_dmk_broyden.set_nks(1); + mix_dmk_broyden.set_mixing(mixing); + mixing->coef = { 1.1, -0.1 }; + std::vector>>> dm_broyden(3); + for (int istep = 0; istep < 3; ++istep) + { + dm_broyden[istep] = std::vector>>(1); + dm_broyden[istep][0] = mix_complexdata_vector[istep]; + mix_dmk_broyden.mix(dm_broyden[istep], (istep == 0)); + } + std::vector>*> dm_broyden_out = mix_dmk_broyden.get_DMk_out(); for (int i = 0; i < nr; ++i) { for (int j = 0; j < nc; ++j) { std::complex first_step_result - = (1 - mixing_beta) * mix_complexdata_vector[0][i * nc + j] + mixing_beta * mix_complexdata_vector[1][i * nc + j]; - EXPECT_DOUBLE_EQ(dm_out[0][0][i * nc + j].real(), first_step_result.real()); - EXPECT_DOUBLE_EQ(dm_out[0][0][i * nc + j].imag(), first_step_result.imag()); + = (1 - mixing_beta) * mix_complexdata_vector[0][i * nc + j] + + mixing_beta * mix_complexdata_vector[1][i * nc + j]; + std::complex second_step_result + = (1 - mixing_beta) * first_step_result + mixing_beta * mix_complexdata_vector[2][i * nc + j]; + std::complex ref = second_step_result * mixing->coef[1] + first_step_result * mixing->coef[0]; + EXPECT_DOUBLE_EQ(dm_broyden_out[0][0][i * nc + j].real(), ref.real()); + EXPECT_DOUBLE_EQ(dm_broyden_out[0][0][i * nc + j].imag(), ref.imag()); } } -} \ No newline at end of file +}