Skip to content

Commit d0708af

Browse files
PeizeLinPeizeLin
committed
1. update LRI_Cal_Aux::add_Ds(fac) and LRI::cal_loop3(fac_add_Ds)
1 parent c5ff3ab commit d0708af

File tree

6 files changed

+82
-57
lines changed

6 files changed

+82
-57
lines changed

‎include/RI/global/Tensor.h‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -55,6 +55,7 @@ class Tensor
5555
bool empty() const { return shape.empty(); }
5656

5757
Tensor & operator += (const Tensor &);
58+
Tensor & operator -= (const Tensor &);
5859
Tensor operator-() const;
5960

6061
template <class Archive> void serialize( Archive & ar ){ ar(shape, data); } // for cereal

‎include/RI/global/Tensor.hpp‎

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -159,6 +159,14 @@ Tensor<T> &Tensor<T>::operator+= (const Tensor &t)
159159
return *this;
160160
}
161161

162+
template<typename T>
163+
Tensor<T> &Tensor<T>::operator-= (const Tensor &t)
164+
{
165+
assert(same_shape(*this,t));
166+
*this->data -= *t.data;
167+
return *this;
168+
}
169+
162170
template<typename T>
163171
Tensor<T> Tensor<T>::transpose() const
164172
{

‎include/RI/physics/Exx.hpp‎

Lines changed: 24 additions & 36 deletions
Original file line numberDiff line numberDiff line change
@@ -136,8 +136,6 @@ void Exx<TA,Tcell,Ndim,Tdata>::cal_Hs(
136136
this->lri.save_load.save("Ds_"+save_names_suffix[2], {Label::ab::a1b1, Label::ab::a1b2, Label::ab::a2b1, Label::ab::a2b2});
137137
}
138138

139-
140-
141139
template<typename TA, typename Tcell, std::size_t Ndim, typename Tdata>
142140
void Exx<TA,Tcell,Ndim,Tdata>::cal_force(
143141
const std::array<std::string,5> &save_names_suffix) // "Cs","Vs","Ds","dCs","dVs"
@@ -161,34 +159,29 @@ void Exx<TA,Tcell,Ndim,Tdata>::cal_force(
161159
this->lri.save_load.load("Vs_"+save_names_suffix[1], Label::ab::a0b0);
162160
this->lri.save_load.load("Cs_"+save_names_suffix[0], Label::ab::b);
163161

164-
this->lri.coefficients = {[](const Label::ab_ab &label, const TA &Aa01, const TAC &Aa2, const TAC &Ab01, const TAC &Ab2) -> Tdata
165-
{
166-
switch(label)
167-
{
168-
case Label::ab_ab::a0b0_a1b1: case Label::ab_ab::a0b0_a1b2: return -1;
169-
case Label::ab_ab::a0b0_a2b1: case Label::ab_ab::a0b0_a2b2: return 1;
170-
default: throw std::invalid_argument(std::string(__FILE__)+" line "+std::to_string(__LINE__));
171-
}
172-
}};
173-
174-
this->lri.cal(
162+
this->lri.cal_loop3(
175163
{Label::ab_ab::a0b0_a1b1,
176-
Label::ab_ab::a0b0_a1b2,
177-
Label::ab_ab::a0b0_a2b1,
164+
Label::ab_ab::a0b0_a1b2,},
165+
dHs_vec,
166+
-1.0);
167+
168+
this->lri.cal_loop3(
169+
{Label::ab_ab::a0b0_a2b1,
178170
Label::ab_ab::a0b0_a2b2},
179-
dHs_vec);
171+
dHs_vec,
172+
1.0);
180173

181174
this->lri.save_load.save("dCs_"+std::to_string(ipos)+"_"+save_names_suffix[3], Label::ab::a);
182175
this->lri.save_load.save("Vs_"+save_names_suffix[1], Label::ab::a0b0);
183176

184177
this->lri.save_load.load("Cs_"+save_names_suffix[0], Label::ab::a);
185178
this->lri.save_load.load("dVs_"+std::to_string(ipos)+"_"+save_names_suffix[4], Label::ab::a0b0);
186179

187-
this->lri.coefficients = {nullptr};
188-
this->lri.cal(
180+
this->lri.cal_loop3(
189181
{Label::ab_ab::a0b0_a2b2,
190182
Label::ab_ab::a0b0_a2b1},
191-
dHs_vec);
183+
dHs_vec,
184+
1.0);
192185

193186
this->post_2D.cal_force(
194187
this->post_2D.saves["Ds_"+save_names_suffix[2]],
@@ -205,34 +198,29 @@ void Exx<TA,Tcell,Ndim,Tdata>::cal_force(
205198
{
206199
std::vector<std::map<TA,std::map<TAC,Tensor<Tdata>>>> dHs_vec(1);
207200

208-
this->lri.coefficients = {nullptr};
209-
this->lri.cal(
201+
this->lri.cal_loop3(
210202
{Label::ab_ab::a0b0_a2b2,
211203
Label::ab_ab::a0b0_a1b2},
212-
dHs_vec);
204+
dHs_vec,
205+
1.0);
213206

214207
this->lri.save_load.save("dVs_"+std::to_string(ipos)+"_"+save_names_suffix[4], Label::ab::a0b0);
215208
this->lri.save_load.save("Cs_"+save_names_suffix[0], Label::ab::b);
216209

217210
this->lri.save_load.load("Vs_"+save_names_suffix[1], Label::ab::a0b0);
218211
this->lri.save_load.load("dCs_"+std::to_string(ipos)+"_"+save_names_suffix[3], Label::ab::b);
219212

220-
this->lri.coefficients = {[](const Label::ab_ab &label, const TA &Aa01, const TAC &Aa2, const TAC &Ab01, const TAC &Ab2) -> Tdata
221-
{
222-
switch(label)
223-
{
224-
case Label::ab_ab::a0b0_a1b1: case Label::ab_ab::a0b0_a2b1: return 1;
225-
case Label::ab_ab::a0b0_a1b2: case Label::ab_ab::a0b0_a2b2: return -1;
226-
default: throw std::invalid_argument(std::string(__FILE__)+" line "+std::to_string(__LINE__));
227-
}
228-
}};
229-
230-
this->lri.cal(
213+
this->lri.cal_loop3(
231214
{Label::ab_ab::a0b0_a1b1,
232-
Label::ab_ab::a0b0_a1b2,
233-
Label::ab_ab::a0b0_a2b1,
215+
Label::ab_ab::a0b0_a2b1},
216+
dHs_vec,
217+
1.0);
218+
219+
this->lri.cal_loop3(
220+
{Label::ab_ab::a0b0_a1b2,
234221
Label::ab_ab::a0b0_a2b2},
235-
dHs_vec);
222+
dHs_vec,
223+
-1.0);
236224

237225
this->lri.save_load.save("Cs_"+save_names_suffix[0], Label::ab::a);
238226
this->lri.save_load.save("dCs_"+std::to_string(ipos)+"_"+save_names_suffix[3], Label::ab::b);

‎include/RI/ri/LRI-cal_loop3.hpp‎

Lines changed: 7 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -21,7 +21,8 @@ namespace RI
2121
template<typename TA, typename Tcell, std::size_t Ndim, typename Tdata>
2222
void LRI<TA,Tcell,Ndim,Tdata>::cal_loop3(
2323
const std::vector<Label::ab_ab> &labels,
24-
std::vector<std::map<TA, std::map<TAC, Tensor<Tdata>>>> &Ds_result)
24+
std::vector<std::map<TA, std::map<TAC, Tensor<Tdata>>>> &Ds_result,
25+
const double fac_add_Ds)
2526
{
2627
using namespace Array_Operator;
2728

@@ -118,7 +119,7 @@ void LRI<TA,Tcell,Ndim,Tdata>::cal_loop3(
118119
if(!Ds_result_fixed.empty())
119120
LRI_Cal_Aux::add_Ds( LRI_Cal_Aux::Ds_translate(std::move(Ds_result_fixed), Aa2.second, this->period),
120121
Ds_result_thread[0][Aa2.first]);
121-
LRI_Cal_Aux::add_Ds_omp_try(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add);
122+
LRI_Cal_Aux::add_Ds_omp_try(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add, fac_add_Ds);
122123
} // end for Ab2
123124
break;
124125
} // end case a0b0_a1b1
@@ -167,7 +168,7 @@ void LRI<TA,Tcell,Ndim,Tdata>::cal_loop3(
167168
if(!Ds_result_fixed.empty())
168169
LRI_Cal_Aux::add_Ds( LRI_Cal_Aux::Ds_exchange(std::move(Ds_result_fixed), Ab01, this->period),
169170
Ds_result_thread[0]);
170-
LRI_Cal_Aux::add_Ds_omp_try(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add);
171+
LRI_Cal_Aux::add_Ds_omp_try(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add, fac_add_Ds);
171172
} // end for Ab01
172173
break;
173174
} // end case a0b0_a1b2
@@ -217,7 +218,7 @@ void LRI<TA,Tcell,Ndim,Tdata>::cal_loop3(
217218
if(!Ds_result_fixed.empty())
218219
LRI_Cal_Aux::add_Ds( std::move(Ds_result_fixed),
219220
Ds_result_thread[0][Aa01]);
220-
LRI_Cal_Aux::add_Ds_omp_try(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add);
221+
LRI_Cal_Aux::add_Ds_omp_try(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add, fac_add_Ds);
221222
} // end for Aa01
222223
break;
223224
} // end case a0b0_a2b1
@@ -267,7 +268,7 @@ void LRI<TA,Tcell,Ndim,Tdata>::cal_loop3(
267268
if(!Ds_result_fixed.empty())
268269
LRI_Cal_Aux::add_Ds( std::move(Ds_result_fixed),
269270
Ds_result_thread[0][Aa01]);
270-
LRI_Cal_Aux::add_Ds_omp_try(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add);
271+
LRI_Cal_Aux::add_Ds_omp_try(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add, fac_add_Ds);
271272
} // end for Aa01
272273
break;
273274
} // end case a0b0_a2b2
@@ -277,7 +278,7 @@ void LRI<TA,Tcell,Ndim,Tdata>::cal_loop3(
277278
} // end switch(label)
278279
} // end for label
279280

280-
LRI_Cal_Aux::add_Ds_omp_wait(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add);
281+
LRI_Cal_Aux::add_Ds_omp_wait(std::move(Ds_result_thread), Ds_result, lock_Ds_result_add, fac_add_Ds);
281282
} // end #pragma omp parallel
282283

283284
omp_destroy_lock(&lock_Ds_result_add);

‎include/RI/ri/LRI.h‎

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -58,7 +58,8 @@ class LRI
5858
std::vector<std::map<TA, std::map<TAC, Tensor<Tdata>>>> &Ds_result);
5959
void cal_loop3(
6060
const std::vector<Label::ab_ab> &labels,
61-
std::vector<std::map<TA, std::map<TAC, Tensor<Tdata>>>> &Ds_result);
61+
std::vector<std::map<TA, std::map<TAC, Tensor<Tdata>>>> &Ds_result,
62+
const double fac_add_Ds = 1.0);
6263

6364
public:
6465
std::shared_ptr<Parallel_LRI<TA,Tcell,Ndim,Tdata>>

‎include/RI/ri/LRI_Cal_Aux.h‎

Lines changed: 40 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,7 @@
1313
#include <memory.h>
1414
#include <cassert>
1515
#include <stdexcept>
16+
#include <omp.h>
1617

1718
namespace RI
1819
{
@@ -66,43 +67,66 @@ namespace LRI_Cal_Aux
6667
*/
6768

6869
template<typename Tdata>
69-
inline void add_Ds(Tensor<Tdata> &&D_add, Tensor<Tdata> &D_result)
70+
inline void add_Ds(
71+
Tensor<Tdata> &&D_add,
72+
Tensor<Tdata> &D_result,
73+
const double fac = 1.0)
7074
{
7175
if(D_result.empty())
72-
D_result = std::move(D_add);
76+
{
77+
if(1.0==fac)
78+
D_result = std::move(D_add);
79+
else if(-1.0==fac)
80+
D_result = -D_add;
81+
else
82+
D_result = Tdata(fac) * D_add;
83+
}
7384
else
74-
D_result += D_add;
85+
{
86+
if(1.0==fac)
87+
D_result += D_add;
88+
else if(-1.0==fac)
89+
D_result -= D_add;
90+
else
91+
D_result += Tdata(fac) * D_add;
92+
}
7593
}
7694

7795
template<typename Tkey, typename Tvalue>
7896
void add_Ds(
7997
std::map<Tkey, Tvalue> &&Ds_add,
80-
std::map<Tkey, Tvalue> &Ds_result)
98+
std::map<Tkey, Tvalue> &Ds_result,
99+
const double fac = 1.0)
81100
{
82-
if(Ds_result.empty())
101+
if(Ds_result.empty() && 1.0==fac)
83102
Ds_result = std::move(Ds_add);
84103
else
85104
{
86105
for(auto &&Ds_add_A : Ds_add)
87-
add_Ds(std::move(Ds_add_A.second), Ds_result[Ds_add_A.first]);
106+
add_Ds(std::move(Ds_add_A.second), Ds_result[Ds_add_A.first], fac);
88107
Ds_add.clear();
89108
}
90109
}
110+
91111
template<typename Tvalue>
92112
void add_Ds(
93113
std::vector<Tvalue> &&Ds_add,
94-
std::vector<Tvalue> &Ds_result)
114+
std::vector<Tvalue> &Ds_result,
115+
const double fac = 1.0)
95116
{
96-
if(Ds_result.empty())
117+
if(Ds_result.empty() && 1.0==fac)
97118
{
98119
Ds_result = std::move(Ds_add);
99120
Ds_add.resize(Ds_result.size());
100121
}
101122
else
102123
{
103-
assert(Ds_add.size()==Ds_result.size());
124+
if(Ds_result.empty())
125+
Ds_result.resize(Ds_add.size());
126+
else
127+
assert(Ds_add.size()==Ds_result.size());
104128
for(std::size_t i=0; i<Ds_result.size(); ++i)
105-
add_Ds(std::move(Ds_add[i]), Ds_result[i]);
129+
add_Ds(std::move(Ds_add[i]), Ds_result[i], fac);
106130
Ds_add.clear();
107131
Ds_add.resize(Ds_result.size());
108132
}
@@ -199,11 +223,12 @@ namespace LRI_Cal_Aux
199223
void add_Ds_omp_try(
200224
std::vector<std::map<TA, std::map<TAC, Tensor<Tdata>>>> &&Ds_result_thread,
201225
std::vector<std::map<TA, std::map<TAC, Tensor<Tdata>>>> &Ds_result,
202-
omp_lock_t &lock_Ds_result_add)
226+
omp_lock_t &lock_Ds_result_add,
227+
const double &fac)
203228
{
204229
if( !LRI_Cal_Aux::judge_Ds_empty(Ds_result_thread) && omp_test_lock(&lock_Ds_result_add) )
205230
{
206-
LRI_Cal_Aux::add_Ds(std::move(Ds_result_thread), Ds_result);
231+
LRI_Cal_Aux::add_Ds(std::move(Ds_result_thread), Ds_result, fac);
207232
omp_unset_lock(&lock_Ds_result_add);
208233
Ds_result_thread.clear();
209234
Ds_result_thread.resize(Ds_result.size());
@@ -214,12 +239,13 @@ namespace LRI_Cal_Aux
214239
void add_Ds_omp_wait(
215240
std::vector<std::map<TA, std::map<TAC, Tensor<Tdata>>>> &&Ds_result_thread,
216241
std::vector<std::map<TA, std::map<TAC, Tensor<Tdata>>>> &Ds_result,
217-
omp_lock_t &lock_Ds_result_add)
242+
omp_lock_t &lock_Ds_result_add,
243+
const double &fac)
218244
{
219245
if(!LRI_Cal_Aux::judge_Ds_empty(Ds_result_thread))
220246
{
221247
omp_set_lock(&lock_Ds_result_add);
222-
LRI_Cal_Aux::add_Ds(std::move(Ds_result_thread), Ds_result);
248+
LRI_Cal_Aux::add_Ds(std::move(Ds_result_thread), Ds_result, fac);
223249
omp_unset_lock(&lock_Ds_result_add);
224250
Ds_result_thread.clear();
225251
Ds_result_thread.resize(Ds_result.size());

0 commit comments

Comments
 (0)