Skip to content

Commit 2366837

Browse files
QianruipkuQianruipku
committed
merge
2 parents 3495f9f + 53d83b2 commit 2366837

File tree

4 files changed

+104
-188
lines changed

4 files changed

+104
-188
lines changed

‎source/module_esolver/esolver_fp.cpp‎

Lines changed: 74 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,11 @@
11
#include "esolver_fp.h"
22

33
#include "module_base/global_variable.h"
4+
#include "module_elecstate/cal_ux.h"
45
#include "module_elecstate/module_charge/symmetry_rho.h"
56
#include "module_elecstate/read_pseudo.h"
7+
#include "module_hamilt_general/module_ewald/H_Ewald_pw.h"
8+
#include "module_hamilt_general/module_vdw/vdw.h"
69
#include "module_hamilt_pw/hamilt_pwdft/global.h"
710
#include "module_io/cif_io.h"
811
#include "module_io/cube_io.h"
@@ -283,6 +286,77 @@ void ESolver_FP::before_scf(UnitCell& ucell, const int istep)
283286
ModuleBase::GlobalFunc::DONE(GlobalV::ofs_running, "INIT K-POINTS");
284287
}
285288

289+
// charge extrapolation
290+
if (ucell.ionic_position_updated)
291+
{
292+
this->CE.update_all_dis(ucell);
293+
this->CE.extrapolate_charge(&(this->Pgrid),
294+
ucell,
295+
this->pelec->charge,
296+
&(this->sf),
297+
GlobalV::ofs_running,
298+
GlobalV::ofs_warning);
299+
}
300+
301+
//----------------------------------------------------------
302+
// about vdw, jiyy add vdwd3 and linpz add vdwd2
303+
//----------------------------------------------------------
304+
auto vdw_solver = vdw::make_vdw(ucell, PARAM.inp, &(GlobalV::ofs_running));
305+
if (vdw_solver != nullptr)
306+
{
307+
this->pelec->f_en.evdw = vdw_solver->get_energy();
308+
}
309+
310+
// calculate ewald energy
311+
if (!PARAM.inp.test_skip_ewald)
312+
{
313+
this->pelec->f_en.ewald_energy = H_Ewald_pw::compute_ewald(ucell, this->pw_rhod, this->sf.strucFac);
314+
}
315+
316+
//----------------------------------------------------------
317+
//! cal_ux should be called before init_scf because
318+
//! the direction of ux is used in noncoline_rho
319+
//----------------------------------------------------------
320+
elecstate::cal_ux(ucell);
321+
322+
//! output the initial charge density
323+
if (PARAM.inp.out_chg[0] == 2)
324+
{
325+
for (int is = 0; is < PARAM.inp.nspin; is++)
326+
{
327+
std::stringstream ss;
328+
ss << PARAM.globalv.global_out_dir << "SPIN" << is + 1 << "_CHG_INI.cube";
329+
ModuleIO::write_vdata_palgrid(this->Pgrid,
330+
this->pelec->charge->rho[is],
331+
is,
332+
PARAM.inp.nspin,
333+
istep,
334+
ss.str(),
335+
this->pelec->eferm.ef,
336+
&(ucell));
337+
}
338+
}
339+
340+
//! output total local potential of the initial charge density
341+
if (PARAM.inp.out_pot == 3)
342+
{
343+
for (int is = 0; is < PARAM.inp.nspin; is++)
344+
{
345+
std::stringstream ss;
346+
ss << PARAM.globalv.global_out_dir << "SPIN" << is + 1 << "_POT_INI.cube";
347+
ModuleIO::write_vdata_palgrid(this->Pgrid,
348+
this->pelec->pot->get_effective_v(is),
349+
is,
350+
PARAM.inp.nspin,
351+
istep,
352+
ss.str(),
353+
0.0, // efermi
354+
&(ucell),
355+
11, // precsion
356+
0); // out_fermi
357+
}
358+
}
359+
286360
return;
287361
}
288362

‎source/module_esolver/esolver_ks_pw.cpp‎

Lines changed: 30 additions & 107 deletions
Original file line numberDiff line numberDiff line change
@@ -2,6 +2,7 @@
22

33
#include "module_base/formatter.h"
44
#include "module_base/global_variable.h"
5+
#include "module_base/kernels/math_kernel_op.h"
56
#include "module_base/memory.h"
67
#include "module_elecstate/cal_ux.h"
78
#include "module_elecstate/elecstate_pw.h"
@@ -99,7 +100,7 @@ ESolver_KS_PW<T, Device>::~ESolver_KS_PW()
99100
std::cout << " ** Closing DSP Hardware..." << std::endl;
100101
mtfunc::dspDestoryHandle(GlobalV::MY_RANK);
101102
#endif
102-
if(PARAM.inp.device == "gpu" || PARAM.inp.precision == "single")
103+
if (PARAM.inp.device == "gpu" || PARAM.inp.precision == "single")
103104
{
104105
delete this->kspw_psi;
105106
}
@@ -181,7 +182,6 @@ void ESolver_KS_PW<T, Device>::before_all_runners(UnitCell& ucell, const Input_p
181182
&(this->pelec->f_en.vtxc));
182183
}
183184

184-
185185
//! initalize local pseudopotential
186186
this->locpp.init_vloc(ucell, this->pw_rhod);
187187
ModuleBase::GlobalFunc::DONE(GlobalV::ofs_running, "LOCAL POTENTIAL");
@@ -243,16 +243,6 @@ void ESolver_KS_PW<T, Device>::before_scf(UnitCell& ucell, const int istep)
243243

244244
this->p_psi_init->prepare_init(PARAM.inp.pw_seed);
245245
}
246-
if (ucell.ionic_position_updated)
247-
{
248-
this->CE.update_all_dis(ucell);
249-
this->CE.extrapolate_charge(&this->Pgrid,
250-
ucell,
251-
this->pelec->charge,
252-
&this->sf,
253-
GlobalV::ofs_running,
254-
GlobalV::ofs_warning);
255-
}
256246

257247
// init Hamilt, this should be allocated before each scf loop
258248
// Operators in HamiltPW should be reallocated once cell changed
@@ -262,77 +252,12 @@ void ESolver_KS_PW<T, Device>::before_scf(UnitCell& ucell, const int istep)
262252
// allocate HamiltPW
263253
this->allocate_hamilt(ucell);
264254

265-
//----------------------------------------------------------
266-
// about vdw, jiyy add vdwd3 and linpz add vdwd2
267-
//----------------------------------------------------------
268-
auto vdw_solver = vdw::make_vdw(ucell, PARAM.inp, &(GlobalV::ofs_running));
269-
if (vdw_solver != nullptr)
270-
{
271-
this->pelec->f_en.evdw = vdw_solver->get_energy();
272-
}
273-
274-
//----------------------------------------------------------
275-
// calculate ewald energy
276-
//----------------------------------------------------------
277-
if (!PARAM.inp.test_skip_ewald)
278-
{
279-
this->pelec->f_en.ewald_energy = H_Ewald_pw::compute_ewald(ucell, this->pw_rhod, this->sf.strucFac);
280-
}
281-
282-
//----------------------------------------------------------
283-
//! cal_ux should be called before init_scf because
284-
//! the direction of ux is used in noncoline_rho
285-
//----------------------------------------------------------
286-
elecstate::cal_ux(ucell);
287-
288255
//----------------------------------------------------------
289256
//! calculate the total local pseudopotential in real space
290257
//----------------------------------------------------------
291258
this->pelec
292259
->init_scf(istep, ucell, this->Pgrid, this->sf.strucFac, this->locpp.numeric, ucell.symm, (void*)this->pw_wfc);
293260

294-
//----------------------------------------------------------
295-
//! output the initial charge density
296-
//----------------------------------------------------------
297-
if (PARAM.inp.out_chg[0] == 2)
298-
{
299-
for (int is = 0; is < PARAM.inp.nspin; is++)
300-
{
301-
std::stringstream ss;
302-
ss << PARAM.globalv.global_out_dir << "SPIN" << is + 1 << "_CHG_INI.cube";
303-
ModuleIO::write_vdata_palgrid(this->Pgrid,
304-
this->pelec->charge->rho[is],
305-
is,
306-
PARAM.inp.nspin,
307-
istep,
308-
ss.str(),
309-
this->pelec->eferm.ef,
310-
&(ucell));
311-
}
312-
}
313-
314-
//----------------------------------------------------------
315-
//! output total local potential of the initial charge density
316-
//----------------------------------------------------------
317-
if (PARAM.inp.out_pot == 3)
318-
{
319-
for (int is = 0; is < PARAM.inp.nspin; is++)
320-
{
321-
std::stringstream ss;
322-
ss << PARAM.globalv.global_out_dir << "SPIN" << is + 1 << "_POT_INI.cube";
323-
ModuleIO::write_vdata_palgrid(this->Pgrid,
324-
this->pelec->pot->get_effective_v(is),
325-
is,
326-
PARAM.inp.nspin,
327-
istep,
328-
ss.str(),
329-
0.0, // efermi
330-
&(ucell),
331-
11, // precsion
332-
0); // out_fermi
333-
}
334-
}
335-
336261
//----------------------------------------------------------
337262
//! Symmetry_rho should behind init_scf, because charge should be
338263
//! initialized first. liuyu comment: Symmetry_rho should be located between
@@ -583,7 +508,7 @@ void ESolver_KS_PW<T, Device>::iter_finish(UnitCell& ucell, const int istep, int
583508

584509
// 2) Update USPP-related quantities
585510
// D in USPP needs vloc, thus needs update when veff updated
586-
// calculate the effective coefficient matrix for non-local pp projectors
511+
// calculate the effective coefficient matrix for non-local pp projectors
587512
// liuyu 2023-10-24
588513
if (PARAM.globalv.use_uspp)
589514
{
@@ -594,9 +519,7 @@ void ESolver_KS_PW<T, Device>::iter_finish(UnitCell& ucell, const int istep, int
594519
// 3) Print out electronic wavefunctions in pw basis
595520
if (PARAM.inp.out_wfc_pw == 1 || PARAM.inp.out_wfc_pw == 2)
596521
{
597-
if (iter % PARAM.inp.out_freq_elec == 0 ||
598-
iter == PARAM.inp.scf_nmax ||
599-
conv_esolver)
522+
if (iter % PARAM.inp.out_freq_elec == 0 || iter == PARAM.inp.scf_nmax || conv_esolver)
600523
{
601524
std::stringstream ssw;
602525
ssw << PARAM.globalv.global_out_dir << "WAVEFUNC";
@@ -931,35 +854,35 @@ void ESolver_KS_PW<T, Device>::after_all_runners(UnitCell& ucell)
931854

932855
#ifdef __MLKEDF
933856
// generate training data for ML-KEDF
934-
if(PARAM.inp.of_ml_gene_data == 1)
857+
if (PARAM.inp.of_ml_gene_data == 1)
935858
{
936859
this->pelec->pot->update_from_charge(this->pelec->charge, &ucell);
937860

938-
ML_data ml_data;
939-
ml_data.set_para(this->pelec->charge->nrxx,
940-
PARAM.inp.nelec,
941-
PARAM.inp.of_tf_weight,
942-
PARAM.inp.of_vw_weight,
943-
PARAM.inp.of_ml_chi_p,
944-
PARAM.inp.of_ml_chi_q,
945-
PARAM.inp.of_ml_chi_xi,
946-
PARAM.inp.of_ml_chi_pnl,
947-
PARAM.inp.of_ml_chi_qnl,
948-
PARAM.inp.of_ml_nkernel,
949-
PARAM.inp.of_ml_kernel,
950-
PARAM.inp.of_ml_kernel_scaling,
951-
PARAM.inp.of_ml_yukawa_alpha,
952-
PARAM.inp.of_ml_kernel_file,
953-
ucell.omega,
954-
this->pw_rho);
955-
956-
ml_data.generateTrainData_KS(this->kspw_psi,
957-
this->pelec,
958-
this->pw_wfc,
959-
this->pw_rho,
960-
ucell,
961-
this->pelec->pot->get_effective_v(0));
962-
}
861+
ML_data ml_data;
862+
ml_data.set_para(this->pelec->charge->nrxx,
863+
PARAM.inp.nelec,
864+
PARAM.inp.of_tf_weight,
865+
PARAM.inp.of_vw_weight,
866+
PARAM.inp.of_ml_chi_p,
867+
PARAM.inp.of_ml_chi_q,
868+
PARAM.inp.of_ml_chi_xi,
869+
PARAM.inp.of_ml_chi_pnl,
870+
PARAM.inp.of_ml_chi_qnl,
871+
PARAM.inp.of_ml_nkernel,
872+
PARAM.inp.of_ml_kernel,
873+
PARAM.inp.of_ml_kernel_scaling,
874+
PARAM.inp.of_ml_yukawa_alpha,
875+
PARAM.inp.of_ml_kernel_file,
876+
ucell.omega,
877+
this->pw_rho);
878+
879+
ml_data.generateTrainData_KS(this->kspw_psi,
880+
this->pelec,
881+
this->pw_wfc,
882+
this->pw_rho,
883+
ucell,
884+
this->pelec->pot->get_effective_v(0));
885+
}
963886
#endif
964887
}
965888

‎source/module_esolver/esolver_of.cpp‎

Lines changed: 0 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -258,17 +258,9 @@ void ESolver_OF::before_opt(const int istep, UnitCell& ucell)
258258
this->precip_dir_[is] = new std::complex<double>[pw_rho->npw];
259259
}
260260
}
261-
if (ucell.ionic_position_updated)
262-
{
263-
CE.update_all_dis(ucell);
264-
CE.extrapolate_charge(&Pgrid, ucell, pelec->charge, &sf, GlobalV::ofs_running, GlobalV::ofs_warning);
265-
}
266261

267262
this->pelec->init_scf(istep, ucell, Pgrid, sf.strucFac, locpp.numeric, ucell.symm);
268263

269-
// calculate ewald energy
270-
this->pelec->f_en.ewald_energy = H_Ewald_pw::compute_ewald(ucell, this->pw_rho, sf.strucFac);
271-
272264
Symmetry_rho srho;
273265
for (int is = 0; is < PARAM.inp.nspin; is++)
274266
{

0 commit comments

Comments
 (0)