From 242184cdbc2fa0b59f46fdacd7c1701704ae3222 Mon Sep 17 00:00:00 2001 From: mohanchen Date: Thu, 16 Oct 2025 14:47:46 +0800 Subject: [PATCH 1/4] add update_pot in source_estate --- source/source_esolver/esolver_ks.cpp | 3 ++- source/source_esolver/esolver_ks_lcao.cpp | 16 ---------------- source/source_esolver/esolver_ks_pw.cpp | 21 ++------------------- source/source_estate/update_pot.cpp | 20 ++++++++++++++++++++ source/source_estate/update_pot.h | 16 ++++++++++++++++ 5 files changed, 40 insertions(+), 36 deletions(-) create mode 100644 source/source_estate/update_pot.cpp create mode 100644 source/source_estate/update_pot.h diff --git a/source/source_esolver/esolver_ks.cpp b/source/source_esolver/esolver_ks.cpp index a4e8126ad0..d98571d541 100644 --- a/source/source_esolver/esolver_ks.cpp +++ b/source/source_esolver/esolver_ks.cpp @@ -23,6 +23,7 @@ #include "source_io/json_output/init_info.h" #include "source_io/json_output/output_info.h" +#include "source_estate/update_pot.h" // mohan add 20251016 namespace ModuleESolver { @@ -374,7 +375,7 @@ void ESolver_KS::iter_finish(UnitCell& ucell, const int istep, int& i if (this->scf_ene_thr > 0.0) { // calculate energy of output charge density - this->update_pot(ucell, istep, iter, conv_esolver); + elecstate::update_pot(ucell, istep, iter, conv_esolver); this->pelec->cal_energies(2); // 2 means Kohn-Sham functional // now, etot_old is the energy of input density, while etot is the energy of output density this->pelec->f_en.etot_delta = this->pelec->f_en.etot - this->pelec->f_en.etot_old; diff --git a/source/source_esolver/esolver_ks_lcao.cpp b/source/source_esolver/esolver_ks_lcao.cpp index b7601da811..c7171d0de7 100644 --- a/source/source_esolver/esolver_ks_lcao.cpp +++ b/source/source_esolver/esolver_ks_lcao.cpp @@ -446,22 +446,6 @@ void ESolver_KS_LCAO::hamilt2rho_single(UnitCell& ucell, int istep, int this->pelec->f_en.deband = this->pelec->cal_delta_eband(ucell); } -template -void ESolver_KS_LCAO::update_pot(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver) -{ - ModuleBase::TITLE("ESolver_KS_LCAO", "update_pot"); - - if (!conv_esolver) - { - elecstate::cal_ux(ucell); - this->pelec->pot->update_from_charge(&this->chr, &ucell); - this->pelec->f_en.descf = this->pelec->cal_delta_escf(); - } - else - { - this->pelec->cal_converged(); - } -} template void ESolver_KS_LCAO::iter_finish(UnitCell& ucell, const int istep, int& iter, bool& conv_esolver) diff --git a/source/source_esolver/esolver_ks_pw.cpp b/source/source_esolver/esolver_ks_pw.cpp index d2e2612070..3fd406d086 100644 --- a/source/source_esolver/esolver_ks_pw.cpp +++ b/source/source_esolver/esolver_ks_pw.cpp @@ -28,6 +28,7 @@ #include "source_estate/setup_estate_pw.h" // mohan add 20251005 #include "source_io/ctrl_output_pw.h" // mohan add 20250927 #include "source_estate/module_charge/chgmixing.h" // use charge mixing, mohan add 20251006 +#include "source_estate/update_pot.h" // mohan add 20251016 namespace ModuleESolver { @@ -268,24 +269,6 @@ void ESolver_KS_PW::hamilt2rho_single(UnitCell& ucell, const int iste ModuleBase::timer::tick("ESolver_KS_PW", "hamilt2rho_single"); } -// Temporary, it should be rewritten with Hamilt class. -template -void ESolver_KS_PW::update_pot(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver) -{ - if (!conv_esolver) - { - elecstate::cal_ux(ucell); - this->pelec->pot->update_from_charge(&this->chr, &ucell); - this->pelec->f_en.descf = this->pelec->cal_delta_escf(); -#ifdef __MPI - MPI_Bcast(&(this->pelec->f_en.descf), 1, MPI_DOUBLE, 0, BP_WORLD); -#endif - } - else - { - this->pelec->cal_converged(); - } -} template void ESolver_KS_PW::iter_finish(UnitCell& ucell, const int istep, int& iter, bool& conv_esolver) @@ -343,7 +326,7 @@ void ESolver_KS_PW::iter_finish(UnitCell& ucell, const int istep, int << std::endl; exx_helper.op_exx->first_iter = false; XC_Functional::set_xc_type(ucell.atoms[0].ncpp.xc_func); - update_pot(ucell, istep, iter, conv_esolver); + elecstate::update_pot(ucell, istep, iter, conv_esolver); exx_helper.iter_inc(); } } diff --git a/source/source_estate/update_pot.cpp b/source/source_estate/update_pot.cpp new file mode 100644 index 0000000000..4f4c6791f7 --- /dev/null +++ b/source/source_estate/update_pot.cpp @@ -0,0 +1,20 @@ +#include "source_estate/update_pot.h" + +void elecstate::update_pot(UnitCell& ucell, // unitcell + elecstate::ElecState* &pelec, // pointer of electrons + Charge &chr) // charge density +{ + if (!conv_esolver) + { + elecstate::cal_ux(ucell); + pelec->pot->update_from_charge(&chr, &ucell); + pelec->f_en.descf = pelec->cal_delta_escf(); +#ifdef __MPI + MPI_Bcast(&(pelec->f_en.descf), 1, MPI_DOUBLE, 0, BP_WORLD); +#endif + } + else + { + pelec->cal_converged(); + } +} diff --git a/source/source_estate/update_pot.h b/source/source_estate/update_pot.h new file mode 100644 index 0000000000..33eca12f7c --- /dev/null +++ b/source/source_estate/update_pot.h @@ -0,0 +1,16 @@ +#ifndef UPDATE_POT_H +#define UPDATE_POT_H + +#include "source_cell/unitcell.h" +#include "source_estate/elecstate.h" + +namespace elecstate +{ + +void update_pot(UnitCell& ucell, // unitcell + elecstate::ElecState* &pelec, // pointer of electrons + Charge &chr); // charge density +} + + +#endif From 9441396323cd0405e9d84d63e7dda55879e8645b Mon Sep 17 00:00:00 2001 From: mohanchen Date: Fri, 17 Oct 2025 16:59:50 +0800 Subject: [PATCH 2/4] split update_pot in rt-TDDFT into two functions, one is the original update_pot, which is replaced by update_pot function in estate, the other is save2 function used in after_iter --- source/Makefile.Objects | 1 + source/source_esolver/esolver_ks.cpp | 4 ++-- source/source_esolver/esolver_ks.h | 3 --- source/source_esolver/esolver_ks_lcao.h | 2 -- source/source_esolver/esolver_ks_lcao_tddft.cpp | 6 +++++- source/source_esolver/esolver_ks_lcao_tddft.h | 3 ++- source/source_esolver/esolver_ks_pw.cpp | 4 ++-- source/source_esolver/esolver_ks_pw.h | 2 -- source/source_esolver/rho_restart.cpp | 2 ++ source/source_estate/CMakeLists.txt | 1 + source/source_estate/update_pot.cpp | 6 ++++-- source/source_estate/update_pot.h | 5 +++-- 12 files changed, 22 insertions(+), 17 deletions(-) diff --git a/source/Makefile.Objects b/source/Makefile.Objects index b05e69f8f9..e612261c6e 100644 --- a/source/Makefile.Objects +++ b/source/Makefile.Objects @@ -250,6 +250,7 @@ OBJS_ELECSTAT=elecstate.o\ read_pseudo.o\ cal_wfc.o\ setup_estate_pw.o\ + update_pot.o\ OBJS_ELECSTAT_LCAO=elecstate_lcao.o\ elecstate_lcao_cal_tau.o\ diff --git a/source/source_esolver/esolver_ks.cpp b/source/source_esolver/esolver_ks.cpp index d98571d541..d9544209d3 100644 --- a/source/source_esolver/esolver_ks.cpp +++ b/source/source_esolver/esolver_ks.cpp @@ -375,7 +375,7 @@ void ESolver_KS::iter_finish(UnitCell& ucell, const int istep, int& i if (this->scf_ene_thr > 0.0) { // calculate energy of output charge density - elecstate::update_pot(ucell, istep, iter, conv_esolver); + elecstate::update_pot(ucell, this->pelec, this->chr, conv_esolver); this->pelec->cal_energies(2); // 2 means Kohn-Sham functional // now, etot_old is the energy of input density, while etot is the energy of output density this->pelec->f_en.etot_delta = this->pelec->f_en.etot - this->pelec->f_en.etot_old; @@ -433,7 +433,7 @@ void ESolver_KS::iter_finish(UnitCell& ucell, const int istep, int& i #endif // 4) Update potentials (should be done every SF iter) - this->update_pot(ucell, istep, iter, conv_esolver); + elecstate::update_pot(ucell, this->pelec, this->chr, conv_esolver); // 5) calculate energies // 1 means Harris-Foulkes functional diff --git a/source/source_esolver/esolver_ks.h b/source/source_esolver/esolver_ks.h index d54794dd89..3f271dc4cc 100644 --- a/source/source_esolver/esolver_ks.h +++ b/source/source_esolver/esolver_ks.h @@ -46,9 +46,6 @@ class ESolver_KS : public ESolver_FP //! Something to do after SCF iterations when SCF is converged or comes to the max iter step. virtual void after_scf(UnitCell& ucell, const int istep, const bool conv_esolver) override; - //! It should be replaced by a function in Hamilt Class - virtual void update_pot(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver){}; - //! Hamiltonian hamilt::Hamilt* p_hamilt = nullptr; diff --git a/source/source_esolver/esolver_ks_lcao.h b/source/source_esolver/esolver_ks_lcao.h index 13359f7c8e..67ca7a0482 100644 --- a/source/source_esolver/esolver_ks_lcao.h +++ b/source/source_esolver/esolver_ks_lcao.h @@ -52,8 +52,6 @@ class ESolver_KS_LCAO : public ESolver_KS virtual void hamilt2rho_single(UnitCell& ucell, const int istep, const int iter, const double ethr) override; - virtual void update_pot(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver) override; - virtual void iter_finish(UnitCell& ucell, const int istep, int& iter, bool& conv_esolver) override; virtual void after_scf(UnitCell& ucell, const int istep, const bool conv_esolver) override; diff --git a/source/source_esolver/esolver_ks_lcao_tddft.cpp b/source/source_esolver/esolver_ks_lcao_tddft.cpp index c68ff0955b..20f1feddce 100644 --- a/source/source_esolver/esolver_ks_lcao_tddft.cpp +++ b/source/source_esolver/esolver_ks_lcao_tddft.cpp @@ -336,16 +336,19 @@ void ESolver_KS_LCAO_TDDFT::iter_finish( GlobalV::ofs_running << std::endl; } + this->save2(ucell, istep, iter, conv_esolver); + ESolver_KS_LCAO, TR>::iter_finish(ucell, istep, iter, conv_esolver); } template -void ESolver_KS_LCAO_TDDFT::update_pot(UnitCell& ucell, +void ESolver_KS_LCAO_TDDFT::save2(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver) { // Calculate new potential according to new Charge Density +/* if (!conv_esolver) { elecstate::cal_ux(ucell); @@ -356,6 +359,7 @@ void ESolver_KS_LCAO_TDDFT::update_pot(UnitCell& ucell, { this->pelec->cal_converged(); } +*/ const int nloc = this->pv.nloc; const int ncol_nbands = this->pv.ncol_bands; diff --git a/source/source_esolver/esolver_ks_lcao_tddft.h b/source/source_esolver/esolver_ks_lcao_tddft.h index fc0a62e41d..fe3f55cf29 100644 --- a/source/source_esolver/esolver_ks_lcao_tddft.h +++ b/source/source_esolver/esolver_ks_lcao_tddft.h @@ -64,7 +64,8 @@ class ESolver_KS_LCAO_TDDFT : public ESolver_KS_LCAO, TR> virtual void hamilt2rho_single(UnitCell& ucell, const int istep, const int iter, const double ethr) override; - virtual void update_pot(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver) override; + // mohan change update_pot to save2, 2025-10-17 + virtual void save2(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver); virtual void iter_finish(UnitCell& ucell, const int istep, int& iter, bool& conv_esolver) override; diff --git a/source/source_esolver/esolver_ks_pw.cpp b/source/source_esolver/esolver_ks_pw.cpp index 3fd406d086..a8e5d553b7 100644 --- a/source/source_esolver/esolver_ks_pw.cpp +++ b/source/source_esolver/esolver_ks_pw.cpp @@ -326,8 +326,8 @@ void ESolver_KS_PW::iter_finish(UnitCell& ucell, const int istep, int << std::endl; exx_helper.op_exx->first_iter = false; XC_Functional::set_xc_type(ucell.atoms[0].ncpp.xc_func); - elecstate::update_pot(ucell, istep, iter, conv_esolver); - exx_helper.iter_inc(); + elecstate::update_pot(ucell, this->pelec, this->chr, conv_esolver); + exx_helper.iter_inc(); } } } diff --git a/source/source_esolver/esolver_ks_pw.h b/source/source_esolver/esolver_ks_pw.h index acdd7083ee..09fbefc175 100644 --- a/source/source_esolver/esolver_ks_pw.h +++ b/source/source_esolver/esolver_ks_pw.h @@ -41,8 +41,6 @@ class ESolver_KS_PW : public ESolver_KS virtual void iter_init(UnitCell& ucell, const int istep, const int iter) override; - virtual void update_pot(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver) override; - virtual void iter_finish(UnitCell& ucell, const int istep, int& iter, bool& conv_esolver) override; virtual void after_scf(UnitCell& ucell, const int istep, const bool conv_esolver) override; diff --git a/source/source_esolver/rho_restart.cpp b/source/source_esolver/rho_restart.cpp index bd7497f34b..3a8c4b40b5 100644 --- a/source/source_esolver/rho_restart.cpp +++ b/source/source_esolver/rho_restart.cpp @@ -1,5 +1,6 @@ #include "rho_restart.h" +/* void ModuleESolver::rho_restart(const Input_para& inp, const UnitCell& ucell, const elecstate::ElecState& elec, @@ -145,3 +146,4 @@ void ModuleESolver::rho_restart(const Input_para& inp, #endif } +*/ diff --git a/source/source_estate/CMakeLists.txt b/source/source_estate/CMakeLists.txt index a8380ade4e..8a54110e55 100644 --- a/source/source_estate/CMakeLists.txt +++ b/source/source_estate/CMakeLists.txt @@ -39,6 +39,7 @@ list(APPEND objects read_pseudo.cpp cal_wfc.cpp setup_estate_pw.cpp + update_pot.cpp ) if(ENABLE_LCAO) diff --git a/source/source_estate/update_pot.cpp b/source/source_estate/update_pot.cpp index 4f4c6791f7..1cc51da836 100644 --- a/source/source_estate/update_pot.cpp +++ b/source/source_estate/update_pot.cpp @@ -1,8 +1,10 @@ #include "source_estate/update_pot.h" -void elecstate::update_pot(UnitCell& ucell, // unitcell +void elecstate::update_pot(const UnitCell& ucell, // unitcell elecstate::ElecState* &pelec, // pointer of electrons - Charge &chr) // charge density + const Charge &chr, + const bool conv_esolver + ) // charge density { if (!conv_esolver) { diff --git a/source/source_estate/update_pot.h b/source/source_estate/update_pot.h index 33eca12f7c..27c55bf582 100644 --- a/source/source_estate/update_pot.h +++ b/source/source_estate/update_pot.h @@ -7,9 +7,10 @@ namespace elecstate { -void update_pot(UnitCell& ucell, // unitcell +void update_pot(const UnitCell& ucell, // unitcell elecstate::ElecState* &pelec, // pointer of electrons - Charge &chr); // charge density + const Charge &chr, + const bool conv_esolver); // charge density } From dc200841cf7166f985e1e61c7109005f08cb04fe Mon Sep 17 00:00:00 2001 From: mohanchen Date: Fri, 17 Oct 2025 20:56:48 +0800 Subject: [PATCH 3/4] fix bugs --- source/source_estate/update_pot.cpp | 3 ++- source/source_estate/update_pot.h | 2 +- 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/source/source_estate/update_pot.cpp b/source/source_estate/update_pot.cpp index 1cc51da836..8a49f66e0c 100644 --- a/source/source_estate/update_pot.cpp +++ b/source/source_estate/update_pot.cpp @@ -1,6 +1,7 @@ #include "source_estate/update_pot.h" +#include "source_estate/cal_ux.h" -void elecstate::update_pot(const UnitCell& ucell, // unitcell +void elecstate::update_pot(UnitCell& ucell, // unitcell elecstate::ElecState* &pelec, // pointer of electrons const Charge &chr, const bool conv_esolver diff --git a/source/source_estate/update_pot.h b/source/source_estate/update_pot.h index 27c55bf582..cce4c42428 100644 --- a/source/source_estate/update_pot.h +++ b/source/source_estate/update_pot.h @@ -7,7 +7,7 @@ namespace elecstate { -void update_pot(const UnitCell& ucell, // unitcell +void update_pot(UnitCell& ucell, // unitcell elecstate::ElecState* &pelec, // pointer of electrons const Charge &chr, const bool conv_esolver); // charge density From 7829228d9a1f125351fe7fa466d8d3c81094d678 Mon Sep 17 00:00:00 2001 From: mohanchen Date: Fri, 17 Oct 2025 22:54:53 +0800 Subject: [PATCH 4/4] update --- source/source_esolver/esolver_ks_lcao_tddft.cpp | 3 ++- source/source_esolver/esolver_ks_lcao_tddft.h | 2 +- 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/source/source_esolver/esolver_ks_lcao_tddft.cpp b/source/source_esolver/esolver_ks_lcao_tddft.cpp index 20f1feddce..3072579164 100644 --- a/source/source_esolver/esolver_ks_lcao_tddft.cpp +++ b/source/source_esolver/esolver_ks_lcao_tddft.cpp @@ -336,9 +336,10 @@ void ESolver_KS_LCAO_TDDFT::iter_finish( GlobalV::ofs_running << std::endl; } + ESolver_KS_LCAO, TR>::iter_finish(ucell, istep, iter, conv_esolver); + this->save2(ucell, istep, iter, conv_esolver); - ESolver_KS_LCAO, TR>::iter_finish(ucell, istep, iter, conv_esolver); } template diff --git a/source/source_esolver/esolver_ks_lcao_tddft.h b/source/source_esolver/esolver_ks_lcao_tddft.h index fe3f55cf29..75bd287591 100644 --- a/source/source_esolver/esolver_ks_lcao_tddft.h +++ b/source/source_esolver/esolver_ks_lcao_tddft.h @@ -65,7 +65,7 @@ class ESolver_KS_LCAO_TDDFT : public ESolver_KS_LCAO, TR> virtual void hamilt2rho_single(UnitCell& ucell, const int istep, const int iter, const double ethr) override; // mohan change update_pot to save2, 2025-10-17 - virtual void save2(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver); + void save2(UnitCell& ucell, const int istep, const int iter, const bool conv_esolver); virtual void iter_finish(UnitCell& ucell, const int istep, int& iter, bool& conv_esolver) override;