Refactor:remove cal_tau from ElecStateLCAO (#5802)

A-006 · web-flow · commit 414446ac18aa · 2025-01-04T17:33:07.000+08:00
* modify the cal_tau in lcao

* add template for cal_tau

* updatea func for cal_tau
diff --git a/source/module_elecstate/elecstate_lcao.cpp b/source/module_elecstate/elecstate_lcao.cpp
@@ -8,6 +8,7 @@
 #include "module_hamilt_lcao/module_gint/grid_technique.h"
 #include "module_hamilt_pw/hamilt_pwdft/global.h"
 #include "module_parameter/parameter.h"
+#include "elecstate_lcao_cal_tau.h"
 
 #include <vector>
 
@@ -64,7 +65,7 @@ void ElecStateLCAO<std::complex<double>>::psiToRho(const psi::Psi<std::complex<d
 
     if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
     {
-        this->cal_tau(psi);
+        elecstate::lcao_cal_tau_k(gint_k, this->charge);
     }
 
     this->charge->renormalize_rho();
@@ -99,7 +100,7 @@ void ElecStateLCAO<double>::psiToRho(const psi::Psi<double>& psi)
 
     if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
     {
-        this->cal_tau(psi);
+        elecstate::lcao_cal_tau_gamma(gint_gamma, this->charge);
     }
 
     this->charge->renormalize_rho();
diff --git a/source/module_elecstate/elecstate_lcao.h b/source/module_elecstate/elecstate_lcao.h
@@ -46,7 +46,6 @@ class ElecStateLCAO : public ElecState
     // virtual void psiToRho(const psi::Psi<double>& psi) override;
     //  return current electronic density rho, as a input for constructing Hamiltonian
     //  const double* getRho(int spin) const override;
-    virtual void cal_tau(const psi::Psi<TK>& psi) override;
 
     // update charge density for next scf step
     // void getNewRho() override;
diff --git a/source/module_elecstate/elecstate_lcao_cal_tau.cpp b/source/module_elecstate/elecstate_lcao_cal_tau.cpp
@@ -1,41 +1,56 @@
 #include "elecstate_lcao.h"
-
+#include "elecstate_lcao_cal_tau.h"
 #include "module_base/timer.h"
 
 namespace elecstate
 {
 
 // calculate the kinetic energy density tau, multi-k case
-template <>
-void ElecStateLCAO<std::complex<double>>::cal_tau(const psi::Psi<std::complex<double>>& psi)
+void lcao_cal_tau_k(Gint_k* gint_k, 
+                    Charge* charge)
 {
     ModuleBase::timer::tick("ElecStateLCAO", "cal_tau");
 
     for (int is = 0; is < PARAM.inp.nspin; is++)
     {
-        ModuleBase::GlobalFunc::ZEROS(this->charge->kin_r[is], this->charge->nrxx);
+        ModuleBase::GlobalFunc::ZEROS(charge->kin_r[is], charge->nrxx);
     }
-    Gint_inout inout1(this->charge->kin_r, Gint_Tools::job_type::tau, PARAM.inp.nspin);
-    this->gint_k->cal_gint(&inout1);
+    Gint_inout inout1(charge->kin_r, Gint_Tools::job_type::tau, PARAM.inp.nspin);
+    gint_k->cal_gint(&inout1);
 
     ModuleBase::timer::tick("ElecStateLCAO", "cal_tau");
     return;
 }
 
 // calculate the kinetic energy density tau, gamma-only case
-template <>
-void ElecStateLCAO<double>::cal_tau(const psi::Psi<double>& psi)
+void lcao_cal_tau_gamma(Gint_Gamma* gint_gamma,
+                        Charge* charge)
 {
     ModuleBase::timer::tick("ElecStateLCAO", "cal_tau");
 
     for (int is = 0; is < PARAM.inp.nspin; is++)
     {
-        ModuleBase::GlobalFunc::ZEROS(this->charge->kin_r[is], this->charge->nrxx);
+        ModuleBase::GlobalFunc::ZEROS(charge->kin_r[is], charge->nrxx);
     }
-    Gint_inout inout1(this->charge->kin_r, Gint_Tools::job_type::tau, PARAM.inp.nspin);
-    this->gint_gamma->cal_gint(&inout1);
+    Gint_inout inout1(charge->kin_r, Gint_Tools::job_type::tau, PARAM.inp.nspin);
+    gint_gamma->cal_gint(&inout1);
 
     ModuleBase::timer::tick("ElecStateLCAO", "cal_tau");
     return;
 }
-}
+template <> 
+void lcao_cal_tau<double>(Gint_Gamma* gint_gamma, 
+                  Gint_k* gint_k, 
+                  Charge* charge)
+{
+    lcao_cal_tau_gamma(gint_gamma, charge);
+}
+template <> 
+void lcao_cal_tau<complex<double>>(Gint_Gamma* gint_gamma, 
+                    Gint_k* gint_k, 
+                    Charge* charge)
+{
+    lcao_cal_tau_k(gint_k, charge);
+}
+
+} // namespace elecstate
diff --git a/source/module_elecstate/elecstate_lcao_cal_tau.h b/source/module_elecstate/elecstate_lcao_cal_tau.h
@@ -0,0 +1,21 @@
+#ifndef ELECSTATE_LCAO_CAL_TAU_H
+#define ELECSTATE_LCAO_CAL_TAU_H
+#include "module_elecstate/module_charge/charge.h"
+#include "module_hamilt_lcao/module_gint/gint_gamma.h"
+#include "module_hamilt_lcao/module_gint/gint_k.h"
+namespace elecstate
+{
+
+    void lcao_cal_tau_k(Gint_k* gint_k, 
+                        Charge* charge);
+
+    void lcao_cal_tau_gamma(Gint_Gamma* gint_gamma,
+                            Charge* charge);
+
+    template <typename T>
+    void lcao_cal_tau(Gint_Gamma* gint_gamma, 
+                      Gint_k* gint_k, 
+                      Charge* charge);
+
+}
+#endif
diff --git a/source/module_esolver/esolver_ks_lcao.cpp b/source/module_esolver/esolver_ks_lcao.cpp
@@ -34,6 +34,7 @@
 #include "module_base/global_function.h"
 #include "module_cell/module_neighbor/sltk_grid_driver.h"
 #include "module_elecstate/cal_ux.h"
+#include "module_elecstate/elecstate_lcao_cal_tau.h"
 #include "module_elecstate/module_charge/symmetry_rho.h"
 #include "module_elecstate/occupy.h"
 #include "module_hamilt_lcao/hamilt_lcaodft/LCAO_domain.h" // need divide_HS_in_frag
@@ -927,8 +928,9 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
     // 1) calculate the kinetic energy density tau, sunliang 2024-09-18
     if (PARAM.inp.out_elf[0] > 0)
     {
-        assert(this->psi != nullptr);
-        this->pelec->cal_tau(*(this->psi));
+        elecstate::lcao_cal_tau<TK>(&(this->GG), 
+                                    &(this->GK),
+                                    this->pelec->charge);
     }
 
     //! 2) call after_scf() of ESolver_KS
diff --git a/source/module_rdmft/update_state_rdmft.cpp b/source/module_rdmft/update_state_rdmft.cpp
@@ -8,7 +8,7 @@
 #include "module_elecstate/module_dm/cal_dm_psi.h"
 #include "module_elecstate/module_dm/density_matrix.h"
 #include "module_elecstate/module_charge/symmetry_rho.h"
-
+#include "module_elecstate/elecstate_lcao_cal_tau.h"
 
 namespace rdmft
 {
@@ -118,7 +118,7 @@ void RDMFT<TK, TR>::update_charge(UnitCell& ucell)
             // }
             // Gint_inout inout1(charge->kin_r, Gint_Tools::job_type::tau);
             // GG->cal_gint(&inout1);
-            this->pelec->cal_tau(wfc);
+            elecstate::lcao_cal_tau_gamma(GG, charge);
         }
 
         charge->renormalize_rho();
@@ -148,7 +148,7 @@ void RDMFT<TK, TR>::update_charge(UnitCell& ucell)
             // }
             // Gint_inout inout1(charge->kin_r, Gint_Tools::job_type::tau);
             // GK->cal_gint(&inout1);
-            this->pelec->cal_tau(wfc);
+            elecstate::lcao_cal_tau_k(GK, charge);
         }
 
         charge->renormalize_rho();

Original file line number	Diff line number	Diff line change
`@@ -8,6 +8,7 @@`
`8`	`8`	`#include "module_hamilt_lcao/module_gint/grid_technique.h"`
`9`	`9`	`#include "module_hamilt_pw/hamilt_pwdft/global.h"`
`10`	`10`	`#include "module_parameter/parameter.h"`
	`11`	`+#include "elecstate_lcao_cal_tau.h"`
`11`	`12`
`12`	`13`	`#include <vector>`
`13`	`14`
`@@ -64,7 +65,7 @@ void ElecStateLCAO<std::complex<double>>::psiToRho(const psi::Psi<std::complex<d`
`64`	`65`
`65`	`66`	`if (XC_Functional::get_func_type() == 3 \|\| XC_Functional::get_func_type() == 5)`
`66`	`67`	`{`
`67`		`- this->cal_tau(psi);`
	`68`	`+ elecstate::lcao_cal_tau_k(gint_k, this->charge);`
`68`	`69`	`}`
`69`	`70`
`70`	`71`	`this->charge->renormalize_rho();`
`@@ -99,7 +100,7 @@ void ElecStateLCAO<double>::psiToRho(const psi::Psi<double>& psi)`
`99`	`100`
`100`	`101`	`if (XC_Functional::get_func_type() == 3 \|\| XC_Functional::get_func_type() == 5)`
`101`	`102`	`{`
`102`		`- this->cal_tau(psi);`
	`103`	`+ elecstate::lcao_cal_tau_gamma(gint_gamma, this->charge);`
`103`	`104`	`}`
`104`	`105`
`105`	`106`	`this->charge->renormalize_rho();`