@@ -19,8 +19,7 @@ Get_pchg_lcao::Get_pchg_lcao(const psi::Psi<std::complex<double>>& psi, const Pa
1919}
2020
2121// For gamma_only
22- void Get_pchg_lcao::begin_gamma (double * const * rho,
23- const UnitCell& ucell,
22+ void Get_pchg_lcao::begin_gamma (const UnitCell& ucell,
2423 const Parallel_Grid& pgrid,
2524 const Grid_Driver& grid_driver,
2625 const std::vector<int >& out_pchg,
@@ -34,9 +33,14 @@ void Get_pchg_lcao::begin_gamma(double* const* rho,
3433 prepare_get_pchg (ofs_running);
3534
3635 assert (psi_gamma_ != nullptr );
37- assert (rho != nullptr );
3836 const int nrxx = pgrid.get_nrxx ();
39- std::vector<double *> rho_pointers (rho, rho + nspin_);
37+ const int precision = 11 ;
38+ std::vector<std::vector<double >> rho (nspin_, std::vector<double >(nrxx));
39+ std::vector<double *> rho_pointers (nspin_);
40+ for (int is = 0 ; is < nspin_; ++is)
41+ {
42+ rho_pointers[is] = rho[is].data ();
43+ }
4044 const std::vector<int > bands_picked = select_bands (out_pchg);
4145
4246 for (int ib = 0 ; ib < nbands_; ++ib)
@@ -55,23 +59,13 @@ void Get_pchg_lcao::begin_gamma(double* const* rho,
5559
5660 for (int is = 0 ; is < nspin_; ++is)
5761 {
58- ModuleBase::GlobalFunc::ZEROS (rho[is], nrxx );
62+ std::fill (rho[is]. begin (), rho[is]. end (), 0.0 );
5963 }
6064
6165 DM .init_DMR (&grid_driver, &ucell);
6266 DM .cal_DMR ();
6367 ModuleGint::cal_gint_rho (DM .get_DMR_vector (), nspin_, rho_pointers.data ());
6468
65- // A solution to replace the original implementation of the following code:
66- // pelec->charge->save_rho_before_sum_band();
67- // Using std::vector to replace the original double** rho_save
68- std::vector<std::vector<double >> rho_save (nspin_, std::vector<double >(nrxx));
69-
70- for (int is = 0 ; is < nspin_; ++is)
71- {
72- ModuleBase::GlobalFunc::DCOPY (rho[is], rho_save[is].data (), nrxx); // Copy data
73- }
74-
7569 for (int is = 0 ; is < nspin_; ++is)
7670 {
7771 // ssc should be inside the inner loop to reset the string stream each time
@@ -80,19 +74,7 @@ void Get_pchg_lcao::begin_gamma(double* const* rho,
8074
8175 ofs_running << " Writing cube file " << ssc.str () << std::endl;
8276
83- const int precision = 6 ;
84- ModuleIO::write_vdata_palgrid (pgrid,
85- rho_save[is].data (),
86- is,
87- nspin_,
88- 0 ,
89- ssc.str (),
90- 0.0 ,
91- &ucell,
92- precision,
93- 0 ,
94- false ,
95- false );
77+ ModuleIO::write_vdata_palgrid (pgrid, rho[is].data (), is, nspin_, 0 , ssc.str (), 0.0 , &ucell, precision, 0 , false , false );
9678 }
9779 }
9880 }
@@ -101,9 +83,7 @@ void Get_pchg_lcao::begin_gamma(double* const* rho,
10183}
10284
10385// For multi-k
104- void Get_pchg_lcao::begin_k (double * const * rho,
105- std::complex <double >* const * rhog,
106- const ModulePW::PW_Basis& rho_pw,
86+ void Get_pchg_lcao::begin_k (const ModulePW::PW_Basis& rho_pw,
10787 UnitCell& ucell,
10888 const Parallel_Grid& pgrid,
10989 const Grid_Driver& grid_driver,
@@ -120,12 +100,28 @@ void Get_pchg_lcao::begin_k(double* const* rho,
120100 prepare_get_pchg (ofs_running);
121101
122102 assert (psi_k_ != nullptr );
123- assert (rho != nullptr );
124- assert (rhog != nullptr );
125103 assert (pgrid.get_nrxx () == rho_pw.nrxx );
126104 const int nrxx = pgrid.get_nrxx ();
127- std::vector<double *> rho_pointers (rho, rho + nspin_);
128- std::vector<std::complex <double >*> rhog_pointers (rhog, rhog + nspin_);
105+ const int precision = 11 ;
106+ std::vector<std::vector<double >> rho (nspin_, std::vector<double >(nrxx));
107+ std::vector<double *> rho_pointers (nspin_);
108+ for (int is = 0 ; is < nspin_; ++is)
109+ {
110+ rho_pointers[is] = rho[is].data ();
111+ }
112+
113+ const bool needs_symmetry = !if_separate_k && ModuleSymmetry::Symmetry::symm_flag == 1 ;
114+ std::vector<std::vector<std::complex <double >>> rhog;
115+ std::vector<std::complex <double >*> rhog_pointers;
116+ if (needs_symmetry)
117+ {
118+ rhog.resize (nspin_, std::vector<std::complex <double >>(rho_pw.npw ));
119+ rhog_pointers.resize (nspin_);
120+ for (int is = 0 ; is < nspin_; ++is)
121+ {
122+ rhog_pointers[is] = rhog[is].data ();
123+ }
124+ }
129125 const std::vector<int > bands_picked = select_bands (out_pchg);
130126
131127 for (int ib = 0 ; ib < nbands_; ++ib)
@@ -152,21 +148,13 @@ void Get_pchg_lcao::begin_k(double* const* rho,
152148 {
153149 for (int is = 0 ; is < nspin_; ++is)
154150 {
155- ModuleBase::GlobalFunc::ZEROS (rho[is], nrxx );
151+ std::fill (rho[is]. begin (), rho[is]. end (), 0.0 );
156152 }
157153
158154 DM .init_DMR (&grid_driver, &ucell);
159155 DM .cal_DMR (ik);
160156 ModuleGint::cal_gint_rho (DM .get_DMR_vector (), nspin_, rho_pointers.data ());
161157
162- // Using std::vector to replace the original double** rho_save
163- std::vector<std::vector<double >> rho_save (nspin_, std::vector<double >(nrxx));
164-
165- for (int is = 0 ; is < nspin_; ++is)
166- {
167- ModuleBase::GlobalFunc::DCOPY (rho[is], rho_save[is].data (), nrxx); // Copy data
168- }
169-
170158 for (int is = 0 ; is < nspin_; ++is)
171159 {
172160 // ssc should be inside the inner loop to reset the string stream each time
@@ -175,9 +163,8 @@ void Get_pchg_lcao::begin_k(double* const* rho,
175163
176164 ofs_running << " Writing cube file " << ssc.str () << std::endl;
177165
178- const int precision = 6 ;
179166 ModuleIO::write_vdata_palgrid (pgrid,
180- rho_save [is].data (),
167+ rho [is].data (),
181168 is,
182169 nspin_,
183170 0 ,
@@ -195,30 +182,29 @@ void Get_pchg_lcao::begin_k(double* const* rho,
195182 {
196183 for (int is = 0 ; is < nspin_; ++is)
197184 {
198- ModuleBase::GlobalFunc::ZEROS (rho[is], nrxx );
185+ std::fill (rho[is]. begin (), rho[is]. end (), 0.0 );
199186 }
200187
201188 DM .init_DMR (&grid_driver, &ucell);
202189 DM .cal_DMR ();
203190 ModuleGint::cal_gint_rho (DM .get_DMR_vector (), nspin_, rho_pointers.data ());
204- // Using std::vector to replace the original double** rho_save
205- std::vector<std::vector<double >> rho_save (nspin_, std::vector<double >(nrxx));
206-
207- for (int is = 0 ; is < nspin_; ++is)
208- {
209- ModuleBase::GlobalFunc::DCOPY (rho[is], rho_save[is].data (), nrxx); // Copy data
210- }
211191
212192 // Symmetrize the charge density, otherwise the results are incorrect if the symmetry is on
213- Symmetry_rho srho;
214- for (int is = 0 ; is < nspin_; ++is)
193+ if (needs_symmetry)
215194 {
216- std::vector<double *> rho_save_pointers (nspin_);
217- for (int i = 0 ; i < nspin_; ++i)
195+ Symmetry_rho srho;
196+ if (nspin_ == 4 )
197+ {
198+ srho.begin (0 , rho_pointers.data (), rhog_pointers.data (), rho_pw.npw , nullptr , &rho_pw, ucell.symm );
199+ srho.begin_soc (rho_pointers.data (), rhog_pointers.data (), &rho_pw, ucell.symm );
200+ }
201+ else
218202 {
219- rho_save_pointers[i] = rho_save[i].data ();
203+ for (int is = 0 ; is < nspin_; ++is)
204+ {
205+ srho.begin (is, rho_pointers.data (), rhog_pointers.data (), rho_pw.npw , nullptr , &rho_pw, ucell.symm );
206+ }
220207 }
221- srho.begin (is, rho_save_pointers.data (), rhog_pointers.data (), rho_pw.npw , nullptr , &rho_pw, ucell.symm );
222208 }
223209
224210 for (int is = 0 ; is < nspin_; ++is)
@@ -229,19 +215,7 @@ void Get_pchg_lcao::begin_k(double* const* rho,
229215
230216 ofs_running << " Writing cube file " << ssc.str () << std::endl;
231217
232- const int precision = 6 ;
233- ModuleIO::write_vdata_palgrid (pgrid,
234- rho_save[is].data (),
235- is,
236- nspin_,
237- 0 ,
238- ssc.str (),
239- 0.0 ,
240- &ucell,
241- precision,
242- 0 ,
243- false ,
244- false );
218+ ModuleIO::write_vdata_palgrid (pgrid, rho[is].data (), is, nspin_, 0 , ssc.str (), 0.0 , &ucell, precision, 0 , false , false );
245219 }
246220 }
247221 }
0 commit comments