Skip to content

Commit 23e677e

Browse files
authored
Fix: Prevent missing EXX in RT-TDDFT evolution when using estep_per_md > 1 evolution strategy (#7360)
* Remove unnecessary cout in TDDFT current file * Fix RT-TDDFT EXX bug when using estep_per_md * Modify cout format * Fix a compiling issue with respect to std::vector * Update test 08_EXX/14_NO_TDDFT_PBE0
1 parent 8e50659 commit 23e677e

7 files changed

Lines changed: 204 additions & 137 deletions

File tree

source/source_estate/module_pot/H_TDDFT_pw.cpp

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -16,7 +16,7 @@ bool H_TDDFT_pw::is_initialized = false;
1616
double H_TDDFT_pw::amp;
1717

1818
// Used for calculating electric field force on ions, summing over directions
19-
vector<double> H_TDDFT_pw::global_vext_time = {0.0, 0.0, 0.0};
19+
std::vector<double> H_TDDFT_pw::global_vext_time = {0.0, 0.0, 0.0};
2020

2121
int H_TDDFT_pw::stype; // 0 : length gauge 1: velocity gauge
2222

source/source_estate/module_pot/H_TDDFT_pw.h

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -4,6 +4,8 @@
44
#include "pot_base.h"
55
#include "source_io/module_parameter/parameter.h"
66

7+
#include <vector>
8+
79
namespace elecstate
810
{
911

@@ -120,7 +122,7 @@ class H_TDDFT_pw : public PotBase
120122
static bool is_initialized; // static flag variable, used to ensure initialization only once
121123

122124
static double amp;
123-
static vector<double> global_vext_time;
125+
static std::vector<double> global_vext_time;
124126

125127
const UnitCell* ucell_ = nullptr;
126128

source/source_io/module_current/td_current_io_comm.cpp

Lines changed: 134 additions & 86 deletions
Large diffs are not rendered by default.

source/source_lcao/module_operator_lcao/op_exx_lcao.h

Lines changed: 49 additions & 40 deletions
Original file line numberDiff line numberDiff line change
@@ -3,12 +3,13 @@
33

44
#ifdef __EXX
55

6-
#include <RI/global/Tensor.h>
7-
#include <RI/ri/Cell_Nearest.h>
86
#include "operator_lcao.h"
97
#include "source_cell/klist.h"
108
#include "source_hamilt/module_xc/exx_info.h"
119

10+
#include <RI/global/Tensor.h>
11+
#include <RI/ri/Cell_Nearest.h>
12+
1213
namespace hamilt
1314
{
1415

@@ -21,13 +22,19 @@ class OperatorEXX : public T
2122
};
2223

2324
#endif
24-
enum Add_Hexx_Type { R, k };
25+
enum Add_Hexx_Type
26+
{
27+
R,
28+
k
29+
};
2530
template <typename TK, typename TR>
2631
class OperatorEXX<OperatorLCAO<TK, TR>> : public OperatorLCAO<TK, TR>
2732
{
2833
using TAC = std::pair<int, std::array<int, 3>>;
29-
public:
30-
OperatorEXX<OperatorLCAO<TK, TR>>(HS_Matrix_K<TK>* hsk_in,
34+
35+
public:
36+
OperatorEXX<OperatorLCAO<TK, TR>>(
37+
HS_Matrix_K<TK>* hsk_in,
3138
hamilt::HContainer<TR>* hR_in,
3239
const UnitCell& ucell,
3340
const K_Vectors& kv_in,
@@ -42,36 +49,37 @@ class OperatorEXX<OperatorLCAO<TK, TR>> : public OperatorLCAO<TK, TR>
4249
virtual void contributeHR() override;
4350

4451
private:
45-
Add_Hexx_Type add_hexx_type = Add_Hexx_Type::R;
46-
int current_spin = 0;
47-
bool HR_fixed_done = false;
48-
49-
std::vector<std::map<int, std::map<TAC, RI::Tensor<double>>>>* Hexxd = nullptr;
50-
std::vector<std::map<int, std::map<TAC, RI::Tensor<std::complex<double>>>>>* Hexxc = nullptr;
51-
52-
/// @brief the step of the outer loop.
53-
/// nullptr: no dependence on the number of two_level_step, contributeHk will do enerything normally.
54-
/// 0: the first outer loop. If restart, contributeHk will directly add Hexx to Hloc. else, do nothing.
55-
/// >0: not the first outer loop. contributeHk will do enerything normally.
56-
int* two_level_step = nullptr;
57-
/// @brief if restart, read and save Hexx, and directly use it during the first outer loop.
58-
bool restart = false;
59-
60-
const int istep = 0; // the ion step
61-
62-
void add_loaded_Hexx(const int ik);
63-
64-
const UnitCell& ucell;
65-
66-
const K_Vectors& kv;
67-
68-
// if k points has no shift, use cell_nearest to reduce the memory cost
69-
RI::Cell_Nearest<int, int, 3, double, 3> cell_nearest;
70-
bool use_cell_nearest = true;
71-
72-
/// @brief Hexxk for all k-points, only for the 1st scf loop ofrestart load
73-
std::vector<std::vector<double>> Hexxd_k_load;
74-
std::vector<std::vector<std::complex<double>>> Hexxc_k_load;
52+
Add_Hexx_Type add_hexx_type = Add_Hexx_Type::R;
53+
int current_spin = 0;
54+
bool HR_fixed_done = false;
55+
bool initial_gga_done = false; // Taoni Bao add 2026-05-18, to fix RT-TDDFT EXX missing problem in the evolution
56+
57+
std::vector<std::map<int, std::map<TAC, RI::Tensor<double>>>>* Hexxd = nullptr;
58+
std::vector<std::map<int, std::map<TAC, RI::Tensor<std::complex<double>>>>>* Hexxc = nullptr;
59+
60+
/// @brief the step of the outer loop.
61+
/// nullptr: no dependence on the number of two_level_step, contributeHk will do enerything normally.
62+
/// 0: the first outer loop. If restart, contributeHk will directly add Hexx to Hloc. else, do nothing.
63+
/// >0: not the first outer loop. contributeHk will do enerything normally.
64+
int* two_level_step = nullptr;
65+
/// @brief if restart, read and save Hexx, and directly use it during the first outer loop.
66+
bool restart = false;
67+
68+
const int istep = 0; // the ion step
69+
70+
void add_loaded_Hexx(const int ik);
71+
72+
const UnitCell& ucell;
73+
74+
const K_Vectors& kv;
75+
76+
// if k points has no shift, use cell_nearest to reduce the memory cost
77+
RI::Cell_Nearest<int, int, 3, double, 3> cell_nearest;
78+
bool use_cell_nearest = true;
79+
80+
/// @brief Hexxk for all k-points, only for the 1st scf loop ofrestart load
81+
std::vector<std::vector<double>> Hexxd_k_load;
82+
std::vector<std::vector<std::complex<double>>> Hexxc_k_load;
7583
};
7684

7785
using TAC = std::pair<int, std::array<int, 3>>;
@@ -81,14 +89,15 @@ RI::Cell_Nearest<int, int, 3, double, 3> init_cell_nearest(const UnitCell& ucell
8189
// allocate according to the read-in HexxR, used in nscf
8290
template <typename Tdata, typename TR>
8391
void reallocate_hcontainer(const std::vector<std::map<int, std::map<TAC, RI::Tensor<Tdata>>>>& Hexxs,
84-
HContainer<TR>* hR,
85-
const RI::Cell_Nearest<int, int, 3, double, 3>* const cell_nearest = nullptr);
92+
HContainer<TR>* hR,
93+
const RI::Cell_Nearest<int, int, 3, double, 3>* const cell_nearest = nullptr);
8694

8795
/// allocate according to BvK cells, used in scf
8896
template <typename TR>
89-
void reallocate_hcontainer(const int nat, HContainer<TR>* hR,
90-
const std::array<int, 3>& Rs_period,
91-
const RI::Cell_Nearest<int, int, 3, double, 3>* const cell_nearest = nullptr);
97+
void reallocate_hcontainer(const int nat,
98+
HContainer<TR>* hR,
99+
const std::array<int, 3>& Rs_period,
100+
const RI::Cell_Nearest<int, int, 3, double, 3>* const cell_nearest = nullptr);
92101

93102
} // namespace hamilt
94103
#endif // __EXX

source/source_lcao/module_operator_lcao/op_exx_lcao.hpp

Lines changed: 9 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -417,14 +417,22 @@ void OperatorEXX<OperatorLCAO<TK, TR>>::contributeHk(int ik)
417417
// 2. For the first ionic step:
418418
else if (this->istep == 0)
419419
{
420+
// If EXX is once turned on (two_level_step > 0), let OperatorEXX remember this
421+
if (this->two_level_step != nullptr && *this->two_level_step > 0)
422+
{
423+
this->initial_gga_done = true;
424+
}
425+
420426
// Check if we are in the pre-convergence stage of the two-level SCF (i.e., the pure GGA loop)
421427
bool in_gga_pre_loop = (this->two_level_step != nullptr && *this->two_level_step == 0);
422428

423429
// Check if a high-quality initial guess is missing
424430
bool lacks_good_guess = (!this->restart);
425431

426432
// If in the pre-convergence loop and lacking a good initial guess, skip adding the EXX contribution
427-
if (in_gga_pre_loop && lacks_good_guess)
433+
// Taoni Bao add 2026-05-18, only skip EXX if initial GGA loop is not done
434+
// Fix RT-TDDFT EXX missing problem in the evolution
435+
if (in_gga_pre_loop && lacks_good_guess && !this->initial_gga_done)
428436
{
429437
return; // In the non-EXX loop, skip adding EXX contribution
430438
}
Lines changed: 6 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
1-
1 -2.4064298950993433e-19 2.4067857486819294e-19 2.4048227372973307e-19
2-
2 1.3583754132510189e-08 -1.3583754049585592e-08 2.4306592013239020e-06
3-
3 3.8454029971142809e-08 -3.8454029713369349e-08 7.2734897464115248e-06
4-
4 4.3394109439952714e-08 -4.3394109112344272e-08 9.6350822407378473e-06
5-
5 2.6144935011100098e-08 -2.6144935077756982e-08 9.4987415747287647e-06
6-
6 6.1114518547755749e-10 -6.1114527641090651e-10 9.2969143096824206e-06
1+
1 -2.4064812749552146e-19 2.4058918971736318e-19 2.4056235004085640e-19
2+
2 1.3401891378677289e-08 -1.3401891246959429e-08 2.4369369671292067e-06
3+
3 3.9314244407525817e-08 -3.9314246669204259e-08 7.3211566892522424e-06
4+
4 4.7934170184878929e-08 -4.7934175958037426e-08 9.7956029317437335e-06
5+
5 3.4376283041277918e-08 -3.4376284895741181e-08 9.8560550504479175e-06
6+
6 5.9843794461236954e-09 -5.9843813295237610e-09 9.9145073758045465e-06
Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
etotref -170.67031862
2-
etotperatomref -85.33515931
1+
etotref -203.7710904354790
2+
etotperatomref -101.8855452177
33
CompareCurrent_pass 0
44
totaltimeref 16.38

0 commit comments

Comments
 (0)