Skip to content

Commit

Permalink
remove LM
Browse files Browse the repository at this point in the history
  • Loading branch information
@maki49
maki49 committed Jul 12, 2024
1 parent e7a76bd commit 6ffbea4
Show file tree
Hide file tree
Showing 20 changed files with 306 additions and 340 deletions.
132 changes: 63 additions & 69 deletions source/module_esolver/esolver_ks_lcao.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -63,22 +63,6 @@ ESolver_KS_LCAO<TK, TR>::ESolver_KS_LCAO()
{
this->classname = "ESolver_KS_LCAO";
this->basisname = "LCAO";

// the following EXX code should be removed to other places, mohan 2024/05/11
#ifdef __EXX
if (GlobalC::exx_info.info_ri.real_number)
{
this->exx_lri_double = std::make_shared<Exx_LRI<double>>(GlobalC::exx_info.info_ri);
this->exd = std::make_shared<Exx_LRI_Interface<TK, double>>(this->exx_lri_double);
this->LM.Hexxd = &this->exd->get_Hexxs();
}
else
{
this->exx_lri_complex = std::make_shared<Exx_LRI<std::complex<double>>>(GlobalC::exx_info.info_ri);
this->exc = std::make_shared<Exx_LRI_Interface<TK, std::complex<double>>>(this->exx_lri_complex);
this->LM.Hexxc = &this->exc->get_Hexxs();
}
#endif
}

//------------------------------------------------------------------------------
Expand Down Expand Up @@ -162,9 +146,6 @@ void ESolver_KS_LCAO<TK, TR>::before_all_runners(Input& inp, UnitCell& ucell)
this->init_basis_lcao(inp, ucell);
//------------------init Basis_lcao----------------------

//! pass basis-pointer to EState and Psi
this->LM.ParaV = &(this->ParaV);

// 5) initialize density matrix
dynamic_cast<elecstate::ElecStateLCAO<TK>*>(this->pelec)
->init_DM(&this->kv, &(this->ParaV), GlobalV::NSPIN);
Expand All @@ -186,14 +167,20 @@ void ESolver_KS_LCAO<TK, TR>::before_all_runners(Input& inp, UnitCell& ucell)
// PLEASE simplify the Exx_Global interface
if (GlobalV::CALCULATION == "scf" || GlobalV::CALCULATION == "relax"
|| GlobalV::CALCULATION == "cell-relax"
|| GlobalV::CALCULATION == "md") {
if (GlobalC::exx_info.info_global.cal_exx) {
|| GlobalV::CALCULATION == "md")
{
if (GlobalC::exx_info.info_global.cal_exx)
{
XC_Functional::set_xc_first_loop(ucell);
if (GlobalC::exx_info.info_ri.real_number) {
// initialize 2-center radial tables for EXX-LRI
if (GlobalC::exx_info.info_ri.real_number)
{
this->exx_lri_double = std::make_shared<Exx_LRI<double>>(GlobalC::exx_info.info_ri);
this->exx_lri_double->init(MPI_COMM_WORLD, this->kv);
}
else
{
this->exx_lri_complex = std::make_shared<Exx_LRI<std::complex<double>>>(GlobalC::exx_info.info_ri);
this->exx_lri_complex->init(MPI_COMM_WORLD, this->kv);
}
}
Expand All @@ -202,7 +189,7 @@ void ESolver_KS_LCAO<TK, TR>::before_all_runners(Input& inp, UnitCell& ucell)

// 8) initialize DFT+U
if (GlobalV::dft_plus_u) {
GlobalC::dftu.init(ucell, this->LM.ParaV, this->kv.get_nks());
GlobalC::dftu.init(ucell, &this->ParaV, this->kv.get_nks());
}

// 9) initialize ppcell
Expand Down Expand Up @@ -276,7 +263,7 @@ void ESolver_KS_LCAO<TK, TR>::init_after_vc(Input& inp, UnitCell& ucell)
this->pw_big);

dynamic_cast<elecstate::ElecStateLCAO<TK>*>(this->pelec)
->init_DM(&this->kv, this->LM.ParaV, GlobalV::NSPIN);
->init_DM(&this->kv, &this->ParaV, GlobalV::NSPIN);

GlobalC::ppcell.init_vloc(GlobalC::ppcell.vloc, this->pw_rho);

Expand Down Expand Up @@ -473,21 +460,26 @@ void ESolver_KS_LCAO<TK, TR>::after_all_runners()
if (PARAM.inp.out_mat_xc)
{
ModuleIO::write_Vxc<TK, TR>(GlobalV::NSPIN,
GlobalV::NLOCAL,
GlobalV::DRANK,
*this->psi,
GlobalC::ucell,
this->sf,
*this->pw_rho,
*this->pw_rhod,
GlobalC::ppcell.vloc,
*this->pelec->charge,
this->GG,
this->GK,
this->LM,
this->kv,
this->pelec->wg,
GlobalC::GridD);
GlobalV::NLOCAL,
GlobalV::DRANK,
&this->ParaV,
*this->psi,
GlobalC::ucell,
this->sf,
*this->pw_rho,
*this->pw_rhod,
GlobalC::ppcell.vloc,
*this->pelec->charge,
this->GG,
this->GK,
this->kv,
this->pelec->wg,
GlobalC::GridD
#ifdef __EXX
, this->exx_lri_double ? &this->exx_lri_double->Hexxs : nullptr
, this->exx_lri_complex ? &this->exx_lri_complex->Hexxs : nullptr
#endif
);
}

ModuleBase::timer::tick("ESolver_KS_LCAO", "after_all_runners");
Expand Down Expand Up @@ -848,7 +840,7 @@ void ESolver_KS_LCAO<TK, TR>::hamilt2density(int istep, int iter, double ethr)
if (GlobalV::sc_mag_switch)
{
SpinConstrain<TK, base_device::DEVICE_CPU>& sc = SpinConstrain<TK, base_device::DEVICE_CPU>::getScInstance();
sc.cal_MW(iter, &(this->LM));
sc.cal_MW(iter, this->p_hamilt);
}

// 9) use new charge density to calculate energy
Expand Down Expand Up @@ -1015,8 +1007,7 @@ void ESolver_KS_LCAO<TK, TR>::iter_finish(int iter)
for (int ik = 0; ik < this->kv.get_nks(); ++ik) {
Hexxk_save.set_zero_hk();
hamilt::OperatorEXX<hamilt::OperatorLCAO<TK, TR>> opexx_save(&Hexxk_save,
&this->LM,
hamilt::OperatorEXX<hamilt::OperatorLCAO<TK, TR>> opexx_save(&Hexxk_save,
nullptr,
this->kv);
Expand Down Expand Up @@ -1069,7 +1060,7 @@ void ESolver_KS_LCAO<TK, TR>::iter_finish(int iter)
if (GlobalV::sc_mag_switch)
{
SpinConstrain<TK, base_device::DEVICE_CPU>& sc = SpinConstrain<TK, base_device::DEVICE_CPU>::getScInstance();
sc.cal_MW(iter, &(this->LM));
sc.cal_MW(iter, this->p_hamilt);
}

// 6) calculate the total energy.
Expand Down Expand Up @@ -1315,20 +1306,24 @@ bool ESolver_KS_LCAO<TK, TR>::do_after_converge(int& iter)
}
#endif
#ifdef __EXX
if (GlobalC::exx_info.info_ri.real_number) {
return this->exd->exx_after_converge(
*this->p_hamilt,
*dynamic_cast<const elecstate::ElecStateLCAO<TK>*>(this->pelec)
->get_DM(),
this->kv,
iter);
} else {
return this->exc->exx_after_converge(
*this->p_hamilt,
*dynamic_cast<const elecstate::ElecStateLCAO<TK>*>(this->pelec)
->get_DM(),
this->kv,
iter);
if (GlobalC::exx_info.info_global.cal_exx)
{
if (GlobalC::exx_info.info_ri.real_number) {
return this->exd->exx_after_converge(
*this->p_hamilt,
*dynamic_cast<const elecstate::ElecStateLCAO<TK>*>(this->pelec)
->get_DM(),
this->kv,
iter);
}
else {
return this->exc->exx_after_converge(
*this->p_hamilt,
*dynamic_cast<const elecstate::ElecStateLCAO<TK>*>(this->pelec)
->get_DM(),
this->kv,
iter);
}
}
#endif // __EXX

Expand All @@ -1353,18 +1348,17 @@ template <typename TK, typename TR>
ModuleIO::Output_Mat_Sparse<TK> ESolver_KS_LCAO<TK, TR>::create_Output_Mat_Sparse(int istep)
{
return ModuleIO::Output_Mat_Sparse<TK>(hsolver::HSolverLCAO<TK>::out_mat_hsR,
hsolver::HSolverLCAO<TK>::out_mat_dh,
hsolver::HSolverLCAO<TK>::out_mat_t,
INPUT.out_mat_r,
istep,
this->pelec->pot->get_effective_v(),
this->ParaV,
this->GK, // mohan add 2024-04-01
two_center_bundle_,
this->LM,
GlobalC::GridD, // mohan add 2024-04-06
this->kv,
this->p_hamilt);
hsolver::HSolverLCAO<TK>::out_mat_dh,
hsolver::HSolverLCAO<TK>::out_mat_t,
INPUT.out_mat_r,
istep,
this->pelec->pot->get_effective_v(),
this->ParaV,
this->GK, // mohan add 2024-04-01
two_center_bundle_,
GlobalC::GridD, // mohan add 2024-04-06
this->kv,
this->p_hamilt);
}

//------------------------------------------------------------------------------
Expand Down
3 changes: 0 additions & 3 deletions source/module_esolver/esolver_ks_lcao.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -76,9 +76,6 @@ class ESolver_KS_LCAO : public ESolver_KS<TK> {
// used for gamma only algorithms.
Gint_Gamma GG;

// we will get rid of this class soon, don't use it, mohan 2024-03-28
LCAO_Matrix LM;

Grid_Technique GridT;

TwoCenterBundle two_center_bundle_;
Expand Down
17 changes: 9 additions & 8 deletions source/module_esolver/esolver_ks_lcao_elec.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -100,7 +100,7 @@ void ESolver_KS_LCAO<TK, TR>::set_matrix_grid(Record_adj& ra)
// (2)For each atom, calculate the adjacent atoms in different cells
// and allocate the space for H(R) and S(R).
// If k point is used here, allocate HlocR after atom_arrange.
Parallel_Orbitals* pv = this->LM.ParaV;
Parallel_Orbitals* pv = &this->ParaV;
ra.for_2d(*pv, GlobalV::GAMMA_ONLY_LOCAL);
if (!GlobalV::GAMMA_ONLY_LOCAL)
{
Expand Down Expand Up @@ -183,25 +183,25 @@ void ESolver_KS_LCAO<TK, TR>::beforesolver(const int istep)
this->p_hamilt = new hamilt::HamiltLCAO<TK, TR>(
GlobalV::GAMMA_ONLY_LOCAL ? &(this->GG) : nullptr,
GlobalV::GAMMA_ONLY_LOCAL ? nullptr : &(this->GK),
&(this->LM),
&this->ParaV,
this->pelec->pot,
this->kv,
two_center_bundle_,
DM
#ifdef __EXX
DM,
GlobalC::exx_info.info_ri.real_number ? &this->exd->two_level_step : &this->exc->two_level_step);
#else
DM);
, GlobalC::exx_info.info_ri.real_number ? &this->exd->two_level_step : &this->exc->two_level_step
, GlobalC::exx_info.info_ri.real_number ? &exx_lri_double->Hexxs : nullptr
, GlobalC::exx_info.info_ri.real_number ? nullptr : &exx_lri_complex->Hexxs
#endif
);
}

#ifdef __DEEPKS
// for each ionic step, the overlap <psi|alpha> must be rebuilt
// since it depends on ionic positions
if (GlobalV::deepks_setorb)
{
const Parallel_Orbitals* pv = this->LM.ParaV;
const Parallel_Orbitals* pv = &this->ParaV;
// build and save <psi(0)|alpha(R)> at beginning
GlobalC::ld.build_psialpha(GlobalV::CAL_FORCE,
GlobalC::ucell,
Expand Down Expand Up @@ -284,10 +284,12 @@ void ESolver_KS_LCAO<TK, TR>::before_scf(const int istep)
#ifdef __EXX // set xc type before the first cal of xc in pelec->init_scf
if (GlobalC::exx_info.info_ri.real_number)
{
this->exd = std::make_shared<Exx_LRI_Interface<TK, double>>(exx_lri_double);
this->exd->exx_beforescf(this->kv, *this->p_chgmix);
}
else
{
this->exc = std::make_shared<Exx_LRI_Interface<TK, std::complex<double>>>(exx_lri_complex);
this->exc->exx_beforescf(this->kv, *this->p_chgmix);
}
#endif // __EXX
Expand Down Expand Up @@ -583,7 +585,6 @@ void ESolver_KS_LCAO<std::complex<double>, std::complex<double>>::get_S(void)
GlobalV::test_atom_input);

this->RA.for_2d(this->ParaV, GlobalV::GAMMA_ONLY_LOCAL);
this->LM.ParaV = &this->ParaV;
if (this->p_hamilt == nullptr) {
this->p_hamilt = new hamilt::HamiltLCAO<std::complex<double>,
std::complex<double>>(
Expand Down
52 changes: 24 additions & 28 deletions source/module_esolver/esolver_ks_lcao_tddft.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -89,18 +89,14 @@ void ESolver_KS_LCAO_TDDFT::before_all_runners(Input& inp, UnitCell& ucell)
// 5) allocate H and S matrices according to computational resources
LCAO_domain::divide_HS_in_frag(GlobalV::GAMMA_ONLY_LOCAL, ParaV, kv.get_nks());

// this part will be updated soon
// pass Hamilt-pointer to Operator
this->LM.ParaV = &(this->ParaV);

// 6) initialize Density Matrix
dynamic_cast<elecstate::ElecStateLCAO<std::complex<double>>*>(this->pelec)
->init_DM(&kv, this->LM.ParaV, GlobalV::NSPIN);
->init_DM(&kv, &this->ParaV, GlobalV::NSPIN);

// 7) initialize Hsolver
if (this->phsol == nullptr)
{
this->phsol = new hsolver::HSolverLCAO<std::complex<double>>(this->LM.ParaV);
this->phsol = new hsolver::HSolverLCAO<std::complex<double>>(&this->ParaV);
this->phsol->method = GlobalV::KS_SOLVER;
}

Expand Down Expand Up @@ -256,28 +252,28 @@ void ESolver_KS_LCAO_TDDFT::update_pot(const int istep, const int iter)
if (hsolver::HSolverLCAO<std::complex<double>>::out_mat_hs[0])
{
ModuleIO::save_mat(istep,
h_mat.p,
GlobalV::NLOCAL,
bit,
hsolver::HSolverLCAO<std::complex<double>>::out_mat_hs[1],
1,
GlobalV::out_app_flag,
"H",
"data-" + std::to_string(ik),
*this->LM.ParaV,
GlobalV::DRANK);
h_mat.p,
GlobalV::NLOCAL,
bit,
hsolver::HSolverLCAO<std::complex<double>>::out_mat_hs[1],
1,
GlobalV::out_app_flag,
"H",
"data-" + std::to_string(ik),
this->ParaV,
GlobalV::DRANK);

ModuleIO::save_mat(istep,
s_mat.p,
GlobalV::NLOCAL,
bit,
hsolver::HSolverLCAO<std::complex<double>>::out_mat_hs[1],
1,
GlobalV::out_app_flag,
"S",
"data-" + std::to_string(ik),
*this->LM.ParaV,
GlobalV::DRANK);
s_mat.p,
GlobalV::NLOCAL,
bit,
hsolver::HSolverLCAO<std::complex<double>>::out_mat_hs[1],
1,
GlobalV::out_app_flag,
"S",
"data-" + std::to_string(ik),
this->ParaV,
GlobalV::DRANK);
}
}
}
Expand Down Expand Up @@ -311,8 +307,8 @@ void ESolver_KS_LCAO_TDDFT::update_pot(const int istep, const int iter)
}

const int nloc = this->ParaV.nloc;
const int ncol_nbands = this->LM.ParaV->ncol_bands;
const int nrow = this->LM.ParaV->nrow;
const int ncol_nbands = this->ParaV.ncol_bands;
const int nrow = this->ParaV.nrow;
const int nbands = GlobalV::NBANDS;
const int nlocal = GlobalV::NLOCAL;

Expand Down
13 changes: 6 additions & 7 deletions source/module_esolver/io_npz.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -458,17 +458,16 @@ void ESolver_KS_LCAO<TK, TR>::output_mat_npz(std::string& zipname, const hamilt:

}
#else
const Parallel_Orbitals* paraV = this->LM->ParaV;
auto row_indexes = paraV->get_indexes_row();
auto col_indexes = paraV->get_indexes_col();
auto row_indexes = paraV.get_indexes_row();
auto col_indexes = paraV.get_indexes_col();
for(int iap=0;iap<hR.size_atom_pairs();++iap)
{
int atom_i = hR.get_atom_pair(iap).get_atom_i();
int atom_j = hR.get_atom_pair(iap).get_atom_j();
int start_i = paraV->atom_begin_row[atom_i];
int start_j = paraV->atom_begin_col[atom_j];
int row_size = paraV->get_row_size(atom_i);
int col_size = paraV->get_col_size(atom_j);
int start_i = paraV.atom_begin_row[atom_i];
int start_j = paraV.atom_begin_col[atom_j];
int row_size = paraV.get_row_size(atom_i);
int col_size = paraV.get_col_size(atom_j);
for(int iR=0;iR<hR.get_atom_pair(iap).get_R_size();++iR)
{
auto& matrix = hR.get_atom_pair(iap).get_HR_values(iR);
Expand Down
Loading

0 comments on commit 6ffbea4

Please sign in to comment.