Skip to content

Commit 80eb6a1

Browse files
committed
update reference of lat0 in ri
1 parent 704f0b8 commit 80eb6a1

File tree

9 files changed

+115
-75
lines changed

9 files changed

+115
-75
lines changed

source/module_ri/Matrix_Orbs11.cpp

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -19,7 +19,7 @@ void Matrix_Orbs11::init(const int mode,
1919
ModuleBase::timer::tick("Matrix_Orbs11", "init");
2020

2121
int Lmax_used, Lmax;
22-
this->lat0 = ucell.lat0;
22+
this->lat0 = &ucell.lat0;
2323
const int ntype = orb.get_ntype();
2424
int lmax_orb = -1, lmax_beta = -1;
2525
for (int it = 0; it < ntype; it++)
@@ -127,6 +127,7 @@ void Matrix_Orbs11::init_radial_table(const std::map<size_t, std::map<size_t, st
127127
{
128128
ModuleBase::TITLE("Matrix_Orbs11", "init_radial_table_Rs");
129129
ModuleBase::timer::tick("Matrix_Orbs11", "init_radial_table");
130+
const double lat0 = *this->lat0;
130131
for (const auto& RsA: Rs) {
131132
for (const auto& RsB: RsA.second)
132133
{

source/module_ri/Matrix_Orbs11.h

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -71,7 +71,7 @@ class Matrix_Orbs11
7171
ModuleBase::Sph_Bessel_Recursive::D2* psb_ = nullptr;
7272
ORB_gaunt_table MGT;
7373
const double lcao_dr_ = 0.01;
74-
double lat0 = 0.0; // restore ucell.lat0
74+
double* lat0=nullptr; // restore ucell.lat0
7575
std::map<size_t, // TA
7676
std::map<size_t, // TB
7777
std::map<int, // LA

source/module_ri/Matrix_Orbs11.hpp

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -21,6 +21,7 @@ RI::Tensor<Tdata> Matrix_Orbs11::cal_overlap_matrix(
2121
const Matrix_Order &matrix_order) const
2222
{
2323
RI::Tensor<Tdata> m;
24+
const double lat0 = *this->lat0;
2425
const size_t sizeA = index_A[TA].count_size;
2526
const size_t sizeB = index_B[TB].count_size;
2627
switch(matrix_order)
@@ -75,6 +76,7 @@ std::array<RI::Tensor<Tdata>,3> Matrix_Orbs11::cal_grad_overlap_matrix(
7576
const Matrix_Order &matrix_order) const
7677
{
7778
std::array<RI::Tensor<Tdata>,3> m;
79+
const double lat0 = *this->lat0;
7880
const size_t sizeA = index_A[TA].count_size;
7981
const size_t sizeB = index_B[TB].count_size;
8082
for(int i=0; i<m.size(); ++i)

source/module_ri/Matrix_Orbs21.cpp

Lines changed: 102 additions & 70 deletions
Original file line numberDiff line numberDiff line change
@@ -21,7 +21,7 @@ void Matrix_Orbs21::init(const int mode,
2121

2222
const int ntype = orb.get_ntype();
2323
int lmax_orb = -1, lmax_beta = -1;
24-
this->lat0 = ucell.lat0;
24+
this->lat0 = &ucell.lat0;
2525
for (int it = 0; it < ntype; it++)
2626
{
2727
lmax_orb = std::max(lmax_orb, orb.Phi[it].getLmax());
@@ -62,29 +62,37 @@ void Matrix_Orbs21::init_radial(const std::vector<std::vector<std::vector<Numeri
6262
ModuleBase::TITLE("Matrix_Orbs21", "init_radial");
6363
ModuleBase::timer::tick("Matrix_Orbs21", "init_radial");
6464
assert(orb_A1.size() == orb_A2.size());
65-
for (size_t TA = 0; TA != orb_A1.size(); ++TA) {
66-
for (size_t TB = 0; TB != orb_B.size(); ++TB) {
67-
for (int LA1 = 0; LA1 != orb_A1[TA].size(); ++LA1) {
68-
for (size_t NA1 = 0; NA1 != orb_A1[TA][LA1].size(); ++NA1) {
69-
for (int LA2 = 0; LA2 != orb_A2[TA].size(); ++LA2) {
70-
for (size_t NA2 = 0; NA2 != orb_A2[TA][LA2].size(); ++NA2) {
71-
for (int LB = 0; LB != orb_B[TB].size(); ++LB) {
72-
for (size_t NB = 0; NB != orb_B[TB][LB].size(); ++NB) {
65+
for (size_t TA = 0; TA != orb_A1.size(); ++TA)
66+
{
67+
for (size_t TB = 0; TB != orb_B.size(); ++TB)
68+
{
69+
for (int LA1 = 0; LA1 != orb_A1[TA].size(); ++LA1)
70+
{
71+
for (size_t NA1 = 0; NA1 != orb_A1[TA][LA1].size(); ++NA1)
72+
{
73+
for (int LA2 = 0; LA2 != orb_A2[TA].size(); ++LA2)
74+
{
75+
for (size_t NA2 = 0; NA2 != orb_A2[TA][LA2].size(); ++NA2)
76+
{
77+
for (int LB = 0; LB != orb_B[TB].size(); ++LB)
78+
{
79+
for (size_t NB = 0; NB != orb_B[TB][LB].size(); ++NB)
80+
{
7381
center2_orb21_s[TA][TB][LA1][NA1][LA2][NA2][LB].insert(
7482
std::make_pair(NB,
7583
Center2_Orb::Orb21(orb_A1[TA][LA1][NA1],
7684
orb_A2[TA][LA2][NA2],
7785
orb_B[TB][LB][NB],
7886
psb_,
79-
this->MGT)));
80-
}
81-
}
82-
}
83-
}
84-
}
85-
}
86-
}
87-
}
87+
this->MGT)));
88+
}
89+
}
90+
}
91+
}
92+
}
93+
}
94+
}
95+
}
8896
ModuleBase::timer::tick("Matrix_Orbs21", "init_radial");
8997
}
9098

@@ -95,60 +103,77 @@ void Matrix_Orbs21::init_radial(const std::vector<std::vector<std::vector<Numeri
95103
ModuleBase::TITLE("Matrix_Orbs21", "init_radial");
96104
ModuleBase::timer::tick("Matrix_Orbs21", "init_radial");
97105
assert(orb_A1.size() == orb_A2.get_ntype());
98-
for (size_t TA = 0; TA != orb_A1.size(); ++TA) {
99-
for (size_t TB = 0; TB != orb_B.get_ntype(); ++TB) {
100-
for (int LA1 = 0; LA1 != orb_A1[TA].size(); ++LA1) {
101-
for (size_t NA1 = 0; NA1 != orb_A1[TA][LA1].size(); ++NA1) {
102-
for (int LA2 = 0; LA2 <= orb_A2.Phi[TA].getLmax(); ++LA2) {
103-
for (size_t NA2 = 0; NA2 != orb_A2.Phi[TA].getNchi(LA2); ++NA2) {
104-
for (int LB = 0; LB <= orb_B.Phi[TB].getLmax(); ++LB) {
105-
for (size_t NB = 0; NB != orb_B.Phi[TB].getNchi(LB); ++NB) {
106+
for (size_t TA = 0; TA != orb_A1.size(); ++TA)
107+
{
108+
for (size_t TB = 0; TB != orb_B.get_ntype(); ++TB)
109+
{
110+
for (int LA1 = 0; LA1 != orb_A1[TA].size(); ++LA1)
111+
{
112+
for (size_t NA1 = 0; NA1 != orb_A1[TA][LA1].size(); ++NA1)
113+
{
114+
for (int LA2 = 0; LA2 <= orb_A2.Phi[TA].getLmax(); ++LA2)
115+
{
116+
for (size_t NA2 = 0; NA2 != orb_A2.Phi[TA].getNchi(LA2); ++NA2)
117+
{
118+
for (int LB = 0; LB <= orb_B.Phi[TB].getLmax(); ++LB)
119+
{
120+
for (size_t NB = 0; NB != orb_B.Phi[TB].getNchi(LB); ++NB)
121+
{
106122
center2_orb21_s[TA][TB][LA1][NA1][LA2][NA2][LB].insert(
107123
std::make_pair(NB,
108124
Center2_Orb::Orb21(orb_A1[TA][LA1][NA1],
109125
orb_A2.Phi[TA].PhiLN(LA2, NA2),
110126
orb_B.Phi[TB].PhiLN(LB, NB),
111127
psb_,
112-
this->MGT)));
113-
}
114-
}
115-
}
116-
}
117-
}
118-
}
119-
}
120-
}
128+
this->MGT)));
129+
}
130+
}
131+
}
132+
}
133+
}
134+
}
135+
}
136+
}
121137
ModuleBase::timer::tick("Matrix_Orbs21", "init_radial");
122138
}
123139

124140
void Matrix_Orbs21::init_radial_table()
125141
{
126142
ModuleBase::TITLE("Matrix_Orbs21", "init_radial");
127143
ModuleBase::timer::tick("Matrix_Orbs21", "init_radial_table");
128-
for (auto& coA: center2_orb21_s) {
129-
for (auto& coB: coA.second) {
130-
for (auto& coC: coB.second) {
131-
for (auto& coD: coC.second) {
132-
for (auto& coE: coD.second) {
133-
for (auto& coF: coE.second) {
134-
for (auto& coG: coF.second) {
135-
for (auto& coH: coG.second) {
136-
coH.second.init_radial_table();
137-
}
138-
}
139-
}
140-
}
141-
}
142-
}
143-
}
144-
}
144+
for (auto& coA: center2_orb21_s)
145+
{
146+
for (auto& coB: coA.second)
147+
{
148+
for (auto& coC: coB.second)
149+
{
150+
for (auto& coD: coC.second)
151+
{
152+
for (auto& coE: coD.second)
153+
{
154+
for (auto& coF: coE.second)
155+
{
156+
for (auto& coG: coF.second)
157+
{
158+
for (auto& coH: coG.second)
159+
{
160+
coH.second.init_radial_table();
161+
}
162+
}
163+
}
164+
}
165+
}
166+
}
167+
}
168+
}
145169
ModuleBase::timer::tick("Matrix_Orbs21", "init_radial_table");
146170
}
147171

148172
void Matrix_Orbs21::init_radial_table(const std::map<size_t, std::map<size_t, std::set<double>>>& Rs)
149173
{
150174
ModuleBase::TITLE("Matrix_Orbs21", "init_radial_table_Rs");
151175
ModuleBase::timer::tick("Matrix_Orbs21", "init_radial_table");
176+
const double lat0 = *this->lat0;
152177
for (const auto& RsA: Rs) {
153178
for (const auto& RsB: RsA.second)
154179
{
@@ -164,25 +189,32 @@ void Matrix_Orbs21::init_radial_table(const std::map<size_t, std::map<size_t, st
164189
{
165190
const double position = R * lat0 / lcao_dr_;
166191
const size_t iq = static_cast<size_t>(position);
167-
for (size_t i = 0; i != 4; ++i) {
168-
radials.insert(iq + i);
169-
}
192+
for (size_t i = 0; i != 4; ++i)
193+
{
194+
radials.insert(iq + i);
195+
}
196+
}
197+
for (auto& coC: *center2_orb21_sAB)
198+
{
199+
for (auto& coD: coC.second)
200+
{
201+
for (auto& coE: coD.second)
202+
{
203+
for (auto& coF: coE.second)
204+
{
205+
for (auto& coG: coF.second)
206+
{
207+
for (auto& coH: coG.second)
208+
{
209+
coH.second.init_radial_table(radials);
210+
}
211+
}
212+
}
213+
}
214+
}
170215
}
171-
for (auto& coC: *center2_orb21_sAB) {
172-
for (auto& coD: coC.second) {
173-
for (auto& coE: coD.second) {
174-
for (auto& coF: coE.second) {
175-
for (auto& coG: coF.second) {
176-
for (auto& coH: coG.second) {
177-
coH.second.init_radial_table(radials);
178-
}
179-
}
180-
}
181-
}
182-
}
183-
}
184216
}
185-
}
186-
}
217+
}
218+
}
187219
ModuleBase::timer::tick("Matrix_Orbs21", "init_radial_table");
188220
}

source/module_ri/Matrix_Orbs21.h

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -79,7 +79,7 @@ class Matrix_Orbs21
7979
ModuleBase::Sph_Bessel_Recursive::D2* psb_ = nullptr;
8080
ORB_gaunt_table MGT;
8181
const double lcao_dr_ = 0.01;
82-
double lat0 = 0; // restore ucell.lat0
82+
double* lat0 = nullptr; // restore ucell.lat0
8383
std::map<size_t, // TA
8484
std::map<size_t, // TB
8585
std::map<int, // LA1

source/module_ri/Matrix_Orbs21.hpp

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -22,6 +22,7 @@ RI::Tensor<Tdata> Matrix_Orbs21::cal_overlap_matrix(
2222
const Matrix_Order &matrix_order) const
2323
{
2424
RI::Tensor<Tdata> m;
25+
const double lat0 = *this->lat0;
2526
const size_t sizeA1 = index_A1[TA].count_size;
2627
const size_t sizeA2 = index_A2[TA].count_size;
2728
const size_t sizeB = index_B[TB].count_size;
@@ -98,6 +99,7 @@ std::array<RI::Tensor<Tdata>,3> Matrix_Orbs21::cal_grad_overlap_matrix(
9899
const Matrix_Order &matrix_order) const
99100
{
100101
std::array<RI::Tensor<Tdata>,3> m;
102+
const double lat0 = *this->lat0;
101103
const size_t sizeA1 = index_A1[TA].count_size;
102104
const size_t sizeA2 = index_A2[TA].count_size;
103105
const size_t sizeB = index_B[TB].count_size;

source/module_ri/Matrix_Orbs22.cpp

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -19,6 +19,7 @@ void Matrix_Orbs22::init(const int mode,
1919
ModuleBase::timer::tick("Matrix_Orbs22", "init");
2020
int Lmax_used, Lmax;
2121

22+
this->lat0 = &ucell.lat0;
2223
const int ntype = orb.get_ntype();
2324
int lmax_orb = -1, lmax_beta = -1;
2425
for (int it = 0; it < ntype; it++)
@@ -137,6 +138,7 @@ void Matrix_Orbs22::init_radial_table(const std::map<size_t, std::map<size_t, st
137138
{
138139
ModuleBase::TITLE("Matrix_Orbs22", "init_radial_table_Rs");
139140
ModuleBase::timer::tick("Matrix_Orbs22", "init_radial_table");
141+
const double lat0 = *this->lat0;
140142
for (const auto& RsA: Rs)
141143
for (const auto& RsB: RsA.second)
142144
{

source/module_ri/Matrix_Orbs22.h

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -103,7 +103,7 @@ class Matrix_Orbs22
103103
ModuleBase::Sph_Bessel_Recursive::D2* psb_ = nullptr;
104104
ORB_gaunt_table MGT;
105105
const double lcao_dr_ = 0.01;
106-
double lat0 = 0; // restore ucell.lat0
106+
double* lat0 = nullptr; // restore ucell.lat0
107107
std::map<
108108
size_t, // TA
109109
std::map<size_t, // TB

source/module_ri/Matrix_Orbs22.hpp

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -22,6 +22,7 @@ RI::Tensor<Tdata> Matrix_Orbs22::cal_overlap_matrix(
2222
const ModuleBase::Element_Basis_Index::IndexLNM &index_B2,
2323
const Matrix_Order &matrix_order) const
2424
{
25+
const double lat0 = *this->lat0;
2526
RI::Tensor<Tdata> m;
2627
const size_t sizeA1 = index_A1[TA].count_size;
2728
const size_t sizeA2 = index_A2[TA].count_size;
@@ -184,7 +185,7 @@ std::array<RI::Tensor<Tdata>,3> Matrix_Orbs22::cal_grad_overlap_matrix(
184185
default: throw std::invalid_argument(std::string(__FILE__)+" line "+std::to_string(__LINE__));
185186
}
186187
}
187-
188+
const double lat0 = *this->lat0;
188189
for( const auto &co3 : center2_orb22_s.at(TA).at(TB) )
189190
{
190191
const int LA1 = co3.first;

0 commit comments

Comments
 (0)