Begin EXX KPAR

Author: Flying-dragon-boxing

source/source_pw/module_pwdft/operator_pw/exx_pw_pot.cpp

@@ -32,16 +32,10 @@ void get_exx_potential(const K_Vectors* kv,
 
     // calculate Fock pot
     auto param_fock = GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock];
-    for (auto param: param_fock)
+    for (int i = 0; i < param_fock.size(); i++)
     {
-        double exx_div = exx_divergence(Conv_Coulomb_Pot_K::Coulomb_Type::Fock,
-                                        0.0,
-                                        kv,
-                                        wfcpw,
-                                        rhopw_dev,
-                                        tpiba,
-                                        gamma_extrapolation,
-                                        ucell_omega);
+        auto param = param_fock[i];
+        double exx_div = OperatorEXXPW<std::complex<Real>, Device>::fock_div[i];
         double alpha = std::stod(param["alpha"]);
         const ModuleBase::Vector3<double> k_c = wfcpw->kvec_c[ik];
         const ModuleBase::Vector3<double> k_d = wfcpw->kvec_d[ik];
@@ -98,20 +92,21 @@ void get_exx_potential(const K_Vectors* kv,
 
     // calculate erfc pot
     auto param_erfc = GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc];
-    for (auto param: param_erfc)
+    for (int i = 0; i < param_erfc.size(); i++)
     {
+        auto param = param_erfc[i];
         double erfc_omega = std::stod(param["omega"]);
         double erfc_omega2 = erfc_omega * erfc_omega;
         double alpha = std::stod(param["alpha"]);
+        // double exx_div = OperatorEXXPW<std::complex<Real>, Device>::erfc_div[i];
         double exx_div = exx_divergence(Conv_Coulomb_Pot_K::Coulomb_Type::Erfc,
-                                        erfc_omega,
-                                        kv,
-                                        wfcpw,
-                                        rhopw_dev,
-                                        tpiba,
-                                        gamma_extrapolation,
-                                        ucell_omega);
-
+                                          erfc_omega,
+                                          kv,
+                                          wfcpw,
+                                          rhopw_dev,
+                                          tpiba,
+                                          gamma_extrapolation,
+                                          ucell_omega);
         const ModuleBase::Vector3<double> k_c = wfcpw->kvec_c[ik];
         const ModuleBase::Vector3<double> k_d = wfcpw->kvec_d[ik];
         const ModuleBase::Vector3<double> q_c = wfcpw->kvec_c[iq];
@@ -212,14 +207,14 @@ void get_exx_stress_potential(const K_Vectors* kv,
     auto param_fock = GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock];
     for (auto param: param_fock)
     {
-        double exx_div = exx_divergence(Conv_Coulomb_Pot_K::Coulomb_Type::Fock,
-                                        0.0,
-                                        kv,
-                                        wfcpw,
-                                        rhopw_dev,
-                                        tpiba,
-                                        gamma_extrapolation,
-                                        ucell_omega);
+        // double exx_div = exx_divergence(Conv_Coulomb_Pot_K::Coulomb_Type::Fock,
+        //                                 0.0,
+        //                                 kv,
+        //                                 wfcpw,
+        //                                 rhopw_dev,
+        //                                 tpiba,
+        //                                 gamma_extrapolation,
+        //                                 ucell_omega);
         double alpha = std::stod(param["alpha"]);
 
         const ModuleBase::Vector3<double> k_c = wfcpw->kvec_c[ik];

source/source_pw/module_pwdft/operator_pw/op_exx_pw.cpp

@@ -23,6 +23,11 @@
 
 namespace hamilt
 {
+template <typename T, typename Device>
+std::vector<typename GetTypeReal<T>::type> OperatorEXXPW<T, Device>::fock_div = {};
+
+template <typename T, typename Device>
+std::vector<typename GetTypeReal<T>::type> OperatorEXXPW<T, Device>::erfc_div = {};
 
 template <typename T, typename Device>
 OperatorEXXPW<T, Device>::OperatorEXXPW(const int* isk_in,
@@ -83,6 +88,31 @@ OperatorEXXPW<T, Device>::OperatorEXXPW(const int* isk_in,
     rhopw_dev->setuptransform();
     rhopw_dev->collect_local_pw();
 
+    auto param_fock = GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Fock];
+    for (auto param: param_fock)
+    {
+        fock_div.push_back(exx_divergence(Conv_Coulomb_Pot_K::Coulomb_Type::Fock,
+                                          0.0,
+                                          kv,
+                                          wfcpw,
+                                          rhopw_dev,
+                                          tpiba,
+                                          gamma_extrapolation,
+                                          ucell->omega));
+    }
+    auto param_erfc = GlobalC::exx_info.info_global.coulomb_param[Conv_Coulomb_Pot_K::Coulomb_Type::Erfc];
+    for (auto param: param_erfc)
+    {
+        erfc_div.push_back(exx_divergence(Conv_Coulomb_Pot_K::Coulomb_Type::Erfc,
+                                          std::stod(param["omega"]),
+                                          kv,
+                                          wfcpw,
+                                          rhopw_dev,
+                                          tpiba,
+                                          gamma_extrapolation,
+                                          ucell->omega));
+    }
+
 }   // end of constructor
 
 template <typename T, typename Device>

source/source_pw/module_pwdft/operator_pw/op_exx_pw.h

@@ -55,13 +55,14 @@ class OperatorEXXPW : public OperatorPW<T, Device>
 
     bool first_iter = true;
 
+    static std::vector<Real> fock_div, erfc_div;
+
   private:
     const int* isk = nullptr;
     const ModulePW::PW_Basis_K* wfcpw = nullptr;
     const ModulePW::PW_Basis* rhopw = nullptr;
     ModulePW::PW_Basis* rhopw_dev = nullptr; // for device
     const UnitCell *ucell = nullptr;
-//    Real exx_div = 0;
     Real tpiba = 0;
 
     std::vector<int> get_q_points(const int ik) const;