change lat0 in Matrix.cpp

A-006 · A-006 · commit 372de4c069b4 · 2024-12-08T22:16:51.000+08:00
diff --git a/source/module_ri/ABFs_Construct-PCA.cpp b/source/module_ri/ABFs_Construct-PCA.cpp
@@ -104,7 +104,7 @@ namespace PCA
 		for( std::size_t it=0; it!=abfs.size(); ++it ) {
 			delta_R[it][it] = {0.0};
 }
-		m_abfslcaos_lcaos.init_radial_table(ucell.lat0,delta_R);
+		m_abfslcaos_lcaos.init_radial_table(delta_R);
 
 		GlobalC::exx_info.info_ri.abfs_Lmax = Lmax_bak;
 		
diff --git a/source/module_ri/Matrix_Orbs11.cpp b/source/module_ri/Matrix_Orbs11.cpp
@@ -19,7 +19,7 @@ void Matrix_Orbs11::init(const int mode,
     ModuleBase::timer::tick("Matrix_Orbs11", "init");
 
     int Lmax_used, Lmax;
-
+    this->lat0 = ucell.lat0;
     const int ntype = orb.get_ntype();
     int lmax_orb = -1, lmax_beta = -1;
     for (int it = 0; it < ntype; it++)
@@ -123,8 +123,7 @@ void Matrix_Orbs11::init_radial_table()
     ModuleBase::timer::tick("Matrix_Orbs11", "init_radial_table");
 }
 
-void Matrix_Orbs11::init_radial_table(const double lat0,
-                                      const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)
+void Matrix_Orbs11::init_radial_table(const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)
 {
     ModuleBase::TITLE("Matrix_Orbs11", "init_radial_table_Rs");
     ModuleBase::timer::tick("Matrix_Orbs11", "init_radial_table");
diff --git a/source/module_ri/Matrix_Orbs11.h b/source/module_ri/Matrix_Orbs11.h
@@ -35,8 +35,7 @@ class Matrix_Orbs11
     void init_radial(const LCAO_Orbitals& orb_A, const LCAO_Orbitals& orb_B);
 
     void init_radial_table();
-    void init_radial_table(const double lat0,
-                           const std::map<size_t, std::map<size_t, std::set<double>>>& Rs); // unit: ucell.lat0
+    void init_radial_table(const std::map<size_t, std::map<size_t, std::set<double>>>& Rs); // unit: ucell.lat0
 
     enum class Matrix_Order
     {
@@ -71,7 +70,7 @@ class Matrix_Orbs11
     ModuleBase::Sph_Bessel_Recursive::D2* psb_ = nullptr;
     ORB_gaunt_table MGT;
     const double lcao_dr_ = 0.01;
-
+    double lat0 = 0.0;                                            // restore ucell.lat0
     std::map<size_t,                                              // TA
              std::map<size_t,                                     // TB
                       std::map<int,                               // LA
diff --git a/source/module_ri/Matrix_Orbs11.hpp b/source/module_ri/Matrix_Orbs11.hpp
@@ -46,7 +46,7 @@ RI::Tensor<Tdata> Matrix_Orbs11::cal_overlap_matrix(
 						const size_t NB = co6.first;
 						for( size_t MB=0; MB!=2*LB+1; ++MB )
 						{
-							const Tdata overlap = co6.second.cal_overlap( tauA*GlobalC::ucell.lat0, tauB*GlobalC::ucell.lat0, MA, MB );
+							const Tdata overlap = co6.second.cal_overlap( tauA*lat0, tauB*lat0, MA, MB );
 							const size_t iA = index_A[TA][LA][NA][MA];
 							const size_t iB = index_B[TB][LB][NB][MB];
 							switch(matrix_order)
@@ -103,7 +103,7 @@ std::array<RI::Tensor<Tdata>,3> Matrix_Orbs11::cal_grad_overlap_matrix(
 						const size_t NB = co6.first;
 						for( size_t MB=0; MB!=2*LB+1; ++MB )
 						{
-							const std::array<double,3> grad_overlap = RI_Util::Vector3_to_array3(co6.second.cal_grad_overlap( tauA*GlobalC::ucell.lat0, tauB*GlobalC::ucell.lat0, MA, MB ));
+							const std::array<double,3> grad_overlap = RI_Util::Vector3_to_array3(co6.second.cal_grad_overlap( tauA*lat0, tauB*lat0, MA, MB ));
 							const size_t iA = index_A[TA][LA][NA][MA];
 							const size_t iB = index_B[TB][LB][NB][MB];
 							for(size_t i=0; i<m.size(); ++i)
diff --git a/source/module_ri/Matrix_Orbs21.cpp b/source/module_ri/Matrix_Orbs21.cpp
@@ -21,6 +21,7 @@ void Matrix_Orbs21::init(const int mode,
 
     const int ntype = orb.get_ntype();
     int lmax_orb = -1, lmax_beta = -1;
+    this->lat0   = ucell.lat0;
     for (int it = 0; it < ntype; it++)
     {
         lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax());
@@ -120,8 +121,7 @@ void Matrix_Orbs21::init_radial_table()
     ModuleBase::timer::tick("Matrix_Orbs21", "init_radial_table");
 }
 
-void Matrix_Orbs21::init_radial_table(const double lat0,
-                                      const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)
+void Matrix_Orbs21::init_radial_table(const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)
 {
     ModuleBase::TITLE("Matrix_Orbs21", "init_radial_table_Rs");
     ModuleBase::timer::tick("Matrix_Orbs21", "init_radial_table");
diff --git a/source/module_ri/Matrix_Orbs21.h b/source/module_ri/Matrix_Orbs21.h
@@ -36,8 +36,7 @@ class Matrix_Orbs21
                      const LCAO_Orbitals& orb_B);
 
     void init_radial_table();
-    void init_radial_table(const double lat0,
-                           const std::map<size_t, std::map<size_t, std::set<double>>>& Rs); // unit: ucell.lat0
+    void init_radial_table(const std::map<size_t, std::map<size_t, std::set<double>>>& Rs); // unit: ucell.lat0
 
     enum class Matrix_Order
     {
@@ -79,7 +78,7 @@ class Matrix_Orbs21
     ModuleBase::Sph_Bessel_Recursive::D2* psb_ = nullptr;
     ORB_gaunt_table MGT;
     const double lcao_dr_ = 0.01;
-
+    double lat0 = 0;                                                                // restore ucell.lat0
     std::map<size_t,                                                                // TA
              std::map<size_t,                                                       // TB
                       std::map<int,                                                 // LA1
diff --git a/source/module_ri/Matrix_Orbs21.hpp b/source/module_ri/Matrix_Orbs21.hpp
@@ -60,7 +60,7 @@ RI::Tensor<Tdata> Matrix_Orbs21::cal_overlap_matrix(
 									const size_t NB = co8.first;
 									for( size_t MB=0; MB!=2*LB+1; ++MB )
 									{
-										const Tdata overlap = co8.second.cal_overlap( tauA*GlobalC::ucell.lat0, tauB*GlobalC::ucell.lat0, MA1, MA2, MB );
+										const Tdata overlap = co8.second.cal_overlap( tauA*lat0, tauB*lat0, MA1, MA2, MB );
 										const size_t iA1 = index_A1[TA][LA1][NA1][MA1];
 										const size_t iA2 = index_A2[TA][LA2][NA2][MA2];
 										const size_t iB = index_B[TB][LB][NB][MB];
@@ -139,7 +139,7 @@ std::array<RI::Tensor<Tdata>,3> Matrix_Orbs21::cal_grad_overlap_matrix(
 									const size_t NB = co8.first;
 									for( size_t MB=0; MB!=2*LB+1; ++MB )
 									{
-										const std::array<double,3> grad_overlap = RI_Util::Vector3_to_array3(co8.second.cal_grad_overlap( tauA*GlobalC::ucell.lat0, tauB*GlobalC::ucell.lat0, MA1, MA2, MB ));
+										const std::array<double,3> grad_overlap = RI_Util::Vector3_to_array3(co8.second.cal_grad_overlap( tauA*lat0, tauB*lat0, MA1, MA2, MB ));
 										const size_t iA1 = index_A1[TA][LA1][NA1][MA1];
 										const size_t iA2 = index_A2[TA][LA2][NA2][MA2];
 										const size_t iB = index_B[TB][LB][NB][MB];
diff --git a/source/module_ri/Matrix_Orbs22.cpp b/source/module_ri/Matrix_Orbs22.cpp
@@ -133,8 +133,7 @@ void Matrix_Orbs22::init_radial_table()
     ModuleBase::timer::tick("Matrix_Orbs22", "init_radial_table");
 }
 
-void Matrix_Orbs22::init_radial_table(const double lat0,
-                                      const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)
+void Matrix_Orbs22::init_radial_table(const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)
 {
     ModuleBase::TITLE("Matrix_Orbs22", "init_radial_table_Rs");
     ModuleBase::timer::tick("Matrix_Orbs22", "init_radial_table");
diff --git a/source/module_ri/Matrix_Orbs22.h b/source/module_ri/Matrix_Orbs22.h
@@ -39,8 +39,7 @@ class Matrix_Orbs22
                      const LCAO_Orbitals& orb_B2);
 
     void init_radial_table();
-    void init_radial_table(const double lat0,
-                           const std::map<size_t, std::map<size_t, std::set<double>>>& Rs); // unit: ucell.lat0
+    void init_radial_table(const std::map<size_t, std::map<size_t, std::set<double>>>& Rs); // unit: ucell.lat0
 
     enum class Matrix_Order
     {
@@ -103,7 +102,7 @@ class Matrix_Orbs22
     ModuleBase::Sph_Bessel_Recursive::D2* psb_ = nullptr;
     ORB_gaunt_table MGT;
     const double lcao_dr_ = 0.01;
-
+    double lat0 = 0;                                                                             // restore ucell.lat0
     std::map<
         size_t,                                                                                  // TA
         std::map<size_t,                                                                         // TB
diff --git a/source/module_ri/Matrix_Orbs22.hpp b/source/module_ri/Matrix_Orbs22.hpp
@@ -88,7 +88,7 @@ RI::Tensor<Tdata> Matrix_Orbs22::cal_overlap_matrix(
 												const size_t NB2 = co10.first;
 												for( size_t MB2=0; MB2!=2*LB2+1; ++MB2 )
 												{
-													const Tdata overlap = co10.second.cal_overlap( tauA*GlobalC::ucell.lat0, tauB*GlobalC::ucell.lat0, MA1, MA2, MB1, MB2 );
+													const Tdata overlap = co10.second.cal_overlap( tauA*lat0, tauB*lat0, MA1, MA2, MB1, MB2 );
 													const size_t iA1 = index_A1[TA][LA1][NA1][MA1];
 													const size_t iA2 = index_A2[TA][LA2][NA2][MA2];
 													const size_t iB1 = index_B1[TB][LB1][NB1][MB1];
@@ -217,10 +217,10 @@ std::array<RI::Tensor<Tdata>,3> Matrix_Orbs22::cal_grad_overlap_matrix(
 												const size_t NB2 = co10.first;
 												for( size_t MB2=0; MB2!=2*LB2+1; ++MB2 )
 												{
-													const Tdata overlap = co10.second.cal_overlap( tauA*GlobalC::ucell.lat0, tauB*GlobalC::ucell.lat0, MA1, MA2, MB1, MB2 );
+													const Tdata overlap = co10.second.cal_overlap( tauA*lat0, tauB*lat0, MA1, MA2, MB1, MB2 );
 													switch(matrix_order)
 													{										
-														const std::array<double,3> grad_overlap = RI_Util::Vector3_to_array3(co10.second.cal_grad_overlap( tauA*GlobalC::ucell.lat0, tauB*GlobalC::ucell.lat0, MA1, MA2, MB1, MB2 ));
+														const std::array<double,3> grad_overlap = RI_Util::Vector3_to_array3(co10.second.cal_grad_overlap( tauA*lat0, tauB*lat0, MA1, MA2, MB1, MB2 ));
 														const size_t iA1 = index_A1[TA][LA1][NA1][MA1];
 														const size_t iA2 = index_A2[TA][LA2][NA2][MA2];
 														const size_t iB1 = index_B1[TB][LB1][NB1][MB1];
diff --git a/source/module_ri/exx_opt_orb.cpp b/source/module_ri/exx_opt_orb.cpp
@@ -76,7 +76,7 @@ void Exx_Opt_Orb::generate_matrix(const K_Vectors &kv, const UnitCell& ucell, co
 		m_lcaoslcaos_lcaoslcaos.init( 1, ucell,orb, this->kmesh_times, 1 );
 		m_lcaoslcaos_lcaoslcaos.init_radial( lcaos, lcaos, lcaos, lcaos );
 		#if TEST_EXX_RADIAL>=1
-		m_lcaoslcaos_lcaoslcaos.init_radial_table(ucell.lat0,radial_R);
+		m_lcaoslcaos_lcaoslcaos.init_radial_table(radial_R);
 		#else
 		m_lcaoslcaos_lcaoslcaos.init_radial_table();
 		#endif
@@ -92,7 +92,7 @@ void Exx_Opt_Orb::generate_matrix(const K_Vectors &kv, const UnitCell& ucell, co
 		m_jyslcaos_lcaos.init( 1, ucell , orb, this->kmesh_times, 1 );
 		m_jyslcaos_lcaos.init_radial( jle.jle, lcaos, lcaos );
 		#if TEST_EXX_RADIAL>=1
-		m_jyslcaos_lcaos.init_radial_table(ucell.lat0, radial_R);
+		m_jyslcaos_lcaos.init_radial_table( radial_R);
 		#else
 		m_jyslcaos_lcaos.init_radial_table();
 		#endif
@@ -108,7 +108,7 @@ void Exx_Opt_Orb::generate_matrix(const K_Vectors &kv, const UnitCell& ucell, co
 		m_jys_jys.init( 2,ucell,orb, this->kmesh_times, 1 );
 		m_jys_jys.init_radial( jle.jle, jle.jle );
 		#if TEST_EXX_RADIAL>=1
-		m_jys_jys.init_radial_table(ucell.lat0,radial_R);
+		m_jys_jys.init_radial_table(radial_R);
 		#else
 		m_jys_jys.init_radial_table();
 		#endif
@@ -124,7 +124,7 @@ void Exx_Opt_Orb::generate_matrix(const K_Vectors &kv, const UnitCell& ucell, co
 		m_abfs_abfs.init( 2, ucell, orb, this->kmesh_times, 1 );
 		m_abfs_abfs.init_radial( abfs, abfs );
 		#if TEST_EXX_RADIAL>=1
-		m_abfs_abfs.init_radial_table(ucell.lat0,radial_R);
+		m_abfs_abfs.init_radial_table(radial_R);
 		#else
 		m_abfs_abfs.init_radial_table();
 		#endif
@@ -140,7 +140,7 @@ void Exx_Opt_Orb::generate_matrix(const K_Vectors &kv, const UnitCell& ucell, co
 		m_abfslcaos_lcaos.init( 1, ucell , orb, this->kmesh_times, 1 );
 		m_abfslcaos_lcaos.init_radial( abfs, lcaos, lcaos );
 		#if TEST_EXX_RADIAL>=1
-		m_abfslcaos_lcaos.init_radial_table(ucell.lat0,radial_R);
+		m_abfslcaos_lcaos.init_radial_table(radial_R);
 		#else
 		m_abfslcaos_lcaos.init_radial_table();
 		#endif
@@ -156,7 +156,7 @@ void Exx_Opt_Orb::generate_matrix(const K_Vectors &kv, const UnitCell& ucell, co
 		m_jys_abfs.init( 2, ucell,orb, this->kmesh_times, 1 );
 		m_jys_abfs.init_radial( jle.jle, abfs );
 		#if TEST_EXX_RADIAL>=1
-		m_jys_abfs.init_radial_table(ucell.lat0,radial_R);
+		m_jys_abfs.init_radial_table(radial_R);
 		#else
 		m_jys_abfs.init_radial_table();
 		#endif

Original file line number	Diff line number	Diff line change
`@@ -104,7 +104,7 @@ namespace PCA`
`104`	`104`	`for( std::size_t it=0; it!=abfs.size(); ++it ) {`
`105`	`105`	`delta_R[it][it] = {0.0};`
`106`	`106`	`}`
`107`		`- m_abfslcaos_lcaos.init_radial_table(ucell.lat0,delta_R);`
	`107`	`+ m_abfslcaos_lcaos.init_radial_table(delta_R);`
`108`	`108`
`109`	`109`	`GlobalC::exx_info.info_ri.abfs_Lmax = Lmax_bak;`
`110`	`110`
Original file line number	Diff line number	Diff line change
`@@ -21,6 +21,7 @@ void Matrix_Orbs21::init(const int mode,`
`21`	`21`
`22`	`22`	`const int ntype = orb.get_ntype();`
`23`	`23`	`int lmax_orb = -1, lmax_beta = -1;`
	`24`	`+ this->lat0 = ucell.lat0;`
`24`	`25`	`for (int it = 0; it < ntype; it++)`
`25`	`26`	`{`
`26`	`27`	`lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax());`
`@@ -120,8 +121,7 @@ void Matrix_Orbs21::init_radial_table()`
`120`	`121`	`ModuleBase::timer::tick("Matrix_Orbs21", "init_radial_table");`
`121`	`122`	`}`
`122`	`123`
`123`		`-void Matrix_Orbs21::init_radial_table(const double lat0,`
`124`		`- const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)`
	`124`	`+void Matrix_Orbs21::init_radial_table(const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)`
`125`	`125`	`{`
`126`	`126`	`ModuleBase::TITLE("Matrix_Orbs21", "init_radial_table_Rs");`
`127`	`127`	`ModuleBase::timer::tick("Matrix_Orbs21", "init_radial_table");`
Original file line number	Diff line number	Diff line change
`@@ -133,8 +133,7 @@ void Matrix_Orbs22::init_radial_table()`
`133`	`133`	`ModuleBase::timer::tick("Matrix_Orbs22", "init_radial_table");`
`134`	`134`	`}`
`135`	`135`
`136`		`-void Matrix_Orbs22::init_radial_table(const double lat0,`
`137`		`- const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)`
	`136`	`+void Matrix_Orbs22::init_radial_table(const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)`
`138`	`137`	`{`
`139`	`138`	`ModuleBase::TITLE("Matrix_Orbs22", "init_radial_table_Rs");`
`140`	`139`	`ModuleBase::timer::tick("Matrix_Orbs22", "init_radial_table");`