change ucell in module_ri/exx_lri.hpp

Author: A-006

source/module_hamilt_lcao/hamilt_lcaodft/FORCE_STRESS.cpp

@@ -375,25 +375,25 @@ void Force_Stress_LCAO<T>::getForceStress(const bool isforce,
         {
             if (GlobalC::exx_info.info_ri.real_number)
             {
-                exx_lri_double.cal_exx_force();
+                exx_lri_double.cal_exx_force(ucell.nat);
                 force_exx = GlobalC::exx_info.info_global.hybrid_alpha * exx_lri_double.force_exx;
             }
             else
             {
-                exx_lri_complex.cal_exx_force();
+                exx_lri_complex.cal_exx_force(ucell.nat);
                 force_exx = GlobalC::exx_info.info_global.hybrid_alpha * exx_lri_complex.force_exx;
             }
         }
         if (isstress)
         {
             if (GlobalC::exx_info.info_ri.real_number)
             {
-                exx_lri_double.cal_exx_stress();
+                exx_lri_double.cal_exx_stress(ucell.omega,ucell.lat0);
                 stress_exx = GlobalC::exx_info.info_global.hybrid_alpha * exx_lri_double.stress_exx;
             }
             else
             {
-                exx_lri_complex.cal_exx_stress();
+                exx_lri_complex.cal_exx_stress(ucell.omega,ucell.lat0);
                 stress_exx = GlobalC::exx_info.info_global.hybrid_alpha * exx_lri_complex.stress_exx;
             }
         }

source/module_ri/Exx_LRI.h

@@ -60,8 +60,8 @@ class Exx_LRI
 			  const UnitCell &ucell,
 			  const K_Vectors &kv_in, 
 			  const LCAO_Orbitals& orb);
-	void cal_exx_force();
-    void cal_exx_stress();
+	void cal_exx_force(const int& nat);
+    void cal_exx_stress(const double& omega, const double& lat0);
 	void cal_exx_ions(const UnitCell& ucell, const bool write_cv = false);
     void cal_exx_elec(const std::vector<std::map<TA, std::map<TAC, RI::Tensor<Tdata>>>>& Ds,
 		const UnitCell& ucell,

source/module_ri/Exx_LRI.hpp

@@ -219,7 +219,7 @@ void Exx_LRI<Tdata>::cal_exx_elec(const std::vector<std::map<TA, std::map<TAC, R
             // reduce but not repeat
             auto Hs_a2D = this->exx_lri.post_2D.set_tensors_map2(this->exx_lri.Hs);
             // rotate locally without repeat
-            Hs_a2D = p_symrot->restore_HR(GlobalC::ucell.symm, GlobalC::ucell.atoms, GlobalC::ucell.st, 'H', Hs_a2D);
+            Hs_a2D = p_symrot->restore_HR(ucell.symm, ucell.atoms, ucell.st, 'H', Hs_a2D);
             // cal energy using full Hs without repeat
             this->exx_lri.energy = this->exx_lri.post_2D.cal_energy(
                 this->exx_lri.post_2D.saves["Ds_" + suffix],
@@ -274,20 +274,20 @@ post_process_old
 */
 
 template<typename Tdata>
-void Exx_LRI<Tdata>::cal_exx_force()
+void Exx_LRI<Tdata>::cal_exx_force(const int& nat)
 {
 	ModuleBase::TITLE("Exx_LRI","cal_exx_force");
 	ModuleBase::timer::tick("Exx_LRI", "cal_exx_force");
 
-	this->force_exx.create(GlobalC::ucell.nat, Ndim);
+	this->force_exx.create(nat, Ndim);
 	for(int is=0; is<PARAM.inp.nspin; ++is)
 	{
 		this->exx_lri.cal_force({"","",std::to_string(is),"",""});
 		for(std::size_t idim=0; idim<Ndim; ++idim) {
 			for(const auto &force_item : this->exx_lri.force[idim]) {
 				this->force_exx(force_item.first, idim) += std::real(force_item.second);
-}
-}
+			}
+		}
 	}
 
 	const double SPIN_multiple = std::map<int,double>{{1,2}, {2,1}, {4,1}}.at(PARAM.inp.nspin);				// why?
@@ -298,7 +298,7 @@ void Exx_LRI<Tdata>::cal_exx_force()
 
 
 template<typename Tdata>
-void Exx_LRI<Tdata>::cal_exx_stress()
+void Exx_LRI<Tdata>::cal_exx_stress(const double& omega, const double& lat0)
 {
 	ModuleBase::TITLE("Exx_LRI","cal_exx_stress");
 	ModuleBase::timer::tick("Exx_LRI", "cal_exx_stress");
@@ -310,12 +310,12 @@ void Exx_LRI<Tdata>::cal_exx_stress()
 		for(std::size_t idim0=0; idim0<Ndim; ++idim0) {
 			for(std::size_t idim1=0; idim1<Ndim; ++idim1) {
 				this->stress_exx(idim0,idim1) += std::real(this->exx_lri.stress(idim0,idim1));
-}
-}
+		}
+	}
 	}
 
 	const double SPIN_multiple = std::map<int,double>{{1,2}, {2,1}, {4,1}}.at(PARAM.inp.nspin);				// why?
-	const double frac = 2 * SPIN_multiple / GlobalC::ucell.omega * GlobalC::ucell.lat0;		// why?
+	const double frac = 2 * SPIN_multiple / omega * lat0;		// why?
 	this->stress_exx *= frac;
 
 	ModuleBase::timer::tick("Exx_LRI", "cal_exx_stress");