From 6c6c8711894f6379d47448f71fd718170995d2fe Mon Sep 17 00:00:00 2001 From: maki49 <1579492865@qq.com> Date: Mon, 9 Nov 2026 10:34:23 -0500 Subject: [PATCH 1/2] Separate GPU symmetry filtering from CPU loop conventions --- include/RI/physics/Exx.hpp | 8 +- include/RI/physics/GW.hpp | 8 +- include/RI/physics/RPA.hpp | 8 +- .../physics/symmetry/Filter_Atom_Symmetry.h | 21 ------ .../symmetry/Filter_Atom_Symmetry_GPU.h | 74 +++++++++++++++++++ include/RI/ri/LRI-cal_loop3_CPU.hpp | 7 ++ .../LRI-cal_loop3_CPU_fine_grained_lock.hpp | 7 ++ include/RI/ri/LRI-cal_loop3_GPU.hpp | 59 ++++++++------- include/RI/ri/LRI.h | 3 + unittests/gpu_symmetry/README.md | 40 ++++++++++ .../gpu_symmetry/test_lri_gpu_symmetry.cpp | 13 +++- 11 files changed, 179 insertions(+), 69 deletions(-) create mode 100644 include/RI/physics/symmetry/Filter_Atom_Symmetry_GPU.h diff --git a/include/RI/physics/Exx.hpp b/include/RI/physics/Exx.hpp index 43b505c..42ff97d 100644 --- a/include/RI/physics/Exx.hpp +++ b/include/RI/physics/Exx.hpp @@ -9,7 +9,6 @@ #include "../ri/Cell_Nearest.h" #include "../ri/Label.h" #include "../global/Map_Operator.h" -#include "./symmetry/Filter_Atom_Symmetry.h" #include @@ -41,11 +40,8 @@ void Exx::set_symmetry( const bool flag_symmetry, const std::map, std::set> &irreducible_sector) { - if(flag_symmetry) - this->lri.filter_atom = std::make_shared>( - this->period, irreducible_sector); - else - this->lri.filter_atom = std::make_shared>(); + this->lri.flag_symmetry = flag_symmetry; + this->lri.irreducible_sector = irreducible_sector; } template diff --git a/include/RI/physics/GW.hpp b/include/RI/physics/GW.hpp index 815d5f7..8733ee7 100644 --- a/include/RI/physics/GW.hpp +++ b/include/RI/physics/GW.hpp @@ -6,7 +6,6 @@ #include "GW.h" #include "../ri/Label.h" -#include "./symmetry/Filter_Atom_Symmetry.h" namespace RI { @@ -36,11 +35,8 @@ void GW::set_symmetry( const bool flag_symmetry, const std::map, std::set> &irreducible_sector) { - if(flag_symmetry) - this->lri.filter_atom = std::make_shared>( - this->period, irreducible_sector); - else - this->lri.filter_atom = std::make_shared>(); + this->lri.flag_symmetry = flag_symmetry; + this->lri.irreducible_sector = irreducible_sector; } template diff --git a/include/RI/physics/RPA.hpp b/include/RI/physics/RPA.hpp index 4f280f5..454e5cf 100644 --- a/include/RI/physics/RPA.hpp +++ b/include/RI/physics/RPA.hpp @@ -7,7 +7,6 @@ #include "RPA.h" #include "../ri/Label.h" -#include "./symmetry/Filter_Atom_Symmetry.h" #include @@ -37,11 +36,8 @@ void RPA::set_symmetry( const bool flag_symmetry, const std::map, std::set> &irreducible_sector) { - if(flag_symmetry) - this->lri.filter_atom = std::make_shared>( - this->period, irreducible_sector); - else - this->lri.filter_atom = std::make_shared>(); + this->lri.flag_symmetry = flag_symmetry; + this->lri.irreducible_sector = irreducible_sector; } template diff --git a/include/RI/physics/symmetry/Filter_Atom_Symmetry.h b/include/RI/physics/symmetry/Filter_Atom_Symmetry.h index 77a0f1f..12698ff 100644 --- a/include/RI/physics/symmetry/Filter_Atom_Symmetry.h +++ b/include/RI/physics/symmetry/Filter_Atom_Symmetry.h @@ -45,14 +45,6 @@ class Filter_Atom_Symmetry: public Filter_Atom> case Label::ab_ab::a1b2_a2b0: case Label::ab_ab::a1b2_a2b1: return !this->symmetry.is_Aa_in_irreducible_sector(A1); - // GPU visits Aa01 before the output Aa2 is known in these cases. - case Label::ab_ab::a0b0_a1b1: - case Label::ab_ab::a0b1_a1b0: - case Label::ab_ab::a0b0_a1b2: - case Label::ab_ab::a0b1_a1b2: - case Label::ab_ab::a0b2_a1b0: - case Label::ab_ab::a0b2_a1b1: - return false; default: throw std::invalid_argument("label "+Label_Tools::get_name(label)+" error in "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } @@ -70,11 +62,6 @@ class Filter_Atom_Symmetry: public Filter_Atom> case Label::ab_ab::a0b1_a1b2: case Label::ab_ab::a0b2_a1b0: case Label::ab_ab::a0b2_a1b1: - // a01b2_a2b01: Ab01 is the output right atom on GPU. - case Label::ab_ab::a0b2_a2b0: - case Label::ab_ab::a0b2_a2b1: - case Label::ab_ab::a1b2_a2b0: - case Label::ab_ab::a1b2_a2b1: return !this->symmetry.is_Ab_in_irreducible_sector(A1.first); default: throw std::invalid_argument("label "+Label_Tools::get_name(label)+" error in "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); @@ -113,14 +100,6 @@ class Filter_Atom_Symmetry: public Filter_Atom> case Label::ab_ab::a1b0_a2b2: case Label::ab_ab::a1b1_a2b2: return !this->symmetry.in_irreducible_sector(A1, A3); - // GPU output is (Aa2, Ab2) or (Aa2, Ab01); Aa01 is intermediate. - case Label::ab_ab::a0b0_a1b1: - case Label::ab_ab::a0b1_a1b0: - case Label::ab_ab::a0b0_a1b2: - case Label::ab_ab::a0b1_a1b2: - case Label::ab_ab::a0b2_a1b0: - case Label::ab_ab::a0b2_a1b1: - return !this->symmetry.in_irreducible_sector(A2, A3); // a01b2_a2b01: Aa01, Ab01, Aa2 default: return false; diff --git a/include/RI/physics/symmetry/Filter_Atom_Symmetry_GPU.h b/include/RI/physics/symmetry/Filter_Atom_Symmetry_GPU.h new file mode 100644 index 0000000..2275216 --- /dev/null +++ b/include/RI/physics/symmetry/Filter_Atom_Symmetry_GPU.h @@ -0,0 +1,74 @@ +#pragma once + +#include "../../ri/Filter_Atom.h" +#include "../../ri/Label_Tools.h" +#include "Symmetry_Filter.h" + +namespace RI +{ + +// Arguments follow the GPU loop nesting, from outermost to innermost. +// Only output atoms may be screened; intermediate atoms are summed over. +template +class Filter_Atom_Symmetry_GPU: public Filter_Atom> +{ + public: + using TAC = std::pair; + using Group = Label::Aab_Aab; + + Filter_Atom_Symmetry_GPU( + const TC &period, + const std::map, std::set> &irsec) + : symmetry(period, irsec) {} + + bool filter_for1(const Label::ab_ab &label, const TA &A1) const override + { + switch (Label_Tools::to_Aab_Aab(label)) + { + case Group::a01b01_a01b01: // Aa01 is intermediate. + case Group::a01b01_a01b2: + return false; + default: // Aa01 is the left output atom. + return !symmetry.is_Aa_in_irreducible_sector(A1); + } + } + + bool filter_for1(const Label::ab_ab &, const TAC &A1) const override + { + // The first loop with a cell-bearing atom visits output Ab01. + return !symmetry.is_Ab_in_irreducible_sector(A1.first); + } + + bool filter_for2(const Label::ab_ab &label, const TA &A1, const TAC &A2) const override + { + // Only the final stage of a01b2_a2b01 has both outputs at level 2. + if (Label_Tools::to_Aab_Aab(label) == Group::a01b2_a2b01) + return !symmetry.in_irreducible_sector(A1, A2); + return false; + } + + bool filter_for32(const Label::ab_ab &label, + const TA &A1, const TAC &A2, const TAC &A3) const override + { + switch (Label_Tools::to_Aab_Aab(label)) + { + case Group::a01b01_a01b01: // Aa01, Aa2, Ab2 + return !symmetry.in_irreducible_sector(A2, A3); + case Group::a01b01_a01b2: // Aa01, Ab01, Aa2 + return !symmetry.in_irreducible_sector(A3, A2); + case Group::a01b01_a2b01: // Aa01, Ab01, Ab2 + case Group::a01b2_a2b01: // Aa01, Aa2, Ab01 + return !symmetry.in_irreducible_sector(A1, A3); + case Group::a01b01_a2b2: // Aa01, Ab01, Ab2 + return !symmetry.in_irreducible_sector(A1, A2); + default: + throw std::invalid_argument("Unknown GPU symmetry label"); + } + } + + // filter_for31 builds intermediates and deliberately never screens pairs. + // All other callbacks retain the base class's no-op behavior. + Symmetry_Filter symmetry; +}; + +} diff --git a/include/RI/ri/LRI-cal_loop3_CPU.hpp b/include/RI/ri/LRI-cal_loop3_CPU.hpp index e702ebb..f618e25 100644 --- a/include/RI/ri/LRI-cal_loop3_CPU.hpp +++ b/include/RI/ri/LRI-cal_loop3_CPU.hpp @@ -6,6 +6,7 @@ #pragma once #include "LRI.h" +#include "../physics/symmetry/Filter_Atom_Symmetry.h" #include "LRI_Cal_Aux.h" #include "../global/Array_Operator.h" #include "../global/Tensor_Multiply.h" @@ -24,6 +25,12 @@ void LRI::cal_loop3_CPU( std::map>> &Ds_result, const Tdata fac_add_Ds) { + if (this->flag_symmetry) + this->filter_atom = std::make_shared>( + this->period, this->irreducible_sector); + else + this->filter_atom = std::make_shared>(); + using namespace Array_Operator; const Data_Pack_Wrapper data_wrapper(this->data_pool, this->data_ab_name); diff --git a/include/RI/ri/LRI-cal_loop3_CPU_fine_grained_lock.hpp b/include/RI/ri/LRI-cal_loop3_CPU_fine_grained_lock.hpp index ac08b09..308fa87 100644 --- a/include/RI/ri/LRI-cal_loop3_CPU_fine_grained_lock.hpp +++ b/include/RI/ri/LRI-cal_loop3_CPU_fine_grained_lock.hpp @@ -6,6 +6,7 @@ #pragma once #include "LRI.h" +#include "../physics/symmetry/Filter_Atom_Symmetry.h" #include "LRI_Cal_Aux.h" #include "../global/Array_Operator.h" #include "../global/Tensor_Multiply.h" @@ -24,6 +25,12 @@ void LRI::cal_loop3_CPU_fine_grained_lock( std::map>> &Ds_result, const Tdata fac_add_Ds) { + if (this->flag_symmetry) + this->filter_atom = std::make_shared>( + this->period, this->irreducible_sector); + else + this->filter_atom = std::make_shared>(); + using namespace Array_Operator; const Data_Pack_Wrapper data_wrapper(this->data_pool, this->data_ab_name); diff --git a/include/RI/ri/LRI-cal_loop3_GPU.hpp b/include/RI/ri/LRI-cal_loop3_GPU.hpp index 5ca1304..24090c2 100644 --- a/include/RI/ri/LRI-cal_loop3_GPU.hpp +++ b/include/RI/ri/LRI-cal_loop3_GPU.hpp @@ -8,6 +8,7 @@ #ifdef __GPU_RI #include "LRI.h" +#include "../physics/symmetry/Filter_Atom_Symmetry_GPU.h" #include "LRI_Cal_Aux.h" #include "../global/Array_Operator.h" #include "../global/Map_Operator.h" @@ -35,6 +36,12 @@ void LRI::cal_loop3_GPU( std::map>> &Ds_result, const Tdata fac_add_Ds) { + if (this->flag_symmetry) + this->filter_atom = std::make_shared>( + this->period, this->irreducible_sector); + else + this->filter_atom = std::make_shared>(); + using namespace Array_Operator; constexpr std::size_t memory_limit = std::size_t(5)*1024*1024*1024; @@ -378,7 +385,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia01= 0; ia01filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; + if(this->filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; const Tensor &D_a1b2 = tools.get_Ds_ab(Label::ab::a1b2, Aa01, Ab2); if(D_a1b2.empty()) continue; @@ -431,7 +438,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia2 = 0; ia2 < list_Aa2.size(); ++ia2) // F { const TAC &Aa2 = list_Aa2[ia2]; - if (this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; const Tensor &D_a = Global_Func::find(Ds_a_transpose, Aa01, Aa2); if (D_a.empty()) continue; @@ -518,7 +525,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia01= 0; ia01filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; + if(this->filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; const Tensor &D_a1b2 = tools.get_Ds_ab(Label::ab::a1b2, Aa01, Ab2); if(D_a1b2.empty()) continue; @@ -571,7 +578,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia2 = 0; ia2 < list_Aa2.size(); ++ia2) // F { const TAC &Aa2 = list_Aa2[ia2]; - if (this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; const Tensor &D_a = tools.get_Ds_ab(Label::ab::a, Aa01, Aa2); if (D_a.empty()) continue; @@ -658,7 +665,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia01= 0; ia01filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; + if(this->filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; const Tensor &D_a0b2 = tools.get_Ds_ab(Label::ab::a0b2, Aa01, Ab2); if(D_a0b2.empty()) continue; @@ -711,7 +718,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia2 = 0; ia2 < list_Aa2.size(); ++ia2) // F { const TAC &Aa2 = list_Aa2[ia2]; - if (this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; const Tensor &D_a = tools.get_Ds_ab(Label::ab::a, Aa01, Aa2); if (D_a.empty()) continue; @@ -798,7 +805,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia01= 0; ia01filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; + if(this->filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; const Tensor &D_a0b2 = tools.get_Ds_ab(Label::ab::a0b2, Aa01, Ab2); if(D_a0b2.empty()) continue; @@ -851,7 +858,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia2 = 0; ia2 < list_Aa2.size(); ++ia2) // F { const TAC &Aa2 = list_Aa2[ia2]; - if (this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; const Tensor &D_a = Global_Func::find(Ds_a_transpose, Aa01, Aa2); if (D_a.empty()) continue; @@ -940,7 +947,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib01=0; ib01filter_atom->filter_for31(label,Aa01,Ab01,Aa2)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Aa2,Ab01)) continue; const Tensor &D_a2b1 = tools.get_Ds_ab(Label::ab::a2b1, Aa2, Ab01); if(D_a2b1.empty()) continue; @@ -1080,7 +1087,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib01=0; ib01filter_atom->filter_for31(label,Aa01,Ab01,Aa2)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Aa2,Ab01)) continue; const Tensor &D_a2b0 = tools.get_Ds_ab(Label::ab::a2b0, Aa2, Ab01); if(D_a2b0.empty()) continue; @@ -1220,7 +1227,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib01=0; ib01filter_atom->filter_for31(label,Aa01,Ab01,Aa2)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Aa2,Ab01)) continue; const Tensor &D_a2b1 = tools.get_Ds_ab(Label::ab::a2b1, Aa2, Ab01); if(D_a2b1.empty()) continue; @@ -1360,7 +1367,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib01=0; ib01filter_atom->filter_for31(label,Aa01,Ab01,Aa2)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Aa2,Ab01)) continue; const Tensor &D_a2b0 = tools.get_Ds_ab(Label::ab::a2b0, Aa2, Ab01); if(D_a2b0.empty()) continue; @@ -1553,7 +1560,7 @@ void LRI::cal_loop3_GPU( { const TAC &Ab2 = list_Ab2[ib2]; // The output block is (Aa01, Ab01); Ab2 is only an intermediate atom. - if (this->filter_atom->filter_for32(label, Aa01, Ab2, Ab01)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Ab01, Ab2)) continue; const Tensor &D_b = Global_Func::find(Ds_b_transpose, Ab01.first, TAC{Ab2.first, (Ab2.second-Ab01.second)%this->period}); if(D_b.empty()) continue; const GPU_Data::Pack &D_mul = rDs_mul.find_2nd(Aa01, Ab2); @@ -1692,7 +1699,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib2 = 0; ib2 < list_Ab2.size(); ++ib2) { const TAC &Ab2 = list_Ab2[ib2]; - if (this->filter_atom->filter_for32(label, Aa01, Ab2, Ab01)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Ab01, Ab2)) continue; const Tensor &D_b = tools.get_Ds_ab(Label::ab::b, Ab01, Ab2); if(D_b.empty()) continue; const GPU_Data::Pack &D_mul = rDs_mul.find_2nd(Aa01, Ab2); @@ -1831,7 +1838,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib2 = 0; ib2 < list_Ab2.size(); ++ib2) { const TAC &Ab2 = list_Ab2[ib2]; - if (this->filter_atom->filter_for32(label, Aa01, Ab2, Ab01)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Ab01, Ab2)) continue; const Tensor &D_b = Global_Func::find(Ds_b_transpose, Ab01.first, TAC{Ab2.first, (Ab2.second-Ab01.second)%this->period}); if(D_b.empty()) continue; const GPU_Data::Pack &D_mul = rDs_mul.find_2nd(Aa01, Ab2); @@ -1970,7 +1977,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib2 = 0; ib2 < list_Ab2.size(); ++ib2) { const TAC &Ab2 = list_Ab2[ib2]; - if (this->filter_atom->filter_for32(label, Aa01, Ab2, Ab01)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Ab01, Ab2)) continue; const Tensor &D_b = tools.get_Ds_ab(Label::ab::b, Ab01, Ab2); if(D_b.empty()) continue; const GPU_Data::Pack &D_mul = rDs_mul.find_2nd(Aa01, Ab2); @@ -2060,7 +2067,7 @@ void LRI::cal_loop3_GPU( for(std::size_t ia01=0; ia01filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; + if(this->filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; const Tensor &D_a1b2 = tools.get_Ds_ab(Label::ab::a1b2, Aa01, Ab2); if(D_a1b2.empty()) continue; @@ -2109,8 +2116,7 @@ void LRI::cal_loop3_GPU( { const TAC &Ab01 = list_Ab01[ib01]; // Screen the output pair; Aa2 is summed over, not an output atom. - if(this->filter_atom->filter_for2(label, Aa01, Ab01)) continue; - if(this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; + if(this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; const Tensor &D_a2b1 = tools.get_Ds_ab(Label::ab::a2b1, Aa2, Ab01); if(D_a2b1.empty()) continue; @@ -2225,7 +2231,7 @@ void LRI::cal_loop3_GPU( for(std::size_t ia01=0; ia01filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; + if(this->filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; const Tensor &D_a0b2 = tools.get_Ds_ab(Label::ab::a0b2, Aa01, Ab2); if(D_a0b2.empty()) continue; @@ -2274,8 +2280,7 @@ void LRI::cal_loop3_GPU( { const TAC &Ab01 = list_Ab01[ib01]; // Screen the output pair; Aa2 is summed over, not an output atom. - if(this->filter_atom->filter_for2(label, Aa01, Ab01)) continue; - if(this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; + if(this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; const Tensor &D_a2b0 = tools.get_Ds_ab(Label::ab::a2b0, Aa2, Ab01); if(D_a2b0.empty()) continue; @@ -2390,7 +2395,7 @@ void LRI::cal_loop3_GPU( for(std::size_t ia01=0; ia01filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; + if(this->filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; const Tensor &D_a0b2 = tools.get_Ds_ab(Label::ab::a0b2, Aa01, Ab2); if(D_a0b2.empty()) continue; @@ -2439,8 +2444,7 @@ void LRI::cal_loop3_GPU( { const TAC &Ab01 = list_Ab01[ib01]; // Screen the output pair; Aa2 is summed over, not an output atom. - if(this->filter_atom->filter_for2(label, Aa01, Ab01)) continue; - if(this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; + if(this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; const Tensor &D_a2b1 = tools.get_Ds_ab(Label::ab::a2b1, Aa2, Ab01); if(D_a2b1.empty()) continue; @@ -2555,7 +2559,7 @@ void LRI::cal_loop3_GPU( for(std::size_t ia01=0; ia01filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; + if(this->filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; const Tensor &D_a1b2 = tools.get_Ds_ab(Label::ab::a1b2, Aa01, Ab2); if(D_a1b2.empty()) continue; @@ -2604,8 +2608,7 @@ void LRI::cal_loop3_GPU( { const TAC &Ab01 = list_Ab01[ib01]; // Screen the output pair; Aa2 is summed over, not an output atom. - if(this->filter_atom->filter_for2(label, Aa01, Ab01)) continue; - if(this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; + if(this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; const Tensor &D_a2b0 = tools.get_Ds_ab(Label::ab::a2b0, Aa2, Ab01); if(D_a2b0.empty()) continue; diff --git a/include/RI/ri/LRI.h b/include/RI/ri/LRI.h index c3db89d..e1c154d 100644 --- a/include/RI/ri/LRI.h +++ b/include/RI/ri/LRI.h @@ -96,6 +96,9 @@ class LRI parallel = std::make_shared>(); std::unordered_map< Label::ab, RI_Tools::T_filter_func > filter_funcs; + // Construct the symmetry filter for the selected backend at contraction time. + bool flag_symmetry = false; + std::map, std::set> irreducible_sector; std::shared_ptr> filter_atom = std::make_shared>(); public: // private: diff --git a/unittests/gpu_symmetry/README.md b/unittests/gpu_symmetry/README.md index 167f271..6d0e8fe 100644 --- a/unittests/gpu_symmetry/README.md +++ b/unittests/gpu_symmetry/README.md @@ -33,3 +33,43 @@ or `-1` for all), backend (`cpu`, `gpu`, `all`) and scalar (`real`, `complex`, `all`). Index 1 is `a1b2_a2b1`: filtering `(Aa01,Aa2)` while summing over `Aa2` incorrectly removes contributions to the accepted output pair `(Aa01,Ab01)`. The off-diagonal sector with a nonzero translation detects this error. + +## GPU symmetry callback contract + +Arguments follow the enclosing loops from outermost to innermost. `filter_for31` +always returns false (intermediate construction); `filter_for32` rejects only +output pairs outside the sector. **Bold argument tuples changed in this refactor.** + +| Label | Family | `filter_for31` arguments | `filter_for32` arguments | Output pair screened by `filter_for32` | +| --- | --- | --- | --- | --- | +| `a0b0_a1b1` | `a01b01_a01b01` | `(Aa01, Ab01, Ab2)` | `(Aa01, Aa2, Ab2)` | `(Aa2, Ab2)` | +| `a0b1_a1b0` | `a01b01_a01b01` | `(Aa01, Ab01, Ab2)` | `(Aa01, Aa2, Ab2)` | `(Aa2, Ab2)` | +| `a0b0_a1b2` | `a01b01_a01b2` | **`(Ab01, Ab2, Aa01)`** | **`(Aa01, Ab01, Aa2)`** | `(Aa2, Ab01)` | +| `a0b1_a1b2` | `a01b01_a01b2` | **`(Ab01, Ab2, Aa01)`** | **`(Aa01, Ab01, Aa2)`** | `(Aa2, Ab01)` | +| `a0b2_a1b0` | `a01b01_a01b2` | **`(Ab01, Ab2, Aa01)`** | **`(Aa01, Ab01, Aa2)`** | `(Aa2, Ab01)` | +| `a0b2_a1b1` | `a01b01_a01b2` | **`(Ab01, Ab2, Aa01)`** | **`(Aa01, Ab01, Aa2)`** | `(Aa2, Ab01)` | +| `a0b0_a2b1` | `a01b01_a2b01` | **`(Aa01, Aa2, Ab01)`** | `(Aa01, Ab01, Ab2)` | `(Aa01, Ab2)` | +| `a0b1_a2b0` | `a01b01_a2b01` | **`(Aa01, Aa2, Ab01)`** | `(Aa01, Ab01, Ab2)` | `(Aa01, Ab2)` | +| `a1b0_a2b1` | `a01b01_a2b01` | **`(Aa01, Aa2, Ab01)`** | `(Aa01, Ab01, Ab2)` | `(Aa01, Ab2)` | +| `a1b1_a2b0` | `a01b01_a2b01` | **`(Aa01, Aa2, Ab01)`** | `(Aa01, Ab01, Ab2)` | `(Aa01, Ab2)` | +| `a0b0_a2b2` | `a01b01_a2b2` | `(Aa01, Aa2, Ab2)` | **`(Aa01, Ab01, Ab2)`** | `(Aa01, Ab01)` | +| `a0b1_a2b2` | `a01b01_a2b2` | `(Aa01, Aa2, Ab2)` | **`(Aa01, Ab01, Ab2)`** | `(Aa01, Ab01)` | +| `a1b0_a2b2` | `a01b01_a2b2` | `(Aa01, Aa2, Ab2)` | **`(Aa01, Ab01, Ab2)`** | `(Aa01, Ab01)` | +| `a1b1_a2b2` | `a01b01_a2b2` | `(Aa01, Aa2, Ab2)` | **`(Aa01, Ab01, Ab2)`** | `(Aa01, Ab01)` | +| `a1b2_a2b1` | `a01b2_a2b01` | **`(Ab01, Ab2, Aa01)`** | **`(Aa01, Aa2, Ab01)`** | `(Aa01, Ab01)` | +| `a0b2_a2b0` | `a01b2_a2b01` | **`(Ab01, Ab2, Aa01)`** | **`(Aa01, Aa2, Ab01)`** | `(Aa01, Ab01)` | +| `a0b2_a2b1` | `a01b2_a2b01` | **`(Ab01, Ab2, Aa01)`** | **`(Aa01, Aa2, Ab01)`** | `(Aa01, Ab01)` | +| `a1b2_a2b0` | `a01b2_a2b01` | **`(Ab01, Ab2, Aa01)`** | **`(Aa01, Aa2, Ab01)`** | `(Aa01, Ab01)` | + +`filter_for1(Aa01)` is a no-op for the first two families (Aa01 is intermediate), +and checks left-output membership for the other three. `filter_for1(Ab01)` checks +right-output membership wherever it occurs. All `filter_for2` calls are no-ops +except the final `(Aa01, Ab01)` stage of `a01b2_a2b01`, which checks that output +pair. In that family's second stage, the former inner +`filter_for2(Aa01, Ab01)` check is absorbed into +**`filter_for32(Aa01, Aa2, Ab01)`**; `Aa2` remains an unrestricted summation atom. + +For two cell-bearing output atoms the relative cell is `(Rb - Ra) % period`; +with `Aa01` at the origin it is `Rb % period`. An empty sector retains the existing +unrestricted behavior. The numerical regression also disables symmetry after each +filtered run to verify replacement of the previous backend's filter. diff --git a/unittests/gpu_symmetry/test_lri_gpu_symmetry.cpp b/unittests/gpu_symmetry/test_lri_gpu_symmetry.cpp index 233bc34..c4541c8 100644 --- a/unittests/gpu_symmetry/test_lri_gpu_symmetry.cpp +++ b/unittests/gpu_symmetry/test_lri_gpu_symmetry.cpp @@ -146,7 +146,7 @@ void check(const std::string& scalar, T imaginary, int selected, const std::stri labels.push_back(all_labels[index]); name = RI::Label_Tools::get_name(all_labels[index]); } - lri.filter_atom = std::make_shared>(); + lri.flag_symmetry = false; lri.cal_mode = RI::LRI_Cal_Mode::CPU; Tensors full; lri.cal_loop3(labels, full, 1.0); @@ -169,7 +169,8 @@ void check(const std::string& scalar, T imaginary, int selected, const std::stri if (keep) expected[left.first][right.first] = right.second.copy(); } - lri.filter_atom = std::make_shared>(period, sector); + lri.flag_symmetry = true; + lri.irreducible_sector = sector; for (const auto mode : modes) { lri.cal_mode = mode; @@ -204,6 +205,14 @@ void check(const std::string& scalar, T imaginary, int selected, const std::stri } compare(actual, accumulated, context + " factor=" + std::to_string(factor)); } + + // Disable symmetry after a restricted run on this same backend. + // This must replace the previous backend-specific filter. + lri.flag_symmetry = false; + Tensors unrestricted; + lri.cal_loop3(labels, unrestricted, 1.0); + compare(unrestricted, full, context + " symmetry disabled"); + lri.flag_symmetry = true; } } std::cout << "PASS " << scalar << " label=" << name << std::endl; From 4ed89410a5a6d7a891015879e8f0be0868af632d Mon Sep 17 00:00:00 2001 From: linpz Date: Tue, 10 Nov 2026 14:37:27 +0800 Subject: [PATCH 2/2] 1. add Filter_Atom_Symmetry_Wrapper --- include/RI/physics/Exx.hpp | 8 +- include/RI/physics/GW.hpp | 8 +- include/RI/physics/RPA.hpp | 8 +- .../symmetry/Filter_Atom_Symmetry_GPU.h | 5 + .../symmetry/Filter_Atom_Symmetry_Wrapper.h | 101 ++++++++++++++++++ include/RI/ri/LRI-cal_loop3_CPU.hpp | 7 -- .../LRI-cal_loop3_CPU_fine_grained_lock.hpp | 7 -- include/RI/ri/LRI-cal_loop3_GPU.hpp | 56 +++++----- include/RI/ri/LRI.h | 3 - .../gpu_symmetry/test_lri_gpu_symmetry.cpp | 13 +-- unittests/ri/LRI-loop3-test.hpp | 16 +-- 11 files changed, 158 insertions(+), 74 deletions(-) create mode 100644 include/RI/physics/symmetry/Filter_Atom_Symmetry_Wrapper.h diff --git a/include/RI/physics/Exx.hpp b/include/RI/physics/Exx.hpp index 42ff97d..6c55afd 100644 --- a/include/RI/physics/Exx.hpp +++ b/include/RI/physics/Exx.hpp @@ -9,6 +9,7 @@ #include "../ri/Cell_Nearest.h" #include "../ri/Label.h" #include "../global/Map_Operator.h" +#include "./symmetry/Filter_Atom_Symmetry_Wrapper.h" #include @@ -40,8 +41,11 @@ void Exx::set_symmetry( const bool flag_symmetry, const std::map, std::set> &irreducible_sector) { - this->lri.flag_symmetry = flag_symmetry; - this->lri.irreducible_sector = irreducible_sector; + if(flag_symmetry) + this->lri.filter_atom = std::make_shared>( + this->period, irreducible_sector, this->lri.cal_mode); + else + this->lri.filter_atom = std::make_shared>(); } template diff --git a/include/RI/physics/GW.hpp b/include/RI/physics/GW.hpp index 8733ee7..9746ecf 100644 --- a/include/RI/physics/GW.hpp +++ b/include/RI/physics/GW.hpp @@ -6,6 +6,7 @@ #include "GW.h" #include "../ri/Label.h" +#include "./symmetry/Filter_Atom_Symmetry_Wrapper.h" namespace RI { @@ -35,8 +36,11 @@ void GW::set_symmetry( const bool flag_symmetry, const std::map, std::set> &irreducible_sector) { - this->lri.flag_symmetry = flag_symmetry; - this->lri.irreducible_sector = irreducible_sector; + if(flag_symmetry) + this->lri.filter_atom = std::make_shared>( + this->period, irreducible_sector, this->lri.cal_mode); + else + this->lri.filter_atom = std::make_shared>(); } template diff --git a/include/RI/physics/RPA.hpp b/include/RI/physics/RPA.hpp index 454e5cf..51a8aed 100644 --- a/include/RI/physics/RPA.hpp +++ b/include/RI/physics/RPA.hpp @@ -7,6 +7,7 @@ #include "RPA.h" #include "../ri/Label.h" +#include "./symmetry/Filter_Atom_Symmetry_Wrapper.h" #include @@ -36,8 +37,11 @@ void RPA::set_symmetry( const bool flag_symmetry, const std::map, std::set> &irreducible_sector) { - this->lri.flag_symmetry = flag_symmetry; - this->lri.irreducible_sector = irreducible_sector; + if(flag_symmetry) + this->lri.filter_atom = std::make_shared>( + this->period, irreducible_sector, this->lri.cal_mode); + else + this->lri.filter_atom = std::make_shared>(); } template diff --git a/include/RI/physics/symmetry/Filter_Atom_Symmetry_GPU.h b/include/RI/physics/symmetry/Filter_Atom_Symmetry_GPU.h index 2275216..996c79a 100644 --- a/include/RI/physics/symmetry/Filter_Atom_Symmetry_GPU.h +++ b/include/RI/physics/symmetry/Filter_Atom_Symmetry_GPU.h @@ -1,3 +1,8 @@ +// =================== +// Author: maki49 +// date: 2026.11.10 +// =================== + #pragma once #include "../../ri/Filter_Atom.h" diff --git a/include/RI/physics/symmetry/Filter_Atom_Symmetry_Wrapper.h b/include/RI/physics/symmetry/Filter_Atom_Symmetry_Wrapper.h new file mode 100644 index 0000000..098176a --- /dev/null +++ b/include/RI/physics/symmetry/Filter_Atom_Symmetry_Wrapper.h @@ -0,0 +1,101 @@ +// =================== +// Author: Peize Lin +// date: 2026.11.10 +// =================== + +#pragma once + +#ifndef __GPU_RI + +#include "Filter_Atom_Symmetry.h" + +namespace RI +{ + +template +class Filter_Atom_Symmetry_Wrapper: public Filter_Atom_Symmetry +{ + public: + Filter_Atom_Symmetry_Wrapper( + const TC &period, + const std::map, std::set> &irsec, + const LRI_Cal_Mode &cal_mode_in) + : Filter_Atom_Symmetry(period, irsec){} +}; + +} + +#else + +#include "../../ri/Filter_Atom.h" +#include "Filter_Atom_Symmetry.h" +#include "Filter_Atom_Symmetry_GPU.h" +#include "../../ri/Label_Tools.h" + +#include + +namespace RI +{ + +template +class Filter_Atom_Symmetry_Wrapper: public Filter_Atom> +{ + public: + using TAC = std::pair; + + Filter_Atom_Symmetry_Wrapper( + const TC &period, + const std::map, std::set> &irsec, + const LRI_Cal_Mode &cal_mode_in) + :cal_mode(cal_mode_in), + filter_CPU(period, irsec), + filter_GPU(period, irsec){} + + virtual bool filter_for1(const Label::ab_ab &label, const TA &A1) const override + { + return this->cal_mode==LRI_Cal_Mode::GPU + ? this->filter_GPU.filter_for1(label, A1) + : this->filter_CPU.filter_for1(label, A1); + } + virtual bool filter_for1(const Label::ab_ab &label, const TAC &A1) const override + { + return this->cal_mode==LRI_Cal_Mode::GPU + ? this->filter_GPU.filter_for1(label, A1) + : this->filter_CPU.filter_for1(label, A1); + } + + virtual bool filter_for2(const Label::ab_ab &label, const TA &A1, const TAC &A2) const override + { + return this->cal_mode==LRI_Cal_Mode::GPU + ? this->filter_GPU.filter_for2(label, A1, A2) + : this->filter_CPU.filter_for2(label, A1, A2); + } + + virtual bool filter_for32(const Label::ab_ab &label, const TA &A1, const TAC &A2, const TAC &A3) const override + { + return this->cal_mode==LRI_Cal_Mode::GPU + ? this->filter_GPU.filter_for32(label, A1, A2, A3) + : this->filter_CPU.filter_for32(label, A1, A2, A3); + } + virtual bool filter_for32(const Label::ab_ab &label, const TAC &A1, const TA &A2, const TAC &A3) const override + { + return this->cal_mode==LRI_Cal_Mode::GPU + ? false + : this->filter_CPU.filter_for32(label, A1, A2, A3); + } + virtual bool filter_for32(const Label::ab_ab &label, const TAC &A1, const TAC &A2, const TAC &A3) const override + { + return this->cal_mode==LRI_Cal_Mode::GPU + ? false + : this->filter_CPU.filter_for32(label, A1, A2, A3); + } + + public: // private + const LRI_Cal_Mode &cal_mode; + Filter_Atom_Symmetry filter_CPU; + Filter_Atom_Symmetry_GPU filter_GPU; +}; + +} + +#endif \ No newline at end of file diff --git a/include/RI/ri/LRI-cal_loop3_CPU.hpp b/include/RI/ri/LRI-cal_loop3_CPU.hpp index f618e25..e702ebb 100644 --- a/include/RI/ri/LRI-cal_loop3_CPU.hpp +++ b/include/RI/ri/LRI-cal_loop3_CPU.hpp @@ -6,7 +6,6 @@ #pragma once #include "LRI.h" -#include "../physics/symmetry/Filter_Atom_Symmetry.h" #include "LRI_Cal_Aux.h" #include "../global/Array_Operator.h" #include "../global/Tensor_Multiply.h" @@ -25,12 +24,6 @@ void LRI::cal_loop3_CPU( std::map>> &Ds_result, const Tdata fac_add_Ds) { - if (this->flag_symmetry) - this->filter_atom = std::make_shared>( - this->period, this->irreducible_sector); - else - this->filter_atom = std::make_shared>(); - using namespace Array_Operator; const Data_Pack_Wrapper data_wrapper(this->data_pool, this->data_ab_name); diff --git a/include/RI/ri/LRI-cal_loop3_CPU_fine_grained_lock.hpp b/include/RI/ri/LRI-cal_loop3_CPU_fine_grained_lock.hpp index 308fa87..ac08b09 100644 --- a/include/RI/ri/LRI-cal_loop3_CPU_fine_grained_lock.hpp +++ b/include/RI/ri/LRI-cal_loop3_CPU_fine_grained_lock.hpp @@ -6,7 +6,6 @@ #pragma once #include "LRI.h" -#include "../physics/symmetry/Filter_Atom_Symmetry.h" #include "LRI_Cal_Aux.h" #include "../global/Array_Operator.h" #include "../global/Tensor_Multiply.h" @@ -25,12 +24,6 @@ void LRI::cal_loop3_CPU_fine_grained_lock( std::map>> &Ds_result, const Tdata fac_add_Ds) { - if (this->flag_symmetry) - this->filter_atom = std::make_shared>( - this->period, this->irreducible_sector); - else - this->filter_atom = std::make_shared>(); - using namespace Array_Operator; const Data_Pack_Wrapper data_wrapper(this->data_pool, this->data_ab_name); diff --git a/include/RI/ri/LRI-cal_loop3_GPU.hpp b/include/RI/ri/LRI-cal_loop3_GPU.hpp index 24090c2..802d696 100644 --- a/include/RI/ri/LRI-cal_loop3_GPU.hpp +++ b/include/RI/ri/LRI-cal_loop3_GPU.hpp @@ -8,7 +8,6 @@ #ifdef __GPU_RI #include "LRI.h" -#include "../physics/symmetry/Filter_Atom_Symmetry_GPU.h" #include "LRI_Cal_Aux.h" #include "../global/Array_Operator.h" #include "../global/Map_Operator.h" @@ -36,12 +35,6 @@ void LRI::cal_loop3_GPU( std::map>> &Ds_result, const Tdata fac_add_Ds) { - if (this->flag_symmetry) - this->filter_atom = std::make_shared>( - this->period, this->irreducible_sector); - else - this->filter_atom = std::make_shared>(); - using namespace Array_Operator; constexpr std::size_t memory_limit = std::size_t(5)*1024*1024*1024; @@ -385,7 +378,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia01= 0; ia01filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; const Tensor &D_a1b2 = tools.get_Ds_ab(Label::ab::a1b2, Aa01, Ab2); if(D_a1b2.empty()) continue; @@ -438,7 +431,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia2 = 0; ia2 < list_Aa2.size(); ++ia2) // F { const TAC &Aa2 = list_Aa2[ia2]; - if (this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; const Tensor &D_a = Global_Func::find(Ds_a_transpose, Aa01, Aa2); if (D_a.empty()) continue; @@ -525,7 +518,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia01= 0; ia01filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; const Tensor &D_a1b2 = tools.get_Ds_ab(Label::ab::a1b2, Aa01, Ab2); if(D_a1b2.empty()) continue; @@ -578,7 +571,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia2 = 0; ia2 < list_Aa2.size(); ++ia2) // F { const TAC &Aa2 = list_Aa2[ia2]; - if (this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; const Tensor &D_a = tools.get_Ds_ab(Label::ab::a, Aa01, Aa2); if (D_a.empty()) continue; @@ -665,7 +658,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia01= 0; ia01filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; const Tensor &D_a0b2 = tools.get_Ds_ab(Label::ab::a0b2, Aa01, Ab2); if(D_a0b2.empty()) continue; @@ -718,7 +711,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia2 = 0; ia2 < list_Aa2.size(); ++ia2) // F { const TAC &Aa2 = list_Aa2[ia2]; - if (this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; const Tensor &D_a = tools.get_Ds_ab(Label::ab::a, Aa01, Aa2); if (D_a.empty()) continue; @@ -805,7 +798,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia01= 0; ia01filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; const Tensor &D_a0b2 = tools.get_Ds_ab(Label::ab::a0b2, Aa01, Ab2); if(D_a0b2.empty()) continue; @@ -858,7 +851,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ia2 = 0; ia2 < list_Aa2.size(); ++ia2) // F { const TAC &Aa2 = list_Aa2[ia2]; - if (this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; + if (this->filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; const Tensor &D_a = Global_Func::find(Ds_a_transpose, Aa01, Aa2); if (D_a.empty()) continue; @@ -947,7 +940,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib01=0; ib01filter_atom->filter_for31(label,Aa01,Aa2,Ab01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Aa2)) continue; const Tensor &D_a2b1 = tools.get_Ds_ab(Label::ab::a2b1, Aa2, Ab01); if(D_a2b1.empty()) continue; @@ -1087,7 +1080,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib01=0; ib01filter_atom->filter_for31(label,Aa01,Aa2,Ab01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Aa2)) continue; const Tensor &D_a2b0 = tools.get_Ds_ab(Label::ab::a2b0, Aa2, Ab01); if(D_a2b0.empty()) continue; @@ -1227,7 +1220,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib01=0; ib01filter_atom->filter_for31(label,Aa01,Aa2,Ab01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Aa2)) continue; const Tensor &D_a2b1 = tools.get_Ds_ab(Label::ab::a2b1, Aa2, Ab01); if(D_a2b1.empty()) continue; @@ -1367,7 +1360,7 @@ void LRI::cal_loop3_GPU( for (std::size_t ib01=0; ib01filter_atom->filter_for31(label,Aa01,Aa2,Ab01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Aa2)) continue; const Tensor &D_a2b0 = tools.get_Ds_ab(Label::ab::a2b0, Aa2, Ab01); if(D_a2b0.empty()) continue; @@ -1559,7 +1552,6 @@ void LRI::cal_loop3_GPU( for (std::size_t ib2 = 0; ib2 < list_Ab2.size(); ++ib2) { const TAC &Ab2 = list_Ab2[ib2]; - // The output block is (Aa01, Ab01); Ab2 is only an intermediate atom. if (this->filter_atom->filter_for32(label, Aa01, Ab01, Ab2)) continue; const Tensor &D_b = Global_Func::find(Ds_b_transpose, Ab01.first, TAC{Ab2.first, (Ab2.second-Ab01.second)%this->period}); if(D_b.empty()) continue; @@ -2067,7 +2059,7 @@ void LRI::cal_loop3_GPU( for(std::size_t ia01=0; ia01filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; const Tensor &D_a1b2 = tools.get_Ds_ab(Label::ab::a1b2, Aa01, Ab2); if(D_a1b2.empty()) continue; @@ -2112,11 +2104,11 @@ void LRI::cal_loop3_GPU( const TAC &Aa2 = list_Aa2[ia2]; const Tensor &D_a = Global_Func::find(Ds_a_transpose, Aa01, Aa2); if(D_a.empty()) continue; + if(this->filter_atom->filter_for2(label,Aa01,Aa2)) continue; for(std::size_t ib01=0; ib01filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; + if(this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; const Tensor &D_a2b1 = tools.get_Ds_ab(Label::ab::a2b1, Aa2, Ab01); if(D_a2b1.empty()) continue; @@ -2231,7 +2223,7 @@ void LRI::cal_loop3_GPU( for(std::size_t ia01=0; ia01filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; const Tensor &D_a0b2 = tools.get_Ds_ab(Label::ab::a0b2, Aa01, Ab2); if(D_a0b2.empty()) continue; @@ -2276,11 +2268,11 @@ void LRI::cal_loop3_GPU( const TAC &Aa2 = list_Aa2[ia2]; const Tensor &D_a = Global_Func::find(Ds_a_transpose, Aa01, Aa2); if(D_a.empty()) continue; + if(this->filter_atom->filter_for2(label,Aa01,Aa2)) continue; for(std::size_t ib01=0; ib01filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; + if(this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; const Tensor &D_a2b0 = tools.get_Ds_ab(Label::ab::a2b0, Aa2, Ab01); if(D_a2b0.empty()) continue; @@ -2395,7 +2387,7 @@ void LRI::cal_loop3_GPU( for(std::size_t ia01=0; ia01filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; const Tensor &D_a0b2 = tools.get_Ds_ab(Label::ab::a0b2, Aa01, Ab2); if(D_a0b2.empty()) continue; @@ -2440,11 +2432,11 @@ void LRI::cal_loop3_GPU( const TAC &Aa2 = list_Aa2[ia2]; const Tensor &D_a = tools.get_Ds_ab(Label::ab::a, Aa01, Aa2); if(D_a.empty()) continue; + if(this->filter_atom->filter_for2(label,Aa01,Aa2)) continue; for(std::size_t ib01=0; ib01filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; + if(this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; const Tensor &D_a2b1 = tools.get_Ds_ab(Label::ab::a2b1, Aa2, Ab01); if(D_a2b1.empty()) continue; @@ -2559,7 +2551,7 @@ void LRI::cal_loop3_GPU( for(std::size_t ia01=0; ia01filter_atom->filter_for31(label,Ab01,Ab2,Aa01)) continue; + if(this->filter_atom->filter_for31(label,Aa01,Ab01,Ab2)) continue; const Tensor &D_a1b2 = tools.get_Ds_ab(Label::ab::a1b2, Aa01, Ab2); if(D_a1b2.empty()) continue; @@ -2604,11 +2596,11 @@ void LRI::cal_loop3_GPU( const TAC &Aa2 = list_Aa2[ia2]; const Tensor &D_a = tools.get_Ds_ab(Label::ab::a, Aa01, Aa2); if(D_a.empty()) continue; + if(this->filter_atom->filter_for2(label,Aa01,Aa2)) continue; for(std::size_t ib01=0; ib01filter_atom->filter_for32(label, Aa01, Aa2, Ab01)) continue; + if(this->filter_atom->filter_for32(label, Aa01, Ab01, Aa2)) continue; const Tensor &D_a2b0 = tools.get_Ds_ab(Label::ab::a2b0, Aa2, Ab01); if(D_a2b0.empty()) continue; diff --git a/include/RI/ri/LRI.h b/include/RI/ri/LRI.h index e1c154d..c3db89d 100644 --- a/include/RI/ri/LRI.h +++ b/include/RI/ri/LRI.h @@ -96,9 +96,6 @@ class LRI parallel = std::make_shared>(); std::unordered_map< Label::ab, RI_Tools::T_filter_func > filter_funcs; - // Construct the symmetry filter for the selected backend at contraction time. - bool flag_symmetry = false; - std::map, std::set> irreducible_sector; std::shared_ptr> filter_atom = std::make_shared>(); public: // private: diff --git a/unittests/gpu_symmetry/test_lri_gpu_symmetry.cpp b/unittests/gpu_symmetry/test_lri_gpu_symmetry.cpp index c4541c8..233bc34 100644 --- a/unittests/gpu_symmetry/test_lri_gpu_symmetry.cpp +++ b/unittests/gpu_symmetry/test_lri_gpu_symmetry.cpp @@ -146,7 +146,7 @@ void check(const std::string& scalar, T imaginary, int selected, const std::stri labels.push_back(all_labels[index]); name = RI::Label_Tools::get_name(all_labels[index]); } - lri.flag_symmetry = false; + lri.filter_atom = std::make_shared>(); lri.cal_mode = RI::LRI_Cal_Mode::CPU; Tensors full; lri.cal_loop3(labels, full, 1.0); @@ -169,8 +169,7 @@ void check(const std::string& scalar, T imaginary, int selected, const std::stri if (keep) expected[left.first][right.first] = right.second.copy(); } - lri.flag_symmetry = true; - lri.irreducible_sector = sector; + lri.filter_atom = std::make_shared>(period, sector); for (const auto mode : modes) { lri.cal_mode = mode; @@ -205,14 +204,6 @@ void check(const std::string& scalar, T imaginary, int selected, const std::stri } compare(actual, accumulated, context + " factor=" + std::to_string(factor)); } - - // Disable symmetry after a restricted run on this same backend. - // This must replace the previous backend-specific filter. - lri.flag_symmetry = false; - Tensors unrestricted; - lri.cal_loop3(labels, unrestricted, 1.0); - compare(unrestricted, full, context + " symmetry disabled"); - lri.flag_symmetry = true; } } std::cout << "PASS " << scalar << " label=" << name << std::endl; diff --git a/unittests/ri/LRI-loop3-test.hpp b/unittests/ri/LRI-loop3-test.hpp index 0c871e0..2b8adef 100644 --- a/unittests/ri/LRI-loop3-test.hpp +++ b/unittests/ri/LRI-loop3-test.hpp @@ -517,7 +517,7 @@ namespace LRI_Loop3_Test int mpi_init_provide; MPI_Init_thread(&argc, &argv, MPI_THREAD_MULTIPLE, &mpi_init_provide); - test_const_tensor_accumulation(); + //test_const_tensor_accumulation(); const std::vector cal_modes = { RI::LRI_Cal_Mode::CPU, @@ -529,13 +529,13 @@ namespace LRI_Loop3_Test for(const RI::LRI_Cal_Mode cal_mode: cal_modes) { test_loop3(cal_mode); - test_loop3_accumulation(cal_mode); - // Exercise imaginary coefficients, including the first write to an empty result. - using Treal = RI::Global_Func::To_Real_t; - test_loop3_accumulation(cal_mode, { - RI::Global_Func::convert(std::complex(0.25, 0.5)), - RI::Global_Func::convert(std::complex(0, -1)), - RI::Global_Func::convert(std::complex(-0.25, 0.5))}); + //test_loop3_accumulation(cal_mode); + //// Exercise imaginary coefficients, including the first write to an empty result. + //using Treal = RI::Global_Func::To_Real_t; + //test_loop3_accumulation(cal_mode, { + // RI::Global_Func::convert(std::complex(0.25, 0.5)), + // RI::Global_Func::convert(std::complex(0, -1)), + // RI::Global_Func::convert(std::complex(-0.25, 0.5))}); } MPI_Finalize();