Skip to content

Commit 29fc792

Browse files
linpeizePeizeLinmohanchen
authored
Refactor: simplify Mixing DM (#7462)
* Refactor class Mix_DMk_2D and Mix_Matrix * fix plain mixing for seperate_loop * fix mixing in RPA_LRI * update Mix_DMk_2D --------- Co-authored-by: linpz <linpz@mail.ustc.edu.cn> Co-authored-by: Mohan Chen <mohanchen@pku.edu.cn>
1 parent b2e3a37 commit 29fc792

11 files changed

Lines changed: 188 additions & 500 deletions

File tree

source/Makefile.Objects

Lines changed: 0 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -679,7 +679,6 @@ OBJS_MODULE_RI=conv_coulomb_pot_k.o\
679679
Matrix_Orbs22.o\
680680
RI_2D_Comm.o\
681681
Mix_DMk_2D.o\
682-
Mix_Matrix.o\
683682
symmetry_rotation.o\
684683
symmetry_rotation_output.o\
685684
symmetry_irreducible_sector.o\

source/source_lcao/module_ri/CMakeLists.txt

Lines changed: 0 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -8,7 +8,6 @@ if (ENABLE_LIBRI)
88
Matrix_Orbs22.cpp
99
RI_2D_Comm.cpp
1010
Mix_DMk_2D.cpp
11-
Mix_Matrix.cpp
1211
)
1312

1413
if(ENABLE_LCAO)

source/source_lcao/module_ri/Exx_LRI_interface.h

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -124,7 +124,8 @@ class Exx_LRI_Interface
124124
std::shared_ptr<Exx_LRI<Tdata>> exx_ptr;
125125

126126
private:
127-
Mix_DMk_2D mix_DMk_2D;
127+
128+
Mix_DMk_2D<T> mix_DMk_2D;
128129

129130
bool exx_spacegroup_symmetry = false;
130131
ModuleSymmetry::Symmetry_rotation symrot_;

source/source_lcao/module_ri/Exx_LRI_interface.hpp

Lines changed: 41 additions & 65 deletions
Original file line numberDiff line numberDiff line change
@@ -35,9 +35,7 @@ void Exx_LRI_Interface<T, Tdata>::cal_exx_ions(const UnitCell& ucell, const bool
3535
ModuleBase::TITLE("Exx_LRI_Interface","cal_exx_ions");
3636
if(!this->flag_finish.init)
3737
{ throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); }
38-
3938
this->exx_ptr->cal_exx_ions(ucell, write_cv);
40-
4139
this->flag_finish.ions = true;
4240
}
4341

@@ -49,13 +47,9 @@ void Exx_LRI_Interface<T, Tdata>::cal_exx_elec(const std::vector<std::map<TA, st
4947
{
5048
ModuleBase::TITLE("Exx_LRI_Interface","cal_exx_elec");
5149
if(!this->flag_finish.init || !this->flag_finish.ions)
52-
{
53-
throw std::runtime_error("Exx init unfinished when "
54-
+std::string(__FILE__)+" line "+std::to_string(__LINE__));
55-
}
50+
{ throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); }
5651

5752
this->exx_ptr->cal_exx_elec(Ds, ucell, pv, p_symrot);
58-
5953
this->flag_finish.elec = true;
6054
}
6155

@@ -64,17 +58,11 @@ void Exx_LRI_Interface<T, Tdata>::cal_exx_force(const int& nat)
6458
{
6559
ModuleBase::TITLE("Exx_LRI_Interface","cal_exx_force");
6660
if(!this->flag_finish.init || !this->flag_finish.ions)
67-
{
68-
throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__));
69-
}
61+
{ throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); }
7062
if(!this->flag_finish.elec)
71-
{
72-
throw std::runtime_error("Exx Hamiltonian unfinished when "+std::string(__FILE__)
73-
+" line "+std::to_string(__LINE__));
74-
}
63+
{ throw std::runtime_error("Exx Hamiltonian unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); }
7564

7665
this->exx_ptr->cal_exx_force(nat);
77-
7866
this->flag_finish.force = true;
7967
}
8068

@@ -83,24 +71,17 @@ void Exx_LRI_Interface<T, Tdata>::cal_exx_stress(const double& omega, const doub
8371
{
8472
ModuleBase::TITLE("Exx_LRI_Interface","cal_exx_stress");
8573
if(!this->flag_finish.init || !this->flag_finish.ions)
86-
{
87-
throw std::runtime_error("Exx init unfinished when "
88-
+std::string(__FILE__)+" line "+std::to_string(__LINE__));
89-
}
74+
{ throw std::runtime_error("Exx init unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); }
9075
if(!this->flag_finish.elec)
91-
{
92-
throw std::runtime_error("Exx Hamiltonian unfinished when "
93-
+std::string(__FILE__)+" line "+std::to_string(__LINE__));
94-
}
76+
{ throw std::runtime_error("Exx Hamiltonian unfinished when "+std::string(__FILE__)+" line "+std::to_string(__LINE__)); }
9577

9678
this->exx_ptr->cal_exx_stress(omega, lat0);
97-
9879
this->flag_finish.stress = true;
9980
}
10081

10182
template<typename T, typename Tdata>
10283
void Exx_LRI_Interface<T, Tdata>::exx_before_all_runners(
103-
const K_Vectors& kv,
84+
const K_Vectors& kv,
10485
const UnitCell& ucell,
10586
const Parallel_2D& pv)
10687
{
@@ -147,14 +128,15 @@ void Exx_LRI_Interface<T, Tdata>::exx_beforescf(const int istep,
147128
if(GlobalC::exx_info.info_global.cal_exx)
148129
{
149130
if (this->exx_spacegroup_symmetry)
150-
{this->mix_DMk_2D.set_nks(kv.get_nkstot_full() * (PARAM.inp.nspin == 2 ? 2 : 1), PARAM.globalv.gamma_only_local);}
131+
{ this->mix_DMk_2D.set_nks(kv.get_nkstot_full() * (PARAM.inp.nspin == 2 ? 2 : 1)); }
151132
else
152-
{this->mix_DMk_2D.set_nks(kv.get_nks(), PARAM.globalv.gamma_only_local);}
133+
{ this->mix_DMk_2D.set_nks(kv.get_nks()); }
153134

154-
if(GlobalC::exx_info.info_global.separate_loop)
155-
{ this->mix_DMk_2D.set_mixing(nullptr); }
135+
if (GlobalC::exx_info.info_global.separate_loop)
136+
{ this->mix_DMk_2D.set_mixing_plain(GlobalC::exx_info.info_global.mixing_beta_for_loop1); }
156137
else
157138
{ this->mix_DMk_2D.set_mixing(chgmix.get_mixing()); }
139+
158140
// for exx two_level scf
159141
this->two_level_step = 0;
160142
}
@@ -171,13 +153,13 @@ void Exx_LRI_Interface<T, Tdata>::exx_eachiterinit(const int istep,
171153
ModuleBase::TITLE("Exx_LRI_Interface","exx_eachiterinit");
172154
if (GlobalC::exx_info.info_global.cal_exx)
173155
{
174-
if (!GlobalC::exx_info.info_global.separate_loop
175-
&& (this->two_level_step
176-
|| istep > 0
156+
if (!GlobalC::exx_info.info_global.separate_loop
157+
&& (this->two_level_step
158+
|| istep > 0
177159
|| PARAM.inp.init_wfc == "file") // non separate loop case
178-
|| (GlobalC::exx_info.info_global.separate_loop
179-
&& PARAM.inp.init_wfc == "file"
180-
&& this->two_level_step == 0
160+
|| (GlobalC::exx_info.info_global.separate_loop
161+
&& PARAM.inp.init_wfc == "file"
162+
&& this->two_level_step == 0
181163
&& iter == 1)
182164
) // the first iter in separate loop case
183165
{
@@ -186,33 +168,25 @@ void Exx_LRI_Interface<T, Tdata>::exx_eachiterinit(const int istep,
186168
auto cal = [this, &ucell,&kv, &flag_restart](const elecstate::DensityMatrix<T, double>& dm_in)
187169
{
188170
if (this->exx_spacegroup_symmetry)
189-
{ this->mix_DMk_2D.mix(symrot_.restore_dm(kv,dm_in.get_DMK_vector(), *dm_in.get_paraV_pointer()), flag_restart); }
171+
{ this->mix_DMk_2D.mix(symrot_.restore_dm(kv, dm_in.get_DMK_vector(), *dm_in.get_paraV_pointer()), flag_restart); }
190172
else
191173
{ this->mix_DMk_2D.mix(dm_in.get_DMK_vector(), flag_restart); }
192-
const std::vector<std::map<TA, std::map<TAC, RI::Tensor<Tdata>>>>
193-
Ds = PARAM.globalv.gamma_only_local
194-
? RI_2D_Comm::split_m2D_ktoR<Tdata>(
195-
ucell,
196-
*this->exx_ptr->p_kv,
197-
this->mix_DMk_2D.get_DMk_gamma_out(),
198-
*dm_in.get_paraV_pointer(),
199-
PARAM.inp.nspin)
200-
: RI_2D_Comm::split_m2D_ktoR<Tdata>(
201-
ucell,
202-
*this->exx_ptr->p_kv,
203-
this->mix_DMk_2D.get_DMk_k_out(),
204-
*dm_in.get_paraV_pointer(),
205-
PARAM.inp.nspin,
206-
this->exx_spacegroup_symmetry);
207-
208-
if (this->exx_spacegroup_symmetry && GlobalC::exx_info.info_ri.exx_symmetry_realspace)
174+
const std::vector<std::map<TA, std::map<TAC, RI::Tensor<Tdata>>>> Ds =
175+
RI_2D_Comm::split_m2D_ktoR<Tdata>(
176+
ucell,
177+
*this->exx_ptr->p_kv,
178+
this->mix_DMk_2D.get_DMk_out(),
179+
*dm_in.get_paraV_pointer(),
180+
PARAM.inp.nspin,
181+
this->exx_spacegroup_symmetry);
182+
if(this->exx_spacegroup_symmetry && GlobalC::exx_info.info_ri.exx_symmetry_realspace)
209183
{ this->cal_exx_elec(Ds, ucell,*dm_in.get_paraV_pointer(), &this->symrot_); }
210184
else
211185
{ this->cal_exx_elec(Ds, ucell,*dm_in.get_paraV_pointer()); }
212186
};
213187

214188
if(istep > 0 && flag_restart)
215-
{ cal(*dm_last_step); }
189+
{ cal(*this->dm_last_step); }
216190
else
217191
{ cal(dm); }
218192
}
@@ -387,21 +361,23 @@ bool Exx_LRI_Interface<T, Tdata>::exx_after_converge(
387361
// if init_wfc == "file", DM is calculated in the 1st iter of the 1st two-level step, so we mix it here
388362
const bool flag_restart = (this->two_level_step == 0 && PARAM.inp.init_wfc != "file") ? true : false;
389363

390-
if (this->exx_spacegroup_symmetry)
391-
{this->mix_DMk_2D.mix(symrot_.restore_dm(kv, dm.get_DMK_vector(), *dm.get_paraV_pointer()), flag_restart);}
364+
if(this->exx_spacegroup_symmetry)
365+
{ this->mix_DMk_2D.mix(symrot_.restore_dm(kv, dm.get_DMK_vector(), *dm.get_paraV_pointer()), flag_restart); }
392366
else
393-
{this->mix_DMk_2D.mix(dm.get_DMK_vector(), flag_restart);}
394-
395-
// GlobalC::exx_lcao.cal_exx_elec(p_esolver->LOC, p_esolver->LOWF.wfc_k_grid);
396-
const std::vector<std::map<int, std::map<std::pair<int, std::array<int, 3>>, RI::Tensor<Tdata>>>>
397-
Ds = std::is_same<T, double>::value //gamma_only_local
398-
? RI_2D_Comm::split_m2D_ktoR<Tdata>(ucell,*this->exx_ptr->p_kv, this->mix_DMk_2D.get_DMk_gamma_out(), *dm.get_paraV_pointer(), nspin)
399-
: RI_2D_Comm::split_m2D_ktoR<Tdata>(ucell,*this->exx_ptr->p_kv, this->mix_DMk_2D.get_DMk_k_out(), *dm.get_paraV_pointer(), nspin, this->exx_spacegroup_symmetry);
400-
401-
if (this->exx_spacegroup_symmetry && GlobalC::exx_info.info_ri.exx_symmetry_realspace)
367+
{ this->mix_DMk_2D.mix(dm.get_DMK_vector(), flag_restart); }
368+
const std::vector<std::map<TA, std::map<TAC, RI::Tensor<Tdata>>>> Ds =
369+
RI_2D_Comm::split_m2D_ktoR<Tdata>(
370+
ucell,
371+
*this->exx_ptr->p_kv,
372+
this->mix_DMk_2D.get_DMk_out(),
373+
*dm.get_paraV_pointer(),
374+
nspin,
375+
this->exx_spacegroup_symmetry);
376+
if(this->exx_spacegroup_symmetry && GlobalC::exx_info.info_ri.exx_symmetry_realspace)
402377
{ this->cal_exx_elec(Ds, ucell, *dm.get_paraV_pointer(), &this->symrot_); }
403378
else
404379
{ this->cal_exx_elec(Ds, ucell, *dm.get_paraV_pointer()); } // restore DM but not Hexx
380+
405381
iter = 0;
406382
this->two_level_step++;
407383

source/source_lcao/module_ri/Mix_DMk_2D.cpp

Lines changed: 63 additions & 47 deletions
Original file line numberDiff line numberDiff line change
@@ -4,69 +4,85 @@
44
//=======================
55

66
#include "Mix_DMk_2D.h"
7+
#include "source_base/module_mixing/plain_mixing.h"
78
#include "source_base/tool_title.h"
89

9-
Mix_DMk_2D &Mix_DMk_2D::set_nks(const int nks, const bool gamma_only_in)
10+
#include <cassert>
11+
12+
template <typename Tdata>
13+
Mix_DMk_2D<Tdata>::~Mix_DMk_2D<Tdata>()
14+
{
15+
if(this->flag_del_mixing)
16+
delete this->mixing;
17+
}
18+
19+
template <typename Tdata>
20+
void Mix_DMk_2D<Tdata>::set_nks(const int nks)
1021
{
11-
ModuleBase::TITLE("Mix_DMk_2D", "set_nks");
12-
this->gamma_only = gamma_only_in;
13-
if (this->gamma_only)
14-
this->mix_DMk_gamma.resize(nks);
15-
else
16-
this->mix_DMk_k.resize(nks);
17-
return *this;
22+
this->mix_DMk.clear();
23+
this->mix_DMk.resize(nks);
1824
}
1925

20-
Mix_DMk_2D &Mix_DMk_2D::set_mixing(Base_Mixing::Mixing* mixing_in)
26+
template <typename Tdata>
27+
void Mix_DMk_2D<Tdata>::set_mixing(Base_Mixing::Mixing* mixing_in)
2128
{
22-
ModuleBase::TITLE("Mix_DMk_2D","set_mixing");
23-
if(this->gamma_only)
24-
for (Mix_Matrix<std::vector<double>>& mix_one : this->mix_DMk_gamma)
25-
mix_one.init(mixing_in);
26-
else
27-
for (Mix_Matrix<std::vector<std::complex<double>>>& mix_one : this->mix_DMk_k)
28-
mix_one.init(mixing_in);
29-
return *this;
29+
if(this->flag_del_mixing)
30+
delete this->mixing;
31+
this->mixing = mixing_in;
32+
this->flag_del_mixing = false;
3033
}
3134

32-
Mix_DMk_2D &Mix_DMk_2D::set_mixing_beta(const double mixing_beta)
35+
template <typename Tdata>
36+
void Mix_DMk_2D<Tdata>::set_mixing_plain(const double& mixing_beta)
3337
{
34-
ModuleBase::TITLE("Mix_DMk_2D","set_mixing_beta");
35-
if(this->gamma_only)
36-
for (Mix_Matrix<std::vector<double>>& mix_one : this->mix_DMk_gamma)
37-
mix_one.mixing_beta = mixing_beta;
38-
else
39-
for (Mix_Matrix<std::vector<std::complex<double>>>& mix_one : this->mix_DMk_k)
40-
mix_one.mixing_beta = mixing_beta;
41-
return *this;
38+
if(this->flag_del_mixing)
39+
delete this->mixing;
40+
this->mixing = new Base_Mixing::Plain_Mixing(mixing_beta);
41+
this->flag_del_mixing = true;
4242
}
4343

44-
void Mix_DMk_2D::mix(const std::vector<std::vector<double>>& dm, const bool flag_restart)
44+
template <typename Tdata>
45+
void Mix_DMk_2D<Tdata>::mix(const std::vector<std::vector<Tdata>>& dm, const bool flag_restart)
4546
{
46-
ModuleBase::TITLE("Mix_DMk_2D","mix");
47-
assert(this->mix_DMk_gamma.size() == dm.size());
48-
for(int ik=0; ik<dm.size(); ++ik)
49-
this->mix_DMk_gamma[ik].mix(dm[ik], flag_restart);
47+
ModuleBase::TITLE("Mix_DMk_2D", "mix");
48+
if (flag_restart)
49+
{ this->restart_all(dm); }
50+
else
51+
{ this->mix_all(dm); }
5052
}
51-
void Mix_DMk_2D::mix(const std::vector<std::vector<std::complex<double>>>& dm, const bool flag_restart)
53+
54+
template <typename Tdata>
55+
std::vector<const std::vector<Tdata>*> Mix_DMk_2D<Tdata>::get_DMk_out() const
5256
{
53-
ModuleBase::TITLE("Mix_DMk_2D","mix");
54-
assert(this->mix_DMk_k.size() == dm.size());
55-
for(int ik=0; ik<dm.size(); ++ik)
56-
this->mix_DMk_k[ik].mix(dm[ik], flag_restart);
57+
std::vector<const std::vector<Tdata>*> DMk_out(this->mix_DMk.size());
58+
for (int ik = 0; ik < this->mix_DMk.size(); ++ik)
59+
{ DMk_out[ik] = &this->mix_DMk[ik].data_out; }
60+
return DMk_out;
5761
}
5862

59-
std::vector<const std::vector<double>*> Mix_DMk_2D::get_DMk_gamma_out() const
63+
template <typename Tdata>
64+
void Mix_DMk_2D<Tdata>::restart_all(const std::vector<std::vector<Tdata>>& data_in)
6065
{
61-
std::vector<const std::vector<double>*> DMk_out(this->mix_DMk_gamma.size());
62-
for(int ik=0; ik<this->mix_DMk_gamma.size(); ++ik)
63-
DMk_out[ik] = &this->mix_DMk_gamma[ik].get_data_out();
64-
return DMk_out;
66+
assert(this->mix_DMk.size() == data_in.size());
67+
assert(this->mixing != nullptr);
68+
for (int ik = 0; ik < data_in.size(); ++ik)
69+
{
70+
this->mix_DMk[ik].data_out = data_in[ik];
71+
this->mixing->init_mixing_data(this->mix_DMk[ik].mixing_data, data_in[ik].size(), sizeof(Tdata));
72+
}
6573
}
66-
std::vector<const std::vector<std::complex<double>>*> Mix_DMk_2D::get_DMk_k_out() const
74+
75+
template <typename Tdata>
76+
void Mix_DMk_2D<Tdata>::mix_all(const std::vector<std::vector<Tdata>>& data_in)
6777
{
68-
std::vector<const std::vector<std::complex<double>>*> DMk_out(this->mix_DMk_k.size());
69-
for(int ik=0; ik<this->mix_DMk_k.size(); ++ik)
70-
DMk_out[ik] = &this->mix_DMk_k[ik].get_data_out();
71-
return DMk_out;
72-
}
78+
assert(this->mix_DMk.size() == data_in.size());
79+
assert(this->mixing != nullptr);
80+
for (int ik = 0; ik < data_in.size(); ++ik)
81+
{
82+
this->mixing->push_data(this->mix_DMk[ik].mixing_data, this->mix_DMk[ik].data_out.data(), data_in[ik].data(), nullptr, false);
83+
this->mixing->mix_data(this->mix_DMk[ik].mixing_data, this->mix_DMk[ik].data_out.data());
84+
}
85+
}
86+
87+
template class Mix_DMk_2D<double>;
88+
template class Mix_DMk_2D<std::complex<double>>;

0 commit comments

Comments
 (0)