Merge pull request #161 from PeizeLin/develop

1. fix bug in LCAO_Hamilt::calculate_HR_exx_sparse()
abacusmodeling · Jul 11, 2022 · b3287b2 · b3287b2
2 parents fb188b3 + e227fd1
commit b3287b2
Show file tree

Hide file tree

Showing 3 changed files with 27 additions and 8 deletions.
diff --git a/source/Makefile.vars b/source/Makefile.vars
@@ -4,12 +4,12 @@ CPLUSPLUS_MPI = mpiicpc
 
 LAPACK_DIR    = $(MKLROOT)
 
-FFTW_DIR      = /home/qianrui/intelcompile/fftw_3.3.8
+FFTW_DIR      = /public/software/fftw_3.3.8
 
-ELPA_DIR      = /home/qianrui/intelcompile/elpa_21.05.002
+ELPA_DIR      = /public/software/elpa_21.05.002
 ELPA_INCLUDE_DIR = ${ELPA_DIR}/include/elpa-2021.05.002
 
-CEREAL_DIR    = /home/qianrui/headfile/cereal
+CEREAL_DIR    = /public/software/cereal
 
 # LIBXC_DIR     = /public/software/libxc-5.0.0
 

diff --git a/source/src_external/src_test/print_tmp.h b/source/src_external/src_test/print_tmp.h
@@ -23,7 +23,20 @@ static void print_tmp(const T &S)
 	}
 }
 
-void print_matrixess(const std::map<size_t,std::map<size_t,std::map<Abfs::Vector3_Order<double>,std::weak_ptr<ModuleBase::matrix>>>> &ms, std::ofstream &ofs)
+template<typename T1,typename T2,typename T3>
+void print_matrixess(const std::map<T1,std::map<T2,std::map<Abfs::Vector3_Order<T3>,ModuleBase::matrix>>> &ms, std::ofstream &ofs)
+{
+	for( const auto & c1 : ms )
+		for( const auto & c2 : c1.second )
+			for( const auto & c3 : c2.second )
+			{
+				ofs<<c1.first<<"\t"<<c2.first<<"\t"<<c3.first<<std::endl;
+				c3.second.print(ofs, 1E-10)<<std::endl;
+			}
+}
+
+template<typename T1,typename T2,typename T3>
+void print_matrixess(const std::map<T1,std::map<T2,std::map<Abfs::Vector3_Order<T3>,std::weak_ptr<ModuleBase::matrix>>>> &ms, std::ofstream &ofs)
 {
 	for( const auto & c1 : ms )
 		for( const auto & c2 : c1.second )
@@ -34,7 +47,8 @@ void print_matrixess(const std::map<size_t,std::map<size_t,std::map<Abfs::Vector
 			}
 }
 
-void print_matrixess(const std::map<size_t,std::map<size_t,std::map<Abfs::Vector3_Order<int>,std::shared_ptr<ModuleBase::matrix>>>> &ms, std::ofstream &ofs)
+template<typename T1,typename T2,typename T3>
+void print_matrixess(const std::map<T1,std::map<T2,std::map<Abfs::Vector3_Order<T3>,std::shared_ptr<ModuleBase::matrix>>>> &ms, std::ofstream &ofs)
 {
 	for( const auto & c1 : ms )
 		for( const auto & c2 : c1.second )

diff --git a/source/src_lcao/LCAO_hamilt.cpp b/source/src_lcao/LCAO_hamilt.cpp
@@ -1027,6 +1027,8 @@ void LCAO_Hamilt::calculat_HR_dftu_soc_sparse(const int &current_spin, const dou
 
 }
 
+#include "src_external/src_test/src_global/matrix-test.h"
+
 #ifdef __MPI
 // Peize Lin add 2021.11.16
 void LCAO_Hamilt::calculate_HR_exx_sparse(const int &current_spin, const double &sparse_threshold)
@@ -1040,9 +1042,10 @@ void LCAO_Hamilt::calculate_HR_exx_sparse(const int &current_spin, const double
 			+ModuleBase::GlobalFunc::TO_STRING(__FILE__)+" line "+ModuleBase::GlobalFunc::TO_STRING(__LINE__));
 	const std::vector<Abfs::Vector3_Order<int>> Rs = Abfs::get_Born_von_Karmen_boxes( Rs_period );
 
-
 	const int ik_begin = (GlobalV::NSPIN==2) ? (current_spin*GlobalC::kv.nks/2) : 0;
 	const int ik_end = (GlobalV::NSPIN==2) ? ((current_spin+1)*GlobalC::kv.nks/2) : GlobalC::kv.nks;
+	const double frac = (GlobalV::NSPIN==1) ? 0.5 : 1.0;                        // Peize Lin add 2022.07.09
+
 	for(const Abfs::Vector3_Order<int> &R : Rs)
 	{
 		ModuleBase::matrix HexxR;
@@ -1051,11 +1054,13 @@ void LCAO_Hamilt::calculate_HR_exx_sparse(const int &current_spin, const double
 			ModuleBase::matrix HexxR_tmp;
 			if(GlobalV::GAMMA_ONLY_LOCAL)
 				HexxR_tmp = GlobalC::exx_global.info.hybrid_alpha
-					* GlobalC::exx_lcao.Hexx_para.HK_Gamma_m2D[ik];
+					* GlobalC::exx_lcao.Hexx_para.HK_Gamma_m2D[ik]
+					* (GlobalC::kv.wk[ik] * frac);
 			else
 				HexxR_tmp = GlobalC::exx_global.info.hybrid_alpha
 					* (GlobalC::exx_lcao.Hexx_para.HK_K_m2D[ik]
-					* std::exp( ModuleBase::TWO_PI*ModuleBase::IMAG_UNIT * (GlobalC::kv.kvec_c[ik] * (R*GlobalC::ucell.latvec)) )).real();
+					* std::exp( ModuleBase::TWO_PI*ModuleBase::IMAG_UNIT * (GlobalC::kv.kvec_c[ik] * (R*GlobalC::ucell.latvec)) )).real()
+					* (GlobalC::kv.wk[ik] * frac);
 
 			if(HexxR.c)
 				HexxR += HexxR_tmp;