Feature: update LibRI v0.2.0

source/module_hamilt_general/module_xc/exx_info.h

@@ -47,6 +47,8 @@ struct Exx_Info
         double cauchy_threshold = 0;
         double C_grad_threshold = 0;
         double V_grad_threshold = 0;
+        double C_grad_R_threshold = 0;
+        double V_grad_R_threshold = 0;
         double cauchy_force_threshold = 0;
         double cauchy_stress_threshold = 0;
         double ccp_rmesh_times = 10;

source/module_io/input_conv.cpp

@@ -497,6 +497,8 @@ void Input_Conv::Convert()
         GlobalC::exx_info.info_ri.cauchy_threshold = PARAM.inp.exx_cauchy_threshold;
         GlobalC::exx_info.info_ri.C_grad_threshold = PARAM.inp.exx_c_grad_threshold;
         GlobalC::exx_info.info_ri.V_grad_threshold = PARAM.inp.exx_v_grad_threshold;
+        GlobalC::exx_info.info_ri.C_grad_R_threshold = PARAM.inp.exx_c_grad_r_threshold;
+        GlobalC::exx_info.info_ri.V_grad_R_threshold = PARAM.inp.exx_v_grad_r_threshold;
         GlobalC::exx_info.info_ri.cauchy_force_threshold = PARAM.inp.exx_cauchy_force_threshold;
         GlobalC::exx_info.info_ri.cauchy_stress_threshold = PARAM.inp.exx_cauchy_stress_threshold;
         GlobalC::exx_info.info_ri.ccp_rmesh_times = std::stod(PARAM.inp.exx_ccp_rmesh_times);

source/module_io/read_input_item_exx_dftu.cpp

@@ -149,6 +149,18 @@ void ReadInput::item_exx()
         read_sync_double(input.exx_v_grad_threshold);
         this->add_item(item);
     }
+    {
+        Input_Item item("exx_c_grad_r_threshold");
+        item.annotation = "threshold to screen nabla C matrix in exx";
+        read_sync_double(input.exx_c_grad_r_threshold);
+        this->add_item(item);
+    }
+    {
+        Input_Item item("exx_v_grad_r_threshold");
+        item.annotation = "threshold to screen nabla V matrix in exx";
+        read_sync_double(input.exx_v_grad_r_threshold);
+        this->add_item(item);
+    }
     {
         Input_Item item("exx_cauchy_force_threshold");
         item.annotation = "threshold to screen exx force using Cauchy-Schwartz inequality";

source/module_parameter/input_parameter.h

@@ -471,6 +471,8 @@ struct Input_para
     double exx_cauchy_threshold = 1e-07;        ///< threshold to screen exx using Cauchy-Schwartz inequality
     double exx_c_grad_threshold = 0.0001;       ///< threshold to screen nabla C matrix in exx
     double exx_v_grad_threshold = 0.1;          ///< threshold to screen nabla V matrix in exx
+    double exx_c_grad_r_threshold = 0.0001;     ///< threshold to screen nabla C matrix in exx
+    double exx_v_grad_r_threshold = 0.1;        ///< threshold to screen nabla V matrix in exx
     double exx_cauchy_force_threshold = 1e-07;  ///< threshold to screen exx force using Cauchy-Schwartz
                                                 ///< inequality
     double exx_cauchy_stress_threshold = 1e-07; ///< threshold to screen exx stress using Cauchy-Schwartz

source/module_ri/Exx_LRI.hpp

@@ -140,6 +140,11 @@ void Exx_LRI<Tdata>::cal_exx_ions()
 				{{"writable_dVws",true}});
 		this->cv.dVws = LRI_CV_Tools::get_dCVws(dVs);
 		this->exx_lri.set_dVs(std::move(dVs), this->info.V_grad_threshold);
+		if(GlobalV::CAL_STRESS)
+		{
+			std::array<std::array<std::map<TA,std::map<TAC,RI::Tensor<Tdata>>>,3>,3> dVRs = LRI_CV_Tools::cal_dMRs(dVs);
+			this->exx_lri.set_dVRs(std::move(dVRs), this->info.V_grad_R_threshold);
+		}
 	}
 
 	const std::array<Tcell,Ndim> period_Cs = LRI_CV_Tools::cal_latvec_range<Tcell>(2);
@@ -160,6 +165,11 @@ void Exx_LRI<Tdata>::cal_exx_ions()
 		std::array<std::map<TA,std::map<TAC,RI::Tensor<Tdata>>>,3> &dCs = std::get<1>(Cs_dCs);
 		this->cv.dCws = LRI_CV_Tools::get_dCVws(dCs);
 		this->exx_lri.set_dCs(std::move(dCs), this->info.C_grad_threshold);
+		if(GlobalV::CAL_STRESS)
+		{
+			std::array<std::array<std::map<TA,std::map<TAC,RI::Tensor<Tdata>>>,3>,3> dCRs = LRI_CV_Tools::cal_dMRs(dCs);
+			this->exx_lri.set_dCRs(std::move(dCRs), this->info.C_grad_R_threshold);
+		}
 	}
 	ModuleBase::timer::tick("Exx_LRI", "cal_exx_ions");
 }
@@ -172,8 +182,6 @@ void Exx_LRI<Tdata>::cal_exx_elec(const std::vector<std::map<TA,std::map<TAC,RI:
 
 	const std::vector<std::tuple<std::set<TA>, std::set<TA>>> judge = RI_2D_Comm::get_2D_judge(pv);
 
-	this->exx_lri.set_csm_threshold(this->info.cauchy_threshold);
-
 	this->Hexxs.resize(GlobalV::NSPIN);
 	this->Eexx = 0;
 	for(int is=0; is<GlobalV::NSPIN; ++is)
@@ -209,7 +217,7 @@ void Exx_LRI<Tdata>::post_process_Hexx( std::map<TA, std::map<TAC, RI::Tensor<Td
 }
 
 template<typename Tdata>
-double Exx_LRI<Tdata>::post_process_Eexx(const double& Eexx_in) const
+double Exx_LRI<Tdata>::post_process_Eexx(const double& Eexx_in) const
 {
 	ModuleBase::TITLE("Exx_LRI","post_process_Eexx");
     const double SPIN_multiple = std::map<int, double>{ {1,2}, {2,1}, {4,1} }.at(GlobalV::NSPIN);				// why?
@@ -239,8 +247,6 @@ void Exx_LRI<Tdata>::cal_exx_force()
 {
 	ModuleBase::TITLE("Exx_LRI","cal_exx_force");
 	ModuleBase::timer::tick("Exx_LRI", "cal_exx_force");
-
-	this->exx_lri.set_csm_threshold(this->info.cauchy_force_threshold);
 
 	this->force_exx.create(GlobalC::ucell.nat, Ndim);
 	for(int is=0; is<GlobalV::NSPIN; ++is)
@@ -263,8 +269,6 @@ void Exx_LRI<Tdata>::cal_exx_stress()
 {
 	ModuleBase::TITLE("Exx_LRI","cal_exx_stress");
 	ModuleBase::timer::tick("Exx_LRI", "cal_exx_stress");
-
-	this->exx_lri.set_csm_threshold(this->info.cauchy_stress_threshold);
 
 	this->stress_exx.create(Ndim, Ndim);
 	for(int is=0; is<GlobalV::NSPIN; ++is)

source/module_ri/LRI_CV.hpp

@@ -8,10 +8,11 @@
 
 #include "LRI_CV.h"
 #include "LRI_CV_Tools.h"
-#include "module_ri/exx_abfs-abfs_index.h"
+#include "exx_abfs-abfs_index.h"
 #include "RI_Util.h"
-#include "module_base/tool_title.h"
-#include "module_base/timer.h"
+#include "../module_base/tool_title.h"
+#include "../module_base/timer.h"
+#include "../module_hamilt_pw/hamilt_pwdft/global.h"
 #include <RI/global/Global_Func-1.h>
 #include <omp.h>
 

source/module_ri/LRI_CV_Tools.h

@@ -6,11 +6,14 @@
 #ifndef LRI_CV_TOOLS_H
 #define LRI_CV_TOOLS_H
 
+#include "module_base/abfs-vector3_order.h"
+
 #include <RI/global/Tensor.h>
 
 #include <cstddef>
 #include <array>
 #include <vector>
+#include <map>
 
 namespace LRI_CV_Tools
 {
@@ -84,7 +87,12 @@ namespace LRI_CV_Tools
 	template<typename TA, typename Tcell, typename Tdata>
 	extern std::map<int,std::map<int,std::map<Abfs::Vector3_Order<double>,std::array<RI::Tensor<Tdata>,3>>>>
 	get_dCVws(
-		const std::array<std::map<TA,std::map<std::pair<TA,std::array<Tcell,3>>,RI::Tensor<Tdata>>>,3> &dCVs);	
+		const std::array<std::map<TA,std::map<std::pair<TA,std::array<Tcell,3>>,RI::Tensor<Tdata>>>,3> &dCVs);
+
+	template<typename TA, typename TC, typename Tdata>
+	extern std::array<std::array<std::map<TA,std::map<std::pair<TA,TC>,RI::Tensor<Tdata>>>,3>,3>
+	cal_dMRs(
+		const std::array<std::map<TA,std::map<std::pair<TA,TC>,RI::Tensor<Tdata>>>,3> &dMs);
 }
 
 #include "LRI_CV_Tools.hpp"

source/module_ri/LRI_CV_Tools.hpp

@@ -9,6 +9,7 @@
 #include "LRI_CV_Tools.h"
 #include "Inverse_Matrix.h"
 #include "../module_base/mathzone.h"
+#include "../module_hamilt_pw/hamilt_pwdft/global.h"
 
 template<typename Tdata>
 RI::Tensor<Tdata>
@@ -54,8 +55,10 @@ template<typename T, std::size_t N>
 bool LRI_CV_Tools::exist(const std::array<T,N> &dV)
 {
 	for(size_t i=0; i<3; ++i)
+	{
 		if(!dV[i].empty())
 			return true;
+	}
 	return false;
 }
 
@@ -312,4 +315,45 @@ LRI_CV_Tools::get_dCVws(
 	return dCVws;
 }
 
+
+// dMRs[ipos0][ipos1] = \nabla_{ipos0} M R_{ipos1}
+template<typename TA, typename TC, typename Tdata>
+std::array<std::array<std::map<TA,std::map<std::pair<TA,TC>,RI::Tensor<Tdata>>>,3>,3>
+LRI_CV_Tools::cal_dMRs(
+	const std::array<std::map<TA,std::map<std::pair<TA,TC>,RI::Tensor<Tdata>>>,3> &dMs)
+{
+	auto get_R_delta = [&](const TA &iat0, const std::pair<TA,TC> &A1) -> std::array<Tdata,3>
+	{
+		const TA iat1 = A1.first;
+		const TC &cell1 = A1.second;
+		const int it0 = GlobalC::ucell.iat2it[iat0];
+		const int ia0 = GlobalC::ucell.iat2ia[iat0];
+		const int it1 = GlobalC::ucell.iat2it[iat1];
+		const int ia1 = GlobalC::ucell.iat2ia[iat1];
+		const ModuleBase::Vector3<double> tau0 = GlobalC::ucell.atoms[it0].tau[ia0];
+		const ModuleBase::Vector3<double> tau1 = GlobalC::ucell.atoms[it1].tau[ia1];
+		const Abfs::Vector3_Order<double> R_delta = -tau0+tau1+(RI_Util::array3_to_Vector3(cell1)*GlobalC::ucell.latvec);
+		return std::array<Tdata,3>{R_delta.x, R_delta.y, R_delta.z};
+	};
+	constexpr int Npos = 3;
+	std::array<std::array<std::map<TA,std::map<std::pair<TA,TC>,RI::Tensor<Tdata>>>,Npos>,Npos> dMRs;
+	for(int ipos0=0; ipos0<Npos; ++ipos0)
+	{
+		for(int ipos1=0; ipos1<Npos; ++ipos1)
+		{
+			for(const auto &dMs_A : dMs[ipos0])
+			{
+				const TA iat0 = dMs_A.first;
+				for(const auto &dMs_B : dMs_A.second)
+				{
+					const std::pair<TA,TC> A1 = dMs_B.first;
+					const RI::Tensor<Tdata> &dM = dMs_B.second;
+					const std::array<Tdata,3> R_delta = get_R_delta(iat0, A1);
+					dMRs[ipos0][ipos1][iat0][A1] = dM * R_delta[ipos1];
+				}
+			}
+		}
+	}
+	return dMRs;
+}
 #endif

tests/integrate/207_NO_KP_OXC/k-1-Vxc.ref

@@ -1,6 +1,6 @@
-6 (-1.039229,-3.6984385e-16) (-1.4337629e-08,8.1948535e-16) (-1.6549262e-15,-2.4275964e-17) (-7.3477682e-14,-2.8317688e-17) (-2.2743308e-08,6.0642134e-16) (-4.7365417e-10,2.7036238e-14)
- (-0.96491479,1.2723875e-15) (-6.2469742e-09,-1.8654043e-12) (-8.1614632e-11,8.9132756e-13) (0.12279288,7.8499298e-15) (9.3436543e-14,2.364213e-16)
- (-1.064337,2.3693393e-16) (1.0547119e-15,1.5451914e-16) (-6.6444545e-09,2.0880435e-12) (-5.7610897e-10,-3.0200736e-12)
- (-1.064337,2.3700845e-16) (6.3717573e-10,-9.9770901e-13) (4.9778198e-10,-3.4757565e-13)
- (-0.65449032,1.2140463e-16) (1.8561136e-14,4.810951e-17)
- (-0.65138796,-7.4288639e-18)
+6 (-1.0394425,-3.7363449e-16) (3.9703113e-10,8.3455762e-16) (4.4417161e-14,5.5568749e-17) (-2.4341466e-13,4.3883506e-17) (-2.0039369e-09,5.8451215e-16) (1.9596449e-11,6.8294051e-12)
+ (-0.96326621,1.2640312e-15) (-1.4831543e-09,-4.0086346e-12) (8.370459e-10,-5.1125035e-13) (0.12202161,-5.3763309e-15) (1.7561473e-14,3.428681e-15)
+ (-1.0644129,2.3743998e-16) (-7.1678774e-15,-2.6869918e-16) (-1.4937099e-09,4.0465037e-12) (2.5345497e-09,8.0997135e-10)
+ (-1.0644129,2.3735151e-16) (8.1785897e-10,5.1608006e-13) (-7.4156745e-10,-2.3462101e-10)
+ (-0.65468695,1.1556487e-16) (4.4158931e-15,1.0408529e-15)
+ (-0.65149008,-4.4004079e-17)

tests/integrate/207_NO_KP_OXC/kinetic_out.ref

@@ -1,15 +1,15 @@
 2
 1
 6
-7.5338982275919106e-02	2.0500788809470163e+00
-6.9290936005296899e-01	1.8855030988507920e+01
-6.9290936005302095e-01	1.8855030988509334e+01
-6.9295771489195779e-01	1.8856346791180162e+01
-6.9099816910638534e-01	1.8803024814828817e+01
-6.8997237978555814e-01	1.8775111655407219e+01
-2.1196370472037387e-01	5.7678283067731630e+00
-2.4311138823535741e-01	6.6154002573820518e+00
-6.0355581862793817e-01	1.6423596388789001e+01
-6.0355581862806040e-01	1.6423596388792330e+01
-6.4467362246896842e-01	1.7542469231757597e+01
-7.5610037427593779e-01	2.0574546700170757e+01
+7.5360628562083129e-02	2.0506679066117548e+00
+6.9287220698058305e-01	1.8854020001542612e+01
+6.9287220698057750e-01	1.8854020001542459e+01
+6.9314985795349449e-01	1.8861575272116287e+01
+6.9028764642341500e-01	1.8783690500735531e+01
+6.8994928867395522e-01	1.8774483314025311e+01
+2.1202140024261817e-01	5.7693982824763790e+00
+2.4283946374398946e-01	6.6080008123653400e+00
+6.0350742927314394e-01	1.6422279646893511e+01
+6.0350742927330381e-01	1.6422279646897863e+01
+6.4498868732620107e-01	1.7551042586353439e+01
+7.5610285478687722e-01	2.0574614198336214e+01

tests/integrate/207_NO_KP_OXC/result.ref

@@ -1,9 +1,9 @@
-etotref -212.6300856413655
-etotperatomref -106.3150428207
+etotref -212.6220609473505
+etotperatomref -106.3110304737
 CompareVXC_pass 0
 CompareOrbXC_pass 0
 CompareOrbKinetic_pass 0
 CompareOrbVL_pass 0
 CompareOrbVNL_pass 0
 CompareOrbVHartree_pass 0
-totaltimeref 22.06
+totaltimeref 13.12

tests/integrate/207_NO_KP_OXC/vhartree_out.ref

@@ -1,15 +1,15 @@
 2
 1
 6
-1.4113685655409342e-01	3.8405308938886313e+00
-2.5168586236213109e-01	6.8487236683374446e+00
-2.5168586236216250e-01	6.8487236683382990e+00
-2.3818785304408960e-01	6.4814239915725276e+00
-9.6061104054610796e-02	2.6139567426272201e+00
-1.0497474887605052e-01	2.8565094616667657e+00
-1.8354710156875595e-01	4.9945728654396397e+00
-1.0537573740554317e-01	2.8674209193342479e+00
-2.1880380947272768e-01	5.9539571058709440e+00
-2.1880380947281258e-01	5.9539571058732541e+00
-9.8617011618821054e-02	2.6835065554963409e+00
-9.3894268939185510e-02	2.5549941342346769e+00
+1.4112466165742804e-01	3.8401990537262907e+00
+2.5163564025411278e-01	6.8473570546682039e+00
+2.5163564025411655e-01	6.8473570546683060e+00
+2.3821003729805079e-01	6.4820276560920300e+00
+9.6583832464817498e-02	2.6281809123978053e+00
+1.0492875802915166e-01	2.8552579865194256e+00
+1.8356635884948339e-01	4.9950968829313966e+00
+1.0498309602344302e-01	2.8567365991999334e+00
+2.1873498376172118e-01	5.9520842621937646e+00
+2.1873498376183687e-01	5.9520842621969123e+00
+9.8983158673885074e-02	2.6934699280059218e+00
+9.3899408260331163e-02	2.5551339823375425e+00

tests/integrate/207_NO_KP_OXC/vpp_local_out.ref

@@ -1,15 +1,15 @@
 2
 1
 6
--3.0098788292896494e-01	-8.1903004735817060e+00
--4.6431550238352493e-01	-1.2634673004297040e+01
--4.6431550238366942e-01	-1.2634673004300973e+01
--4.6265884371492372e-01	-1.2589593009228901e+01
--2.5904047415029896e-01	-7.0488529221315490e+00
--2.5489046490776013e-01	-6.9359253772291645e+00
--3.9897177020111796e-01	-1.0856578831763620e+01
--1.9304303478937274e-01	-5.2529704646953981e+00
--3.9905222628974663e-01	-1.0858768154251907e+01
--3.9905222629003950e-01	-1.0858768154259876e+01
--2.9572559706291357e-01	-8.0471063290153779e+00
--2.2291264124345445e-01	-6.0657641542815712e+00
+-3.0101836258858156e-01	-8.1911298676694777e+00
+-4.6425055402880600e-01	-1.2632905668897235e+01
+-4.6425055402883653e-01	-1.2632905668898067e+01
+-4.6307169914994006e-01	-1.2600827381961883e+01
+-2.6013937714472535e-01	-7.0787556066784711e+00
+-2.5486302990834309e-01	-6.9351788325957680e+00
+-3.9905184265146898e-01	-1.0858757714918813e+01
+-1.9159273894286299e-01	-5.2135058900988662e+00
+-3.9895346713212027e-01	-1.0856080779705110e+01
+-3.9895346713250973e-01	-1.0856080779715708e+01
+-2.9660581100861438e-01	-8.0710581792565659e+00
+-2.2304717215673603e-01	-6.0694249282371180e+00

tests/integrate/207_NO_KP_OXC/vpp_nonlocal_out.ref

@@ -1,15 +1,15 @@
 2
 1
 6
-2.8613596948111919e-01	7.7861591753946486e+00
-2.2849804473444624e-01	6.2177507804947316e+00
-2.2849804473451799e-01	6.2177507804966838e+00
-2.2695283203968664e-01	6.1757033859534012e+00
-1.9157221661016460e-01	5.2129474487769665e+00
-1.8886548044839890e-01	5.1392933792116402e+00
-3.1324706939237390e-01	8.5238900510753659e+00
-2.1372431018924415e-01	5.8157368393863580e+00
-2.2455517390069851e-01	6.1104597608607722e+00
-2.2455517390081420e-01	6.1104597608639200e+00
-2.3364137951596947e-01	6.3577080999953335e+00
-1.3061189823674019e-01	3.5541320852316391e+00
+2.8613536205720358e-01	7.7861426465419417e+00
+2.2848104040381553e-01	6.2172880689202241e+00
+2.2848104040385026e-01	6.2172880689211691e+00
+2.2710998981794001e-01	6.1799798684919338e+00
+1.9348330177195094e-01	5.2649507439040590e+00
+1.8886683703805932e-01	5.1393302939100991e+00
+3.1324386694132939e-01	8.5238029079118220e+00
+2.1271941813565287e-01	5.7883923237788322e+00
+2.2452501763844440e-01	6.1096391668666952e+00
+2.2452501763860333e-01	6.1096391668710206e+00
+2.3398759639736544e-01	6.3671291446568841e+00
+1.3072279135213710e-01	3.5571496417083783e+00

tests/integrate/207_NO_KP_OXC/vxc_out.ref

@@ -1,15 +1,15 @@
 2
 1
 6
--5.1230914735474076e-01	-1.3940647083092204e+01
--5.4586414912548942e-01	-1.4853725524056747e+01
--5.4586414912542969e-01	-1.4853725524055122e+01
--5.2184281554508472e-01	-1.4200071503552257e+01
--3.3292737504158443e-01	-9.0594186414970697e+00
--3.3960159355658315e-01	-9.2410334444992319e+00
--5.1961448164761326e-01	-1.4139435427447937e+01
--4.8245729881713090e-01	-1.3128336611203281e+01
--5.3216848113640292e-01	-1.4481047278921190e+01
--5.3216848113629889e-01	-1.4481047278918359e+01
--3.2724521859846717e-01	-8.9047992323894558e+00
--3.2569397421272883e-01	-8.8625877071163526e+00
+-5.1241966263460559e-01	-1.3943654358136657e+01
+-5.4586098339243305e-01	-1.4853639380040919e+01
+-5.4586098339241196e-01	-1.4853639380040345e+01
+-5.2166040137040381e-01	-1.4195107759209002e+01
+-3.3409473287713004e-01	-9.0911840778338053e+00
+-3.3961143930774479e-01	-9.2413013611330097e+00
+-5.1972127113784583e-01	-1.4142341318555294e+01
+-4.8163310485213451e-01	-1.3105909142840954e+01
+-5.3220647068527449e-01	-1.4482081027579396e+01
+-5.3220647068497240e-01	-1.4482081027571175e+01
+-3.2734347697981392e-01	-8.9074729801146013e+00
+-3.2574503984024938e-01	-8.8639772741288034e+00

tests/integrate/281_NO_KP_HSE/result.ref

@@ -1,5 +1,5 @@
-etotref -205.1553398987854848
-etotperatomref -102.5776699494
-totalforceref 6.064356
-totalstressref 4166.698320
-totaltimeref 29.87
+etotref -205.1574474096462097
+etotperatomref -102.5787237048
+totalforceref 6.063234
+totalstressref 2214.189951
+totaltimeref 17.99