diff --git a/source/module_esolver/esolver_ks.cpp b/source/module_esolver/esolver_ks.cpp index 215182775d..bbbe1bb4ad 100644 --- a/source/module_esolver/esolver_ks.cpp +++ b/source/module_esolver/esolver_ks.cpp @@ -547,6 +547,24 @@ void ESolver_KS::runner(const int istep, UnitCell& ucell) this->conv_esolver = (drho < this->scf_thr && not_restart_step && is_U_converged); + // add energy threshold for SCF convergence + if (this->scf_ene_thr > 0.0) + { + // calculate energy of output charge density + this->update_pot(istep, iter); + 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; + // output etot_delta + GlobalV::ofs_running << " DeltaE_womix = " << this->pelec->f_en.etot_delta * ModuleBase::Ry_to_eV << " eV" << std::endl; + if (iter > 1 && this->conv_esolver == 1) // only check when density is converged + { + // update the convergence flag + this->conv_esolver + = (std::abs(this->pelec->f_en.etot_delta * ModuleBase::Ry_to_eV) < this->scf_ene_thr); + } + } + // If drho < hsolver_error in the first iter or drho < scf_thr, we // do not change rho. if (drho < hsolver_error || this->conv_esolver) @@ -669,13 +687,6 @@ void ESolver_KS::iter_finish(int& iter) } this->pelec->f_en.etot_delta = this->pelec->f_en.etot - this->pelec->f_en.etot_old; this->pelec->f_en.etot_old = this->pelec->f_en.etot; - - // add a energy threshold for SCF convergence - if (this->scf_ene_thr > 0.0 && this->conv_esolver == 1) // only check when density is converged - { - this->conv_esolver - = (std::abs(this->pelec->f_en.etot_delta * ModuleBase::Ry_to_eV) < this->scf_ene_thr); - } } //! Something to do after SCF iterations when SCF is converged or comes to the max iter step.