From 6730529ca008fa6192c3ee3a5fa7aba47a14070a Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Thu, 4 Jan 2024 11:20:47 +0800 Subject: [PATCH 1/8] move the code about recipical magnetic density from Charge to Charge_Mixing --- .../module_charge/charge_mixing.cpp | 42 ++++++++++++++++--- 1 file changed, 36 insertions(+), 6 deletions(-) diff --git a/source/module_elecstate/module_charge/charge_mixing.cpp b/source/module_elecstate/module_charge/charge_mixing.cpp index acece06a04..67372b4079 100755 --- a/source/module_elecstate/module_charge/charge_mixing.cpp +++ b/source/module_elecstate/module_charge/charge_mixing.cpp @@ -319,10 +319,30 @@ void Charge_Mixing::mix_rho_recip_new(Charge* chr) } else if (GlobalV::NSPIN == 2) { - chr->allocate_rhog_mag(); - chr->get_rhog_mag(); - rhog_in = chr->rhog_mag_save; - rhog_out = chr->rhog_mag; + int _ngmc = chr->ngmc; + // magnetic density + std::complex *rhog_mag = nullptr; + std::complex *rhog_mag_save = nullptr; + // allocate rhog_mag[is*ngmc] and rhog_mag_save[is*ngmc] + rhog_mag = new std::complex[_ngmc * GlobalV::NSPIN]; + rhog_mag_save = new std::complex[_ngmc * GlobalV::NSPIN]; + ModuleBase::GlobalFunc::ZEROS(rhog_mag, _ngmc * GlobalV::NSPIN); + ModuleBase::GlobalFunc::ZEROS(rhog_mag_save, _ngmc * GlobalV::NSPIN); + // get rhog_mag[is*ngmc] and rhog_mag_save[is*ngmc] + for (int ig = 0; ig < _ngmc; ig++) + { + rhog_mag[ig] = chr->rhog[0][ig] + chr->rhog[1][ig]; + rhog_mag_save[ig] = chr->rhog_save[0][ig] + chr->rhog_save[1][ig]; + } + for (int ig = 0; ig < _ngmc; ig++) + { + rhog_mag[ig + _ngmc] = chr->rhog[0][ig] - chr->rhog[1][ig]; + rhog_mag_save[ig + _ngmc] = chr->rhog_save[0][ig] - chr->rhog_save[1][ig]; + } + // + rhog_in = rhog_mag_save; + rhog_out = rhog_mag; + // const int npw = this->rhopw->npw; auto screen = std::bind(&Charge_Mixing::Kerker_screen_recip_new, this, std::placeholders::_1); auto twobeta_mix @@ -346,9 +366,19 @@ void Charge_Mixing::mix_rho_recip_new(Charge* chr) this->mixing->push_data(this->rho_mdata, rhog_in, rhog_out, screen, twobeta_mix, true); this->mixing->cal_coef(this->rho_mdata, inner_product_new); this->mixing->mix_data(this->rho_mdata, rhog_out); + // get rhog[is][nnrx] from rhog_mag[is*ngmc] + for (int is = 0; is < GlobalV::NSPIN; is++) + { + ModuleBase::GlobalFunc::ZEROS(chr->rhog[is], _ngmc); + } + for (int ig = 0; ig < _ngmc; ig++) + { + chr->rhog[0][ig] = 0.5 * (rhog_mag[ig] + rhog_mag[ig+_ngmc]); + chr->rhog[1][ig] = 0.5 * (rhog_mag[ig] - rhog_mag[ig+_ngmc]); + } // delete - chr->get_rhog_from_mag(); - chr->destroy_rhog_mag(); + delete[] rhog_mag; + delete[] rhog_mag_save; } else if (GlobalV::NSPIN == 4 && GlobalV::MIXING_ANGLE <= 0) { From e9a7a31225ce1192554a0abaa336c2cb87c7fbbc Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Thu, 4 Jan 2024 13:11:18 +0800 Subject: [PATCH 2/8] move the code about real magnetic density from Charge to Charge_Mixing --- .../module_charge/charge_mixing.cpp | 41 ++++++++++++++++--- 1 file changed, 35 insertions(+), 6 deletions(-) diff --git a/source/module_elecstate/module_charge/charge_mixing.cpp b/source/module_elecstate/module_charge/charge_mixing.cpp index 67372b4079..4726822ebf 100755 --- a/source/module_elecstate/module_charge/charge_mixing.cpp +++ b/source/module_elecstate/module_charge/charge_mixing.cpp @@ -366,7 +366,7 @@ void Charge_Mixing::mix_rho_recip_new(Charge* chr) this->mixing->push_data(this->rho_mdata, rhog_in, rhog_out, screen, twobeta_mix, true); this->mixing->cal_coef(this->rho_mdata, inner_product_new); this->mixing->mix_data(this->rho_mdata, rhog_out); - // get rhog[is][nnrx] from rhog_mag[is*ngmc] + // get rhog[is][ngmc] from rhog_mag[is*ngmc] for (int is = 0; is < GlobalV::NSPIN; is++) { ModuleBase::GlobalFunc::ZEROS(chr->rhog[is], _ngmc); @@ -572,11 +572,29 @@ void Charge_Mixing::mix_rho_real(Charge* chr) } else if (GlobalV::NSPIN == 2) { - chr->allocate_rho_mag(); - chr->get_rho_mag(); + // magnetic density + double *rho_mag = nullptr; + double *rho_mag_save = nullptr; + const int nrxx = this->rhopw->nrxx; + // allocate rho_mag[is*nnrx] and rho_mag_save[is*nnrx] + rho_mag = new double[nrxx * GlobalV::NSPIN]; + rho_mag_save = new double[nrxx * GlobalV::NSPIN]; + ModuleBase::GlobalFunc::ZEROS(rho_mag, nrxx * GlobalV::NSPIN); + ModuleBase::GlobalFunc::ZEROS(rho_mag_save, nrxx * GlobalV::NSPIN); + // get rho_mag[is*nnrx] and rho_mag_save[is*nnrx] + for (int ir = 0; ir < nrxx; ir++) + { + rho_mag[ir] = chr->rho[0][ir] + chr->rho[1][ir]; + rho_mag_save[ir] = chr->rho_save[0][ir] + chr->rho_save[1][ir]; + } + for (int ir = 0; ir < nrxx; ir++) + { + rho_mag[ir + nrxx] = chr->rho[0][ir] - chr->rho[1][ir]; + rho_mag_save[ir + nrxx] = chr->rho_save[0][ir] - chr->rho_save[1][ir]; + } + // rhor_in = chr->rho_mag_save; rhor_out = chr->rho_mag; - const int nrxx = this->rhopw->nrxx; auto screen = std::bind(&Charge_Mixing::Kerker_screen_real, this, std::placeholders::_1); auto twobeta_mix = [this, nrxx](double* out, const double* in, const double* sres) { @@ -601,9 +619,20 @@ void Charge_Mixing::mix_rho_real(Charge* chr) = std::bind(&Charge_Mixing::inner_product_real, this, std::placeholders::_1, std::placeholders::_2); this->mixing->cal_coef(this->rho_mdata, inner_product); this->mixing->mix_data(this->rho_mdata, rhor_out); + // get new rho[is][nrxx] from rho_mag[is*nrxx] + for (int is = 0; is < GlobalV::NSPIN; is++) + { + ModuleBase::GlobalFunc::ZEROS(chr->rho[is], nrxx); + //ModuleBase::GlobalFunc::ZEROS(rho_save[is], nrxx); + } + for (int ir = 0; ir < nrxx; ir++) + { + chr->rho[0][ir] = 0.5 * (rho_mag[ir] + rho_mag[ir+nrxx]); + chr->rho[1][ir] = 0.5 * (rho_mag[ir] - rho_mag[ir+nrxx]); + } // delete - chr->get_rho_from_mag(); - chr->destroy_rho_mag(); + delete[] rho_mag; + delete[] rho_mag_save; } else if (GlobalV::NSPIN == 4 && GlobalV::MIXING_ANGLE <= 0) { From 7603be30c0986b20094b2537d7076d9530604f22 Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Thu, 4 Jan 2024 13:15:01 +0800 Subject: [PATCH 3/8] replace _ngmc by npw --- .../module_charge/charge_mixing.cpp | 26 +++++++++---------- 1 file changed, 13 insertions(+), 13 deletions(-) diff --git a/source/module_elecstate/module_charge/charge_mixing.cpp b/source/module_elecstate/module_charge/charge_mixing.cpp index 4726822ebf..3cac601d32 100755 --- a/source/module_elecstate/module_charge/charge_mixing.cpp +++ b/source/module_elecstate/module_charge/charge_mixing.cpp @@ -319,25 +319,25 @@ void Charge_Mixing::mix_rho_recip_new(Charge* chr) } else if (GlobalV::NSPIN == 2) { - int _ngmc = chr->ngmc; // magnetic density std::complex *rhog_mag = nullptr; std::complex *rhog_mag_save = nullptr; + const int npw = this->rhopw->npw; // allocate rhog_mag[is*ngmc] and rhog_mag_save[is*ngmc] - rhog_mag = new std::complex[_ngmc * GlobalV::NSPIN]; - rhog_mag_save = new std::complex[_ngmc * GlobalV::NSPIN]; - ModuleBase::GlobalFunc::ZEROS(rhog_mag, _ngmc * GlobalV::NSPIN); - ModuleBase::GlobalFunc::ZEROS(rhog_mag_save, _ngmc * GlobalV::NSPIN); + rhog_mag = new std::complex[npw * GlobalV::NSPIN]; + rhog_mag_save = new std::complex[npw * GlobalV::NSPIN]; + ModuleBase::GlobalFunc::ZEROS(rhog_mag, npw * GlobalV::NSPIN); + ModuleBase::GlobalFunc::ZEROS(rhog_mag_save, npw * GlobalV::NSPIN); // get rhog_mag[is*ngmc] and rhog_mag_save[is*ngmc] - for (int ig = 0; ig < _ngmc; ig++) + for (int ig = 0; ig < npw; ig++) { rhog_mag[ig] = chr->rhog[0][ig] + chr->rhog[1][ig]; rhog_mag_save[ig] = chr->rhog_save[0][ig] + chr->rhog_save[1][ig]; } - for (int ig = 0; ig < _ngmc; ig++) + for (int ig = 0; ig < npw; ig++) { - rhog_mag[ig + _ngmc] = chr->rhog[0][ig] - chr->rhog[1][ig]; - rhog_mag_save[ig + _ngmc] = chr->rhog_save[0][ig] - chr->rhog_save[1][ig]; + rhog_mag[ig + npw] = chr->rhog[0][ig] - chr->rhog[1][ig]; + rhog_mag_save[ig + npw] = chr->rhog_save[0][ig] - chr->rhog_save[1][ig]; } // rhog_in = rhog_mag_save; @@ -369,12 +369,12 @@ void Charge_Mixing::mix_rho_recip_new(Charge* chr) // get rhog[is][ngmc] from rhog_mag[is*ngmc] for (int is = 0; is < GlobalV::NSPIN; is++) { - ModuleBase::GlobalFunc::ZEROS(chr->rhog[is], _ngmc); + ModuleBase::GlobalFunc::ZEROS(chr->rhog[is], npw); } - for (int ig = 0; ig < _ngmc; ig++) + for (int ig = 0; ig < npw; ig++) { - chr->rhog[0][ig] = 0.5 * (rhog_mag[ig] + rhog_mag[ig+_ngmc]); - chr->rhog[1][ig] = 0.5 * (rhog_mag[ig] - rhog_mag[ig+_ngmc]); + chr->rhog[0][ig] = 0.5 * (rhog_mag[ig] + rhog_mag[ig+npw]); + chr->rhog[1][ig] = 0.5 * (rhog_mag[ig] - rhog_mag[ig+npw]); } // delete delete[] rhog_mag; From 577459434b892c2c30b03e011434a3e04c9ad967 Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Thu, 4 Jan 2024 13:18:05 +0800 Subject: [PATCH 4/8] delete rho_mag and rhog_mag in Charge --- .../module_elecstate/module_charge/charge.cpp | 107 ------------------ .../module_elecstate/module_charge/charge.h | 47 -------- 2 files changed, 154 deletions(-) diff --git a/source/module_elecstate/module_charge/charge.cpp b/source/module_elecstate/module_charge/charge.cpp index d94a9abb9f..709c49de47 100644 --- a/source/module_elecstate/module_charge/charge.cpp +++ b/source/module_elecstate/module_charge/charge.cpp @@ -258,113 +258,6 @@ void Charge::renormalize_rho(void) return; } -// for real space magnetic density -void Charge::get_rho_mag(void) -{ - ModuleBase::TITLE("Charge", "get_rho_tot_mag"); - - for (int ir = 0; ir < nrxx; ir++) - { - rho_mag[ir] = rho[0][ir] + rho[1][ir]; - rho_mag_save[ir] = rho_save[0][ir] + rho_save[1][ir]; - } - for (int ir = 0; ir < nrxx; ir++) - { - rho_mag[ir + nrxx] = rho[0][ir] - rho[1][ir]; - rho_mag_save[ir + nrxx] = rho_save[0][ir] - rho_save[1][ir]; - } - return; -} - -void Charge::get_rho_from_mag(void) -{ - ModuleBase::TITLE("Charge", "get_rho_from_mag"); - for (int is = 0; is < nspin; is++) - { - ModuleBase::GlobalFunc::ZEROS(rho[is], nrxx); - //ModuleBase::GlobalFunc::ZEROS(rho_save[is], nrxx); - } - for (int ir = 0; ir < nrxx; ir++) - { - rho[0][ir] = 0.5 * (rho_mag[ir] + rho_mag[ir+nrxx]); - rho[1][ir] = 0.5 * (rho_mag[ir] - rho_mag[ir+nrxx]); - } - - return; -} - -void Charge::allocate_rho_mag(void) -{ - rho_mag = new double[nrxx * nspin]; - rho_mag_save = new double[nrxx * nspin]; - - ModuleBase::GlobalFunc::ZEROS(rho_mag, nrxx * nspin); - ModuleBase::GlobalFunc::ZEROS(rho_mag_save, nrxx * nspin); - - return; -} - -void Charge::destroy_rho_mag(void) -{ - delete[] rho_mag; - delete[] rho_mag_save; - - return; -} - -// for reciprocal space magnetic density -void Charge::get_rhog_mag(void) -{ - ModuleBase::TITLE("Charge", "get_rhog_tot_mag"); - - for (int ig = 0; ig < ngmc; ig++) - { - rhog_mag[ig] = rhog[0][ig] + rhog[1][ig]; - rhog_mag_save[ig] = rhog_save[0][ig] + rhog_save[1][ig]; - } - for (int ig = 0; ig < ngmc; ig++) - { - rhog_mag[ig + ngmc] = rhog[0][ig] - rhog[1][ig]; - rhog_mag_save[ig + ngmc] = rhog_save[0][ig] - rhog_save[1][ig]; - } - return; -} - -void Charge::get_rhog_from_mag(void) -{ - ModuleBase::TITLE("Charge", "get_rhog_from_mag"); - for (int is = 0; is < nspin; is++) - { - ModuleBase::GlobalFunc::ZEROS(rhog[is], ngmc); - } - for (int ig = 0; ig < ngmc; ig++) - { - rhog[0][ig] = 0.5 * (rhog_mag[ig] + rhog_mag[ig+ngmc]); - rhog[1][ig] = 0.5 * (rhog_mag[ig] - rhog_mag[ig+ngmc]); - } - - return; -} - -void Charge::allocate_rhog_mag(void) -{ - rhog_mag = new std::complex[ngmc * nspin]; - rhog_mag_save = new std::complex[ngmc * nspin]; - - ModuleBase::GlobalFunc::ZEROS(rhog_mag, ngmc * nspin); - ModuleBase::GlobalFunc::ZEROS(rhog_mag_save, ngmc * nspin); - - return; -} - -void Charge::destroy_rhog_mag(void) -{ - delete[] rhog_mag; - delete[] rhog_mag_save; - - return; -} - //------------------------------------------------------- // superposition of atomic charges contained in the array // rho_at (read from pseudopotential files) diff --git a/source/module_elecstate/module_charge/charge.h b/source/module_elecstate/module_charge/charge.h index 55bf2ee369..58bc2ef615 100644 --- a/source/module_elecstate/module_charge/charge.h +++ b/source/module_elecstate/module_charge/charge.h @@ -36,12 +36,6 @@ class Charge double **rho = nullptr; double **rho_save = nullptr; - // for magnetic density, onyl support nspin=2 now - double *rho_mag = nullptr; - double *rho_mag_save = nullptr; - std::complex *rhog_mag = nullptr; - std::complex *rhog_mag_save = nullptr; - std::complex **rhog = nullptr; std::complex **rhog_save = nullptr; @@ -89,47 +83,6 @@ class Charge void renormalize_rho(void); - // for magnetic density - /** - * @brief allocate rho_mag[is*nnrx] and rho_mag_save[is*nnrx] - */ - void allocate_rho_mag(void); - - /** - * @brief destroy rho_mag[is*nnrx] and rho_mag_save[is*nnrx] - */ - void destroy_rho_mag(void); - - /** - * @brief get rho_mag[is*nnrx] - */ - void get_rho_mag(void); - - /** - * @brief get rho[is][nnrx] from rho_mag[is*nnrx] - */ - void get_rho_from_mag(void); - - /** - * @brief // allocate rhog_mag[is*ngmc] and rhog_mag_save[is*ngmc] - */ - void allocate_rhog_mag(void); - - /** - * @brief destroy rhog_mag[is*ngmc] and rhog_mag_save[is*ngmc] - */ - void destroy_rhog_mag(void); - - /** - * @brief get rhog_mag[is*ngmc] - */ - void get_rhog_mag(void); - - /** - * @brief get rhog[is][nnrx] from rhog_mag[is*ngmc] - */ - void get_rhog_from_mag(void); - double sum_rho(void) const; void save_rho_before_sum_band(void); From 88b90573dc28849021eadf72f524799e52a8f2a5 Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Thu, 4 Jan 2024 13:22:59 +0800 Subject: [PATCH 5/8] fix some bugs --- source/module_elecstate/module_charge/charge_mixing.cpp | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/source/module_elecstate/module_charge/charge_mixing.cpp b/source/module_elecstate/module_charge/charge_mixing.cpp index 3cac601d32..0b761e55f2 100755 --- a/source/module_elecstate/module_charge/charge_mixing.cpp +++ b/source/module_elecstate/module_charge/charge_mixing.cpp @@ -343,7 +343,6 @@ void Charge_Mixing::mix_rho_recip_new(Charge* chr) rhog_in = rhog_mag_save; rhog_out = rhog_mag; // - const int npw = this->rhopw->npw; auto screen = std::bind(&Charge_Mixing::Kerker_screen_recip_new, this, std::placeholders::_1); auto twobeta_mix = [this, npw](std::complex* out, const std::complex* in, const std::complex* sres) { @@ -593,8 +592,8 @@ void Charge_Mixing::mix_rho_real(Charge* chr) rho_mag_save[ir + nrxx] = chr->rho_save[0][ir] - chr->rho_save[1][ir]; } // - rhor_in = chr->rho_mag_save; - rhor_out = chr->rho_mag; + rhor_in = rho_mag_save; + rhor_out = rho_mag; auto screen = std::bind(&Charge_Mixing::Kerker_screen_real, this, std::placeholders::_1); auto twobeta_mix = [this, nrxx](double* out, const double* in, const double* sres) { From bffce55769f8b7edf1043bbd7aa030e79cf4d918 Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Thu, 4 Jan 2024 13:23:36 +0800 Subject: [PATCH 6/8] delete AutoSetTest in charge_mixing_test --- .../test/charge_mixing_test.cpp | 34 ------------------- 1 file changed, 34 deletions(-) diff --git a/source/module_elecstate/test/charge_mixing_test.cpp b/source/module_elecstate/test/charge_mixing_test.cpp index fb20b2519a..09e25ccc99 100644 --- a/source/module_elecstate/test/charge_mixing_test.cpp +++ b/source/module_elecstate/test/charge_mixing_test.cpp @@ -315,40 +315,6 @@ TEST_F(ChargeMixingTest, SetMixingTest) EXPECT_THAT(output, testing::HasSubstr("This Mixing mode is not implemended yet,coming soon.")); } -/** -TEST_F(ChargeMixingTest, AutoSetTest) -{ - Charge_Mixing CMtest; - CMtest.set_rhopw(&pw_basis, &pw_basis); - GlobalV::SCF_THR_TYPE = 1; - GlobalV::NSPIN = 1; - - CMtest.set_mixing("broyden", 1.0, 1, 0.2, false); - CMtest.auto_set(0.0, GlobalC::ucell); - EXPECT_EQ(CMtest.mixing_beta, 1.0); - EXPECT_EQ(CMtest.mixing_gg0, 0.2); - - CMtest.need_auto_set(); - CMtest.auto_set(0.0, GlobalC::ucell); - EXPECT_EQ(CMtest.mixing_beta, 0.2); - EXPECT_EQ(CMtest.mixing->mixing_beta, 0.2); - EXPECT_EQ(CMtest.mixing_gg0, 1.0); - - CMtest.need_auto_set(); - CMtest.auto_set(1.0, GlobalC::ucell); - EXPECT_EQ(CMtest.mixing_beta, 0.7); - EXPECT_EQ(CMtest.mixing->mixing_beta, 0.7); - EXPECT_EQ(CMtest.mixing_gg0, 1.0); - - GlobalC::ucell.atoms = new Atom[1]; - GlobalC::ucell.ntype = 1; - GlobalC::ucell.atoms[0].ncpp.psd = "Sc"; - CMtest.need_auto_set(); - CMtest.auto_set(1.0, GlobalC::ucell); - EXPECT_EQ(CMtest.mixing_gg0, 1.0); -} -**/ - TEST_F(ChargeMixingTest, KerkerScreenRecipTest) { Charge_Mixing CMtest; From a34cb03cabff8896bf18e669c6cb6e72974e45b6 Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Thu, 4 Jan 2024 15:12:51 +0800 Subject: [PATCH 7/8] move renormalize_rho() outside Charge_Mixing --- source/module_elecstate/module_charge/charge_mixing.cpp | 4 ---- source/module_esolver/esolver_ks.cpp | 1 + 2 files changed, 1 insertion(+), 4 deletions(-) diff --git a/source/module_elecstate/module_charge/charge_mixing.cpp b/source/module_elecstate/module_charge/charge_mixing.cpp index 0b761e55f2..227d9a050a 100755 --- a/source/module_elecstate/module_charge/charge_mixing.cpp +++ b/source/module_elecstate/module_charge/charge_mixing.cpp @@ -501,9 +501,6 @@ void Charge_Mixing::mix_rho_recip_new(Charge* chr) } } - // renormalize rho in R-space would induce a error in K-space - //chr->renormalize_rho(); - // For kinetic energy density if ((XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5) && mixing_tau) { @@ -726,7 +723,6 @@ void Charge_Mixing::mix_rho_real(Charge* chr) delete[] rho_magabs_save; } - chr->renormalize_rho(); double *taur_out, *taur_in; if ((XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5) && mixing_tau) { diff --git a/source/module_esolver/esolver_ks.cpp b/source/module_esolver/esolver_ks.cpp index 97da14c4e7..d9877a0969 100644 --- a/source/module_esolver/esolver_ks.cpp +++ b/source/module_esolver/esolver_ks.cpp @@ -438,6 +438,7 @@ namespace ModuleESolver // } p_chgmix->mix_rho(pelec->charge); + if (GlobalV::SCF_THR_TYPE == 2) pelec->charge->renormalize_rho(); // renormalize rho in R-space would induce a error in K-space //----------charge mixing done----------- } } From 4aad817a37ad2662a282d1bfaed3bd3190865dc6 Mon Sep 17 00:00:00 2001 From: weiqingzhou Date: Thu, 4 Jan 2024 15:28:24 +0800 Subject: [PATCH 8/8] delete the mocks of rho_mag/renormalize_rho/sum_rho in charge_mixing_test --- .../test/charge_mixing_test.cpp | 169 ------------------ 1 file changed, 169 deletions(-) diff --git a/source/module_elecstate/test/charge_mixing_test.cpp b/source/module_elecstate/test/charge_mixing_test.cpp index 09e25ccc99..4233d0d99d 100644 --- a/source/module_elecstate/test/charge_mixing_test.cpp +++ b/source/module_elecstate/test/charge_mixing_test.cpp @@ -23,175 +23,6 @@ Charge::Charge() { } -double Charge::sum_rho(void) const -{ - ModuleBase::TITLE("Charge", "sum_rho"); - - double sum_rho = 0.0; - int nspin0 = (nspin == 2) ? 2 : 1; - - for (int is = 0; is < nspin0; is++) - { - for (int ir = 0; ir < nrxx; ir++) - { - if(GlobalV::use_paw) - { - sum_rho += this->rho[is][ir] + this->nhat[is][ir]; - } - else - { - sum_rho += this->rho[is][ir]; - } - } - } - - // multiply the sum of charge density by a factor - sum_rho *= 500.0 / static_cast(this->rhopw->nxyz); - -#ifdef __MPI - Parallel_Reduce::reduce_pool(sum_rho); -#endif - - // mohan fixed bug 2010-01-18, - // sum_rho may be smaller than 1, like Na bcc. - //if (sum_rho <= 0.1) - //{ - //GlobalV::ofs_warning << " sum_rho=" << sum_rho << std::endl; - //ModuleBase::WARNING_QUIT("Charge::renormalize_rho", "Can't find even an electron!"); - //} - - return sum_rho; -} - -void Charge::renormalize_rho(void) -{ - ModuleBase::TITLE("Charge", "renormalize_rho"); - - const double sr = this->sum_rho(); - GlobalV::ofs_warning << std::setprecision(15); - ModuleBase::GlobalFunc::OUT(GlobalV::ofs_warning, "charge before normalized", sr); - const double normalize_factor = GlobalV::nelec / sr; - - for (int is = 0; is < nspin; is++) - { - for (int ir = 0; ir < nrxx; ir++) - { - rho[is][ir] *= normalize_factor; - } - } - - ModuleBase::GlobalFunc::OUT(GlobalV::ofs_warning, "charge after normalized", this->sum_rho()); - - GlobalV::ofs_running << std::setprecision(6); - return; -} - -// for real space magnetic density -void Charge::get_rho_mag(void) -{ - ModuleBase::TITLE("Charge", "get_rho_tot_mag"); - - for (int ir = 0; ir < nrxx; ir++) - { - rho_mag[ir] = rho[0][ir] + rho[1][ir]; - rho_mag_save[ir] = rho_save[0][ir] + rho_save[1][ir]; - } - for (int ir = 0; ir < nrxx; ir++) - { - rho_mag[ir + nrxx] = rho[0][ir] - rho[1][ir]; - rho_mag_save[ir + nrxx] = rho_save[0][ir] - rho_save[1][ir]; - } - return; -} - -void Charge::get_rho_from_mag(void) -{ - ModuleBase::TITLE("Charge", "get_rho_from_mag"); - for (int is = 0; is < nspin; is++) - { - ModuleBase::GlobalFunc::ZEROS(rho[is], nrxx); - //ModuleBase::GlobalFunc::ZEROS(rho_save[is], nrxx); - } - for (int ir = 0; ir < nrxx; ir++) - { - rho[0][ir] = 0.5 * (rho_mag[ir] + rho_mag[ir+nrxx]); - rho[1][ir] = 0.5 * (rho_mag[ir] - rho_mag[ir+nrxx]); - } - - return; -} - -void Charge::allocate_rho_mag(void) -{ - rho_mag = new double[nrxx * nspin]; - rho_mag_save = new double[nrxx * nspin]; - - ModuleBase::GlobalFunc::ZEROS(rho_mag, nrxx * nspin); - ModuleBase::GlobalFunc::ZEROS(rho_mag_save, nrxx * nspin); - - return; -} - -void Charge::destroy_rho_mag(void) -{ - delete[] rho_mag; - delete[] rho_mag_save; - - return; -} - -// for reciprocal space magnetic density -void Charge::get_rhog_mag(void) -{ - ModuleBase::TITLE("Charge", "get_rhog_tot_mag"); - - for (int ig = 0; ig < ngmc; ig++) - { - rhog_mag[ig] = rhog[0][ig] + rhog[1][ig]; - rhog_mag_save[ig] = rhog_save[0][ig] + rhog_save[1][ig]; - } - for (int ig = 0; ig < ngmc; ig++) - { - rhog_mag[ig + ngmc] = rhog[0][ig] - rhog[1][ig]; - rhog_mag_save[ig + ngmc] = rhog_save[0][ig] - rhog_save[1][ig]; - } - return; -} - -void Charge::get_rhog_from_mag(void) -{ - ModuleBase::TITLE("Charge", "get_rhog_from_mag"); - for (int is = 0; is < nspin; is++) - { - ModuleBase::GlobalFunc::ZEROS(rhog[is], ngmc); - } - for (int ig = 0; ig < ngmc; ig++) - { - rhog[0][ig] = 0.5 * (rhog_mag[ig] + rhog_mag[ig+ngmc]); - rhog[1][ig] = 0.5 * (rhog_mag[ig] - rhog_mag[ig+ngmc]); - } - - return; -} - -void Charge::allocate_rhog_mag(void) -{ - rhog_mag = new std::complex[ngmc * nspin]; - rhog_mag_save = new std::complex[ngmc * nspin]; - - ModuleBase::GlobalFunc::ZEROS(rhog_mag, ngmc * nspin); - ModuleBase::GlobalFunc::ZEROS(rhog_mag_save, ngmc * nspin); - - return; -} - -void Charge::destroy_rhog_mag(void) -{ - delete[] rhog_mag; - delete[] rhog_mag_save; - - return; -} void Charge::set_rhopw(ModulePW::PW_Basis* rhopw_in) { this->rhopw = rhopw_in;