Skip to content

Commit 84dfe33

Browse files
dyzhengLiu-RXA-006YuLiu98maki49
authored
Fix 1. get_pchg with nspin=4 and 2. out_mul with DFT+U and 3. fix updates in develop branch (#5859)
* Fix: get_pchg_lcao with nspin=4 * Fix: out_mul with DFT+U * Enable the support of singular number of electrions in DeepKS orbital label (#5793) * delete test print (#5795) * delete test print * change the boundry condition * Refactor: update logic of init_chg (#5801) * Fix: optimize lr_spectrum (#5805) * update the broadening function in lr_spectrum * optimize transition analysis and fix norm bug * fix some warnings about esolver (#5807) * Fix the Segmentation fault caused by zero atom case under LCAO basis. (#5821) * add const in module_lr (#5817) * Fix: support negative value in parse_expression (#5826) * Use template to reconstruct parse_expression * Feature: output R matrix at each MD step * Modify'matrix_HS' to 'matrix' for R matrix output * Merge branches 'develop' and 'develop' of https://github.com/1041176461/abacus-develop into develop * Fix: modify index in parse_expression * Fix: add regfree for parse_expression * Doc: update phonopy doc * Doc: update phonopy doc * fix tdos plot for nspin=2 * optimize dosplot for nspin=2 * fix legend for dosplot * Add files via upload * Update cal_edm_tddft.cpp * Refactor: modify exx loop for relax/md * Update result.ref * Fix wrong parameters in integrate test * Update result.ref * Update result.ref * Update result.ref * Update result.ref * Update result.ref * Update result.ref * Update result.ref * Update result.ref * add exx when istep>0 * Update dos.py * Update esolver_sdft_pw.cpp * Update lcao_before_scf.cpp * Update Exx_LRI_interface.h * Update Exx_LRI_interface.hpp * Fix: compile error * Fix: compile error * Fix: change HSE relax/md result.ref for new framework * Fix: compile error * compatible with exx_iter_finish * Add files via upload * Update esolver_ks_lcao_tddft.h * Update esolver_ks_lcao_tddft.cpp * Fix: support negative value in parse_expression * [pre-commit.ci lite] apply automatic fixes --------- Co-authored-by: jiyuang <jiyuyang@mail.ustc.com> Co-authored-by: Qianrui <76200646+Qianruipku@users.noreply.github.com> Co-authored-by: HTZhao <104255052+ESROAMER@users.noreply.github.com> Co-authored-by: maki49 <1579492865@qq.com> Co-authored-by: pre-commit-ci-lite[bot] <117423508+pre-commit-ci-lite[bot]@users.noreply.github.com> * Fix: dmin is replaced by dmax (#5829) * fix #5819 * Fix a bug and a magic number in module_exx_symmetry (#5848) * fix a magic number in get_euler_angle * do not allow higher symmetry of bvk supercell than the original cell --------- Co-authored-by: Liu Renxi <75369672+Liu-RX@users.noreply.github.com> Co-authored-by: liiutao <74701833+A-006@users.noreply.github.com> Co-authored-by: Yu Liu <77716030+YuLiu98@users.noreply.github.com> Co-authored-by: LUNASEA <33978601+maki49@users.noreply.github.com> Co-authored-by: Liang Sun <50293369+sunliang98@users.noreply.github.com> Co-authored-by: jiyuyang <1041176461@qq.com> Co-authored-by: jiyuang <jiyuyang@mail.ustc.com> Co-authored-by: Qianrui <76200646+Qianruipku@users.noreply.github.com> Co-authored-by: HTZhao <104255052+ESROAMER@users.noreply.github.com> Co-authored-by: maki49 <1579492865@qq.com> Co-authored-by: pre-commit-ci-lite[bot] <117423508+pre-commit-ci-lite[bot]@users.noreply.github.com> Co-authored-by: Haozhi Han <haozhi.han@stu.pku.edu.cn>
1 parent c53f445 commit 84dfe33

File tree

30 files changed

+248
-156
lines changed

30 files changed

+248
-156
lines changed

‎source/module_cell/module_neighbor/sltk_grid.cpp‎

Lines changed: 13 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -152,12 +152,19 @@ void Grid::setMemberVariables(std::ofstream& ofs_in, // output data to ofs
152152
this->clear_atoms();
153153

154154
// random selection, in order to estimate again.
155-
this->x_min = ucell.atoms[0].tau[0].x;
156-
this->y_min = ucell.atoms[0].tau[0].y;
157-
this->z_min = ucell.atoms[0].tau[0].z;
158-
this->x_max = ucell.atoms[0].tau[0].x;
159-
this->y_max = ucell.atoms[0].tau[0].y;
160-
this->z_max = ucell.atoms[0].tau[0].z;
155+
for (int it = 0; it < ucell.ntype; it++)
156+
{
157+
if (ucell.atoms[it].na > 0)
158+
{
159+
this->x_min = ucell.atoms[it].tau[0].x;
160+
this->y_min = ucell.atoms[it].tau[0].y;
161+
this->z_min = ucell.atoms[it].tau[0].z;
162+
this->x_max = ucell.atoms[it].tau[0].x;
163+
this->y_max = ucell.atoms[it].tau[0].y;
164+
this->z_max = ucell.atoms[it].tau[0].z;
165+
break;
166+
}
167+
}
161168

162169
ModuleBase::Vector3<double> vec1(ucell.latvec.e11, ucell.latvec.e12, ucell.latvec.e13);
163170
ModuleBase::Vector3<double> vec2(ucell.latvec.e21, ucell.latvec.e22, ucell.latvec.e23);

‎source/module_cell/update_cell.cpp‎

Lines changed: 9 additions & 29 deletions
Original file line numberDiff line numberDiff line change
@@ -384,41 +384,21 @@ void periodic_boundary_adjustment(Atom* atoms,
384384
// first adjust direct coordinates,
385385
// then update them into cartesian coordinates,
386386
//----------------------------------------------
387-
for (int i=0;i<ntype;i++) {
388-
Atom* atom = &atoms[i];
389-
for (int j=0;j<atom->na;j++) {
390-
printf("the taud is %f %f %f\n",atom->taud[j].x,atom->taud[j].y,atom->taud[j].z);
391-
}
392-
}
393387
for (int it = 0; it < ntype; it++) {
394388
Atom* atom = &atoms[it];
395389
for (int ia = 0; ia < atom->na; ia++) {
396390
// mohan update 2011-03-21
397-
if (atom->taud[ia].x < 0)
391+
for (int ik = 0; ik < 3; ik++)
398392
{
399-
atom->taud[ia].x += 1.0;
400-
}
401-
if (atom->taud[ia].y < 0)
402-
{
403-
atom->taud[ia].y += 1.0;
404-
}
405-
if (atom->taud[ia].z < 0)
406-
{
407-
atom->taud[ia].z += 1.0;
408-
}
409-
if (atom->taud[ia].x >= 1.0)
410-
{
411-
atom->taud[ia].x -= 1.0;
412-
}
413-
if (atom->taud[ia].y >= 1.0)
414-
{
415-
atom->taud[ia].y -= 1.0;
416-
}
417-
if (atom->taud[ia].z >= 1.0)
418-
{
419-
atom->taud[ia].z -= 1.0;
393+
if (atom->taud[ia][ik] < 0)
394+
{
395+
atom->taud[ia][ik] += 1.0;
396+
}
397+
if (atom->taud[ia][ik] >= 1.0)
398+
{
399+
atom->taud[ia][ik] -= 1.0;
400+
}
420401
}
421-
422402
if (atom->taud[ia].x < 0
423403
|| atom->taud[ia].y < 0
424404
|| atom->taud[ia].z < 0

‎source/module_elecstate/module_charge/charge_init.cpp‎

Lines changed: 74 additions & 45 deletions
Original file line numberDiff line numberDiff line change
@@ -37,6 +37,7 @@ void Charge::init_rho(elecstate::efermi& eferm_iout,
3737
this->pgrid = &pgrid;
3838

3939
bool read_error = false;
40+
bool read_kin_error = false;
4041
if (PARAM.inp.init_chg == "file" || PARAM.inp.init_chg == "auto")
4142
{
4243
GlobalV::ofs_running << " try to read charge from file" << std::endl;
@@ -101,72 +102,100 @@ void Charge::init_rho(elecstate::efermi& eferm_iout,
101102
}
102103
}
103104

104-
if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
105+
if (read_error)
105106
{
106-
GlobalV::ofs_running << " try to read kinetic energy density from file" << std::endl;
107-
// try to read charge from binary file first, which is the same as QE
108-
std::vector<std::complex<double>> kin_g_space(PARAM.inp.nspin * this->ngmc, {0.0, 0.0});
109-
std::vector<std::complex<double>*> kin_g;
110-
for (int is = 0; is < PARAM.inp.nspin; is++)
107+
const std::string warn_msg
108+
= " WARNING: \"init_chg\" is enabled but ABACUS failed to read charge density from file.\n"
109+
" Please check if there is SPINX_CHG.cube (X=1,...) or {suffix}-CHARGE-DENSITY.restart in the "
110+
"directory.\n";
111+
std::cout << std::endl << warn_msg;
112+
if (PARAM.inp.init_chg == "file")
111113
{
112-
kin_g.push_back(kin_g_space.data() + is * this->ngmc);
114+
ModuleBase::WARNING_QUIT("Charge::init_rho",
115+
"Failed to read in charge density from file.\nIf you want to use atomic "
116+
"charge initialization, \nplease set init_chg to atomic in INPUT.");
113117
}
118+
}
114119

115-
std::stringstream binary;
116-
binary << PARAM.globalv.global_readin_dir << PARAM.inp.suffix + "-TAU-DENSITY.restart";
117-
if (ModuleIO::read_rhog(binary.str(), rhopw, kin_g.data()))
120+
if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
121+
{
122+
// If the charge density is not read in, then the kinetic energy density is not read in either
123+
if (!read_error)
118124
{
119-
GlobalV::ofs_running << " Read in the kinetic energy density: " << binary.str() << std::endl;
120-
for (int is = 0; is < PARAM.inp.nspin; ++is)
125+
GlobalV::ofs_running << " try to read kinetic energy density from file" << std::endl;
126+
// try to read charge from binary file first, which is the same as QE
127+
std::vector<std::complex<double>> kin_g_space(PARAM.inp.nspin * this->ngmc, {0.0, 0.0});
128+
std::vector<std::complex<double>*> kin_g;
129+
for (int is = 0; is < PARAM.inp.nspin; is++)
121130
{
122-
rhopw->recip2real(kin_g[is], this->kin_r[is]);
131+
kin_g.push_back(kin_g_space.data() + is * this->ngmc);
123132
}
124-
}
125-
else
126-
{
127-
for (int is = 0; is < PARAM.inp.nspin; is++)
133+
134+
std::stringstream binary;
135+
binary << PARAM.globalv.global_readin_dir << PARAM.inp.suffix + "-TAU-DENSITY.restart";
136+
if (ModuleIO::read_rhog(binary.str(), rhopw, kin_g.data()))
128137
{
129-
std::stringstream ssc;
130-
ssc << PARAM.globalv.global_readin_dir << "SPIN" << is + 1 << "_TAU.cube";
131-
// mohan update 2012-02-10, sunliang update 2023-03-09
132-
if (ModuleIO::read_vdata_palgrid(
133-
pgrid,
134-
(PARAM.inp.esolver_type == "sdft" ? GlobalV::RANK_IN_STOGROUP : GlobalV::MY_RANK),
135-
GlobalV::ofs_running,
136-
ssc.str(),
137-
this->kin_r[is],
138-
ucell.nat))
138+
GlobalV::ofs_running << " Read in the kinetic energy density: " << binary.str() << std::endl;
139+
for (int is = 0; is < PARAM.inp.nspin; ++is)
139140
{
140-
GlobalV::ofs_running << " Read in the kinetic energy density: " << ssc.str() << std::endl;
141+
rhopw->recip2real(kin_g[is], this->kin_r[is]);
141142
}
142143
}
144+
else
145+
{
146+
for (int is = 0; is < PARAM.inp.nspin; is++)
147+
{
148+
std::stringstream ssc;
149+
ssc << PARAM.globalv.global_readin_dir << "SPIN" << is + 1 << "_TAU.cube";
150+
// mohan update 2012-02-10, sunliang update 2023-03-09
151+
if (ModuleIO::read_vdata_palgrid(
152+
pgrid,
153+
(PARAM.inp.esolver_type == "sdft" ? GlobalV::RANK_IN_STOGROUP : GlobalV::MY_RANK),
154+
GlobalV::ofs_running,
155+
ssc.str(),
156+
this->kin_r[is],
157+
ucell.nat))
158+
{
159+
GlobalV::ofs_running << " Read in the kinetic energy density: " << ssc.str() << std::endl;
160+
}
161+
else
162+
{
163+
read_kin_error = true;
164+
std::cout << " WARNING: \"init_chg\" is enabled but ABACUS failed to read kinetic energy "
165+
"density from file.\n"
166+
" Please check if there is SPINX_TAU.cube (X=1,...) or "
167+
"{suffix}-TAU-DENSITY.restart in the directory.\n"
168+
<< std::endl;
169+
break;
170+
}
171+
}
172+
}
173+
}
174+
else
175+
{
176+
read_kin_error = true;
143177
}
144178
}
145179
}
146-
if (read_error)
180+
181+
if (PARAM.inp.init_chg == "atomic" || read_error) // mohan add 2007-10-17
147182
{
148-
const std::string warn_msg = " WARNING: \"init_chg\" is enabled but ABACUS failed to read charge density from file.\n"
149-
" Please check if there is SPINX_CHG.cube (X=1,...) or {suffix}-CHARGE-DENSITY.restart in the directory.\n";
150-
std::cout << std::endl << warn_msg;
151-
if (PARAM.inp.init_chg == "auto")
183+
if (read_error)
152184
{
153-
std::cout << " Charge::init_rho: use atomic initialization instead." << std::endl << std::endl;
154-
}
155-
else if (PARAM.inp.init_chg == "file")
156-
{
157-
ModuleBase::WARNING_QUIT("Charge::init_rho", "Failed to read in charge density from file.\nIf you want to use atomic charge initialization, \nplease set init_chg to atomic in INPUT.");
185+
std::cout << " Charge::init_rho: use atomic initialization instead." << std::endl;
158186
}
187+
this->atomic_rho(PARAM.inp.nspin, ucell.omega, rho, strucFac, ucell);
159188
}
160189

161-
if (PARAM.inp.init_chg == "atomic" ||
162-
(PARAM.inp.init_chg == "auto" && read_error)) // mohan add 2007-10-17
190+
// wenfei 2021-7-29 : initial tau = 3/5 rho^2/3, Thomas-Fermi
191+
if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
163192
{
164-
this->atomic_rho(PARAM.inp.nspin, ucell.omega, rho, strucFac, ucell);
165-
166-
// liuyu 2023-06-29 : move here from atomic_rho(), which will be called several times in charge extrapolation
167-
// wenfei 2021-7-29 : initial tau = 3/5 rho^2/3, Thomas-Fermi
168-
if (XC_Functional::get_func_type() == 3 || XC_Functional::get_func_type() == 5)
193+
if (PARAM.inp.init_chg == "atomic" || read_kin_error)
169194
{
195+
if (read_kin_error)
196+
{
197+
std::cout << " Charge::init_rho: init kinetic energy density from rho." << std::endl;
198+
}
170199
const double pi = 3.141592653589790;
171200
const double fact = (3.0 / 5.0) * pow(3.0 * pi * pi, 2.0 / 3.0);
172201
for (int is = 0; is < PARAM.inp.nspin; ++is)

‎source/module_esolver/esolver_gets.cpp‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -147,7 +147,7 @@ void ESolver_GetS::runner(UnitCell& ucell, const int istep)
147147
}
148148

149149
void ESolver_GetS::after_all_runners(UnitCell& ucell) {};
150-
double ESolver_GetS::cal_energy() {};
150+
double ESolver_GetS::cal_energy() { return 0.0; };
151151
void ESolver_GetS::cal_force(UnitCell& ucell, ModuleBase::matrix& force) {};
152152
void ESolver_GetS::cal_stress(UnitCell& ucell, ModuleBase::matrix& stress) {};
153153

‎source/module_esolver/esolver_ks.h‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -43,7 +43,7 @@ class ESolver_KS : public ESolver_FP
4343
virtual void iter_init(UnitCell& ucell, const int istep, const int iter);
4444

4545
//! Something to do after hamilt2density function in each iter loop.
46-
virtual void iter_finish(UnitCell& ucell, const int istep, int& iter);
46+
virtual void iter_finish(UnitCell& ucell, const int istep, int& iter) override;
4747

4848
// calculate electron density from a specific Hamiltonian with ethr
4949
virtual void hamilt2density_single(UnitCell& ucell, const int istep, const int iter, const double ethr);

‎source/module_hamilt_lcao/module_deepks/LCAO_deepks_interface.cpp‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -69,7 +69,7 @@ void LCAO_Deepks_Interface<TK, TR>::out_deepks_labels(const double& etot,
6969

7070
if (PARAM.inp.deepks_bandgap)
7171
{
72-
const int nocc = PARAM.inp.nelec / 2;
72+
const int nocc = (PARAM.inp.nelec+1) / 2;
7373
std::vector<double> o_tot(nks);
7474
for (int iks = 0; iks < nks; ++iks)
7575
{

‎source/module_hsolver/diago_iter_assist.cpp‎

Lines changed: 5 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -118,6 +118,9 @@ void DiagoIterAssist<T, Device>::diagH_subspace(const hamilt::Hamilt<T, Device>*
118118
// after generation of H and S matrix, diag them
119119
DiagoIterAssist::diagH_LAPACK(nstart, n_band, hcc, scc, nstart, en, vcc);
120120

121+
122+
const int ld_temp = in_place ? dmax : dmin;
123+
121124
{ // code block to calculate evc
122125
gemm_op<T, Device>()(ctx,
123126
'N',
@@ -132,12 +135,12 @@ void DiagoIterAssist<T, Device>::diagH_subspace(const hamilt::Hamilt<T, Device>*
132135
nstart,
133136
&zero,
134137
temp,
135-
dmin);
138+
ld_temp);
136139
}
137140

138141
if (!in_place)
139142
{
140-
matrixSetToAnother<T, Device>()(ctx, n_band, temp, dmin, evc.get_pointer(), dmax);
143+
matrixSetToAnother<T, Device>()(ctx, n_band, temp, ld_temp, evc.get_pointer(), dmax);
141144
delmem_complex_op()(ctx, temp);
142145
}
143146
delmem_complex_op()(ctx, hcc);

‎source/module_io/get_pchg_lcao.cpp‎

Lines changed: 8 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -200,7 +200,9 @@ void IState_Charge::begin(Gint_k& gk,
200200
mode = 3;
201201
}
202202

203-
const int fermi_band = static_cast<int>((nelec + 1) / 2 + 1.0e-8);
203+
int fermi_band = static_cast<int>((nelec + 1) / 2 + 1.0e-8);
204+
// for nspin = 4 case, fermi_band is the number of electrons because of each band has 1 electron
205+
if(PARAM.inp.nspin == 4) fermi_band = nelec;
204206
std::cout << " number of electrons = " << nelec << std::endl;
205207
std::cout << " number of occupied bands = " << fermi_band << std::endl;
206208

@@ -216,15 +218,15 @@ void IState_Charge::begin(Gint_k& gk,
216218
elecstate::DensityMatrix<std::complex<double>, double> DM(this->ParaV, nspin_dm, kv.kvec_d, kv.get_nks() / nspin_dm);
217219

218220
#ifdef __MPI
219-
this->idmatrix(ib, nspin, nelec, nlocal, wg, DM, kv, if_separate_k);
221+
this->idmatrix(ib, nspin_dm, nelec, nlocal, wg, DM, kv, if_separate_k);
220222
#else
221223
ModuleBase::WARNING_QUIT("IState_Charge::begin", "The `pchg` calculation is only available for MPI now!");
222224
#endif
223225
// If contribution from different k-points need to be output separately
224226
if (if_separate_k)
225227
{
226228
// For multi-k, loop over all real k-points
227-
for (int ik = 0; ik < kv.get_nks() / nspin; ++ik)
229+
for (int ik = 0; ik < kv.get_nks() / nspin_dm; ++ik)
228230
{
229231
for (int is = 0; is < nspin; ++is)
230232
{
@@ -506,16 +508,8 @@ void IState_Charge::idmatrix(const int& ib,
506508
ModuleBase::TITLE("IState_Charge", "idmatrix");
507509
assert(wg.nr == kv.get_nks());
508510

509-
const int fermi_band = static_cast<int>((nelec + 1) / 2 + 1.0e-8);
510-
511511
// To ensure the normalization of charge density in multi-k calculation (if if_separate_k is true)
512-
double wg_sum_k = 0;
513-
double wg_sum_k_homo = 0;
514-
for (int ik = 0; ik < kv.get_nks() / nspin; ++ik)
515-
{
516-
wg_sum_k += wg(ik, ib);
517-
wg_sum_k_homo += wg(ik, fermi_band - 1);
518-
}
512+
double wg_sum_k = 1.0;
519513

520514
for (int ik = 0; ik < kv.get_nks(); ++ik)
521515
{
@@ -530,11 +524,11 @@ void IState_Charge::idmatrix(const int& ib,
530524
double wg_value;
531525
if (if_separate_k)
532526
{
533-
wg_value = (ib < fermi_band) ? wg_sum_k : wg_sum_k_homo;
527+
wg_value = wg_sum_k;
534528
}
535529
else
536530
{
537-
wg_value = (ib < fermi_band) ? wg(ik, ib) : wg(ik, fermi_band - 1);
531+
wg_value = wg(ik, 0);
538532
}
539533
wg_local[ib_local] = wg_value;
540534
}

‎source/module_io/input_conv.h‎

Lines changed: 5 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -12,8 +12,8 @@
1212
#include <iomanip>
1313
#include <iostream>
1414
#include <regex.h>
15-
#include <stdio.h>
16-
#include <string.h>
15+
#include <cstdio>
16+
#include <cstring>
1717
#include <string>
1818
#include <vector>
1919
#include <algorithm>
@@ -31,7 +31,7 @@ void tmp_convert();
3131
* @brief Pass the data members from the INPUT instance(defined in
3232
* module_io/input.cpp) to GlobalV and GlobalC.
3333
*/
34-
void Convert(void);
34+
void Convert();
3535

3636
/**
3737
* @brief To parse input parameters as expressions into vectors
@@ -47,7 +47,7 @@ void parse_expression(const std::string& fn, std::vector<T>& vec)
4747
{
4848
ModuleBase::TITLE("Input_Conv", "parse_expression");
4949
int count = 0;
50-
std::string pattern("([0-9]+\\*[0-9.]+|[0-9,.]+)");
50+
std::string pattern("([-+]?[0-9]+\\*[-+]?[0-9.]+|[-+]?[0-9,.]+)");
5151
std::vector<std::string> str;
5252
std::stringstream ss(fn);
5353
std::string section;
@@ -103,6 +103,7 @@ void parse_expression(const std::string& fn, std::vector<T>& vec)
103103
{
104104
size_t pos = sub_str.find("*");
105105
int num = stoi(sub_str.substr(0, pos));
106+
assert(num>=0);
106107
T occ = stof(sub_str.substr(pos + 1, sub_str.size()));
107108
// std::vector<double> ocp_temp(num, occ);
108109
// const std::vector<double>::iterator dest = vec.begin() + count;

‎source/module_io/output_mulliken.h‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -138,7 +138,7 @@ void cal_mag(Parallel_Orbitals* pv,
138138
auto sc_lambda
139139
= new hamilt::DeltaSpin<hamilt::OperatorLCAO<TK, double>>(nullptr,
140140
kv.kvec_d,
141-
nullptr,
141+
dynamic_cast<hamilt::HamiltLCAO<TK, double>*>(p_ham)->getHR(),
142142
ucell,
143143
&gd,
144144
two_center_bundle.overlap_orb_onsite.get(),
@@ -164,7 +164,7 @@ void cal_mag(Parallel_Orbitals* pv,
164164
auto sc_lambda = new hamilt::DeltaSpin<hamilt::OperatorLCAO<std::complex<double>, std::complex<double>>>(
165165
nullptr,
166166
kv.kvec_d,
167-
nullptr,
167+
dynamic_cast<hamilt::HamiltLCAO<std::complex<double>, std::complex<double>>*>(p_ham)->getHR(),
168168
ucell,
169169
&gd,
170170
two_center_bundle.overlap_orb_onsite.get(),

0 commit comments

Comments
 (0)