Skip to content

Commit 09fa28b

Browse files
ErjieWuErjieWu
committed
Refactor: move writing H(k)&S(k) from update_pot to after_scf.
1 parent a363fdf commit 09fa28b

File tree

2 files changed

+87
-153
lines changed

2 files changed

+87
-153
lines changed

‎source/module_esolver/esolver_ks_lcao.cpp‎

Lines changed: 87 additions & 91 deletions
Original file line numberDiff line numberDiff line change
@@ -757,86 +757,13 @@ void ESolver_KS_LCAO<TK, TR>::hamilt2density_single(UnitCell& ucell, int istep,
757757
//------------------------------------------------------------------------------
758758
//! the 12th function of ESolver_KS_LCAO: update_pot
759759
//! mohan add 2024-05-11
760-
//! 1) print Hamiltonian and Overlap matrix (why related to update_pot()?)
761-
//! 2) print wavefunctions (why related to update_pot()?)
762-
//! 3) print potential
760+
//! 1) print potential
763761
//------------------------------------------------------------------------------
764762
template <typename TK, typename TR>
765763
void ESolver_KS_LCAO<TK, TR>::update_pot(UnitCell& ucell, const int istep, const int iter)
766764
{
767765
ModuleBase::TITLE("ESolver_KS_LCAO", "update_pot");
768766

769-
// 1) print Hamiltonian and Overlap matrix
770-
if (this->conv_esolver || iter == PARAM.inp.scf_nmax)
771-
{
772-
if (!PARAM.globalv.gamma_only_local && (PARAM.inp.out_mat_hs[0] || PARAM.inp.deepks_v_delta))
773-
{
774-
this->GK.renew(true);
775-
}
776-
for (int ik = 0; ik < this->kv.get_nks(); ++ik)
777-
{
778-
if (PARAM.inp.out_mat_hs[0] || PARAM.inp.deepks_v_delta)
779-
{
780-
this->p_hamilt->updateHk(ik);
781-
}
782-
bool bit = false; // LiuXh, 2017-03-21
783-
// if set bit = true, there would be error in soc-multi-core
784-
// calculation, noted by zhengdy-soc
785-
if (this->psi != nullptr && (istep % PARAM.inp.out_interval == 0))
786-
{
787-
hamilt::MatrixBlock<TK> h_mat;
788-
hamilt::MatrixBlock<TK> s_mat;
789-
790-
this->p_hamilt->matrix(h_mat, s_mat);
791-
792-
if (PARAM.inp.out_mat_hs[0])
793-
{
794-
ModuleIO::save_mat(istep,
795-
h_mat.p,
796-
PARAM.globalv.nlocal,
797-
bit,
798-
PARAM.inp.out_mat_hs[1],
799-
1,
800-
PARAM.inp.out_app_flag,
801-
"H",
802-
"data-" + std::to_string(ik),
803-
this->pv,
804-
GlobalV::DRANK);
805-
ModuleIO::save_mat(istep,
806-
s_mat.p,
807-
PARAM.globalv.nlocal,
808-
bit,
809-
PARAM.inp.out_mat_hs[1],
810-
1,
811-
PARAM.inp.out_app_flag,
812-
"S",
813-
"data-" + std::to_string(ik),
814-
this->pv,
815-
GlobalV::DRANK);
816-
}
817-
#ifdef __DEEPKS
818-
if (PARAM.inp.deepks_out_labels && PARAM.inp.deepks_v_delta)
819-
{
820-
DeePKS_domain::save_h_mat(h_mat.p, this->pv.nloc);
821-
}
822-
#endif
823-
}
824-
}
825-
}
826-
827-
// 2) print wavefunctions
828-
if (elecstate::ElecStateLCAO<TK>::out_wfc_lcao && (this->conv_esolver || iter == PARAM.inp.scf_nmax)
829-
&& (istep % PARAM.inp.out_interval == 0))
830-
{
831-
ModuleIO::write_wfc_nao(elecstate::ElecStateLCAO<TK>::out_wfc_lcao,
832-
this->psi[0],
833-
this->pelec->ekb,
834-
this->pelec->wg,
835-
this->pelec->klist->kvec_c,
836-
this->pv,
837-
istep);
838-
}
839-
840767
if (!this->conv_esolver)
841768
{
842769
elecstate::cal_ux(ucell);
@@ -962,7 +889,9 @@ void ESolver_KS_LCAO<TK, TR>::iter_finish(UnitCell& ucell, const int istep, int&
962889
//! 1) call after_scf() of ESolver_KS
963890
//! 2) write density matrix for sparse matrix
964891
//! 4) write density matrix
965-
//! 6) write Exx matrix
892+
//! 5) write Exx matrix
893+
//! 6) write Hamiltonian and Overlap matrix
894+
//! 7) write wavefunctions
966895
//! 11) write deepks information
967896
//! 12) write rpa information
968897
//! 13) write HR in npz format
@@ -982,18 +911,18 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
982911
this->pelec->cal_tau(*(this->psi));
983912
}
984913

985-
//! 5) call after_scf() of ESolver_KS
914+
//! 2) call after_scf() of ESolver_KS
986915
ESolver_KS<TK>::after_scf(ucell, istep);
987916

988-
//! 6) write density matrix for sparse matrix
917+
//! 3) write density matrix for sparse matrix
989918
ModuleIO::write_dmr(dynamic_cast<const elecstate::ElecStateLCAO<TK>*>(this->pelec)->get_DM()->get_DMR_vector(),
990919
this->pv,
991920
PARAM.inp.out_dm1,
992921
false,
993922
PARAM.inp.out_app_flag,
994923
istep);
995924

996-
//! 7) write density matrix
925+
//! 4) write density matrix
997926
if (PARAM.inp.out_dm)
998927
{
999928
std::vector<double> efermis(PARAM.inp.nspin == 2 ? 2 : 1);
@@ -1010,7 +939,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
1010939
}
1011940

1012941
#ifdef __EXX
1013-
//! 8) write Hexx matrix for NSCF (see `out_chg` in docs/advanced/input_files/input-main.md)
942+
//! 5) write Hexx matrix for NSCF (see `out_chg` in docs/advanced/input_files/input-main.md)
1014943
if (PARAM.inp.calculation != "nscf")
1015944
{
1016945
if (GlobalC::exx_info.info_global.cal_exx && PARAM.inp.out_chg[0]
@@ -1029,7 +958,74 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
1029958
}
1030959
#endif
1031960

1032-
//! 9) Write DeePKS information
961+
// 6) write Hamiltonian and Overlap matrix
962+
if (!PARAM.globalv.gamma_only_local && (PARAM.inp.out_mat_hs[0] || PARAM.inp.deepks_v_delta || PARAM.inp.esolver_type == 'tddft'))
963+
{
964+
this->GK.renew(true);
965+
}
966+
for (int ik = 0; ik < this->kv.get_nks(); ++ik)
967+
{
968+
if (PARAM.inp.out_mat_hs[0] || PARAM.inp.deepks_v_delta)
969+
{
970+
this->p_hamilt->updateHk(ik);
971+
}
972+
bool bit = false; // LiuXh, 2017-03-21
973+
// if set bit = true, there would be error in soc-multi-core
974+
// calculation, noted by zhengdy-soc
975+
if (this->psi != nullptr && (istep % PARAM.inp.out_interval == 0))
976+
{
977+
hamilt::MatrixBlock<TK> h_mat;
978+
hamilt::MatrixBlock<TK> s_mat;
979+
980+
this->p_hamilt->matrix(h_mat, s_mat);
981+
982+
if (PARAM.inp.out_mat_hs[0])
983+
{
984+
ModuleIO::save_mat(istep,
985+
h_mat.p,
986+
PARAM.globalv.nlocal,
987+
bit,
988+
PARAM.inp.out_mat_hs[1],
989+
1,
990+
PARAM.inp.out_app_flag,
991+
"H",
992+
"data-" + std::to_string(ik),
993+
this->pv,
994+
GlobalV::DRANK);
995+
ModuleIO::save_mat(istep,
996+
s_mat.p,
997+
PARAM.globalv.nlocal,
998+
bit,
999+
PARAM.inp.out_mat_hs[1],
1000+
1,
1001+
PARAM.inp.out_app_flag,
1002+
"S",
1003+
"data-" + std::to_string(ik),
1004+
this->pv,
1005+
GlobalV::DRANK);
1006+
}
1007+
#ifdef __DEEPKS
1008+
if (PARAM.inp.deepks_out_labels && PARAM.inp.deepks_v_delta)
1009+
{
1010+
DeePKS_domain::save_h_mat(h_mat.p, this->pv.nloc);
1011+
}
1012+
#endif
1013+
}
1014+
}
1015+
1016+
// 7) write wavefunctions
1017+
if (elecstate::ElecStateLCAO<TK>::out_wfc_lcao && (istep % PARAM.inp.out_interval == 0))
1018+
{
1019+
ModuleIO::write_wfc_nao(elecstate::ElecStateLCAO<TK>::out_wfc_lcao,
1020+
this->psi[0],
1021+
this->pelec->ekb,
1022+
this->pelec->wg,
1023+
this->pelec->klist->kvec_c,
1024+
this->pv,
1025+
istep);
1026+
}
1027+
1028+
//! 8) Write DeePKS information
10331029
#ifdef __DEEPKS
10341030
std::shared_ptr<LCAO_Deepks> ld_shared_ptr(&GlobalC::ld, [](LCAO_Deepks*) {});
10351031
LCAO_Deepks_Interface LDI = LCAO_Deepks_Interface(ld_shared_ptr);
@@ -1051,7 +1047,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
10511047
ModuleBase::timer::tick("ESolver_KS_LCAO", "out_deepks_labels");
10521048
#endif
10531049

1054-
//! 10) Perform RDMFT calculations
1050+
//! 9) Perform RDMFT calculations
10551051
/******** test RDMFT *********/
10561052
if ( PARAM.inp.rdmft == true ) // rdmft, added by jghan, 2024-10-17
10571053
{
@@ -1077,7 +1073,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
10771073

10781074

10791075
#ifdef __EXX
1080-
// 11) Write RPA information.
1076+
// 10) Write RPA information.
10811077
if (PARAM.inp.rpa)
10821078
{
10831079
// ModuleRPA::DFT_RPA_interface
@@ -1092,7 +1088,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
10921088
}
10931089
#endif
10941090

1095-
// 12) write HR in npz format.
1091+
// 11) write HR in npz format.
10961092
if (PARAM.inp.out_hr_npz)
10971093
{
10981094
this->p_hamilt->updateHk(0); // first k point, up spin
@@ -1111,7 +1107,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
11111107
}
11121108
}
11131109

1114-
// 13) write density matrix in the 'npz' format.
1110+
// 12) write density matrix in the 'npz' format.
11151111
if (PARAM.inp.out_dm_npz)
11161112
{
11171113
const elecstate::DensityMatrix<TK, double>* dm
@@ -1126,7 +1122,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
11261122
}
11271123
}
11281124

1129-
//! 14) Print out information every 'out_interval' steps.
1125+
//! 13) Print out information every 'out_interval' steps.
11301126
if (PARAM.inp.calculation != "md" || istep % PARAM.inp.out_interval == 0)
11311127
{
11321128
//! Print out sparse matrix
@@ -1152,7 +1148,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
11521148
}
11531149
}
11541150

1155-
//! 15) Print out atomic magnetization only when 'spin_constraint' is on.
1151+
//! 14) Print out atomic magnetization only when 'spin_constraint' is on.
11561152
if (PARAM.inp.sc_mag_switch)
11571153
{
11581154
spinconstrain::SpinConstrain<TK>& sc = spinconstrain::SpinConstrain<TK>::getScInstance();
@@ -1161,14 +1157,14 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
11611157
sc.print_Mag_Force(GlobalV::ofs_running);
11621158
}
11631159

1164-
//! 16) Clean up RA.
1160+
//! 15) Clean up RA.
11651161
//! this should be last function and put it in the end, mohan request 2024-11-28
11661162
if (!PARAM.inp.cal_force && !PARAM.inp.cal_stress)
11671163
{
11681164
RA.delete_grid();
11691165
}
11701166

1171-
//! 17) Print out quasi-orbitals.
1167+
//! 16) Print out quasi-orbitals.
11721168
if (PARAM.inp.qo_switch)
11731169
{
11741170
toQO tqo(PARAM.inp.qo_basis, PARAM.inp.qo_strategy, PARAM.inp.qo_thr, PARAM.inp.qo_screening_coeff);
@@ -1183,7 +1179,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
11831179
tqo.calculate();
11841180
}
11851181

1186-
//! 18) Print out kinetic matrix.
1182+
//! 17) Print out kinetic matrix.
11871183
if (PARAM.inp.out_mat_tk[0])
11881184
{
11891185
hamilt::HS_Matrix_K<TK> hsk(&pv, true);
@@ -1218,7 +1214,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
12181214
delete ekinetic;
12191215
}
12201216

1221-
//! 19) Wannier 90 function, added by jingan in 2018.11.7
1217+
//! 18) Wannier 90 function, added by jingan in 2018.11.7
12221218
if (PARAM.inp.calculation == "nscf" && PARAM.inp.towannier90)
12231219
{
12241220
std::cout << FmtCore::format("\n * * * * * *\n << Start %s.\n", "Wave function to Wannier90");
@@ -1256,7 +1252,7 @@ void ESolver_KS_LCAO<TK, TR>::after_scf(UnitCell& ucell, const int istep)
12561252
std::cout << FmtCore::format(" >> Finish %s.\n * * * * * *\n", "Wave function to Wannier90");
12571253
}
12581254

1259-
//! 20) berry phase calculations, added by jingan
1255+
//! 19) berry phase calculations, added by jingan
12601256
if (PARAM.inp.calculation == "nscf" &&
12611257
berryphase::berry_phase_flag &&
12621258
ModuleSymmetry::Symmetry::symm_flag != 1)

‎source/module_esolver/esolver_ks_lcao_tddft.cpp‎

Lines changed: 0 additions & 62 deletions
Original file line numberDiff line numberDiff line change
@@ -173,68 +173,6 @@ void ESolver_KS_LCAO_TDDFT::iter_finish(UnitCell& ucell, const int istep, int& i
173173

174174
void ESolver_KS_LCAO_TDDFT::update_pot(UnitCell& ucell, const int istep, const int iter)
175175
{
176-
// print Hamiltonian and Overlap matrix
177-
if (this->conv_esolver)
178-
{
179-
if (!PARAM.globalv.gamma_only_local)
180-
{
181-
this->GK.renew(true);
182-
}
183-
for (int ik = 0; ik < kv.get_nks(); ++ik)
184-
{
185-
if (PARAM.inp.out_mat_hs[0])
186-
{
187-
this->p_hamilt->updateHk(ik);
188-
}
189-
bool bit = false; // LiuXh, 2017-03-21
190-
// if set bit = true, there would be error in soc-multi-core
191-
// calculation, noted by zhengdy-soc
192-
if (this->psi != nullptr && (istep % PARAM.inp.out_interval == 0))
193-
{
194-
hamilt::MatrixBlock<complex<double>> h_mat, s_mat;
195-
this->p_hamilt->matrix(h_mat, s_mat);
196-
if (PARAM.inp.out_mat_hs[0])
197-
{
198-
ModuleIO::save_mat(istep,
199-
h_mat.p,
200-
PARAM.globalv.nlocal,
201-
bit,
202-
PARAM.inp.out_mat_hs[1],
203-
1,
204-
PARAM.inp.out_app_flag,
205-
"H",
206-
"data-" + std::to_string(ik),
207-
this->pv,
208-
GlobalV::DRANK);
209-
210-
ModuleIO::save_mat(istep,
211-
s_mat.p,
212-
PARAM.globalv.nlocal,
213-
bit,
214-
PARAM.inp.out_mat_hs[1],
215-
1,
216-
PARAM.inp.out_app_flag,
217-
"S",
218-
"data-" + std::to_string(ik),
219-
this->pv,
220-
GlobalV::DRANK);
221-
}
222-
}
223-
}
224-
}
225-
226-
if (elecstate::ElecStateLCAO<std::complex<double>>::out_wfc_lcao
227-
&& (this->conv_esolver || iter == PARAM.inp.scf_nmax) && (istep % PARAM.inp.out_interval == 0))
228-
{
229-
ModuleIO::write_wfc_nao(elecstate::ElecStateLCAO<std::complex<double>>::out_wfc_lcao,
230-
this->psi[0],
231-
this->pelec->ekb,
232-
this->pelec->wg,
233-
this->pelec->klist->kvec_c,
234-
this->pv,
235-
istep);
236-
}
237-
238176
// Calculate new potential according to new Charge Density
239177
if (!this->conv_esolver)
240178
{

0 commit comments

Comments
 (0)