From 9319fe553b1ccc630c06ec654510a1660ef3cc08 Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Mon, 5 Feb 2024 11:04:05 +0800 Subject: [PATCH 1/3] define Charge_Mixing::init_mixing() --- .../module_charge/charge_mixing.cpp | 29 +++++++++++++++++-- .../module_charge/charge_mixing.h | 6 ++++ 2 files changed, 32 insertions(+), 3 deletions(-) diff --git a/source/module_elecstate/module_charge/charge_mixing.cpp b/source/module_elecstate/module_charge/charge_mixing.cpp index 6a489b05b3..03f8523293 100755 --- a/source/module_elecstate/module_charge/charge_mixing.cpp +++ b/source/module_elecstate/module_charge/charge_mixing.cpp @@ -70,6 +70,14 @@ void Charge_Mixing::set_mixing(const std::string& mixing_mode_in, GlobalV::ofs_running<<"mixing_ndim: "<< this->mixing_ndim <mixing_mode == "broyden") { delete this->mixing; @@ -97,6 +105,8 @@ void Charge_Mixing::set_mixing(const std::string& mixing_mode_in, this->mixing_highf = new Base_Mixing::Plain_Mixing(this->mixing_beta); } + // allocate memory for mixing data, if exists, free it first and then allocate new memory + // initailize rho_mdata if (GlobalV::SCF_THR_TYPE == 1) { if (GlobalV::NSPIN == 4 && GlobalV::MIXING_ANGLE > 0 ) @@ -124,13 +134,26 @@ void Charge_Mixing::set_mixing(const std::string& mixing_mode_in, } } + // initailize tau_mdata + if ((XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5) && mixing_tau) + { + if (GlobalV::SCF_THR_TYPE == 1) + { + this->mixing->init_mixing_data(this->tau_mdata, + this->rhopw->npw * GlobalV::NSPIN, + sizeof(std::complex)); + } + else + { + this->mixing->init_mixing_data(this->tau_mdata, this->rhopw->nrxx * GlobalV::NSPIN, sizeof(double)); + } + } + + // initailize nhat_mdata #ifdef USE_PAW if(GlobalV::use_paw) this->mixing->init_mixing_data(this->nhat_mdata, this->rhopw->nrxx * GlobalV::NSPIN, sizeof(double)); #endif - // Note: we can not init tau_mdata here temporarily, since set_xc_type() is after it. - // you can find initalize tau_mdata in mix_reset(); - // this->mixing->init_mixing_data(this->tau_mdata, this->rhopw->nrxx * GlobalV::NSPIN, sizeof(double)); return; } diff --git a/source/module_elecstate/module_charge/charge_mixing.h b/source/module_elecstate/module_charge/charge_mixing.h index 90af6ae685..269ead2a1b 100644 --- a/source/module_elecstate/module_charge/charge_mixing.h +++ b/source/module_elecstate/module_charge/charge_mixing.h @@ -66,6 +66,12 @@ class Charge_Mixing const double& mixing_angle_in, const bool& mixing_dmr_in); + /** + * @brief initialize mixing, including constructing mixing and allocating memory for mixing data + * @brief this function should be called at eachiterinit() + */ + void init_mixing(); + /** * @brief allocate memory of dmr_mdata * @param nnr size of real-space density matrix From 250709027c00ddc5dd93af5c4bbb4df6bc102bfe Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Mon, 5 Feb 2024 11:10:59 +0800 Subject: [PATCH 2/3] use init_mixing() to replace set_mixing() in eachiterinit() --- source/module_esolver/esolver_ks_lcao.cpp | 25 ++++++----------------- source/module_esolver/esolver_ks_pw.cpp | 15 +------------- 2 files changed, 7 insertions(+), 33 deletions(-) diff --git a/source/module_esolver/esolver_ks_lcao.cpp b/source/module_esolver/esolver_ks_lcao.cpp index ea7b3996ea..fde972f1b1 100644 --- a/source/module_esolver/esolver_ks_lcao.cpp +++ b/source/module_esolver/esolver_ks_lcao.cpp @@ -493,28 +493,15 @@ namespace ModuleESolver { if (iter == 1 || iter == GlobalV::MIXING_RESTART) { - if (iter == GlobalV::MIXING_RESTART) // delete mixing and re-construct it to restart + this->p_chgmix->init_mixing(); + if (iter == GlobalV::MIXING_RESTART && GlobalV::MIXING_DMR) // for mixing_dmr { - this->p_chgmix->set_mixing(GlobalV::MIXING_MODE, - GlobalV::MIXING_BETA, - GlobalV::MIXING_NDIM, - GlobalV::MIXING_GG0, - GlobalV::MIXING_TAU, - GlobalV::MIXING_BETA_MAG, - GlobalV::MIXING_GG0_MAG, - GlobalV::MIXING_GG0_MIN, - GlobalV::MIXING_ANGLE, - GlobalV::MIXING_DMR); // allocate memory for dmr_mdata - if (GlobalV::MIXING_DMR) - { - const elecstate::DensityMatrix* dm - = dynamic_cast*>(this->pelec)->get_DM(); - int nnr_tmp = dm->get_DMR_pointer(1)->get_nnr(); - this->p_chgmix->allocate_mixing_dmr(nnr_tmp); - } + const elecstate::DensityMatrix* dm + = dynamic_cast*>(this->pelec)->get_DM(); + int nnr_tmp = dm->get_DMR_pointer(1)->get_nnr(); + this->p_chgmix->allocate_mixing_dmr(nnr_tmp); } - this->p_chgmix->mix_reset(); } // mohan update 2012-06-05 diff --git a/source/module_esolver/esolver_ks_pw.cpp b/source/module_esolver/esolver_ks_pw.cpp index 9861fabb41..89511f410b 100644 --- a/source/module_esolver/esolver_ks_pw.cpp +++ b/source/module_esolver/esolver_ks_pw.cpp @@ -494,20 +494,7 @@ void ESolver_KS_PW::eachiterinit(const int istep, const int iter) { if (iter == 1 || iter == GlobalV::MIXING_RESTART) { - if (iter == GlobalV::MIXING_RESTART) // delete mixing and re-construct it to restart - { - this->p_chgmix->set_mixing(GlobalV::MIXING_MODE, - GlobalV::MIXING_BETA, - GlobalV::MIXING_NDIM, - GlobalV::MIXING_GG0, - GlobalV::MIXING_TAU, - GlobalV::MIXING_BETA_MAG, - GlobalV::MIXING_GG0_MAG, - GlobalV::MIXING_GG0_MIN, - GlobalV::MIXING_ANGLE, - GlobalV::MIXING_DMR); - } - this->p_chgmix->mix_reset(); + this->p_chgmix->init_mixing(); } // mohan move harris functional to here, 2012-06-05 // use 'rho(in)' and 'v_h and v_xc'(in) From 715cf9b12408f985a6f6b2351c8d2b7015319df2 Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Mon, 5 Feb 2024 12:48:58 +0800 Subject: [PATCH 3/3] polish charge_mixing_test to protect init_mixing() --- .../module_charge/charge_mixing.cpp | 24 ++-- .../module_charge/charge_mixing.h | 36 +++--- .../test/charge_mixing_test.cpp | 107 +++++++++++------- source/module_esolver/esolver_ks_lcao.cpp | 16 +-- 4 files changed, 106 insertions(+), 77 deletions(-) diff --git a/source/module_elecstate/module_charge/charge_mixing.cpp b/source/module_elecstate/module_charge/charge_mixing.cpp index 03f8523293..505fac7772 100755 --- a/source/module_elecstate/module_charge/charge_mixing.cpp +++ b/source/module_elecstate/module_charge/charge_mixing.cpp @@ -52,6 +52,11 @@ void Charge_Mixing::set_mixing(const std::string& mixing_mode_in, ModuleBase::WARNING_QUIT("Charge_Mixing", "You'd better set mixing_beta_mag >= 0.0!"); } + if (!(this->mixing_mode == "plain" || this->mixing_mode == "broyden" || this->mixing_mode == "pulay")) + { + ModuleBase::WARNING_QUIT("Charge_Mixing", "This Mixing mode is not implemended yet,coming soon."); + } + // print into running.log GlobalV::ofs_running<<"\n----------- Double Check Mixing Parameters Begin ------------"<mixing_mode <mixing_mode == "broyden") { @@ -97,6 +105,7 @@ void Charge_Mixing::init_mixing() { ModuleBase::WARNING_QUIT("Charge_Mixing", "This Mixing mode is not implemended yet,coming soon."); } + if (GlobalV::double_grid) { // ONLY smooth part of charge density is mixed by specific mixing method @@ -133,7 +142,7 @@ void Charge_Mixing::init_mixing() this->mixing->init_mixing_data(this->rho_mdata, this->rhopw->nrxx * GlobalV::NSPIN, sizeof(double)); } } - + // initailize tau_mdata if ((XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5) && mixing_tau) { @@ -154,6 +163,8 @@ void Charge_Mixing::init_mixing() if(GlobalV::use_paw) this->mixing->init_mixing_data(this->nhat_mdata, this->rhopw->nrxx * GlobalV::NSPIN, sizeof(double)); #endif + ModuleBase::timer::tick("Charge_Mixing", "init_mixing"); + return; } @@ -956,16 +967,7 @@ void Charge_Mixing::mix_reset() // initailize tau_mdata if ((XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5) && mixing_tau) { - if (GlobalV::SCF_THR_TYPE == 1) - { - this->mixing->init_mixing_data(this->tau_mdata, - this->rhopw->npw * GlobalV::NSPIN, - sizeof(std::complex)); - } - else - { - this->mixing->init_mixing_data(this->tau_mdata, this->rhopw->nrxx * GlobalV::NSPIN, sizeof(double)); - } + this->tau_mdata.reset(); } // reset for paw #ifdef USE_PAW diff --git a/source/module_elecstate/module_charge/charge_mixing.h b/source/module_elecstate/module_charge/charge_mixing.h index 269ead2a1b..3c8cade931 100644 --- a/source/module_elecstate/module_charge/charge_mixing.h +++ b/source/module_elecstate/module_charge/charge_mixing.h @@ -24,24 +24,6 @@ class Charge_Mixing Charge_Mixing(); ~Charge_Mixing(); - /** - * @brief reset mixing - */ - void mix_reset(); - - /** - * @brief charge mixing - * @param chr pointer of Charge object - */ - void mix_rho(Charge* chr); - - /** - * @brief density matrix mixing, only for LCAO - * @param DM pointer of DensityMatrix object - */ - void mix_dmr(elecstate::DensityMatrix* DM); - void mix_dmr(elecstate::DensityMatrix, double>* DM); - /** * @brief Set all private mixing paramters * @param mixing_mode_in mixing mode: "plain", "broyden", "pulay" @@ -78,11 +60,29 @@ class Charge_Mixing */ void allocate_mixing_dmr(int nnr); + /** + * @brief charge mixing + * @param chr pointer of Charge object + */ + void mix_rho(Charge* chr); + + /** + * @brief density matrix mixing, only for LCAO + * @param DM pointer of DensityMatrix object + */ + void mix_dmr(elecstate::DensityMatrix* DM); + void mix_dmr(elecstate::DensityMatrix, double>* DM); + /** * @brief Get the drho * */ double get_drho(Charge* chr, const double nelec); + + /** + * @brief reset mixing, actually we only call init_mixing() to reset mixing instead of this function + */ + void mix_reset(); /** * @brief Set the smooth and dense grids diff --git a/source/module_elecstate/test/charge_mixing_test.cpp b/source/module_elecstate/test/charge_mixing_test.cpp index a01a288b01..52eaaeb62f 100644 --- a/source/module_elecstate/test/charge_mixing_test.cpp +++ b/source/module_elecstate/test/charge_mixing_test.cpp @@ -52,15 +52,13 @@ UnitCell ucell; * - Tested Functions: * - SetMixingTest: * Charge_Mixing::set_mixing(mixing_mode_in,mixing_beta_in,mixing_ndim_in,mixing_gg0_in,mixing_tau_in) + * Charge_Mixing::init_mixing() * Charge_Mixing::set_rhopw(rhopw_in) * Charge_Mixing::get_mixing_mode() * Charge_Mixing::get_mixing_beta() * Charge_Mixing::get_mixing_ndim() * Charge_Mixing::get_mixing_gg0() * - set the basic parameters of class charge_mixing - * - AutoSetTest: Charge_Mixing::auto_set(bandgap_in, ucell_) - * Charge_Mixing::need_auto_set() - * - auto-set the parameters of class charge_mixing * - KerkerScreenTest: Charge_Mixing::Kerker_screen_recip(drhog) * Charge_Mixing::Kerker_screen_real(drhog) * - screen drho with Kerker method @@ -120,8 +118,6 @@ TEST_F(ChargeMixingTest, SetMixingTest) GlobalV::MIXING_NDIM = 1; GlobalV::MIXING_GG0 = 1.0; - FUNC_TYPE = 1; - GlobalV::SCF_THR_TYPE = 1; CMtest.set_mixing(GlobalV::MIXING_MODE, GlobalV::MIXING_BETA, GlobalV::MIXING_NDIM, @@ -132,7 +128,6 @@ TEST_F(ChargeMixingTest, SetMixingTest) GlobalV::MIXING_GG0_MIN, GlobalV::MIXING_ANGLE, GlobalV::MIXING_DMR); - EXPECT_EQ(CMtest.rho_mdata.length, pw_basis.npw); EXPECT_EQ(CMtest.get_mixing_mode(), "broyden"); EXPECT_EQ(CMtest.get_mixing_beta(), 1.0); EXPECT_EQ(CMtest.get_mixing_ndim(), 1); @@ -144,23 +139,8 @@ TEST_F(ChargeMixingTest, SetMixingTest) EXPECT_EQ(CMtest.mixing_angle, -10.0); EXPECT_EQ(CMtest.mixing_dmr, false); - GlobalV::SCF_THR_TYPE = 2; - CMtest.set_mixing(GlobalV::MIXING_MODE, - GlobalV::MIXING_BETA, - GlobalV::MIXING_NDIM, - GlobalV::MIXING_GG0, - GlobalV::MIXING_TAU, - GlobalV::MIXING_BETA_MAG, - GlobalV::MIXING_GG0_MAG, - GlobalV::MIXING_GG0_MIN, - GlobalV::MIXING_ANGLE, - GlobalV::MIXING_DMR); - EXPECT_EQ(CMtest.rho_mdata.length, pw_basis.nrxx); - - FUNC_TYPE = 3; GlobalV::MIXING_TAU = true; GlobalV::MIXING_MODE = "plain"; - GlobalV::SCF_THR_TYPE = 1; CMtest.set_mixing(GlobalV::MIXING_MODE, GlobalV::MIXING_BETA, GlobalV::MIXING_NDIM, @@ -171,25 +151,9 @@ TEST_F(ChargeMixingTest, SetMixingTest) GlobalV::MIXING_GG0_MIN, GlobalV::MIXING_ANGLE, GlobalV::MIXING_DMR); - CMtest.mix_reset(); - EXPECT_EQ(CMtest.tau_mdata.length, pw_basis.npw); EXPECT_EQ(CMtest.mixing_mode, "plain"); EXPECT_EQ(CMtest.mixing_tau, true); - GlobalV::SCF_THR_TYPE = 2; - CMtest.set_mixing(GlobalV::MIXING_MODE, - GlobalV::MIXING_BETA, - GlobalV::MIXING_NDIM, - GlobalV::MIXING_GG0, - GlobalV::MIXING_TAU, - GlobalV::MIXING_BETA_MAG, - GlobalV::MIXING_GG0_MAG, - GlobalV::MIXING_GG0_MIN, - GlobalV::MIXING_ANGLE, - GlobalV::MIXING_DMR); - CMtest.mix_reset(); - EXPECT_EQ(CMtest.tau_mdata.length, pw_basis.nrxx); - GlobalV::MIXING_BETA = 1.1; std::string output; testing::internal::CaptureStdout(); @@ -242,6 +206,69 @@ TEST_F(ChargeMixingTest, SetMixingTest) EXPECT_THAT(output, testing::HasSubstr("This Mixing mode is not implemended yet,coming soon.")); } +TEST_F(ChargeMixingTest, InitMixingTest) +{ + omp_set_num_threads(1); + GlobalV::NSPIN = 1; + FUNC_TYPE = 1; + Charge_Mixing CMtest; + CMtest.set_rhopw(&pw_basis, &pw_basis); + + CMtest.set_mixing(GlobalV::MIXING_MODE, + GlobalV::MIXING_BETA, + GlobalV::MIXING_NDIM, + GlobalV::MIXING_GG0, + GlobalV::MIXING_TAU, + GlobalV::MIXING_BETA_MAG, + GlobalV::MIXING_GG0_MAG, + GlobalV::MIXING_GG0_MIN, + GlobalV::MIXING_ANGLE, + GlobalV::MIXING_DMR); + + GlobalV::SCF_THR_TYPE = 1; + CMtest.init_mixing(); + EXPECT_EQ(CMtest.rho_mdata.length, pw_basis.npw); + + GlobalV::SCF_THR_TYPE = 2; + CMtest.init_mixing(); + EXPECT_EQ(CMtest.rho_mdata.length, pw_basis.nrxx); + + GlobalV::NSPIN = 4; + CMtest.init_mixing(); + EXPECT_EQ(CMtest.rho_mdata.length, 4 * pw_basis.nrxx); + + GlobalV::NSPIN = 1; + GlobalV::MIXING_TAU = true; + CMtest.set_mixing(GlobalV::MIXING_MODE, + GlobalV::MIXING_BETA, + GlobalV::MIXING_NDIM, + GlobalV::MIXING_GG0, + GlobalV::MIXING_TAU, + GlobalV::MIXING_BETA_MAG, + GlobalV::MIXING_GG0_MAG, + GlobalV::MIXING_GG0_MIN, + GlobalV::MIXING_ANGLE, + GlobalV::MIXING_DMR); + FUNC_TYPE = 3; + CMtest.init_mixing(); + EXPECT_EQ(CMtest.tau_mdata.length, pw_basis.nrxx); + + GlobalV::NSPIN = 4; + GlobalV::MIXING_ANGLE = 1.0; + CMtest.set_mixing(GlobalV::MIXING_MODE, + GlobalV::MIXING_BETA, + GlobalV::MIXING_NDIM, + GlobalV::MIXING_GG0, + GlobalV::MIXING_TAU, + GlobalV::MIXING_BETA_MAG, + GlobalV::MIXING_GG0_MAG, + GlobalV::MIXING_GG0_MIN, + GlobalV::MIXING_ANGLE, + GlobalV::MIXING_DMR); + CMtest.init_mixing(); + EXPECT_EQ(CMtest.rho_mdata.length, 2 * pw_basis.nrxx); +} + TEST_F(ChargeMixingTest, KerkerScreenRecipTest) { Charge_Mixing CMtest; @@ -535,7 +562,7 @@ TEST_F(ChargeMixingTest, MixRhoTest) GlobalV::MIXING_GG0_MIN, GlobalV::MIXING_ANGLE, GlobalV::MIXING_DMR); - CMtest_recip.mix_reset(); + CMtest_recip.init_mixing(); for(int i = 0 ; i < nspin * npw; ++i) { charge._space_rhog[i] = recip_ref[i]; @@ -574,7 +601,7 @@ TEST_F(ChargeMixingTest, MixRhoTest) GlobalV::MIXING_GG0_MIN, GlobalV::MIXING_ANGLE, GlobalV::MIXING_DMR); - CMtest_real.mix_reset(); + CMtest_real.init_mixing(); for(int i = 0 ; i < nspin * nrxx; ++i) { charge._space_rho[i] = real_ref[i]; @@ -669,7 +696,7 @@ TEST_F(ChargeMixingTest, MixDoubleGridRhoTest) GlobalV::MIXING_GG0_MIN, GlobalV::MIXING_ANGLE, GlobalV::MIXING_DMR); - CMtest_recip.mix_reset(); + CMtest_recip.init_mixing(); for (int i = 0; i < nspin * npw; ++i) { charge._space_rhog[i] = recip_ref[i]; diff --git a/source/module_esolver/esolver_ks_lcao.cpp b/source/module_esolver/esolver_ks_lcao.cpp index fde972f1b1..c09d38b656 100644 --- a/source/module_esolver/esolver_ks_lcao.cpp +++ b/source/module_esolver/esolver_ks_lcao.cpp @@ -599,14 +599,14 @@ namespace ModuleESolver void ESolver_KS_LCAO::hamilt2density(int istep, int iter, double ethr) { // save input rho - this->pelec->charge->save_rho_before_sum_band(); - // save density matrix for mixing - if (GlobalV::MIXING_RESTART > 0 && GlobalV::MIXING_DMR && iter >= GlobalV::MIXING_RESTART) - { - elecstate::DensityMatrix* dm - = dynamic_cast*>(this->pelec)->get_DM(); - dm->save_DMR(); - } + this->pelec->charge->save_rho_before_sum_band(); + // save density matrix for mixing + if (GlobalV::MIXING_RESTART > 0 && GlobalV::MIXING_DMR && iter >= GlobalV::MIXING_RESTART) + { + elecstate::DensityMatrix* dm + = dynamic_cast*>(this->pelec)->get_DM(); + dm->save_DMR(); + } // using HSolverLCAO::solve() if (this->phsol != nullptr)