From 4cbad1fc6088731e4986e0844b079f0738631758 Mon Sep 17 00:00:00 2001 From: linpz Date: Thu, 25 Dec 2025 15:52:30 +0800 Subject: [PATCH 1/9] Refactor: Simplify Lmax in Center2_Orb --- source/source_io/cal_r_overlap_R.cpp | 2 -- source/source_io/to_wannier90_lcao.cpp | 2 -- source/source_io/unk_overlap_lcao.cpp | 3 --- source/source_lcao/center2_orb.cpp | 2 +- source/source_lcao/center2_orb.h | 1 - source/source_lcao/module_ri/Matrix_Orbs11.cpp | 2 -- source/source_lcao/module_ri/Matrix_Orbs21.cpp | 2 -- source/source_lcao/module_ri/Matrix_Orbs22.cpp | 2 -- 8 files changed, 1 insertion(+), 15 deletions(-) diff --git a/source/source_io/cal_r_overlap_R.cpp b/source/source_io/cal_r_overlap_R.cpp index f81bf8026f..8da2dbd02e 100644 --- a/source/source_io/cal_r_overlap_R.cpp +++ b/source/source_io/cal_r_overlap_R.cpp @@ -18,7 +18,6 @@ cal_r_overlap_R::~cal_r_overlap_R() void cal_r_overlap_R::initialize_orb_table(const UnitCell& ucell, const LCAO_Orbitals& orb) { - int Lmax_used = 0; int Lmax = 0; int exx_lmax = 0; #ifdef __EXX @@ -40,7 +39,6 @@ void cal_r_overlap_R::initialize_orb_table(const UnitCell& ucell, Center2_Orb::init_Table_Spherical_Bessel(2, 3, - Lmax_used, Lmax, exx_lmax, lmax_orb, diff --git a/source/source_io/to_wannier90_lcao.cpp b/source/source_io/to_wannier90_lcao.cpp index 1388dc6975..0509c9d25a 100644 --- a/source/source_io/to_wannier90_lcao.cpp +++ b/source/source_io/to_wannier90_lcao.cpp @@ -265,7 +265,6 @@ void toWannier90_LCAO::out_unk(const psi::Psi>& psi) void toWannier90_LCAO::initialize_orb_table(const UnitCell& ucell) { - int Lmax_used = 0; int Lmax = 0; int exx_lmax = 0; #ifdef __EXX @@ -288,7 +287,6 @@ void toWannier90_LCAO::initialize_orb_table(const UnitCell& ucell) Center2_Orb::init_Table_Spherical_Bessel(2, 3, - Lmax_used, Lmax, exx_lmax, lmax_orb, diff --git a/source/source_io/unk_overlap_lcao.cpp b/source/source_io/unk_overlap_lcao.cpp index a352995a79..538656129c 100644 --- a/source/source_io/unk_overlap_lcao.cpp +++ b/source/source_io/unk_overlap_lcao.cpp @@ -28,8 +28,6 @@ void unkOverlap_lcao::init(const UnitCell& ucell, const int nkstot, const LCAO_Orbitals& orb) { - - int Lmax_used = 0; int Lmax = 0; int exx_lmax = 0; #ifdef __EXX @@ -51,7 +49,6 @@ void unkOverlap_lcao::init(const UnitCell& ucell, Center2_Orb::init_Table_Spherical_Bessel(2, 3, - Lmax_used, Lmax, exx_lmax, lmax_orb, diff --git a/source/source_lcao/center2_orb.cpp b/source/source_lcao/center2_orb.cpp index a2abeaf8bd..54f20c6e93 100644 --- a/source/source_lcao/center2_orb.cpp +++ b/source/source_lcao/center2_orb.cpp @@ -101,7 +101,6 @@ void Center2_Orb::init_Lmax(const int orb_num, // Peize Lin update 2016-01-26 void Center2_Orb::init_Table_Spherical_Bessel(const int orb_num, const int mode, - int& Lmax_used, int& Lmax, const int& Lmax_exx, const int lmax_orb, @@ -114,6 +113,7 @@ void Center2_Orb::init_Table_Spherical_Bessel(const int orb_num, { ModuleBase::TITLE("Center2_Orb", "init_Table_Spherical_Bessel"); + int Lmax_used; init_Lmax(orb_num, mode, Lmax_used, Lmax, Lmax_exx, lmax_orb, lmax_beta); // Peize Lin add 2016-01-26 for (auto& sb: ModuleBase::Sph_Bessel_Recursive_Pool::D2::sb_pool) diff --git a/source/source_lcao/center2_orb.h b/source/source_lcao/center2_orb.h index a10e3dffe5..1e765f7737 100644 --- a/source/source_lcao/center2_orb.h +++ b/source/source_lcao/center2_orb.h @@ -35,7 +35,6 @@ class Center2_Orb static void init_Table_Spherical_Bessel(const int orb_num, const int mode, - int& Lmax_used, int& Lmax, const int& Lmax_exx, const int lmax_orb, diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.cpp b/source/source_lcao/module_ri/Matrix_Orbs11.cpp index 021eb6e5b0..c9158bdc4e 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs11.cpp @@ -19,7 +19,6 @@ void Matrix_Orbs11::init(const int mode, ModuleBase::TITLE("Matrix_Orbs11", "init"); ModuleBase::timer::tick("Matrix_Orbs11", "init"); - int Lmax_used; this->lat0 = &ucell.lat0; const int ntype = orb.get_ntype(); int lmax_orb = -1, lmax_beta = -1; @@ -36,7 +35,6 @@ void Matrix_Orbs11::init(const int mode, Center2_Orb::init_Table_Spherical_Bessel(2, mode, - Lmax_used, Lmax, GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.cpp b/source/source_lcao/module_ri/Matrix_Orbs21.cpp index 170dc9b548..4827d8165c 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs21.cpp @@ -18,7 +18,6 @@ void Matrix_Orbs21::init(const int mode, { ModuleBase::TITLE("Matrix_Orbs21", "init"); ModuleBase::timer::tick("Matrix_Orbs21", "init"); - int Lmax_used; const int ntype = orb.get_ntype(); int lmax_orb = -1, lmax_beta = -1; @@ -36,7 +35,6 @@ void Matrix_Orbs21::init(const int mode, Center2_Orb::init_Table_Spherical_Bessel(3, mode, - Lmax_used, Lmax, GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.cpp b/source/source_lcao/module_ri/Matrix_Orbs22.cpp index 01d0330406..be43aafb02 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs22.cpp @@ -18,7 +18,6 @@ void Matrix_Orbs22::init(const int mode, { ModuleBase::TITLE("Matrix_Orbs22", "init"); ModuleBase::timer::tick("Matrix_Orbs22", "init"); - int Lmax_used; this->lat0 = &ucell.lat0; const int ntype = orb.get_ntype(); @@ -36,7 +35,6 @@ void Matrix_Orbs22::init(const int mode, Center2_Orb::init_Table_Spherical_Bessel(4, mode, - Lmax_used, Lmax, GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, From 54b19796a7477620c488901739636cef23234ec2 Mon Sep 17 00:00:00 2001 From: linpz Date: Thu, 25 Dec 2025 16:50:17 +0800 Subject: [PATCH 2/9] Refactor: Simplify Lmax in Center2_Orb --- source/source_io/cal_r_overlap_R.cpp | 9 +++------ source/source_io/to_wannier90_lcao.cpp | 9 +++------ source/source_io/unk_overlap_lcao.cpp | 9 +++------ source/source_lcao/center2_orb.cpp | 10 +--------- source/source_lcao/center2_orb.h | 7 +------ source/source_lcao/module_ri/Matrix_Orbs11.cpp | 9 +++------ source/source_lcao/module_ri/Matrix_Orbs21.cpp | 9 +++------ source/source_lcao/module_ri/Matrix_Orbs22.cpp | 9 +++------ 8 files changed, 20 insertions(+), 51 deletions(-) diff --git a/source/source_io/cal_r_overlap_R.cpp b/source/source_io/cal_r_overlap_R.cpp index 8da2dbd02e..7f8b9920c9 100644 --- a/source/source_io/cal_r_overlap_R.cpp +++ b/source/source_io/cal_r_overlap_R.cpp @@ -37,12 +37,9 @@ void cal_r_overlap_R::initialize_orb_table(const UnitCell& ucell, int Rmesh = static_cast(orb.get_Rmax() / dr) + 4; Rmesh += 1 - Rmesh % 2; - Center2_Orb::init_Table_Spherical_Bessel(2, - 3, - Lmax, - exx_lmax, - lmax_orb, - lmax_beta, + int Lmax_used; + Center2_Orb::init_Lmax(2, 3, Lmax_used, Lmax, exx_lmax, lmax_orb, lmax_beta); + Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, kmesh, diff --git a/source/source_io/to_wannier90_lcao.cpp b/source/source_io/to_wannier90_lcao.cpp index 0509c9d25a..af9f448c54 100644 --- a/source/source_io/to_wannier90_lcao.cpp +++ b/source/source_io/to_wannier90_lcao.cpp @@ -285,12 +285,9 @@ void toWannier90_LCAO::initialize_orb_table(const UnitCell& ucell) int Rmesh = static_cast(orb_.get_Rmax() / dr) + 4; Rmesh += 1 - Rmesh % 2; - Center2_Orb::init_Table_Spherical_Bessel(2, - 3, - Lmax, - exx_lmax, - lmax_orb, - lmax_beta, + int Lmax_used; + Center2_Orb::init_Lmax(2, 3, Lmax_used, Lmax, exx_lmax, lmax_orb, lmax_beta); + Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, kmesh, diff --git a/source/source_io/unk_overlap_lcao.cpp b/source/source_io/unk_overlap_lcao.cpp index 538656129c..a45cae434c 100644 --- a/source/source_io/unk_overlap_lcao.cpp +++ b/source/source_io/unk_overlap_lcao.cpp @@ -47,12 +47,9 @@ void unkOverlap_lcao::init(const UnitCell& ucell, int Rmesh = static_cast(orb.get_Rmax() / dr) + 4; Rmesh += 1 - Rmesh % 2; - Center2_Orb::init_Table_Spherical_Bessel(2, - 3, - Lmax, - exx_lmax, - lmax_orb, - lmax_beta, + int Lmax_used; + Center2_Orb::init_Lmax(2, 3, Lmax_used, Lmax, exx_lmax, lmax_orb, lmax_beta); + Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, kmesh, diff --git a/source/source_lcao/center2_orb.cpp b/source/source_lcao/center2_orb.cpp index 54f20c6e93..33f7dcd980 100644 --- a/source/source_lcao/center2_orb.cpp +++ b/source/source_lcao/center2_orb.cpp @@ -99,12 +99,7 @@ void Center2_Orb::init_Lmax(const int orb_num, } // Peize Lin update 2016-01-26 -void Center2_Orb::init_Table_Spherical_Bessel(const int orb_num, - const int mode, - int& Lmax, - const int& Lmax_exx, - const int lmax_orb, - const int lmax_beta, +void Center2_Orb::init_Table_Spherical_Bessel(const int Lmax_used, const double dr, const double dk, const int kmesh, @@ -113,9 +108,6 @@ void Center2_Orb::init_Table_Spherical_Bessel(const int orb_num, { ModuleBase::TITLE("Center2_Orb", "init_Table_Spherical_Bessel"); - int Lmax_used; - init_Lmax(orb_num, mode, Lmax_used, Lmax, Lmax_exx, lmax_orb, lmax_beta); // Peize Lin add 2016-01-26 - for (auto& sb: ModuleBase::Sph_Bessel_Recursive_Pool::D2::sb_pool) { if (dr * dk == sb.get_dx()) diff --git a/source/source_lcao/center2_orb.h b/source/source_lcao/center2_orb.h index 1e765f7737..9f2a70ad87 100644 --- a/source/source_lcao/center2_orb.h +++ b/source/source_lcao/center2_orb.h @@ -33,12 +33,7 @@ class Center2_Orb const int lmax_orb, const int lmax_beta); - static void init_Table_Spherical_Bessel(const int orb_num, - const int mode, - int& Lmax, - const int& Lmax_exx, - const int lmax_orb, - const int lmax_beta, + static void init_Table_Spherical_Bessel(const int Lmax_used, const double dr, const double dk, const int kmesh, diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.cpp b/source/source_lcao/module_ri/Matrix_Orbs11.cpp index c9158bdc4e..bd06821b4d 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs11.cpp @@ -33,12 +33,9 @@ void Matrix_Orbs11::init(const int mode, int Rmesh = static_cast(rmax / dr) + 4; Rmesh += 1 - Rmesh % 2; - Center2_Orb::init_Table_Spherical_Bessel(2, - mode, - Lmax, - GlobalC::exx_info.info_ri.abfs_Lmax, - lmax_orb, - lmax_beta, + int Lmax_used; + Center2_Orb::init_Lmax(2, mode, Lmax_used, Lmax, GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, lmax_beta); + Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, kmesh, diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.cpp b/source/source_lcao/module_ri/Matrix_Orbs21.cpp index 4827d8165c..88c794c628 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs21.cpp @@ -33,12 +33,9 @@ void Matrix_Orbs21::init(const int mode, int Rmesh = static_cast(rmax / dr) + 4; Rmesh += 1 - Rmesh % 2; - Center2_Orb::init_Table_Spherical_Bessel(3, - mode, - Lmax, - GlobalC::exx_info.info_ri.abfs_Lmax, - lmax_orb, - lmax_beta, + int Lmax_used; + Center2_Orb::init_Lmax(3, mode, Lmax_used, Lmax, GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, lmax_beta); + Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, kmesh, diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.cpp b/source/source_lcao/module_ri/Matrix_Orbs22.cpp index be43aafb02..9d1c452379 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs22.cpp @@ -33,12 +33,9 @@ void Matrix_Orbs22::init(const int mode, int Rmesh = static_cast(rmax / dr) + 4; Rmesh += 1 - Rmesh % 2; - Center2_Orb::init_Table_Spherical_Bessel(4, - mode, - Lmax, - GlobalC::exx_info.info_ri.abfs_Lmax, - lmax_orb, - lmax_beta, + int Lmax_used; + Center2_Orb::init_Lmax(4, mode, Lmax_used, Lmax, GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, lmax_beta); + Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, kmesh, From 2d26f60090c05bbe23f8c53c1ebab11c514aa643 Mon Sep 17 00:00:00 2001 From: linpz Date: Thu, 25 Dec 2025 18:20:26 +0800 Subject: [PATCH 3/9] Refactor: split Center2_Orb::init_Lmax() --- source/source_io/cal_r_overlap_R.cpp | 12 +- source/source_io/to_wannier90_lcao.cpp | 12 +- source/source_io/unk_overlap_lcao.cpp | 12 +- source/source_lcao/center2_orb.cpp | 111 ++++++++---------- source/source_lcao/center2_orb.h | 29 +++-- .../source_lcao/module_ri/Matrix_Orbs11.cpp | 9 +- .../source_lcao/module_ri/Matrix_Orbs21.cpp | 14 +-- .../source_lcao/module_ri/Matrix_Orbs22.cpp | 12 +- 8 files changed, 100 insertions(+), 111 deletions(-) diff --git a/source/source_io/cal_r_overlap_R.cpp b/source/source_io/cal_r_overlap_R.cpp index 7f8b9920c9..6f8bc5d2d9 100644 --- a/source/source_io/cal_r_overlap_R.cpp +++ b/source/source_io/cal_r_overlap_R.cpp @@ -18,18 +18,11 @@ cal_r_overlap_R::~cal_r_overlap_R() void cal_r_overlap_R::initialize_orb_table(const UnitCell& ucell, const LCAO_Orbitals& orb) { - int Lmax = 0; - int exx_lmax = 0; -#ifdef __EXX - exx_lmax = GlobalC::exx_info.info_ri.abfs_Lmax; -#endif - const int ntype = orb.get_ntype(); - int lmax_orb = -1, lmax_beta = -1; + int lmax_orb = -1; for (int it = 0; it < ntype; it++) { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); - lmax_beta = std::max(lmax_beta, ucell.infoNL.Beta[it].getLmax()); } const double dr = orb.get_dR(); const double dk = orb.get_dk(); @@ -37,8 +30,9 @@ void cal_r_overlap_R::initialize_orb_table(const UnitCell& ucell, int Rmesh = static_cast(orb.get_Rmax() / dr) + 4; Rmesh += 1 - Rmesh % 2; + int Lmax = 0; int Lmax_used; - Center2_Orb::init_Lmax(2, 3, Lmax_used, Lmax, exx_lmax, lmax_orb, lmax_beta); + Center2_Orb::init_Lmax_2_3(lmax_orb, Lmax_used, Lmax); Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, diff --git a/source/source_io/to_wannier90_lcao.cpp b/source/source_io/to_wannier90_lcao.cpp index af9f448c54..17b3eb8888 100644 --- a/source/source_io/to_wannier90_lcao.cpp +++ b/source/source_io/to_wannier90_lcao.cpp @@ -265,19 +265,12 @@ void toWannier90_LCAO::out_unk(const psi::Psi>& psi) void toWannier90_LCAO::initialize_orb_table(const UnitCell& ucell) { - int Lmax = 0; - int exx_lmax = 0; -#ifdef __EXX - exx_lmax = GlobalC::exx_info.info_ri.abfs_Lmax; -#endif - #ifdef __LCAO const int ntype = orb_.get_ntype(); - int lmax_orb = -1, lmax_beta = -1; + int lmax_orb = -1; for (int it = 0; it < ntype; it++) { lmax_orb = std::max(lmax_orb, orb_.Phi[it].getLmax()); - lmax_beta = std::max(lmax_beta, ucell.infoNL.Beta[it].getLmax()); } const double dr = orb_.get_dR(); const double dk = orb_.get_dk(); @@ -285,8 +278,9 @@ void toWannier90_LCAO::initialize_orb_table(const UnitCell& ucell) int Rmesh = static_cast(orb_.get_Rmax() / dr) + 4; Rmesh += 1 - Rmesh % 2; + int Lmax = 0; int Lmax_used; - Center2_Orb::init_Lmax(2, 3, Lmax_used, Lmax, exx_lmax, lmax_orb, lmax_beta); + Center2_Orb::init_Lmax_2_3(lmax_orb, Lmax_used, Lmax); Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, diff --git a/source/source_io/unk_overlap_lcao.cpp b/source/source_io/unk_overlap_lcao.cpp index a45cae434c..4ebf98c44f 100644 --- a/source/source_io/unk_overlap_lcao.cpp +++ b/source/source_io/unk_overlap_lcao.cpp @@ -28,18 +28,11 @@ void unkOverlap_lcao::init(const UnitCell& ucell, const int nkstot, const LCAO_Orbitals& orb) { - int Lmax = 0; - int exx_lmax = 0; -#ifdef __EXX - exx_lmax = GlobalC::exx_info.info_ri.abfs_Lmax; -#endif - const int ntype = orb.get_ntype(); - int lmax_orb = -1, lmax_beta = -1; + int lmax_orb = -1; for (int it = 0; it < ntype; it++) { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); - lmax_beta = std::max(lmax_beta, ucell.infoNL.Beta[it].getLmax()); } const double dr = orb.get_dR(); const double dk = orb.get_dk(); @@ -47,8 +40,9 @@ void unkOverlap_lcao::init(const UnitCell& ucell, int Rmesh = static_cast(orb.get_Rmax() / dr) + 4; Rmesh += 1 - Rmesh % 2; + int Lmax = 0; int Lmax_used; - Center2_Orb::init_Lmax(2, 3, Lmax_used, Lmax, exx_lmax, lmax_orb, lmax_beta); + Center2_Orb::init_Lmax_2_3(lmax_orb, Lmax_used, Lmax); Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, diff --git a/source/source_lcao/center2_orb.cpp b/source/source_lcao/center2_orb.cpp index 33f7dcd980..473b45580f 100644 --- a/source/source_lcao/center2_orb.cpp +++ b/source/source_lcao/center2_orb.cpp @@ -28,73 +28,58 @@ int Center2_Orb::get_rmesh(const double& R1, const double& R2, const double dr) return rmesh; } -// Peize Lin update 2016-01-26 -void Center2_Orb::init_Lmax(const int orb_num, - const int mode, - int& Lmax_used, - int& Lmax, - const int& Lmax_exx, - const int lmax_orb, - const int lmax_beta) +// used in or +void Center2_Orb::init_Lmax_2_1(const int lmax_orb, + const int lmax_beta, + int& Lmax_used, + int& Lmax) { + Lmax = std::max({-1, lmax_orb, lmax_beta}); + Lmax_used = 2 * Lmax + 1; + assert(Lmax_used >= 1); +} - Lmax = -1; +// used in or +void Center2_Orb::init_Lmax_2_2(const int& lmax_exx, + int& Lmax_used, + int& Lmax) +{ + Lmax = std::max(-1, lmax_exx); + Lmax_used = 2 * Lmax + 1; + assert(Lmax_used >= 1); +} - switch (orb_num) - { - case 2: - switch (mode) - { - case 1: // used in or - Lmax = std::max({Lmax, lmax_orb, lmax_beta}); - // use 2lmax+1 in dS - Lmax_used = 2 * Lmax + 1; - break; - case 2: // used in or - Lmax = std::max(Lmax, Lmax_exx); - Lmax_used = 2 * Lmax + 1; - break; - case 3: // used in berryphase by jingan - Lmax = std::max(Lmax, lmax_orb); - Lmax++; - Lmax_used = 2 * Lmax + 1; - break; - default: - throw std::invalid_argument("Center2_Orb::init_Lmax orb_num=2, mode error"); - break; - } - break; - case 3: - switch (mode) - { - case 1: // used in or - Lmax = std::max(Lmax, lmax_orb); - Lmax_used = 2 * Lmax + 1; - Lmax = std::max(Lmax, Lmax_exx); - Lmax_used += Lmax_exx; - break; - default: - throw std::invalid_argument("Center2_Orb::init_Lmax orb_num=3, mode error"); - break; - } - break; - case 4: - switch (mode) - { - case 1: // used in - Lmax = std::max(Lmax, lmax_orb); - Lmax_used = 2 * (2 * Lmax + 1); - break; - default: - throw std::invalid_argument("Center2_Orb::init_Lmax orb_num=4, mode error"); - break; - } - break; - default: - throw std::invalid_argument("Center2_Orb::init_Lmax orb_num error"); - break; - } +// used in berryphase by jingan +void Center2_Orb::init_Lmax_2_3(const int lmax_orb, + int& Lmax_used, + int& Lmax) +{ + Lmax = std::max(-1, lmax_orb); + Lmax++; + Lmax_used = 2 * Lmax + 1; + assert(Lmax_used >= 1); +} + +// used in or +void Center2_Orb::init_Lmax_3_1(const int& lmax_exx, + const int lmax_orb, + int& Lmax_used, + int& Lmax) +{ + Lmax = std::max(-1, lmax_orb); + Lmax_used = 2 * Lmax + 1; + Lmax = std::max(Lmax, lmax_exx); + Lmax_used += lmax_exx; + assert(Lmax_used >= 1); +} +// used in +void Center2_Orb::init_Lmax_4_1(const int lmax_orb, + int& Lmax_used, + int& Lmax) +{ + Lmax = std::max(-1, lmax_orb); + Lmax_used = 2 * (2 * Lmax + 1); assert(Lmax_used >= 1); } diff --git a/source/source_lcao/center2_orb.h b/source/source_lcao/center2_orb.h index 9f2a70ad87..55783187b1 100644 --- a/source/source_lcao/center2_orb.h +++ b/source/source_lcao/center2_orb.h @@ -25,13 +25,28 @@ class Center2_Orb static int get_rmesh(const double& R1, const double& R2, const double dr); - static void init_Lmax(const int orb_num, - const int mode, - int& Lmax_used, - int& Lmax, - const int& Lmax_exx, - const int lmax_orb, - const int lmax_beta); + // used in or + static void init_Lmax_2_1(const int lmax_orb, + const int lmax_beta, + int& Lmax_used, + int& Lmax); + // used in or + static void init_Lmax_2_2(const int& lmax_exx, + int& Lmax_used, + int& Lmax); + // used in berryphase by jingan + static void init_Lmax_2_3(const int lmax_orb, + int& Lmax_used, + int& Lmax); + // used in or + static void init_Lmax_3_1(const int& lmax_exx, + const int lmax_orb, + int& Lmax_used, + int& Lmax); + // used in + static void init_Lmax_4_1(const int lmax_orb, + int& Lmax_used, + int& Lmax); static void init_Table_Spherical_Bessel(const int Lmax_used, const double dr, diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.cpp b/source/source_lcao/module_ri/Matrix_Orbs11.cpp index bd06821b4d..49ba2ce796 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs11.cpp @@ -34,7 +34,14 @@ void Matrix_Orbs11::init(const int mode, Rmesh += 1 - Rmesh % 2; int Lmax_used; - Center2_Orb::init_Lmax(2, mode, Lmax_used, Lmax, GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, lmax_beta); + //if(mode==1) + // { Center2_Orb::init_Lmax_2_1(lmax_orb, lmax_beta, Lmax_used, Lmax); } + if(mode==2) + { Center2_Orb::init_Lmax_2_2(GlobalC::exx_info.info_ri.abfs_Lmax, Lmax_used, Lmax); } + //else if(mode==3) + // { Center2_Orb::init_Lmax_2_3(lmax_orb, Lmax_used, Lmax); } + else + { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.cpp b/source/source_lcao/module_ri/Matrix_Orbs21.cpp index 88c794c628..35cc350da2 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs21.cpp @@ -20,13 +20,10 @@ void Matrix_Orbs21::init(const int mode, ModuleBase::timer::tick("Matrix_Orbs21", "init"); const int ntype = orb.get_ntype(); - int lmax_orb = -1, lmax_beta = -1; - this->lat0 = &ucell.lat0; + this->lat0 = &ucell.lat0; + int lmax_orb = -1; for (int it = 0; it < ntype; it++) - { - lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); - lmax_beta = std::max(lmax_beta, ucell.infoNL.Beta[it].getLmax()); - } + { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); } const double dr = orb.get_dR(); const double dk = orb.get_dk(); const int kmesh = orb.get_kmesh() * kmesh_times + 1; @@ -34,7 +31,10 @@ void Matrix_Orbs21::init(const int mode, Rmesh += 1 - Rmesh % 2; int Lmax_used; - Center2_Orb::init_Lmax(3, mode, Lmax_used, Lmax, GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, lmax_beta); + if(mode==1) + { Center2_Orb::init_Lmax_3_1(GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, Lmax_used, Lmax); } + else + { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.cpp b/source/source_lcao/module_ri/Matrix_Orbs22.cpp index 9d1c452379..af0f0d42f1 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs22.cpp @@ -21,12 +21,9 @@ void Matrix_Orbs22::init(const int mode, this->lat0 = &ucell.lat0; const int ntype = orb.get_ntype(); - int lmax_orb = -1, lmax_beta = -1; + int lmax_orb = -1; for (int it = 0; it < ntype; it++) - { - lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); - lmax_beta = std::max(lmax_beta, ucell.infoNL.Beta[it].getLmax()); - } + { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); } const double dr = orb.get_dR(); const double dk = orb.get_dk(); const int kmesh = orb.get_kmesh() * kmesh_times + 1; @@ -34,7 +31,10 @@ void Matrix_Orbs22::init(const int mode, Rmesh += 1 - Rmesh % 2; int Lmax_used; - Center2_Orb::init_Lmax(4, mode, Lmax_used, Lmax, GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, lmax_beta); + if(mode==1) + { Center2_Orb::init_Lmax_4_1(lmax_orb, Lmax_used, Lmax); } + else + { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, From 395d3f2574e89fa195ec8d9d6de7a4a43a9acd57 Mon Sep 17 00:00:00 2001 From: linpz Date: Thu, 25 Dec 2025 18:35:49 +0800 Subject: [PATCH 4/9] Refactor: simplify Matrix_Orbs::init() --- .../module_ri/ABFs_Construct-PCA.cpp | 5 +-- source/source_lcao/module_ri/LRI_CV.hpp | 7 ++- .../source_lcao/module_ri/Matrix_Orbs11.cpp | 44 ++++++++++--------- source/source_lcao/module_ri/Matrix_Orbs11.h | 5 +-- .../source_lcao/module_ri/Matrix_Orbs21.cpp | 33 +++++++------- source/source_lcao/module_ri/Matrix_Orbs21.h | 5 +-- .../source_lcao/module_ri/Matrix_Orbs22.cpp | 32 +++++++------- source/source_lcao/module_ri/Matrix_Orbs22.h | 5 +-- source/source_lcao/module_ri/exx_opt_orb.cpp | 18 +++----- 9 files changed, 74 insertions(+), 80 deletions(-) diff --git a/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp b/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp index c4f5cd7f13..aaa1ef2f41 100644 --- a/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp +++ b/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp @@ -182,10 +182,9 @@ RI::Tensor get_column_mean0_matrix(const RI::Tensor& m) Matrix_Orbs21 m_abfslcaos_lcaos; ORB_gaunt_table MGT; - int Lmax; - m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, orb.get_Rmax(), Lmax); + const int Lmax = m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, orb.get_Rmax()); MGT.init_Gaunt_CH(Lmax); - MGT.init_Gaunt(Lmax); + MGT.init_Gaunt(Lmax); m_abfslcaos_lcaos.init_radial(abfs, lcaos, lcaos, MGT); std::map>> delta_R; diff --git a/source/source_lcao/module_ri/LRI_CV.hpp b/source/source_lcao/module_ri/LRI_CV.hpp index 53f9931fd9..32e01302e7 100644 --- a/source/source_lcao/module_ri/LRI_CV.hpp +++ b/source/source_lcao/module_ri/LRI_CV.hpp @@ -69,12 +69,11 @@ void LRI_CV::set_orbitals( range_abfs = ModuleBase::Element_Basis_Index::construct_range( abfs ); this->index_abfs = ModuleBase::Element_Basis_Index::construct_index( range_abfs ); - int Lmax_v = std::numeric_limits::min(); - this->m_abfs_abfs.init(2, ucell, orb, kmesh_times, lcaos_rmax + abfs_ccp_rmax, Lmax_v); + const int Lmax_v = this->m_abfs_abfs.init(2, ucell, orb, kmesh_times, lcaos_rmax + abfs_ccp_rmax); int Lmax_c = std::numeric_limits::min(); if (init_C) - this->m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, lcaos_rmax, Lmax_c); - int Lmax = std::max(Lmax_v, Lmax_c); + Lmax_c = this->m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, lcaos_rmax); + const int Lmax = std::max(Lmax_v, Lmax_c); if (init_MGT) { MGT.init_Gaunt_CH(Lmax); diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.cpp b/source/source_lcao/module_ri/Matrix_Orbs11.cpp index 49ba2ce796..73c217649b 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs11.cpp @@ -9,30 +9,25 @@ #include "source_base/tool_title.h" #include "source_pw/module_pwdft/global.h" -void Matrix_Orbs11::init(const int mode, +int Matrix_Orbs11::init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, - const double rmax, - int& Lmax) + const double rmax) { ModuleBase::TITLE("Matrix_Orbs11", "init"); ModuleBase::timer::tick("Matrix_Orbs11", "init"); this->lat0 = &ucell.lat0; - const int ntype = orb.get_ntype(); - int lmax_orb = -1, lmax_beta = -1; - for (int it = 0; it < ntype; it++) - { - lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); - lmax_beta = std::max(lmax_beta, ucell.infoNL.Beta[it].getLmax()); - } - const double dr = orb.get_dR(); - const double dk = orb.get_dk(); - const int kmesh = orb.get_kmesh() * kmesh_times + 1; - int Rmesh = static_cast(rmax / dr) + 4; - Rmesh += 1 - Rmesh % 2; + //const int ntype = orb.get_ntype(); + //int lmax_orb = -1, lmax_beta = -1; + //for (int it = 0; it < ntype; it++) + //{ + // lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); + // lmax_beta = std::max(lmax_beta, ucell.infoNL.Beta[it].getLmax()); + //} + int Lmax; int Lmax_used; //if(mode==1) // { Center2_Orb::init_Lmax_2_1(lmax_orb, lmax_beta, Lmax_used, Lmax); } @@ -42,12 +37,6 @@ void Matrix_Orbs11::init(const int mode, // { Center2_Orb::init_Lmax_2_3(lmax_orb, Lmax_used, Lmax); } else { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } - Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, - dr, - dk, - kmesh, - Rmesh, - psb_); //========================================= // (3) make Gaunt coefficients table @@ -55,7 +44,20 @@ void Matrix_Orbs11::init(const int mode, // this->MGT.init_Gaunt_CH(Lmax); // this->MGT.init_Gaunt(Lmax); + const double dr = orb.get_dR(); + const double dk = orb.get_dk(); + const int kmesh = orb.get_kmesh() * kmesh_times + 1; + int Rmesh = static_cast(rmax / dr) + 4; + Rmesh += 1 - Rmesh % 2; + Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, + dr, + dk, + kmesh, + Rmesh, + psb_); + ModuleBase::timer::tick("Matrix_Orbs11", "init"); + return Lmax; } void Matrix_Orbs11::init_radial(const std::vector>>& orb_A, diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.h b/source/source_lcao/module_ri/Matrix_Orbs11.h index 04b1bd8600..a90f0d3fed 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.h +++ b/source/source_lcao/module_ri/Matrix_Orbs11.h @@ -24,12 +24,11 @@ class Matrix_Orbs11 // mode: // 1: // 2: - void init(const int mode, + int init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, // extend Kcut, keep dK - const double rmax, - int& Lmax); + const double rmax); void init_radial(const std::vector>>& orb_A, const std::vector>>& orb_B, diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.cpp b/source/source_lcao/module_ri/Matrix_Orbs21.cpp index 35cc350da2..539f65d004 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs21.cpp @@ -9,38 +9,26 @@ #include "source_base/tool_title.h" #include "source_pw/module_pwdft/global.h" -void Matrix_Orbs21::init(const int mode, +int Matrix_Orbs21::init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, - const double rmax, - int& Lmax) + const double rmax) { ModuleBase::TITLE("Matrix_Orbs21", "init"); ModuleBase::timer::tick("Matrix_Orbs21", "init"); + this->lat0 = &ucell.lat0; const int ntype = orb.get_ntype(); - this->lat0 = &ucell.lat0; int lmax_orb = -1; for (int it = 0; it < ntype; it++) { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); } - const double dr = orb.get_dR(); - const double dk = orb.get_dk(); - const int kmesh = orb.get_kmesh() * kmesh_times + 1; - int Rmesh = static_cast(rmax / dr) + 4; - Rmesh += 1 - Rmesh % 2; - + int Lmax; int Lmax_used; if(mode==1) { Center2_Orb::init_Lmax_3_1(GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, Lmax_used, Lmax); } else { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } - Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, - dr, - dk, - kmesh, - Rmesh, - psb_); //========================================= // (3) make Gaunt coefficients table @@ -49,7 +37,20 @@ void Matrix_Orbs21::init(const int mode, // this->MGT.init_Gaunt(2 * Lmax + 1); Lmax = 2 * Lmax + 1; + const double dr = orb.get_dR(); + const double dk = orb.get_dk(); + const int kmesh = orb.get_kmesh() * kmesh_times + 1; + int Rmesh = static_cast(rmax / dr) + 4; + Rmesh += 1 - Rmesh % 2; + Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, + dr, + dk, + kmesh, + Rmesh, + psb_); + ModuleBase::timer::tick("Matrix_Orbs21", "init"); + return Lmax; } void Matrix_Orbs21::init_radial(const std::vector>>& orb_A1, diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.h b/source/source_lcao/module_ri/Matrix_Orbs21.h index f1241e87b2..afa6e89ef0 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.h +++ b/source/source_lcao/module_ri/Matrix_Orbs21.h @@ -22,12 +22,11 @@ class Matrix_Orbs21 public: // mode: // 1: - void init(const int mode, + int init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, // extend Kcut, keep dK - const double rmax, - int& Lmax); // extend Rcut, keep dR + const double rmax); // extend Rcut, keep dR void init_radial(const std::vector>>& orb_A1, const std::vector>>& orb_A2, diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.cpp b/source/source_lcao/module_ri/Matrix_Orbs22.cpp index af0f0d42f1..6a95713732 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs22.cpp @@ -9,38 +9,27 @@ #include "source_base/tool_title.h" #include "source_pw/module_pwdft/global.h" -void Matrix_Orbs22::init(const int mode, +int Matrix_Orbs22::init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, - const double rmax, - int& Lmax) + const double rmax) { ModuleBase::TITLE("Matrix_Orbs22", "init"); ModuleBase::timer::tick("Matrix_Orbs22", "init"); this->lat0 = &ucell.lat0; + const int ntype = orb.get_ntype(); int lmax_orb = -1; for (int it = 0; it < ntype; it++) { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); } - const double dr = orb.get_dR(); - const double dk = orb.get_dk(); - const int kmesh = orb.get_kmesh() * kmesh_times + 1; - int Rmesh = static_cast(rmax / dr) + 4; - Rmesh += 1 - Rmesh % 2; - + int Lmax; int Lmax_used; if(mode==1) { Center2_Orb::init_Lmax_4_1(lmax_orb, Lmax_used, Lmax); } else { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } - Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, - dr, - dk, - kmesh, - Rmesh, - psb_); //========================================= // (3) make Gaunt coefficients table @@ -49,7 +38,20 @@ void Matrix_Orbs22::init(const int mode, // this->MGT.init_Gaunt(2 * Lmax + 1); Lmax = 2 * Lmax + 1; + const double dr = orb.get_dR(); + const double dk = orb.get_dk(); + const int kmesh = orb.get_kmesh() * kmesh_times + 1; + int Rmesh = static_cast(rmax / dr) + 4; + Rmesh += 1 - Rmesh % 2; + Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, + dr, + dk, + kmesh, + Rmesh, + psb_); + ModuleBase::timer::tick("Matrix_Orbs22", "init"); + return Lmax; } void Matrix_Orbs22::init_radial(const std::vector>>& orb_A1, diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.h b/source/source_lcao/module_ri/Matrix_Orbs22.h index a7db572176..4e9f49f628 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.h +++ b/source/source_lcao/module_ri/Matrix_Orbs22.h @@ -23,12 +23,11 @@ class Matrix_Orbs22 public: // mode: // 1: - void init(const int mode, + int init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, // extend Kcut, keep dK - const double rmax, - int& Lmax); // extend Rcut, keep dR + const double rmax); // extend Rcut, keep dR void init_radial(const std::vector>>& orb_A1, const std::vector>>& orb_A2, diff --git a/source/source_lcao/module_ri/exx_opt_orb.cpp b/source/source_lcao/module_ri/exx_opt_orb.cpp index f36d65249d..4ba6569107 100644 --- a/source/source_lcao/module_ri/exx_opt_orb.cpp +++ b/source/source_lcao/module_ri/exx_opt_orb.cpp @@ -73,8 +73,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(lcaos)) { return {}; } Matrix_Orbs22 m_lcaoslcaos_lcaoslcaos; ORB_gaunt_table MGT; - int Lmax; - m_lcaoslcaos_lcaoslcaos.init( 1, ucell,orb, info.kmesh_times, orb.get_Rmax(), Lmax ); + const int Lmax = m_lcaoslcaos_lcaoslcaos.init( 1, ucell,orb, info.kmesh_times, orb.get_Rmax() ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_lcaoslcaos_lcaoslcaos.init_radial( lcaos, lcaos, lcaos, lcaos, MGT ); @@ -93,8 +92,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(jle)) { return {}; } Matrix_Orbs21 m_jyslcaos_lcaos; ORB_gaunt_table MGT; - int Lmax; - m_jyslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), Lmax ); + const int Lmax = m_jyslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax() ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_jyslcaos_lcaos.init_radial( jle, lcaos, lcaos, MGT); @@ -112,8 +110,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(jle)) { return {}; } Matrix_Orbs11 m_jys_jys; ORB_gaunt_table MGT; - int Lmax; - m_jys_jys.init( 2,ucell,orb, info.kmesh_times, orb.get_Rmax(), Lmax ); + const int Lmax = m_jys_jys.init( 2,ucell,orb, info.kmesh_times, orb.get_Rmax() ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_jys_jys.init_radial( jle, jle, MGT ); @@ -131,8 +128,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs11 m_abfs_abfs; ORB_gaunt_table MGT; - int Lmax; - m_abfs_abfs.init( 2, ucell, orb, info.kmesh_times, orb.get_Rmax(), Lmax ); + const int Lmax = m_abfs_abfs.init( 2, ucell, orb, info.kmesh_times, orb.get_Rmax() ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_abfs_abfs.init_radial( abfs, abfs, MGT ); @@ -151,8 +147,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs21 m_abfslcaos_lcaos; ORB_gaunt_table MGT; - int Lmax; - m_abfslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), Lmax ); + const int Lmax = m_abfslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax() ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_abfslcaos_lcaos.init_radial( abfs, lcaos, lcaos, MGT ); @@ -171,8 +166,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs11 m_jys_abfs; ORB_gaunt_table MGT; - int Lmax; - m_jys_abfs.init( 2, ucell,orb, info.kmesh_times, orb.get_Rmax(), Lmax ); + const int Lmax = m_jys_abfs.init( 2, ucell,orb, info.kmesh_times, orb.get_Rmax() ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_jys_abfs.init_radial( jle, abfs, MGT ); From aed3d0e636516100ede6460ab4af9c92874db716 Mon Sep 17 00:00:00 2001 From: linpz Date: Thu, 25 Dec 2025 19:03:24 +0800 Subject: [PATCH 5/9] Refactor: Simplify Center2_Orb::init_Lmax() --- source/source_io/cal_r_overlap_R.cpp | 5 +- source/source_io/to_wannier90_lcao.cpp | 5 +- source/source_io/unk_overlap_lcao.cpp | 5 +- source/source_lcao/center2_orb.cpp | 47 ++++++++----------- source/source_lcao/center2_orb.h | 22 ++------- .../source_lcao/module_ri/Matrix_Orbs11.cpp | 9 ++-- .../source_lcao/module_ri/Matrix_Orbs21.cpp | 5 +- .../source_lcao/module_ri/Matrix_Orbs22.cpp | 5 +- 8 files changed, 39 insertions(+), 64 deletions(-) diff --git a/source/source_io/cal_r_overlap_R.cpp b/source/source_io/cal_r_overlap_R.cpp index 6f8bc5d2d9..d947d792c8 100644 --- a/source/source_io/cal_r_overlap_R.cpp +++ b/source/source_io/cal_r_overlap_R.cpp @@ -30,9 +30,8 @@ void cal_r_overlap_R::initialize_orb_table(const UnitCell& ucell, int Rmesh = static_cast(orb.get_Rmax() / dr) + 4; Rmesh += 1 - Rmesh % 2; - int Lmax = 0; - int Lmax_used; - Center2_Orb::init_Lmax_2_3(lmax_orb, Lmax_used, Lmax); + int Lmax, Lmax_used; + std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_3(lmax_orb); Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, diff --git a/source/source_io/to_wannier90_lcao.cpp b/source/source_io/to_wannier90_lcao.cpp index 17b3eb8888..45414d15b9 100644 --- a/source/source_io/to_wannier90_lcao.cpp +++ b/source/source_io/to_wannier90_lcao.cpp @@ -278,9 +278,8 @@ void toWannier90_LCAO::initialize_orb_table(const UnitCell& ucell) int Rmesh = static_cast(orb_.get_Rmax() / dr) + 4; Rmesh += 1 - Rmesh % 2; - int Lmax = 0; - int Lmax_used; - Center2_Orb::init_Lmax_2_3(lmax_orb, Lmax_used, Lmax); + int Lmax, Lmax_used; + std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_3(lmax_orb); Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, diff --git a/source/source_io/unk_overlap_lcao.cpp b/source/source_io/unk_overlap_lcao.cpp index 4ebf98c44f..4d79d78a13 100644 --- a/source/source_io/unk_overlap_lcao.cpp +++ b/source/source_io/unk_overlap_lcao.cpp @@ -40,9 +40,8 @@ void unkOverlap_lcao::init(const UnitCell& ucell, int Rmesh = static_cast(orb.get_Rmax() / dr) + 4; Rmesh += 1 - Rmesh % 2; - int Lmax = 0; - int Lmax_used; - Center2_Orb::init_Lmax_2_3(lmax_orb, Lmax_used, Lmax); + int Lmax, Lmax_used; + std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_3(lmax_orb); Center2_Orb::init_Table_Spherical_Bessel(Lmax_used, dr, dk, diff --git a/source/source_lcao/center2_orb.cpp b/source/source_lcao/center2_orb.cpp index 473b45580f..8df630777b 100644 --- a/source/source_lcao/center2_orb.cpp +++ b/source/source_lcao/center2_orb.cpp @@ -29,58 +29,51 @@ int Center2_Orb::get_rmesh(const double& R1, const double& R2, const double dr) } // used in or -void Center2_Orb::init_Lmax_2_1(const int lmax_orb, - const int lmax_beta, - int& Lmax_used, - int& Lmax) +std::pair Center2_Orb::init_Lmax_2_1(const int lmax_orb, const int lmax_beta) { - Lmax = std::max({-1, lmax_orb, lmax_beta}); - Lmax_used = 2 * Lmax + 1; + const int Lmax = std::max({-1, lmax_orb, lmax_beta}); + const int Lmax_used = 2 * Lmax + 1; assert(Lmax_used >= 1); + return {Lmax_used, Lmax}; } // used in or -void Center2_Orb::init_Lmax_2_2(const int& lmax_exx, - int& Lmax_used, - int& Lmax) +std::pair Center2_Orb::init_Lmax_2_2(const int& lmax_exx) { - Lmax = std::max(-1, lmax_exx); - Lmax_used = 2 * Lmax + 1; + const int Lmax = std::max(-1, lmax_exx); + const int Lmax_used = 2 * Lmax + 1; assert(Lmax_used >= 1); + return {Lmax_used, Lmax}; } // used in berryphase by jingan -void Center2_Orb::init_Lmax_2_3(const int lmax_orb, - int& Lmax_used, - int& Lmax) +std::pair Center2_Orb::init_Lmax_2_3(const int lmax_orb) { - Lmax = std::max(-1, lmax_orb); + int Lmax = std::max(-1, lmax_orb); Lmax++; - Lmax_used = 2 * Lmax + 1; + const int Lmax_used = 2 * Lmax + 1; assert(Lmax_used >= 1); + return {Lmax_used, Lmax}; } // used in or -void Center2_Orb::init_Lmax_3_1(const int& lmax_exx, - const int lmax_orb, - int& Lmax_used, - int& Lmax) +std::pair Center2_Orb::init_Lmax_3_1(const int& lmax_exx, const int lmax_orb) { - Lmax = std::max(-1, lmax_orb); - Lmax_used = 2 * Lmax + 1; + int Lmax = std::max(-1, lmax_orb); + int Lmax_used = 2 * Lmax + 1; Lmax = std::max(Lmax, lmax_exx); Lmax_used += lmax_exx; assert(Lmax_used >= 1); + return {Lmax_used, Lmax}; } // used in -void Center2_Orb::init_Lmax_4_1(const int lmax_orb, - int& Lmax_used, - int& Lmax) +std::pair Center2_Orb::init_Lmax_4_1(const int lmax_orb) { - Lmax = std::max(-1, lmax_orb); - Lmax_used = 2 * (2 * Lmax + 1); + const int Lmax = std::max(-1, lmax_orb); + const int Lmax_used = 2 * (2 * Lmax + 1); assert(Lmax_used >= 1); + return {Lmax_used, Lmax}; } // Peize Lin update 2016-01-26 diff --git a/source/source_lcao/center2_orb.h b/source/source_lcao/center2_orb.h index 55783187b1..92288e97a8 100644 --- a/source/source_lcao/center2_orb.h +++ b/source/source_lcao/center2_orb.h @@ -26,27 +26,15 @@ class Center2_Orb static int get_rmesh(const double& R1, const double& R2, const double dr); // used in or - static void init_Lmax_2_1(const int lmax_orb, - const int lmax_beta, - int& Lmax_used, - int& Lmax); + static std::pair init_Lmax_2_1(const int lmax_orb, const int lmax_beta); // used in or - static void init_Lmax_2_2(const int& lmax_exx, - int& Lmax_used, - int& Lmax); + static std::pair init_Lmax_2_2(const int& lmax_exx); // used in berryphase by jingan - static void init_Lmax_2_3(const int lmax_orb, - int& Lmax_used, - int& Lmax); + static std::pair init_Lmax_2_3(const int lmax_orb); // used in or - static void init_Lmax_3_1(const int& lmax_exx, - const int lmax_orb, - int& Lmax_used, - int& Lmax); + static std::pair init_Lmax_3_1(const int& lmax_exx, const int lmax_orb); // used in - static void init_Lmax_4_1(const int lmax_orb, - int& Lmax_used, - int& Lmax); + static std::pair init_Lmax_4_1(const int lmax_orb); static void init_Table_Spherical_Bessel(const int Lmax_used, const double dr, diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.cpp b/source/source_lcao/module_ri/Matrix_Orbs11.cpp index 73c217649b..906a99f008 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs11.cpp @@ -27,14 +27,13 @@ int Matrix_Orbs11::init(const int mode, // lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); // lmax_beta = std::max(lmax_beta, ucell.infoNL.Beta[it].getLmax()); //} - int Lmax; - int Lmax_used; + int Lmax, Lmax_used; //if(mode==1) - // { Center2_Orb::init_Lmax_2_1(lmax_orb, lmax_beta, Lmax_used, Lmax); } + // { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_1(lmax_orb, lmax_beta); } if(mode==2) - { Center2_Orb::init_Lmax_2_2(GlobalC::exx_info.info_ri.abfs_Lmax, Lmax_used, Lmax); } + { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_2(GlobalC::exx_info.info_ri.abfs_Lmax); } //else if(mode==3) - // { Center2_Orb::init_Lmax_2_3(lmax_orb, Lmax_used, Lmax); } + // { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_3(lmax_orb); } else { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.cpp b/source/source_lcao/module_ri/Matrix_Orbs21.cpp index 539f65d004..2814067f72 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs21.cpp @@ -23,10 +23,9 @@ int Matrix_Orbs21::init(const int mode, int lmax_orb = -1; for (int it = 0; it < ntype; it++) { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); } - int Lmax; - int Lmax_used; + int Lmax, Lmax_used; if(mode==1) - { Center2_Orb::init_Lmax_3_1(GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb, Lmax_used, Lmax); } + { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_3_1(GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb); } else { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.cpp b/source/source_lcao/module_ri/Matrix_Orbs22.cpp index 6a95713732..3c1d0a1e64 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs22.cpp @@ -24,10 +24,9 @@ int Matrix_Orbs22::init(const int mode, int lmax_orb = -1; for (int it = 0; it < ntype; it++) { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); } - int Lmax; - int Lmax_used; + int Lmax, Lmax_used; if(mode==1) - { Center2_Orb::init_Lmax_4_1(lmax_orb, Lmax_used, Lmax); } + { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_4_1(lmax_orb); } else { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } From fedd57c50adf8d343725708d0228047ea4c7a494 Mon Sep 17 00:00:00 2001 From: linpz Date: Fri, 26 Dec 2025 23:48:05 +0800 Subject: [PATCH 6/9] Refactor: change some global variable GlobalC::exx_info.info_ri.abfs_Lmax to local variable --- source/source_lcao/LCAO_init_basis.cpp | 5 --- source/source_lcao/center2_orb.cpp | 39 +++++++++++-------- .../module_ri/ABFs_Construct-PCA.cpp | 16 ++------ source/source_lcao/module_ri/Exx_LRI.hpp | 5 --- source/source_lcao/module_ri/LRI_CV.hpp | 4 +- .../source_lcao/module_ri/Matrix_Orbs11.cpp | 5 ++- source/source_lcao/module_ri/Matrix_Orbs11.h | 12 +++--- .../source_lcao/module_ri/Matrix_Orbs21.cpp | 14 ++++--- source/source_lcao/module_ri/Matrix_Orbs21.h | 12 +++--- source/source_lcao/module_ri/exx_opt_orb.cpp | 16 ++++---- 10 files changed, 61 insertions(+), 67 deletions(-) diff --git a/source/source_lcao/LCAO_init_basis.cpp b/source/source_lcao/LCAO_init_basis.cpp index f8b60b6298..6d198a3a32 100644 --- a/source/source_lcao/LCAO_init_basis.cpp +++ b/source/source_lcao/LCAO_init_basis.cpp @@ -62,11 +62,6 @@ void init_basis_lcao(Parallel_Orbitals& pv, two_center_bundle.build_beta(ucell.ntype, ucell.infoNL.Beta); } - int Lmax = 0; -#ifdef __EXX - Lmax = GlobalC::exx_info.info_ri.abfs_Lmax; -#endif - #ifdef USE_NEW_TWO_CENTER two_center_bundle.tabulate(); #else diff --git a/source/source_lcao/center2_orb.cpp b/source/source_lcao/center2_orb.cpp index 8df630777b..b5ece1d3e3 100644 --- a/source/source_lcao/center2_orb.cpp +++ b/source/source_lcao/center2_orb.cpp @@ -49,8 +49,7 @@ std::pair Center2_Orb::init_Lmax_2_2(const int& lmax_exx) // used in berryphase by jingan std::pair Center2_Orb::init_Lmax_2_3(const int lmax_orb) { - int Lmax = std::max(-1, lmax_orb); - Lmax++; + const int Lmax = std::max(-1, lmax_orb) + 1; const int Lmax_used = 2 * Lmax + 1; assert(Lmax_used >= 1); return {Lmax_used, Lmax}; @@ -172,15 +171,14 @@ void Center2_Orb::cal_ST_Phi12_R(const int& job, // double* integrated_func = new double[kmesh]; - int ll = 0; - if (l != 0) - { - ll = l - 1; - } - - const std::vector>& jlm1 = psb->get_jlx()[ll]; - const std::vector>& jl = psb->get_jlx()[l]; - const std::vector>& jlp1 = psb->get_jlx()[l + 1]; + assert(psb->get_jlx().size()>=l+2); + const int lml = (l>0) ? (l-1) : 0; + const std::vector>& jlm1 = psb->get_jlx().at(lml); + const std::vector>& jl = psb->get_jlx().at(l); + const std::vector>& jlp1 = psb->get_jlx().at(l+1); + assert(jlm1.size()>=rmesh); + assert(jl.size()>=rmesh); + assert(jlp1.size()>=rmesh); #ifdef _OPENMP #pragma omp parallel for schedule(static) @@ -189,6 +187,7 @@ void Center2_Orb::cal_ST_Phi12_R(const int& job, { std::vector integrated_func(kmesh); const std::vector& jl_r = jl[ir]; + assert(jl_r.size()>=kmesh); for (int ik = 0; ik < kmesh; ++ik) { integrated_func[ik] = jl_r[ik] * k1_dot_k2[ik]; @@ -202,6 +201,8 @@ void Center2_Orb::cal_ST_Phi12_R(const int& job, // Peize Lin accelerate 2017-10-02 const std::vector& jlm1_r = jlm1[ir]; const std::vector& jlp1_r = jlp1[ir]; + assert(jlm1_r.size()>=kmesh); + assert(jlp1_r.size()>=kmesh); const double fac = l / (l + 1.0); if (l == 0) { @@ -305,10 +306,11 @@ void Center2_Orb::cal_ST_Phi12_R(const int& job, std::vector integrated_func(kmesh); - const int lm1 = (l > 0 ? l - 1 : 0); - const std::vector>& jlm1 = psb->get_jlx()[lm1]; - const std::vector>& jl = psb->get_jlx()[l]; - const std::vector>& jlp1 = psb->get_jlx()[l + 1]; + assert(psb->get_jlx().size()>=l+2); + const int lm1 = (l>0) ? (l-1) : 0; + const std::vector>& jlm1 = psb->get_jlx().at(lm1); + const std::vector>& jl = psb->get_jlx().at(l); + const std::vector>& jlp1 = psb->get_jlx().at(l+1); for (const size_t& ir: radials) { @@ -320,6 +322,7 @@ void Center2_Orb::cal_ST_Phi12_R(const int& job, } const std::vector& jl_r = jl[ir]; + assert(jl_r.size()>=kmesh); for (int ik = 0; ik < kmesh; ++ik) { integrated_func[ik] = jl_r[ik] * k1_dot_k2[ik]; @@ -329,8 +332,10 @@ void Center2_Orb::cal_ST_Phi12_R(const int& job, ModuleBase::Integral::Simpson_Integral(kmesh, ModuleBase::GlobalFunc::VECTOR_TO_PTR(integrated_func), dk, temp); rs[ir] = temp * ModuleBase::FOUR_PI; - const std::vector& jlm1_r = jlm1[ir]; - const std::vector& jlp1_r = jlp1[ir]; + const std::vector& jlm1_r = jlm1.at(ir); + const std::vector& jlp1_r = jlp1.at(ir); + assert(jlm1_r.size()>=kmesh); + assert(jlp1_r.size()>=kmesh); const double fac = l / (l + 1.0); if (l == 0) { diff --git a/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp b/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp index aaa1ef2f41..9360723464 100644 --- a/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp +++ b/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp @@ -172,29 +172,21 @@ RI::Tensor get_column_mean0_matrix(const RI::Tensor& m) const ModuleBase::Element_Basis_Index::IndexLNM index_abfs = ModuleBase::Element_Basis_Index::construct_index(range_abfs); - const int Lmax_bak = GlobalC::exx_info.info_ri.abfs_Lmax; - GlobalC::exx_info.info_ri.abfs_Lmax = std::numeric_limits::min(); + int lmax_abfs = std::numeric_limits::min(); for (std::size_t T = 0; T != abfs.size(); ++T) - { - GlobalC::exx_info.info_ri.abfs_Lmax - = std::max(GlobalC::exx_info.info_ri.abfs_Lmax, static_cast(abfs[T].size()) - 1); -} + { lmax_abfs = std::max(lmax_abfs, static_cast(abfs[T].size()) - 1); } Matrix_Orbs21 m_abfslcaos_lcaos; ORB_gaunt_table MGT; - const int Lmax = m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, orb.get_Rmax()); + const int Lmax = m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, orb.get_Rmax(), lmax_abfs); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_abfslcaos_lcaos.init_radial(abfs, lcaos, lcaos, MGT); std::map>> delta_R; for (std::size_t it = 0; it != abfs.size(); ++it) - { - delta_R[it][it] = {0.0}; -} + { delta_R[it][it] = {0.0}; } m_abfslcaos_lcaos.init_radial_table(delta_R); - - GlobalC::exx_info.info_ri.abfs_Lmax = Lmax_bak; std::vector, RI::Tensor>>> eig(abfs.size()); for (std::size_t T = 0; T != abfs.size(); ++T) diff --git a/source/source_lcao/module_ri/Exx_LRI.hpp b/source/source_lcao/module_ri/Exx_LRI.hpp index 5cf3a0072a..d271a5c8ce 100644 --- a/source/source_lcao/module_ri/Exx_LRI.hpp +++ b/source/source_lcao/module_ri/Exx_LRI.hpp @@ -73,11 +73,6 @@ void Exx_LRI::cal_exx_ions(const UnitCell& ucell, ModuleBase::TITLE("Exx_LRI","cal_exx_ions"); ModuleBase::timer::tick("Exx_LRI", "cal_exx_ions"); - // init_radial_table_ions( cal_atom_centres_core(atom_pairs_core_origin), atom_pairs_core_origin ); - - // this->m_abfsabfs.init_radial_table(Rradial); - // this->m_abfslcaos_lcaos.init_radial_table(Rradial); - std::vector atoms(ucell.nat); for(int iat=0; iat::set_orbitals( range_abfs = ModuleBase::Element_Basis_Index::construct_range( abfs ); this->index_abfs = ModuleBase::Element_Basis_Index::construct_index( range_abfs ); - const int Lmax_v = this->m_abfs_abfs.init(2, ucell, orb, kmesh_times, lcaos_rmax + abfs_ccp_rmax); + const int Lmax_v = this->m_abfs_abfs.init(2, ucell, orb, kmesh_times, lcaos_rmax + abfs_ccp_rmax, GlobalC::exx_info.info_ri.abfs_Lmax); int Lmax_c = std::numeric_limits::min(); if (init_C) - Lmax_c = this->m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, lcaos_rmax); + Lmax_c = this->m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, lcaos_rmax, GlobalC::exx_info.info_ri.abfs_Lmax); const int Lmax = std::max(Lmax_v, Lmax_c); if (init_MGT) { diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.cpp b/source/source_lcao/module_ri/Matrix_Orbs11.cpp index 906a99f008..506a34d4c0 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs11.cpp @@ -13,7 +13,8 @@ int Matrix_Orbs11::init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, - const double rmax) + const double rmax, + const int lmax_abfs) { ModuleBase::TITLE("Matrix_Orbs11", "init"); ModuleBase::timer::tick("Matrix_Orbs11", "init"); @@ -31,7 +32,7 @@ int Matrix_Orbs11::init(const int mode, //if(mode==1) // { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_1(lmax_orb, lmax_beta); } if(mode==2) - { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_2(GlobalC::exx_info.info_ri.abfs_Lmax); } + { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_2(lmax_abfs); } //else if(mode==3) // { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_3(lmax_orb); } else diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.h b/source/source_lcao/module_ri/Matrix_Orbs11.h index a90f0d3fed..a5e3efc053 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.h +++ b/source/source_lcao/module_ri/Matrix_Orbs11.h @@ -24,11 +24,13 @@ class Matrix_Orbs11 // mode: // 1: // 2: - int init(const int mode, - const UnitCell& ucell, - const LCAO_Orbitals& orb, - const double kmesh_times, // extend Kcut, keep dK - const double rmax); + int init( + const int mode, + const UnitCell& ucell, + const LCAO_Orbitals& orb, + const double kmesh_times, // extend Kcut, keep dK + const double rmax, + const int lmax_abfs); void init_radial(const std::vector>>& orb_A, const std::vector>>& orb_B, diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.cpp b/source/source_lcao/module_ri/Matrix_Orbs21.cpp index 2814067f72..6ce200876a 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs21.cpp @@ -9,11 +9,13 @@ #include "source_base/tool_title.h" #include "source_pw/module_pwdft/global.h" -int Matrix_Orbs21::init(const int mode, - const UnitCell& ucell, - const LCAO_Orbitals& orb, - const double kmesh_times, - const double rmax) +int Matrix_Orbs21::init( + const int mode, + const UnitCell& ucell, + const LCAO_Orbitals& orb, + const double kmesh_times, + const double rmax, + const int lmax_abfs) { ModuleBase::TITLE("Matrix_Orbs21", "init"); ModuleBase::timer::tick("Matrix_Orbs21", "init"); @@ -25,7 +27,7 @@ int Matrix_Orbs21::init(const int mode, { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); } int Lmax, Lmax_used; if(mode==1) - { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_3_1(GlobalC::exx_info.info_ri.abfs_Lmax, lmax_orb); } + { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_3_1(lmax_abfs, lmax_orb); } else { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.h b/source/source_lcao/module_ri/Matrix_Orbs21.h index afa6e89ef0..2552b43fca 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.h +++ b/source/source_lcao/module_ri/Matrix_Orbs21.h @@ -22,11 +22,13 @@ class Matrix_Orbs21 public: // mode: // 1: - int init(const int mode, - const UnitCell& ucell, - const LCAO_Orbitals& orb, - const double kmesh_times, // extend Kcut, keep dK - const double rmax); // extend Rcut, keep dR + int init( + const int mode, + const UnitCell& ucell, + const LCAO_Orbitals& orb, + const double kmesh_times, // extend Kcut, keep dK + const double rmax, // extend Rcut, keep dR + const int lmax_abfs); void init_radial(const std::vector>>& orb_A1, const std::vector>>& orb_A2, diff --git a/source/source_lcao/module_ri/exx_opt_orb.cpp b/source/source_lcao/module_ri/exx_opt_orb.cpp index 4ba6569107..884548688b 100644 --- a/source/source_lcao/module_ri/exx_opt_orb.cpp +++ b/source/source_lcao/module_ri/exx_opt_orb.cpp @@ -47,11 +47,11 @@ void Exx_Opt_Orb::generate_matrix( { jle = Exx_Abfs::IO::construct_abfs( jle, orb, info.files_jles, info.kmesh_times ); } Exx_Abfs::Construct_Orbs::filter_empty_orbs(jle); - // GlobalC::exx_info.info_ri.abfs_Lmax is used for Center2 temporarily + int lmax_abfs = -1; for(const auto &orb_T : abfs) - { GlobalC::exx_info.info_ri.abfs_Lmax = std::max( GlobalC::exx_info.info_ri.abfs_Lmax, static_cast(orb_T.size())-1 ); } + { lmax_abfs = std::max( lmax_abfs, static_cast(orb_T.size())-1 ); } for(const auto &orb_T : jle) - { GlobalC::exx_info.info_ri.abfs_Lmax = std::max( GlobalC::exx_info.info_ri.abfs_Lmax, static_cast(orb_T.size())-1 ); } + { lmax_abfs = std::max( lmax_abfs, static_cast(orb_T.size())-1 ); } const ModuleBase::Element_Basis_Index::Range range_lcaos = ModuleBase::Element_Basis_Index::construct_range( lcaos ); const ModuleBase::Element_Basis_Index::IndexLNM index_lcaos = ModuleBase::Element_Basis_Index::construct_index( range_lcaos ); @@ -92,7 +92,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(jle)) { return {}; } Matrix_Orbs21 m_jyslcaos_lcaos; ORB_gaunt_table MGT; - const int Lmax = m_jyslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax() ); + const int Lmax = m_jyslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_jyslcaos_lcaos.init_radial( jle, lcaos, lcaos, MGT); @@ -110,7 +110,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(jle)) { return {}; } Matrix_Orbs11 m_jys_jys; ORB_gaunt_table MGT; - const int Lmax = m_jys_jys.init( 2,ucell,orb, info.kmesh_times, orb.get_Rmax() ); + const int Lmax = m_jys_jys.init( 2,ucell,orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_jys_jys.init_radial( jle, jle, MGT ); @@ -128,7 +128,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs11 m_abfs_abfs; ORB_gaunt_table MGT; - const int Lmax = m_abfs_abfs.init( 2, ucell, orb, info.kmesh_times, orb.get_Rmax() ); + const int Lmax = m_abfs_abfs.init( 2, ucell, orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_abfs_abfs.init_radial( abfs, abfs, MGT ); @@ -147,7 +147,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs21 m_abfslcaos_lcaos; ORB_gaunt_table MGT; - const int Lmax = m_abfslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax() ); + const int Lmax = m_abfslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_abfslcaos_lcaos.init_radial( abfs, lcaos, lcaos, MGT ); @@ -166,7 +166,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs11 m_jys_abfs; ORB_gaunt_table MGT; - const int Lmax = m_jys_abfs.init( 2, ucell,orb, info.kmesh_times, orb.get_Rmax() ); + const int Lmax = m_jys_abfs.init( 2, ucell,orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); MGT.init_Gaunt_CH(Lmax); MGT.init_Gaunt(Lmax); m_jys_abfs.init_radial( jle, abfs, MGT ); From 7ea4eb3f635d19185aef7879aae2cb90445ffee7 Mon Sep 17 00:00:00 2001 From: linpz Date: Sat, 27 Dec 2025 01:00:20 +0800 Subject: [PATCH 7/9] Refactor: Simplify ORB_gaunt_table in exx --- .../module_ao/ORB_gaunt_table.cpp | 2 + .../source_basis/module_ao/ORB_gaunt_table.h | 6 +++ .../module_ri/ABFs_Construct-PCA.cpp | 7 +--- source/source_lcao/module_ri/Exx_LRI.h | 1 - source/source_lcao/module_ri/Exx_LRI.hpp | 5 +-- source/source_lcao/module_ri/LRI_CV.h | 3 +- source/source_lcao/module_ri/LRI_CV.hpp | 20 +++------ .../source_lcao/module_ri/Matrix_Orbs11.cpp | 20 +++++---- source/source_lcao/module_ri/Matrix_Orbs11.h | 9 ++-- .../source_lcao/module_ri/Matrix_Orbs21.cpp | 23 +++++----- source/source_lcao/module_ri/Matrix_Orbs21.h | 10 ++--- .../source_lcao/module_ri/Matrix_Orbs22.cpp | 23 +++++----- source/source_lcao/module_ri/Matrix_Orbs22.h | 10 ++--- source/source_lcao/module_ri/exx_opt_orb.cpp | 42 ++++++------------- 14 files changed, 81 insertions(+), 100 deletions(-) diff --git a/source/source_basis/module_ao/ORB_gaunt_table.cpp b/source/source_basis/module_ao/ORB_gaunt_table.cpp index b930c6bcba..d05f594f7e 100644 --- a/source/source_basis/module_ao/ORB_gaunt_table.cpp +++ b/source/source_basis/module_ao/ORB_gaunt_table.cpp @@ -20,6 +20,7 @@ void ORB_gaunt_table::init_Gaunt(const int &lmax) ModuleBase::TITLE("ORB_gaunt_table", "init_Gaunt"); ModuleBase::timer::tick("ORB_gaunt_table", "init_Gaunt"); + this->Lmax_Gaunt_Coefficients = lmax; const int nlm = (lmax * 2 + 1) * (lmax * 2 + 1); this->Gaunt_Coefficients.create(nlm, nlm, nlm); @@ -184,6 +185,7 @@ void ORB_gaunt_table::init_Gaunt_CH(const int& Lmax) ModuleBase::TITLE("ORB_gaunt_table","init_Gaunt_CH"); ModuleBase::timer::tick("ORB_gaunt_table","init_Gaunt_CH"); + this->Lmax_Gaunt_CH = Lmax; int L = 2*Lmax + 1; int Eff_Np = this->EP_EL(L); diff --git a/source/source_basis/module_ao/ORB_gaunt_table.h b/source/source_basis/module_ao/ORB_gaunt_table.h index 1e295b9125..0e9aefdb70 100644 --- a/source/source_basis/module_ao/ORB_gaunt_table.h +++ b/source/source_basis/module_ao/ORB_gaunt_table.h @@ -93,6 +93,9 @@ class ORB_gaunt_table static int Index_M(const int& m); + int get_Lmax_Gaunt_Coefficients() const { return Lmax_Gaunt_Coefficients; } + int get_Lmax_Gaunt_CH() const { return Lmax_Gaunt_CH; } + private: // Index Function @@ -122,5 +125,8 @@ class ORB_gaunt_table //direct integral ModuleBase::matrix Ylm_Gaunt; + + int Lmax_Gaunt_Coefficients = -1; + int Lmax_Gaunt_CH = -1; }; #endif diff --git a/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp b/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp index 9360723464..67b89200f2 100644 --- a/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp +++ b/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp @@ -177,11 +177,8 @@ RI::Tensor get_column_mean0_matrix(const RI::Tensor& m) { lmax_abfs = std::max(lmax_abfs, static_cast(abfs[T].size()) - 1); } Matrix_Orbs21 m_abfslcaos_lcaos; - ORB_gaunt_table MGT; - const int Lmax = m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, orb.get_Rmax(), lmax_abfs); - MGT.init_Gaunt_CH(Lmax); - MGT.init_Gaunt(Lmax); - m_abfslcaos_lcaos.init_radial(abfs, lcaos, lcaos, MGT); + m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, orb.get_Rmax(), lmax_abfs); + m_abfslcaos_lcaos.init_radial(abfs, lcaos, lcaos); std::map>> delta_R; for (std::size_t it = 0; it != abfs.size(); ++it) diff --git a/source/source_lcao/module_ri/Exx_LRI.h b/source/source_lcao/module_ri/Exx_LRI.h index 0ca4aff9b6..6f94f1a32b 100644 --- a/source/source_lcao/module_ri/Exx_LRI.h +++ b/source/source_lcao/module_ri/Exx_LRI.h @@ -90,7 +90,6 @@ class Exx_LRI const Exx_Info::Exx_Info_RI &info; MPI_Comm mpi_comm; const K_Vectors *p_kv = nullptr; - ORB_gaunt_table MGT; std::vector orb_cutoff_; std::vector>> lcaos; diff --git a/source/source_lcao/module_ri/Exx_LRI.hpp b/source/source_lcao/module_ri/Exx_LRI.hpp index d271a5c8ce..97e00f8528 100644 --- a/source/source_lcao/module_ri/Exx_LRI.hpp +++ b/source/source_lcao/module_ri/Exx_LRI.hpp @@ -53,14 +53,13 @@ void Exx_LRI::init(const MPI_Comm &mpi_comm_in, this->coulomb_settings = RI_Util::update_coulomb_settings(this->info.coulomb_param, ucell, this->p_kv); - bool init_MGT = true; + std::shared_ptr MGT = std::make_shared(); for(const auto &settings_list : this->coulomb_settings) { this->exx_objs[settings_list.first].abfs_ccp = Conv_Coulomb_Pot_K::cal_orbs_ccp(this->abfs, settings_list.second.second, this->info.ccp_rmesh_times); this->exx_objs[settings_list.first].cv.set_orbitals(ucell, orb, this->lcaos, this->abfs, this->exx_objs[settings_list.first].abfs_ccp, - this->info.kmesh_times, this->MGT, init_MGT, settings_list.second.first ); - init_MGT = false; // only init once + this->info.kmesh_times, MGT, settings_list.second.first ); } ModuleBase::timer::tick("Exx_LRI", "init"); diff --git a/source/source_lcao/module_ri/LRI_CV.h b/source/source_lcao/module_ri/LRI_CV.h index 210fea1afa..aaf00e8656 100644 --- a/source/source_lcao/module_ri/LRI_CV.h +++ b/source/source_lcao/module_ri/LRI_CV.h @@ -40,8 +40,7 @@ class LRI_CV const std::vector>> &abfs_in, const std::vector>> &abfs_ccp_in, const double &kmesh_times, - ORB_gaunt_table& MGT, - const bool& init_MGT, + std::shared_ptr MGT, const bool& init_C); inline std::map>> cal_Vs( diff --git a/source/source_lcao/module_ri/LRI_CV.hpp b/source/source_lcao/module_ri/LRI_CV.hpp index 46b05ccb9f..a945cfeef9 100644 --- a/source/source_lcao/module_ri/LRI_CV.hpp +++ b/source/source_lcao/module_ri/LRI_CV.hpp @@ -44,8 +44,7 @@ void LRI_CV::set_orbitals( const std::vector>> &abfs_in, const std::vector>> &abfs_ccp_in, const double &kmesh_times, - ORB_gaunt_table& MGT, - const bool& init_MGT, + std::shared_ptr MGT, const bool& init_C) { ModuleBase::TITLE("LRI_CV", "set_orbitals"); @@ -69,24 +68,17 @@ void LRI_CV::set_orbitals( range_abfs = ModuleBase::Element_Basis_Index::construct_range( abfs ); this->index_abfs = ModuleBase::Element_Basis_Index::construct_index( range_abfs ); - const int Lmax_v = this->m_abfs_abfs.init(2, ucell, orb, kmesh_times, lcaos_rmax + abfs_ccp_rmax, GlobalC::exx_info.info_ri.abfs_Lmax); - int Lmax_c = std::numeric_limits::min(); + this->m_abfs_abfs.MGT = this->m_abfslcaos_lcaos.MGT = MGT; + this->m_abfs_abfs.init(2, ucell, orb, kmesh_times, lcaos_rmax + abfs_ccp_rmax, GlobalC::exx_info.info_ri.abfs_Lmax); if (init_C) - Lmax_c = this->m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, lcaos_rmax, GlobalC::exx_info.info_ri.abfs_Lmax); - const int Lmax = std::max(Lmax_v, Lmax_c); + this->m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, lcaos_rmax, GlobalC::exx_info.info_ri.abfs_Lmax); - if (init_MGT) { - MGT.init_Gaunt_CH(Lmax); - MGT.init_Gaunt(Lmax); - } - - this->m_abfs_abfs.init_radial(this->abfs_ccp, this->abfs, MGT); + this->m_abfs_abfs.init_radial(this->abfs_ccp, this->abfs); this->m_abfs_abfs.init_radial_table(); if (init_C) { this->m_abfslcaos_lcaos.init_radial(this->abfs_ccp, this->lcaos, - this->lcaos, - MGT); + this->lcaos); this->m_abfslcaos_lcaos.init_radial_table(); } diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.cpp b/source/source_lcao/module_ri/Matrix_Orbs11.cpp index 506a34d4c0..7e293d1b6e 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs11.cpp @@ -9,7 +9,7 @@ #include "source_base/tool_title.h" #include "source_pw/module_pwdft/global.h" -int Matrix_Orbs11::init(const int mode, +void Matrix_Orbs11::init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, @@ -41,8 +41,12 @@ int Matrix_Orbs11::init(const int mode, //========================================= // (3) make Gaunt coefficients table //========================================= - // this->MGT.init_Gaunt_CH(Lmax); - // this->MGT.init_Gaunt(Lmax); + if(!this->MGT) + { this->MGT = std::make_shared(); } + if(this->MGT->get_Lmax_Gaunt_CH() < Lmax) + { this->MGT->init_Gaunt_CH(Lmax); } + if(this->MGT->get_Lmax_Gaunt_Coefficients() < Lmax) + { this->MGT->init_Gaunt(Lmax); } const double dr = orb.get_dR(); const double dk = orb.get_dk(); @@ -57,12 +61,10 @@ int Matrix_Orbs11::init(const int mode, psb_); ModuleBase::timer::tick("Matrix_Orbs11", "init"); - return Lmax; } void Matrix_Orbs11::init_radial(const std::vector>>& orb_A, - const std::vector>>& orb_B, - const ORB_gaunt_table& MGT) + const std::vector>>& orb_B) { ModuleBase::TITLE("Matrix_Orbs11", "init_radial"); ModuleBase::timer::tick("Matrix_Orbs11", "init_radial"); @@ -74,7 +76,7 @@ void Matrix_Orbs11::init_radial(const std::vectorMGT))); } } } @@ -84,7 +86,7 @@ void Matrix_Orbs11::init_radial(const std::vectorMGT))); } } } diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.h b/source/source_lcao/module_ri/Matrix_Orbs11.h index a5e3efc053..8821f029aa 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.h +++ b/source/source_lcao/module_ri/Matrix_Orbs11.h @@ -24,7 +24,7 @@ class Matrix_Orbs11 // mode: // 1: // 2: - int init( + void init( const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, @@ -33,9 +33,8 @@ class Matrix_Orbs11 const int lmax_abfs); void init_radial(const std::vector>>& orb_A, - const std::vector>>& orb_B, - const ORB_gaunt_table& MGT); - void init_radial(const LCAO_Orbitals& orb_A, const LCAO_Orbitals& orb_B, const ORB_gaunt_table& MGT); + const std::vector>>& orb_B); + void init_radial(const LCAO_Orbitals& orb_A, const LCAO_Orbitals& orb_B); void init_radial_table(); void init_radial_table(const std::map>>& Rs); // unit: ucell.lat0 @@ -69,6 +68,8 @@ class Matrix_Orbs11 const UnitCell &ucell, const ModuleBase::Element_Basis_Index::IndexLNM& index_r, const ModuleBase::Element_Basis_Index::IndexLNM& index_c) const; + + std::shared_ptr MGT; private: ModuleBase::Sph_Bessel_Recursive::D2* psb_ = nullptr; diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.cpp b/source/source_lcao/module_ri/Matrix_Orbs21.cpp index 6ce200876a..c3c0f6f554 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs21.cpp @@ -9,7 +9,7 @@ #include "source_base/tool_title.h" #include "source_pw/module_pwdft/global.h" -int Matrix_Orbs21::init( +void Matrix_Orbs21::init( const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, @@ -34,9 +34,13 @@ int Matrix_Orbs21::init( //========================================= // (3) make Gaunt coefficients table //========================================= - // this->MGT.init_Gaunt_CH(2 * Lmax + 1); // why +1 - // this->MGT.init_Gaunt(2 * Lmax + 1); - Lmax = 2 * Lmax + 1; + Lmax = 2*Lmax+1; + if(!this->MGT) + { this->MGT = std::make_shared(); } + if(this->MGT->get_Lmax_Gaunt_CH() < Lmax) + { this->MGT->init_Gaunt_CH(Lmax); } + if(this->MGT->get_Lmax_Gaunt_Coefficients() < Lmax) + { this->MGT->init_Gaunt(Lmax); } const double dr = orb.get_dR(); const double dk = orb.get_dk(); @@ -51,13 +55,11 @@ int Matrix_Orbs21::init( psb_); ModuleBase::timer::tick("Matrix_Orbs21", "init"); - return Lmax; } void Matrix_Orbs21::init_radial(const std::vector>>& orb_A1, const std::vector>>& orb_A2, - const std::vector>>& orb_B, - const ORB_gaunt_table& MGT) + const std::vector>>& orb_B) { ModuleBase::TITLE("Matrix_Orbs21", "init_radial"); ModuleBase::timer::tick("Matrix_Orbs21", "init_radial"); @@ -84,7 +86,7 @@ void Matrix_Orbs21::init_radial(const std::vectorMGT))); } } } @@ -98,8 +100,7 @@ void Matrix_Orbs21::init_radial(const std::vector>>& orb_A1, const LCAO_Orbitals& orb_A2, - const LCAO_Orbitals& orb_B, - const ORB_gaunt_table& MGT) + const LCAO_Orbitals& orb_B) { ModuleBase::TITLE("Matrix_Orbs21", "init_radial"); ModuleBase::timer::tick("Matrix_Orbs21", "init_radial"); @@ -126,7 +127,7 @@ void Matrix_Orbs21::init_radial(const std::vectorMGT))); } } } diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.h b/source/source_lcao/module_ri/Matrix_Orbs21.h index 2552b43fca..6a393f598a 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.h +++ b/source/source_lcao/module_ri/Matrix_Orbs21.h @@ -22,7 +22,7 @@ class Matrix_Orbs21 public: // mode: // 1: - int init( + void init( const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, @@ -32,12 +32,10 @@ class Matrix_Orbs21 void init_radial(const std::vector>>& orb_A1, const std::vector>>& orb_A2, - const std::vector>>& orb_B, - const ORB_gaunt_table& MGT); + const std::vector>>& orb_B); void init_radial(const std::vector>>& orb_A1, const LCAO_Orbitals& orb_A2, - const LCAO_Orbitals& orb_B, - const ORB_gaunt_table& MGT); + const LCAO_Orbitals& orb_B); void init_radial_table(); void init_radial_table(const std::map>>& Rs); // unit: ucell.lat0 @@ -78,6 +76,8 @@ class Matrix_Orbs21 const ModuleBase::Element_Basis_Index::IndexLNM& index_A1, const ModuleBase::Element_Basis_Index::IndexLNM& index_A2, const ModuleBase::Element_Basis_Index::IndexLNM& index_B) const; + + std::shared_ptr MGT; private: ModuleBase::Sph_Bessel_Recursive::D2* psb_ = nullptr; diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.cpp b/source/source_lcao/module_ri/Matrix_Orbs22.cpp index 3c1d0a1e64..cff8c29e81 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs22.cpp @@ -9,7 +9,7 @@ #include "source_base/tool_title.h" #include "source_pw/module_pwdft/global.h" -int Matrix_Orbs22::init(const int mode, +void Matrix_Orbs22::init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, @@ -33,9 +33,13 @@ int Matrix_Orbs22::init(const int mode, //========================================= // (3) make Gaunt coefficients table //========================================= - // this->MGT.init_Gaunt_CH(2 * Lmax + 1); // why +1 - // this->MGT.init_Gaunt(2 * Lmax + 1); - Lmax = 2 * Lmax + 1; + Lmax = 2*Lmax+1; + if(!this->MGT) + { this->MGT = std::make_shared(); } + if(this->MGT->get_Lmax_Gaunt_CH() < Lmax) + { this->MGT->init_Gaunt_CH(Lmax); } + if(this->MGT->get_Lmax_Gaunt_Coefficients() < Lmax) + { this->MGT->init_Gaunt(Lmax); } const double dr = orb.get_dR(); const double dk = orb.get_dk(); @@ -50,14 +54,12 @@ int Matrix_Orbs22::init(const int mode, psb_); ModuleBase::timer::tick("Matrix_Orbs22", "init"); - return Lmax; } void Matrix_Orbs22::init_radial(const std::vector>>& orb_A1, const std::vector>>& orb_A2, const std::vector>>& orb_B1, - const std::vector>>& orb_B2, - const ORB_gaunt_table& MGT) + const std::vector>>& orb_B2) { ModuleBase::TITLE("Matrix_Orbs22", "init_radial"); ModuleBase::timer::tick("Matrix_Orbs22", "init_radial"); @@ -80,15 +82,14 @@ void Matrix_Orbs22::init_radial(const std::vectorMGT))); ModuleBase::timer::tick("Matrix_Orbs22", "init_radial"); } void Matrix_Orbs22::init_radial(const LCAO_Orbitals& orb_A1, const LCAO_Orbitals& orb_A2, const LCAO_Orbitals& orb_B1, - const LCAO_Orbitals& orb_B2, - const ORB_gaunt_table& MGT) + const LCAO_Orbitals& orb_B2) { ModuleBase::TITLE("Matrix_Orbs22", "init_radial"); ModuleBase::timer::tick("Matrix_Orbs22", "init_radial"); @@ -111,7 +112,7 @@ void Matrix_Orbs22::init_radial(const LCAO_Orbitals& orb_A1, orb_B1.Phi[TB].PhiLN(LB1, NB1), orb_B2.Phi[TB].PhiLN(LB2, NB2), psb_, - MGT))); + *this->MGT))); ModuleBase::timer::tick("Matrix_Orbs22", "init_radial"); } diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.h b/source/source_lcao/module_ri/Matrix_Orbs22.h index 4e9f49f628..2ad535f4ee 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.h +++ b/source/source_lcao/module_ri/Matrix_Orbs22.h @@ -23,7 +23,7 @@ class Matrix_Orbs22 public: // mode: // 1: - int init(const int mode, + void init(const int mode, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, // extend Kcut, keep dK @@ -32,13 +32,11 @@ class Matrix_Orbs22 void init_radial(const std::vector>>& orb_A1, const std::vector>>& orb_A2, const std::vector>>& orb_B1, - const std::vector>>& , - const ORB_gaunt_table& MGT); + const std::vector>>& orb_B2); void init_radial(const LCAO_Orbitals& orb_A1, const LCAO_Orbitals& orb_A2, const LCAO_Orbitals& orb_B1, - const LCAO_Orbitals& orb_B2, - const ORB_gaunt_table& MGT); + const LCAO_Orbitals& orb_B2); void init_radial_table(); void init_radial_table(const std::map>>& Rs); // unit: ucell.lat0 @@ -100,6 +98,8 @@ class Matrix_Orbs22 const ModuleBase::Element_Basis_Index::IndexLNM& index_A2, const ModuleBase::Element_Basis_Index::IndexLNM& index_B1, const ModuleBase::Element_Basis_Index::IndexLNM& index_B2) const; + + std::shared_ptr MGT; private: ModuleBase::Sph_Bessel_Recursive::D2* psb_ = nullptr; diff --git a/source/source_lcao/module_ri/exx_opt_orb.cpp b/source/source_lcao/module_ri/exx_opt_orb.cpp index 884548688b..ab68e3e144 100644 --- a/source/source_lcao/module_ri/exx_opt_orb.cpp +++ b/source/source_lcao/module_ri/exx_opt_orb.cpp @@ -72,11 +72,8 @@ void Exx_Opt_Orb::generate_matrix( { if(judge_orbs_empty(lcaos)) { return {}; } Matrix_Orbs22 m_lcaoslcaos_lcaoslcaos; - ORB_gaunt_table MGT; - const int Lmax = m_lcaoslcaos_lcaoslcaos.init( 1, ucell,orb, info.kmesh_times, orb.get_Rmax() ); - MGT.init_Gaunt_CH(Lmax); - MGT.init_Gaunt(Lmax); - m_lcaoslcaos_lcaoslcaos.init_radial( lcaos, lcaos, lcaos, lcaos, MGT ); + m_lcaoslcaos_lcaoslcaos.init( 1, ucell,orb, info.kmesh_times, orb.get_Rmax() ); + m_lcaoslcaos_lcaoslcaos.init_radial( lcaos, lcaos, lcaos, lcaos ); #if TEST_EXX_RADIAL>=1 m_lcaoslcaos_lcaoslcaos.init_radial_table(radial_R); #else @@ -91,11 +88,8 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(lcaos)) { return {}; } if(judge_orbs_empty(jle)) { return {}; } Matrix_Orbs21 m_jyslcaos_lcaos; - ORB_gaunt_table MGT; - const int Lmax = m_jyslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - MGT.init_Gaunt_CH(Lmax); - MGT.init_Gaunt(Lmax); - m_jyslcaos_lcaos.init_radial( jle, lcaos, lcaos, MGT); + m_jyslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); + m_jyslcaos_lcaos.init_radial( jle, lcaos, lcaos); #if TEST_EXX_RADIAL>=1 m_jyslcaos_lcaos.init_radial_table( radial_R); #else @@ -109,11 +103,8 @@ void Exx_Opt_Orb::generate_matrix( { if(judge_orbs_empty(jle)) { return {}; } Matrix_Orbs11 m_jys_jys; - ORB_gaunt_table MGT; - const int Lmax = m_jys_jys.init( 2,ucell,orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - MGT.init_Gaunt_CH(Lmax); - MGT.init_Gaunt(Lmax); - m_jys_jys.init_radial( jle, jle, MGT ); + m_jys_jys.init( 2,ucell,orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); + m_jys_jys.init_radial( jle, jle ); #if TEST_EXX_RADIAL>=1 m_jys_jys.init_radial_table(radial_R); #else @@ -127,11 +118,8 @@ void Exx_Opt_Orb::generate_matrix( { if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs11 m_abfs_abfs; - ORB_gaunt_table MGT; - const int Lmax = m_abfs_abfs.init( 2, ucell, orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - MGT.init_Gaunt_CH(Lmax); - MGT.init_Gaunt(Lmax); - m_abfs_abfs.init_radial( abfs, abfs, MGT ); + m_abfs_abfs.init( 2, ucell, orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); + m_abfs_abfs.init_radial( abfs, abfs ); #if TEST_EXX_RADIAL>=1 m_abfs_abfs.init_radial_table(radial_R); #else @@ -146,11 +134,8 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(lcaos)) { return {}; } if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs21 m_abfslcaos_lcaos; - ORB_gaunt_table MGT; - const int Lmax = m_abfslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - MGT.init_Gaunt_CH(Lmax); - MGT.init_Gaunt(Lmax); - m_abfslcaos_lcaos.init_radial( abfs, lcaos, lcaos, MGT ); + m_abfslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); + m_abfslcaos_lcaos.init_radial( abfs, lcaos, lcaos ); #if TEST_EXX_RADIAL>=1 m_abfslcaos_lcaos.init_radial_table(radial_R); #else @@ -165,11 +150,8 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(jle)) { return {}; } if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs11 m_jys_abfs; - ORB_gaunt_table MGT; - const int Lmax = m_jys_abfs.init( 2, ucell,orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - MGT.init_Gaunt_CH(Lmax); - MGT.init_Gaunt(Lmax); - m_jys_abfs.init_radial( jle, abfs, MGT ); + m_jys_abfs.init( 2, ucell,orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); + m_jys_abfs.init_radial( jle, abfs ); #if TEST_EXX_RADIAL>=1 m_jys_abfs.init_radial_table(radial_R); #else From e25e2bd9c5e05bf6acef9b249dadc75cac54c8cc Mon Sep 17 00:00:00 2001 From: linpz Date: Sat, 27 Dec 2025 02:31:37 +0800 Subject: [PATCH 8/9] Refactor: Simplify Lmax in Matrix_Orbs --- .../module_ri/ABFs_Construct-PCA.cpp | 7 +- source/source_lcao/module_ri/Exx_LRI.hpp | 2 + source/source_lcao/module_ri/LRI_CV.hpp | 12 +-- .../source_lcao/module_ri/Matrix_Orbs11.cpp | 51 ++++-------- source/source_lcao/module_ri/Matrix_Orbs11.h | 15 +--- .../source_lcao/module_ri/Matrix_Orbs21.cpp | 81 +++++++------------ source/source_lcao/module_ri/Matrix_Orbs21.h | 16 +--- .../source_lcao/module_ri/Matrix_Orbs22.cpp | 81 +++++++++---------- source/source_lcao/module_ri/Matrix_Orbs22.h | 25 +++--- source/source_lcao/module_ri/exx_abfs.h | 7 ++ source/source_lcao/module_ri/exx_opt_orb.cpp | 24 ++---- 11 files changed, 118 insertions(+), 203 deletions(-) diff --git a/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp b/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp index 67b89200f2..e0bb58391f 100644 --- a/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp +++ b/source/source_lcao/module_ri/ABFs_Construct-PCA.cpp @@ -172,13 +172,8 @@ RI::Tensor get_column_mean0_matrix(const RI::Tensor& m) const ModuleBase::Element_Basis_Index::IndexLNM index_abfs = ModuleBase::Element_Basis_Index::construct_index(range_abfs); - int lmax_abfs = std::numeric_limits::min(); - for (std::size_t T = 0; T != abfs.size(); ++T) - { lmax_abfs = std::max(lmax_abfs, static_cast(abfs[T].size()) - 1); } - Matrix_Orbs21 m_abfslcaos_lcaos; - m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, orb.get_Rmax(), lmax_abfs); - m_abfslcaos_lcaos.init_radial(abfs, lcaos, lcaos); + m_abfslcaos_lcaos.init(abfs, lcaos, lcaos, ucell, orb, kmesh_times, orb.get_Rmax()); std::map>> delta_R; for (std::size_t it = 0; it != abfs.size(); ++it) diff --git a/source/source_lcao/module_ri/Exx_LRI.hpp b/source/source_lcao/module_ri/Exx_LRI.hpp index 97e00f8528..57cd8440d7 100644 --- a/source/source_lcao/module_ri/Exx_LRI.hpp +++ b/source/source_lcao/module_ri/Exx_LRI.hpp @@ -39,6 +39,7 @@ void Exx_LRI::init(const MPI_Comm &mpi_comm_in, this->orb_cutoff_ = orb.cutoffs(); this->lcaos = Exx_Abfs::Construct_Orbs::change_orbs( orb, this->info.kmesh_times ); + Exx_Abfs::Construct_Orbs::filter_empty_orbs(this->lcaos); const std::vector>> abfs_same_atom = Exx_Abfs::Construct_Orbs::abfs_same_atom(ucell, orb, this->lcaos, this->info.kmesh_times, this->info.pca_threshold ); @@ -46,6 +47,7 @@ void Exx_LRI::init(const MPI_Comm &mpi_comm_in, { this->abfs = abfs_same_atom;} else { this->abfs = Exx_Abfs::IO::construct_abfs( abfs_same_atom, orb, this->info.files_abfs, this->info.kmesh_times ); } + Exx_Abfs::Construct_Orbs::filter_empty_orbs(this->abfs); Exx_Abfs::Construct_Orbs::print_orbs_size(ucell, this->abfs, GlobalV::ofs_running); for( size_t T=0; T!=this->abfs.size(); ++T ) diff --git a/source/source_lcao/module_ri/LRI_CV.hpp b/source/source_lcao/module_ri/LRI_CV.hpp index a945cfeef9..60670d5088 100644 --- a/source/source_lcao/module_ri/LRI_CV.hpp +++ b/source/source_lcao/module_ri/LRI_CV.hpp @@ -69,16 +69,16 @@ void LRI_CV::set_orbitals( this->index_abfs = ModuleBase::Element_Basis_Index::construct_index( range_abfs ); this->m_abfs_abfs.MGT = this->m_abfslcaos_lcaos.MGT = MGT; - this->m_abfs_abfs.init(2, ucell, orb, kmesh_times, lcaos_rmax + abfs_ccp_rmax, GlobalC::exx_info.info_ri.abfs_Lmax); + this->m_abfs_abfs.init( + this->abfs_ccp, this->abfs, + ucell, orb, kmesh_times, lcaos_rmax + abfs_ccp_rmax); if (init_C) - this->m_abfslcaos_lcaos.init(1, ucell, orb, kmesh_times, lcaos_rmax, GlobalC::exx_info.info_ri.abfs_Lmax); + this->m_abfslcaos_lcaos.init( + this->abfs_ccp, this->lcaos, this->lcaos, + ucell, orb, kmesh_times, lcaos_rmax); - this->m_abfs_abfs.init_radial(this->abfs_ccp, this->abfs); this->m_abfs_abfs.init_radial_table(); if (init_C) { - this->m_abfslcaos_lcaos.init_radial(this->abfs_ccp, - this->lcaos, - this->lcaos); this->m_abfslcaos_lcaos.init_radial_table(); } diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.cpp b/source/source_lcao/module_ri/Matrix_Orbs11.cpp index 7e293d1b6e..942f5a4cc6 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs11.cpp @@ -9,34 +9,21 @@ #include "source_base/tool_title.h" #include "source_pw/module_pwdft/global.h" -void Matrix_Orbs11::init(const int mode, - const UnitCell& ucell, - const LCAO_Orbitals& orb, - const double kmesh_times, - const double rmax, - const int lmax_abfs) +void Matrix_Orbs11::init( + const std::vector>>& orb_A, + const std::vector>>& orb_B, + const UnitCell& ucell, + const LCAO_Orbitals& orb, + const double kmesh_times, + const double rmax) { ModuleBase::TITLE("Matrix_Orbs11", "init"); ModuleBase::timer::tick("Matrix_Orbs11", "init"); this->lat0 = &ucell.lat0; - //const int ntype = orb.get_ntype(); - //int lmax_orb = -1, lmax_beta = -1; - //for (int it = 0; it < ntype; it++) - //{ - // lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); - // lmax_beta = std::max(lmax_beta, ucell.infoNL.Beta[it].getLmax()); - //} - int Lmax, Lmax_used; - //if(mode==1) - // { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_1(lmax_orb, lmax_beta); } - if(mode==2) - { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_2(lmax_abfs); } - //else if(mode==3) - // { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_2_3(lmax_orb); } - else - { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } + const int Lmax = std::max({ Exx_Abfs::get_Lmax(orb_A), Exx_Abfs::get_Lmax(orb_B) }); + const int Lmax_used = Exx_Abfs::get_Lmax(orb_A) + Exx_Abfs::get_Lmax(orb_B); //========================================= // (3) make Gaunt coefficients table @@ -60,14 +47,6 @@ void Matrix_Orbs11::init(const int mode, Rmesh, psb_); - ModuleBase::timer::tick("Matrix_Orbs11", "init"); -} - -void Matrix_Orbs11::init_radial(const std::vector>>& orb_A, - const std::vector>>& orb_B) -{ - ModuleBase::TITLE("Matrix_Orbs11", "init_radial"); - ModuleBase::timer::tick("Matrix_Orbs11", "init_radial"); for (size_t TA = 0; TA != orb_A.size(); ++TA) { for (size_t TB = 0; TB != orb_B.size(); ++TB) { for (int LA = 0; LA != orb_A[TA].size(); ++LA) { @@ -77,15 +56,12 @@ void Matrix_Orbs11::init_radial(const std::vectorMGT))); - } - } - } - } - } - } - ModuleBase::timer::tick("Matrix_Orbs11", "init_radial"); + }}}}}} + + ModuleBase::timer::tick("Matrix_Orbs11", "init"); } +/* void Matrix_Orbs11::init_radial(const LCAO_Orbitals& orb_A, const LCAO_Orbitals& orb_B) { ModuleBase::TITLE("Matrix_Orbs11", "init_radial"); @@ -110,6 +86,7 @@ void Matrix_Orbs11::init_radial(const LCAO_Orbitals& orb_A, const LCAO_Orbitals& } ModuleBase::timer::tick("Matrix_Orbs11", "init_radial"); } +*/ void Matrix_Orbs11::init_radial_table() { diff --git a/source/source_lcao/module_ri/Matrix_Orbs11.h b/source/source_lcao/module_ri/Matrix_Orbs11.h index 8821f029aa..069127c355 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs11.h +++ b/source/source_lcao/module_ri/Matrix_Orbs11.h @@ -21,20 +21,13 @@ class Matrix_Orbs11 { public: - // mode: - // 1: - // 2: void init( - const int mode, + const std::vector>>& orb_A, + const std::vector>>& orb_B, const UnitCell& ucell, const LCAO_Orbitals& orb, - const double kmesh_times, // extend Kcut, keep dK - const double rmax, - const int lmax_abfs); - - void init_radial(const std::vector>>& orb_A, - const std::vector>>& orb_B); - void init_radial(const LCAO_Orbitals& orb_A, const LCAO_Orbitals& orb_B); + const double kmesh_times, // extend Kcut, keep dK + const double rmax); // extend Rcut, keep dR void init_radial_table(); void init_radial_table(const std::map>>& Rs); // unit: ucell.lat0 diff --git a/source/source_lcao/module_ri/Matrix_Orbs21.cpp b/source/source_lcao/module_ri/Matrix_Orbs21.cpp index c3c0f6f554..a11a1333b5 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs21.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs21.cpp @@ -10,31 +10,26 @@ #include "source_pw/module_pwdft/global.h" void Matrix_Orbs21::init( - const int mode, + const std::vector>>& orb_A1, + const std::vector>>& orb_A2, + const std::vector>>& orb_B, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, - const double rmax, - const int lmax_abfs) + const double rmax) { ModuleBase::TITLE("Matrix_Orbs21", "init"); ModuleBase::timer::tick("Matrix_Orbs21", "init"); this->lat0 = &ucell.lat0; - const int ntype = orb.get_ntype(); - int lmax_orb = -1; - for (int it = 0; it < ntype; it++) - { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); } - int Lmax, Lmax_used; - if(mode==1) - { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_3_1(lmax_abfs, lmax_orb); } - else - { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } + const int Lmax = std::max({ + Exx_Abfs::get_Lmax(orb_A1) + Exx_Abfs::get_Lmax(orb_A2), + Exx_Abfs::get_Lmax(orb_B) }); + const int Lmax_used = Exx_Abfs::get_Lmax(orb_A1) + Exx_Abfs::get_Lmax(orb_A2) + Exx_Abfs::get_Lmax(orb_B) + 1; //========================================= // (3) make Gaunt coefficients table //========================================= - Lmax = 2*Lmax+1; if(!this->MGT) { this->MGT = std::make_shared(); } if(this->MGT->get_Lmax_Gaunt_CH() < Lmax) @@ -54,50 +49,29 @@ void Matrix_Orbs21::init( Rmesh, psb_); - ModuleBase::timer::tick("Matrix_Orbs21", "init"); -} - -void Matrix_Orbs21::init_radial(const std::vector>>& orb_A1, - const std::vector>>& orb_A2, - const std::vector>>& orb_B) -{ - ModuleBase::TITLE("Matrix_Orbs21", "init_radial"); - ModuleBase::timer::tick("Matrix_Orbs21", "init_radial"); assert(orb_A1.size() == orb_A2.size()); - for (size_t TA = 0; TA != orb_A1.size(); ++TA) - { - for (size_t TB = 0; TB != orb_B.size(); ++TB) - { - for (int LA1 = 0; LA1 != orb_A1[TA].size(); ++LA1) - { - for (size_t NA1 = 0; NA1 != orb_A1[TA][LA1].size(); ++NA1) - { - for (int LA2 = 0; LA2 != orb_A2[TA].size(); ++LA2) - { - for (size_t NA2 = 0; NA2 != orb_A2[TA][LA2].size(); ++NA2) - { - for (int LB = 0; LB != orb_B[TB].size(); ++LB) - { - for (size_t NB = 0; NB != orb_B[TB][LB].size(); ++NB) - { + for (size_t TA = 0; TA != orb_A1.size(); ++TA) { + for (size_t TB = 0; TB != orb_B.size(); ++TB) { + for (int LA1 = 0; LA1 != orb_A1[TA].size(); ++LA1) { + for (size_t NA1 = 0; NA1 != orb_A1[TA][LA1].size(); ++NA1) { + for (int LA2 = 0; LA2 != orb_A2[TA].size(); ++LA2) { + for (size_t NA2 = 0; NA2 != orb_A2[TA][LA2].size(); ++NA2) { + for (int LB = 0; LB != orb_B[TB].size(); ++LB) { + for (size_t NB = 0; NB != orb_B[TB][LB].size(); ++NB) { center2_orb21_s[TA][TB][LA1][NA1][LA2][NA2][LB].insert( - std::make_pair(NB, - Center2_Orb::Orb21(orb_A1[TA][LA1][NA1], - orb_A2[TA][LA2][NA2], - orb_B[TB][LB][NB], - psb_, - *this->MGT))); - } - } - } - } - } - } - } - } - ModuleBase::timer::tick("Matrix_Orbs21", "init_radial"); + std::make_pair( + NB, + Center2_Orb::Orb21( + orb_A1[TA][LA1][NA1], + orb_A2[TA][LA2][NA2], + orb_B[TB][LB][NB], + psb_, + *this->MGT))); + }}}}}}}} + ModuleBase::timer::tick("Matrix_Orbs21", "init"); } +/* void Matrix_Orbs21::init_radial(const std::vector>>& orb_A1, const LCAO_Orbitals& orb_A2, const LCAO_Orbitals& orb_B) @@ -138,6 +112,7 @@ void Matrix_Orbs21::init_radial(const std::vector void init( - const int mode, + const std::vector>>& orb_A1, + const std::vector>>& orb_A2, + const std::vector>>& orb_B, const UnitCell& ucell, const LCAO_Orbitals& orb, const double kmesh_times, // extend Kcut, keep dK - const double rmax, // extend Rcut, keep dR - const int lmax_abfs); - - void init_radial(const std::vector>>& orb_A1, - const std::vector>>& orb_A2, - const std::vector>>& orb_B); - void init_radial(const std::vector>>& orb_A1, - const LCAO_Orbitals& orb_A2, - const LCAO_Orbitals& orb_B); + const double rmax); // extend Rcut, keep dR void init_radial_table(); void init_radial_table(const std::map>>& Rs); // unit: ucell.lat0 diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.cpp b/source/source_lcao/module_ri/Matrix_Orbs22.cpp index cff8c29e81..69ea735fbc 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.cpp +++ b/source/source_lcao/module_ri/Matrix_Orbs22.cpp @@ -9,31 +9,29 @@ #include "source_base/tool_title.h" #include "source_pw/module_pwdft/global.h" -void Matrix_Orbs22::init(const int mode, - const UnitCell& ucell, - const LCAO_Orbitals& orb, - const double kmesh_times, - const double rmax) +void Matrix_Orbs22::init( + const std::vector>>& orb_A1, + const std::vector>>& orb_A2, + const std::vector>>& orb_B1, + const std::vector>>& orb_B2, + const UnitCell& ucell, + const LCAO_Orbitals& orb, + const double kmesh_times, + const double rmax) { ModuleBase::TITLE("Matrix_Orbs22", "init"); ModuleBase::timer::tick("Matrix_Orbs22", "init"); - this->lat0 = &ucell.lat0; + this->lat0 = &ucell.lat0; - const int ntype = orb.get_ntype(); - int lmax_orb = -1; - for (int it = 0; it < ntype; it++) - { lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax()); } - int Lmax, Lmax_used; - if(mode==1) - { std::tie(Lmax_used, Lmax) = Center2_Orb::init_Lmax_4_1(lmax_orb); } - else - { throw std::invalid_argument("mode = "+std::to_string(mode)+"in file "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); } + const int Lmax = std::max({ + Exx_Abfs::get_Lmax(orb_A1) + Exx_Abfs::get_Lmax(orb_A2), + Exx_Abfs::get_Lmax(orb_B1) + Exx_Abfs::get_Lmax(orb_B2) }); + const int Lmax_used = Exx_Abfs::get_Lmax(orb_A1) + Exx_Abfs::get_Lmax(orb_A2) + Exx_Abfs::get_Lmax(orb_B1) + Exx_Abfs::get_Lmax(orb_B2) + 2; //========================================= // (3) make Gaunt coefficients table //========================================= - Lmax = 2*Lmax+1; if(!this->MGT) { this->MGT = std::make_shared(); } if(this->MGT->get_Lmax_Gaunt_CH() < Lmax) @@ -53,39 +51,33 @@ void Matrix_Orbs22::init(const int mode, Rmesh, psb_); - ModuleBase::timer::tick("Matrix_Orbs22", "init"); -} - -void Matrix_Orbs22::init_radial(const std::vector>>& orb_A1, - const std::vector>>& orb_A2, - const std::vector>>& orb_B1, - const std::vector>>& orb_B2) -{ - ModuleBase::TITLE("Matrix_Orbs22", "init_radial"); - ModuleBase::timer::tick("Matrix_Orbs22", "init_radial"); assert(orb_A1.size() == orb_A2.size()); assert(orb_B1.size() == orb_B2.size()); - for (size_t TA = 0; TA != orb_A1.size(); ++TA) - for (size_t TB = 0; TB != orb_B1.size(); ++TB) - for (int LA1 = 0; LA1 != orb_A1[TA].size(); ++LA1) - for (size_t NA1 = 0; NA1 != orb_A1[TA][LA1].size(); ++NA1) - for (int LA2 = 0; LA2 != orb_A2[TA].size(); ++LA2) - for (size_t NA2 = 0; NA2 != orb_A2[TA][LA2].size(); ++NA2) - for (int LB1 = 0; LB1 != orb_B1[TB].size(); ++LB1) - for (size_t NB1 = 0; NB1 != orb_B1[TB][LB1].size(); ++NB1) - for (int LB2 = 0; LB2 != orb_B2[TB].size(); ++LB2) - for (size_t NB2 = 0; NB2 != orb_B2[TB][LB2].size(); ++NB2) + for (size_t TA = 0; TA != orb_A1.size(); ++TA) { + for (size_t TB = 0; TB != orb_B1.size(); ++TB) { + for (int LA1 = 0; LA1 != orb_A1[TA].size(); ++LA1) { + for (size_t NA1 = 0; NA1 != orb_A1[TA][LA1].size(); ++NA1) { + for (int LA2 = 0; LA2 != orb_A2[TA].size(); ++LA2) { + for (size_t NA2 = 0; NA2 != orb_A2[TA][LA2].size(); ++NA2) { + for (int LB1 = 0; LB1 != orb_B1[TB].size(); ++LB1) { + for (size_t NB1 = 0; NB1 != orb_B1[TB][LB1].size(); ++NB1) { + for (int LB2 = 0; LB2 != orb_B2[TB].size(); ++LB2) { + for (size_t NB2 = 0; NB2 != orb_B2[TB][LB2].size(); ++NB2) { center2_orb22_s[TA][TB][LA1][NA1][LA2][NA2][LB1][NB1][LB2].insert( - std::make_pair(NB2, - Center2_Orb::Orb22(orb_A1[TA][LA1][NA1], - orb_A2[TA][LA2][NA2], - orb_B1[TB][LB1][NB1], - orb_B2[TB][LB2][NB2], - psb_, - *this->MGT))); - ModuleBase::timer::tick("Matrix_Orbs22", "init_radial"); + std::make_pair( + NB2, + Center2_Orb::Orb22( + orb_A1[TA][LA1][NA1], + orb_A2[TA][LA2][NA2], + orb_B1[TB][LB1][NB1], + orb_B2[TB][LB2][NB2], + psb_, + *this->MGT))); + }}}}}}}}}} + ModuleBase::timer::tick("Matrix_Orbs22", "init"); } +/* void Matrix_Orbs22::init_radial(const LCAO_Orbitals& orb_A1, const LCAO_Orbitals& orb_A2, const LCAO_Orbitals& orb_B1, @@ -115,6 +107,7 @@ void Matrix_Orbs22::init_radial(const LCAO_Orbitals& orb_A1, *this->MGT))); ModuleBase::timer::tick("Matrix_Orbs22", "init_radial"); } +*/ void Matrix_Orbs22::init_radial_table() { diff --git a/source/source_lcao/module_ri/Matrix_Orbs22.h b/source/source_lcao/module_ri/Matrix_Orbs22.h index 2ad535f4ee..e0c3790289 100644 --- a/source/source_lcao/module_ri/Matrix_Orbs22.h +++ b/source/source_lcao/module_ri/Matrix_Orbs22.h @@ -21,22 +21,15 @@ class Matrix_Orbs22 { public: - // mode: - // 1: - void init(const int mode, - const UnitCell& ucell, - const LCAO_Orbitals& orb, - const double kmesh_times, // extend Kcut, keep dK - const double rmax); // extend Rcut, keep dR - - void init_radial(const std::vector>>& orb_A1, - const std::vector>>& orb_A2, - const std::vector>>& orb_B1, - const std::vector>>& orb_B2); - void init_radial(const LCAO_Orbitals& orb_A1, - const LCAO_Orbitals& orb_A2, - const LCAO_Orbitals& orb_B1, - const LCAO_Orbitals& orb_B2); + void init( + const std::vector>>& orb_A1, + const std::vector>>& orb_A2, + const std::vector>>& orb_B1, + const std::vector>>& orb_B2, + const UnitCell& ucell, + const LCAO_Orbitals& orb, + const double kmesh_times, // extend Kcut, keep dK + const double rmax); // extend Rcut, keep dR void init_radial_table(); void init_radial_table(const std::map>>& Rs); // unit: ucell.lat0 diff --git a/source/source_lcao/module_ri/exx_abfs.h b/source/source_lcao/module_ri/exx_abfs.h index f08462c94e..42cc41a01d 100644 --- a/source/source_lcao/module_ri/exx_abfs.h +++ b/source/source_lcao/module_ri/exx_abfs.h @@ -23,6 +23,13 @@ class Exx_Abfs int rmesh_times = 5; // Peize Lin test int kmesh_times = 1; // Peize Lin test + static int get_Lmax(const std::vector>>& orb) + { + int Lmax = -1; + for( const auto &orb_T : orb ) + { Lmax = std::max( Lmax, static_cast(orb_T.size())-1 ); } + return Lmax; + } }; #endif diff --git a/source/source_lcao/module_ri/exx_opt_orb.cpp b/source/source_lcao/module_ri/exx_opt_orb.cpp index ab68e3e144..83ea26c18d 100644 --- a/source/source_lcao/module_ri/exx_opt_orb.cpp +++ b/source/source_lcao/module_ri/exx_opt_orb.cpp @@ -47,12 +47,6 @@ void Exx_Opt_Orb::generate_matrix( { jle = Exx_Abfs::IO::construct_abfs( jle, orb, info.files_jles, info.kmesh_times ); } Exx_Abfs::Construct_Orbs::filter_empty_orbs(jle); - int lmax_abfs = -1; - for(const auto &orb_T : abfs) - { lmax_abfs = std::max( lmax_abfs, static_cast(orb_T.size())-1 ); } - for(const auto &orb_T : jle) - { lmax_abfs = std::max( lmax_abfs, static_cast(orb_T.size())-1 ); } - const ModuleBase::Element_Basis_Index::Range range_lcaos = ModuleBase::Element_Basis_Index::construct_range( lcaos ); const ModuleBase::Element_Basis_Index::IndexLNM index_lcaos = ModuleBase::Element_Basis_Index::construct_index( range_lcaos ); @@ -72,8 +66,7 @@ void Exx_Opt_Orb::generate_matrix( { if(judge_orbs_empty(lcaos)) { return {}; } Matrix_Orbs22 m_lcaoslcaos_lcaoslcaos; - m_lcaoslcaos_lcaoslcaos.init( 1, ucell,orb, info.kmesh_times, orb.get_Rmax() ); - m_lcaoslcaos_lcaoslcaos.init_radial( lcaos, lcaos, lcaos, lcaos ); + m_lcaoslcaos_lcaoslcaos.init( lcaos, lcaos, lcaos, lcaos, ucell,orb, info.kmesh_times, orb.get_Rmax() ); #if TEST_EXX_RADIAL>=1 m_lcaoslcaos_lcaoslcaos.init_radial_table(radial_R); #else @@ -88,8 +81,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(lcaos)) { return {}; } if(judge_orbs_empty(jle)) { return {}; } Matrix_Orbs21 m_jyslcaos_lcaos; - m_jyslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - m_jyslcaos_lcaos.init_radial( jle, lcaos, lcaos); + m_jyslcaos_lcaos.init( jle, lcaos, lcaos, ucell , orb, info.kmesh_times, orb.get_Rmax() ); #if TEST_EXX_RADIAL>=1 m_jyslcaos_lcaos.init_radial_table( radial_R); #else @@ -103,8 +95,7 @@ void Exx_Opt_Orb::generate_matrix( { if(judge_orbs_empty(jle)) { return {}; } Matrix_Orbs11 m_jys_jys; - m_jys_jys.init( 2,ucell,orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - m_jys_jys.init_radial( jle, jle ); + m_jys_jys.init( jle, jle, ucell,orb, info.kmesh_times, orb.get_Rmax() ); #if TEST_EXX_RADIAL>=1 m_jys_jys.init_radial_table(radial_R); #else @@ -118,8 +109,7 @@ void Exx_Opt_Orb::generate_matrix( { if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs11 m_abfs_abfs; - m_abfs_abfs.init( 2, ucell, orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - m_abfs_abfs.init_radial( abfs, abfs ); + m_abfs_abfs.init( abfs, abfs, ucell, orb, info.kmesh_times, orb.get_Rmax() ); #if TEST_EXX_RADIAL>=1 m_abfs_abfs.init_radial_table(radial_R); #else @@ -134,8 +124,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(lcaos)) { return {}; } if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs21 m_abfslcaos_lcaos; - m_abfslcaos_lcaos.init( 1, ucell , orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - m_abfslcaos_lcaos.init_radial( abfs, lcaos, lcaos ); + m_abfslcaos_lcaos.init( abfs, lcaos, lcaos, ucell , orb, info.kmesh_times, orb.get_Rmax() ); #if TEST_EXX_RADIAL>=1 m_abfslcaos_lcaos.init_radial_table(radial_R); #else @@ -150,8 +139,7 @@ void Exx_Opt_Orb::generate_matrix( if(judge_orbs_empty(jle)) { return {}; } if(judge_orbs_empty(abfs)) { return {}; } Matrix_Orbs11 m_jys_abfs; - m_jys_abfs.init( 2, ucell,orb, info.kmesh_times, orb.get_Rmax(), lmax_abfs ); - m_jys_abfs.init_radial( jle, abfs ); + m_jys_abfs.init( jle, abfs, ucell,orb, info.kmesh_times, orb.get_Rmax() ); #if TEST_EXX_RADIAL>=1 m_jys_abfs.init_radial_table(radial_R); #else From b6beea283a2e5d51b043d41ae8b6571757c9e4df Mon Sep 17 00:00:00 2001 From: linpz Date: Sat, 27 Dec 2025 02:52:14 +0800 Subject: [PATCH 9/9] Refactor: remove global variable GlobalC::exx_info.info_ri.abfs_Lmax --- source/source_hamilt/module_xc/exx_info.h | 4 +--- source/source_lcao/module_ri/Exx_LRI.hpp | 3 --- source/source_lcao/module_ri/Exx_LRI_interface.hpp | 2 +- source/source_lcao/module_ri/RPA_LRI.hpp | 2 +- 4 files changed, 3 insertions(+), 8 deletions(-) diff --git a/source/source_hamilt/module_xc/exx_info.h b/source/source_hamilt/module_xc/exx_info.h index 9cbd51c18f..5cca433b94 100644 --- a/source/source_hamilt/module_xc/exx_info.h +++ b/source/source_hamilt/module_xc/exx_info.h @@ -69,8 +69,6 @@ struct Exx_Info double kmesh_times = 4; double Cs_inv_thr = -1; - int abfs_Lmax = 0; // tmp - Exx_Info_RI(const Exx_Info::Exx_Info_Global& info_global) : coulomb_param(info_global.coulomb_param) { @@ -80,7 +78,7 @@ struct Exx_Info struct Exx_Info_Opt_ABFs { - int abfs_Lmax = 0; // tmp + int abfs_Lmax = 0; double ecut_exx = 60; double tolerence = 1E-12; std::vector files_jles; diff --git a/source/source_lcao/module_ri/Exx_LRI.hpp b/source/source_lcao/module_ri/Exx_LRI.hpp index 57cd8440d7..ee8c62a718 100644 --- a/source/source_lcao/module_ri/Exx_LRI.hpp +++ b/source/source_lcao/module_ri/Exx_LRI.hpp @@ -50,9 +50,6 @@ void Exx_LRI::init(const MPI_Comm &mpi_comm_in, Exx_Abfs::Construct_Orbs::filter_empty_orbs(this->abfs); Exx_Abfs::Construct_Orbs::print_orbs_size(ucell, this->abfs, GlobalV::ofs_running); - for( size_t T=0; T!=this->abfs.size(); ++T ) - { GlobalC::exx_info.info_ri.abfs_Lmax = std::max( GlobalC::exx_info.info_ri.abfs_Lmax, static_cast(this->abfs[T].size())-1 ); } - this->coulomb_settings = RI_Util::update_coulomb_settings(this->info.coulomb_param, ucell, this->p_kv); std::shared_ptr MGT = std::make_shared(); diff --git a/source/source_lcao/module_ri/Exx_LRI_interface.hpp b/source/source_lcao/module_ri/Exx_LRI_interface.hpp index 6f133c05df..0494dce96b 100644 --- a/source/source_lcao/module_ri/Exx_LRI_interface.hpp +++ b/source/source_lcao/module_ri/Exx_LRI_interface.hpp @@ -111,7 +111,7 @@ void Exx_LRI_Interface::exx_before_all_runners( this->symrot_.find_irreducible_sector( ucell.symm, ucell.atoms, ucell.st, RI_Util::get_Born_von_Karmen_cells(period), period, ucell.lat); - this->symrot_.set_abfs_Lmax(GlobalC::exx_info.info_ri.abfs_Lmax); + this->symrot_.set_abfs_Lmax(Exx_Abfs::get_Lmax(this->exx_ptr->abfs)); this->symrot_.cal_Ms(kv, ucell, pv); } } diff --git a/source/source_lcao/module_ri/RPA_LRI.hpp b/source/source_lcao/module_ri/RPA_LRI.hpp index d17479b239..641e07a7f4 100644 --- a/source/source_lcao/module_ri/RPA_LRI.hpp +++ b/source/source_lcao/module_ri/RPA_LRI.hpp @@ -95,7 +95,7 @@ void RPA_LRI::cal_postSCF_exx(const elecstate::DensityMatrix symrot.find_irreducible_sector(ucell.symm, ucell.atoms, ucell.st, Rs, period, ucell.lat); // set Lmax of the rotation matrices to max(l_ao, l_abf), to support rotation under ABF - symrot.set_abfs_Lmax(GlobalC::exx_info.info_ri.abfs_Lmax); + symrot.set_abfs_Lmax(Exx_Abfs::get_Lmax(exx_lri_rpa.abfs)); symrot.cal_Ms(kv, ucell, *dm.get_paraV_pointer()); // output Ts (symrot_R.txt) and Ms (symrot_k.txt) ModuleSymmetry::print_symrot_info_R(symrot, ucell.symm, ucell.lmax, Rs);