Skip to content

Commit

Permalink
Refactor: remove GlobalC::ORB in module_hamilt_lcao (#5068)
Browse files Browse the repository at this point in the history
* remove GlobalC::ORB in module_hamilt_lcao

* [pre-commit.ci lite] apply automatic fixes

---------

Co-authored-by: pre-commit-ci-lite[bot] <117423508+pre-commit-ci-lite[bot]@users.noreply.github.com>
  • Loading branch information
jinzx10 and pre-commit-ci-lite[bot] authored Sep 11, 2024
1 parent bd0e6a6 commit 8f86cad
Show file tree
Hide file tree
Showing 42 changed files with 310 additions and 225 deletions.
8 changes: 8 additions & 0 deletions source/module_basis/module_ao/ORB_read.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -46,6 +46,14 @@ const LCAO_Orbitals& LCAO_Orbitals::get_const_instance()
return GlobalC::ORB;
}

std::vector<double> LCAO_Orbitals::cutoffs() const {
std::vector<double> cutoffs(ntype);
for (int it = 0; it < ntype; ++it) {
cutoffs[it] = Phi[it].getRcut();
}
return cutoffs;
}

void LCAO_Orbitals::init(
std::ofstream& ofs_in,
const int& ntype,
Expand Down
2 changes: 2 additions & 0 deletions source/module_basis/module_ao/ORB_read.h
Original file line number Diff line number Diff line change
Expand Up @@ -83,6 +83,8 @@ class LCAO_Orbitals
//caoyu add 2021-05-24
const double& get_rcutmax_Phi() const { return rcutmax_Phi; }

std::vector<double> cutoffs() const;

/// numerical atomic orbitals
Numerical_Orbital* Phi;

Expand Down
8 changes: 5 additions & 3 deletions source/module_esolver/esolver_ks_lcao.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -171,7 +171,8 @@ void ESolver_KS_LCAO<TK, TR>::before_all_runners(const Input_para& inp, UnitCell
inp.lcao_dr,
inp.lcao_rmax,
ucell,
two_center_bundle_);
two_center_bundle_,
orb_);
//------------------init Basis_lcao----------------------

// 5) initialize density matrix
Expand All @@ -190,7 +191,7 @@ void ESolver_KS_LCAO<TK, TR>::before_all_runners(const Input_para& inp, UnitCell
// 6) initialize Hamilt in LCAO
// * allocate H and S matrices according to computational resources
// * set the 'trace' between local H/S and global H/S
LCAO_domain::divide_HS_in_frag(PARAM.globalv.gamma_only_local, pv, this->kv.get_nks());
LCAO_domain::divide_HS_in_frag(PARAM.globalv.gamma_only_local, pv, this->kv.get_nks(), orb_);

#ifdef __EXX
// 7) initialize exx
Expand All @@ -211,7 +212,7 @@ void ESolver_KS_LCAO<TK, TR>::before_all_runners(const Input_para& inp, UnitCell

// 8) initialize DFT+U
if (PARAM.inp.dft_plus_u) {
GlobalC::dftu.init(ucell, &this->pv, this->kv.get_nks());
GlobalC::dftu.init(ucell, &this->pv, this->kv.get_nks(), orb_);
}

// 9) initialize ppcell
Expand Down Expand Up @@ -309,6 +310,7 @@ void ESolver_KS_LCAO<TK, TR>::cal_force(ModuleBase::matrix& force)
this->GG, // mohan add 2024-04-01
this->GK, // mohan add 2024-04-01
two_center_bundle_,
orb_,
force,
this->scs,
this->sf,
Expand Down
6 changes: 4 additions & 2 deletions source/module_esolver/esolver_ks_lcao_tddft.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -92,10 +92,12 @@ void ESolver_KS_LCAO_TDDFT::before_all_runners(const Input_para& inp, UnitCell&
inp.lcao_dr,
inp.lcao_rmax,
ucell,
two_center_bundle_);
two_center_bundle_,
orb_);

// 5) allocate H and S matrices according to computational resources
LCAO_domain::divide_HS_in_frag(PARAM.globalv.gamma_only_local, this->pv, kv.get_nks());
LCAO_domain::divide_HS_in_frag(PARAM.globalv.gamma_only_local, this->pv, kv.get_nks(), orb_);


// 6) initialize Density Matrix
dynamic_cast<elecstate::ElecStateLCAO<std::complex<double>>*>(this->pelec)->init_DM(&kv, &this->pv, GlobalV::NSPIN);
Expand Down
3 changes: 2 additions & 1 deletion source/module_esolver/lcao_before_scf.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -94,7 +94,7 @@ void ESolver_KS_LCAO<TK, TR>::beforesolver(const int istep)
}

// prepare grid in Gint
LCAO_domain::grid_prepare(this->GridT, this->GG, this->GK, *this->pw_rho, *this->pw_big);
LCAO_domain::grid_prepare(this->GridT, this->GG, this->GK, orb_, *this->pw_rho, *this->pw_big);

// init Hamiltonian
if (this->p_hamilt != nullptr)
Expand All @@ -112,6 +112,7 @@ void ESolver_KS_LCAO<TK, TR>::beforesolver(const int istep)
this->pelec->pot,
this->kv,
two_center_bundle_,
orb_,
DM
#ifdef __EXX
, GlobalC::exx_info.info_ri.real_number ? &this->exd->two_level_step : &this->exc->two_level_step
Expand Down
4 changes: 2 additions & 2 deletions source/module_esolver/lcao_gets.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,7 @@ void ESolver_KS_LCAO<std::complex<double>, double>::get_S(void)
GlobalV::SEARCH_RADIUS,
PARAM.inp.test_atom_input);

this->RA.for_2d(this->pv, PARAM.globalv.gamma_only_local);
this->RA.for_2d(this->pv, PARAM.globalv.gamma_only_local, orb_.cutoffs());

if (this->p_hamilt == nullptr) {
this->p_hamilt = new hamilt::HamiltLCAO<std::complex<double>, double>(
Expand Down Expand Up @@ -97,7 +97,7 @@ void ESolver_KS_LCAO<std::complex<double>, std::complex<double>>::get_S(void)
GlobalV::SEARCH_RADIUS,
PARAM.inp.test_atom_input);

this->RA.for_2d(this->pv, PARAM.globalv.gamma_only_local);
this->RA.for_2d(this->pv, PARAM.globalv.gamma_only_local, orb_.cutoffs());
if (this->p_hamilt == nullptr) {
this->p_hamilt = new hamilt::HamiltLCAO<std::complex<double>,
std::complex<double>>(
Expand Down
6 changes: 3 additions & 3 deletions source/module_esolver/set_matrix_grid.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,7 @@ void ESolver_KS_LCAO<TK, TR>::set_matrix_grid(Record_adj& ra)
std::vector<std::vector<double>> dpsi_u;
std::vector<std::vector<double>> d2psi_u;

Gint_Tools::init_orb(dr_uniform, rcuts, GlobalC::ucell, psi_u, dpsi_u, d2psi_u);
Gint_Tools::init_orb(dr_uniform, rcuts, GlobalC::ucell, orb_, psi_u, dpsi_u, d2psi_u);

this->GridT.set_pbc_grid(this->pw_rho->nx,
this->pw_rho->ny,
Expand Down Expand Up @@ -84,13 +84,13 @@ 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.
ra.for_2d(this->pv, PARAM.globalv.gamma_only_local);
ra.for_2d(this->pv, PARAM.globalv.gamma_only_local, orb_.cutoffs());

if (!PARAM.globalv.gamma_only_local)
{
// need to first calculae lgd.
// using GridT.init.
this->GridT.cal_nnrg(&this->pv);
this->GridT.cal_nnrg(&this->pv, orb_.cutoffs());
}

ModuleBase::timer::tick("ESolver_KS_LCAO", "set_matrix_grid");
Expand Down
4 changes: 3 additions & 1 deletion source/module_hamilt_lcao/hamilt_lcaodft/FORCE.h
Original file line number Diff line number Diff line change
Expand Up @@ -66,6 +66,7 @@ class Force_LCAO
#endif
typename TGint<T>::type& gint,
const TwoCenterBundle& two_center_bundle,
const LCAO_Orbitals& orb,
const Parallel_Orbitals& pv,
const K_Vectors* kv = nullptr,
Record_adj* ra = nullptr);
Expand All @@ -74,14 +75,15 @@ class Force_LCAO
void allocate(const Parallel_Orbitals& pv,
ForceStressArrays& fsr, // mohan add 2024-06-15
const TwoCenterBundle& two_center_bundle,
const LCAO_Orbitals& orb,
const int& nks = 0,
const std::vector<ModuleBase::Vector3<double>>& kvec_d = {});

void finish_ftable(ForceStressArrays& fsr);

void average_force(double* fm);

void test(Parallel_Orbitals& pv, double* mm, const std::string& name);
//void test(Parallel_Orbitals& pv, double* mm, const std::string& name);

//-------------------------------------------------------------
// forces reated to overlap matrix
Expand Down
10 changes: 8 additions & 2 deletions source/module_hamilt_lcao/hamilt_lcaodft/FORCE_STRESS.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -38,6 +38,7 @@ void Force_Stress_LCAO<T>::getForceStress(const bool isforce,
Gint_Gamma& gint_gamma, // mohan add 2024-04-01
Gint_k& gint_k, // mohan add 2024-04-01
const TwoCenterBundle& two_center_bundle,
const LCAO_Orbitals& orb,
ModuleBase::matrix& fcs,
ModuleBase::matrix& scs,
const Structure_Factor& sf,
Expand Down Expand Up @@ -164,6 +165,7 @@ void Force_Stress_LCAO<T>::getForceStress(const bool isforce,
gint_gamma,
gint_k,
two_center_bundle,
orb,
pv,
kv);

Expand Down Expand Up @@ -419,7 +421,7 @@ void Force_Stress_LCAO<T>::getForceStress(const bool isforce,
{
const std::vector<std::vector<double>>& dm_gamma
= dynamic_cast<const elecstate::ElecStateLCAO<double>*>(pelec)->get_DM()->get_DMK_vector();
GlobalC::ld.cal_gdmx(dm_gamma[0], GlobalC::ucell, GlobalC::ORB, GlobalC::GridD, isstress);
GlobalC::ld.cal_gdmx(dm_gamma[0], GlobalC::ucell, orb, GlobalC::GridD, isstress);
}
else
{
Expand All @@ -430,7 +432,7 @@ void Force_Stress_LCAO<T>::getForceStress(const bool isforce,

GlobalC::ld.cal_gdmx_k(dm_k,
GlobalC::ucell,
GlobalC::ORB,
orb,
GlobalC::GridD,
kv.get_nks(),
kv.kvec_d,
Expand Down Expand Up @@ -791,6 +793,7 @@ void Force_Stress_LCAO<double>::integral_part(const bool isGammaOnly,
Gint_Gamma& gint_gamma, // mohan add 2024-04-01
Gint_k& gint_k, // mohan add 2024-04-01
const TwoCenterBundle& two_center_bundle,
const LCAO_Orbitals& orb,
const Parallel_Orbitals& pv,
const K_Vectors& kv)
{
Expand All @@ -814,6 +817,7 @@ void Force_Stress_LCAO<double>::integral_part(const bool isGammaOnly,
#endif
gint_gamma,
two_center_bundle,
orb,
pv);
return;
}
Expand All @@ -839,6 +843,7 @@ void Force_Stress_LCAO<std::complex<double>>::integral_part(const bool isGammaOn
Gint_Gamma& gint_gamma,
Gint_k& gint_k,
const TwoCenterBundle& two_center_bundle,
const LCAO_Orbitals& orb,
const Parallel_Orbitals& pv,
const K_Vectors& kv)
{
Expand All @@ -861,6 +866,7 @@ void Force_Stress_LCAO<std::complex<double>>::integral_part(const bool isGammaOn
#endif
gint_k,
two_center_bundle,
orb,
pv,
&kv,
this->RA);
Expand Down
2 changes: 2 additions & 0 deletions source/module_hamilt_lcao/hamilt_lcaodft/FORCE_STRESS.h
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,7 @@ class Force_Stress_LCAO
Gint_Gamma& gint_gamma, // mohan add 2024-04-01
Gint_k& gint_k, // mohan add 2024-04-01
const TwoCenterBundle& two_center_bundle,
const LCAO_Orbitals& orb,
ModuleBase::matrix& fcs,
ModuleBase::matrix& scs,
const Structure_Factor& sf,
Expand Down Expand Up @@ -90,6 +91,7 @@ class Force_Stress_LCAO
Gint_Gamma& gint_gamma,
Gint_k& gint_k,
const TwoCenterBundle& two_center_bundle,
const LCAO_Orbitals& orb,
const Parallel_Orbitals& pv,
const K_Vectors& kv);

Expand Down
60 changes: 31 additions & 29 deletions source/module_hamilt_lcao/hamilt_lcaodft/FORCE_gamma.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@ template <>
void Force_LCAO<double>::allocate(const Parallel_Orbitals& pv,
ForceStressArrays& fsr, // mohan add 2024-06-15
const TwoCenterBundle& two_center_bundle,
const LCAO_Orbitals& orb,
const int& nks,
const std::vector<ModuleBase::Vector3<double>>& kvec_d)
{
Expand Down Expand Up @@ -72,7 +73,7 @@ void Force_LCAO<double>::allocate(const Parallel_Orbitals& pv,
'S',
cal_deri,
GlobalC::ucell,
GlobalC::ORB,
orb,
pv,
two_center_bundle,
&GlobalC::GridD,
Expand All @@ -95,7 +96,7 @@ void Force_LCAO<double>::allocate(const Parallel_Orbitals& pv,
'T',
cal_deri,
GlobalC::ucell,
GlobalC::ORB,
orb,
pv,
two_center_bundle,
&GlobalC::GridD,
Expand All @@ -106,7 +107,7 @@ void Force_LCAO<double>::allocate(const Parallel_Orbitals& pv,
nullptr,
cal_deri,
GlobalC::ucell,
GlobalC::ORB,
orb,
*(two_center_bundle.overlap_orb_beta),
&GlobalC::GridD);

Expand Down Expand Up @@ -149,28 +150,28 @@ void Force_LCAO<double>::finish_ftable(ForceStressArrays& fsr)
return;
}

template <>
void Force_LCAO<double>::test(Parallel_Orbitals& pv, double* mm, const std::string& name)
{
std::cout << "\n PRINT " << name << std::endl;
std::cout << std::setprecision(6) << std::endl;
for (int i = 0; i < GlobalV::NLOCAL; i++)
{
for (int j = 0; j < GlobalV::NLOCAL; j++)
{
if (std::abs(mm[i * GlobalV::NLOCAL + j]) > 1.0e-5)
{
std::cout << std::setw(12) << mm[i * GlobalV::NLOCAL + j];
}
else
{
std::cout << std::setw(12) << "0";
}
}
std::cout << std::endl;
}
return;
}
//template <>
//void Force_LCAO<double>::test(Parallel_Orbitals& pv, double* mm, const std::string& name)
//{
// std::cout << "\n PRINT " << name << std::endl;
// std::cout << std::setprecision(6) << std::endl;
// for (int i = 0; i < GlobalV::NLOCAL; i++)
// {
// for (int j = 0; j < GlobalV::NLOCAL; j++)
// {
// if (std::abs(mm[i * GlobalV::NLOCAL + j]) > 1.0e-5)
// {
// std::cout << std::setw(12) << mm[i * GlobalV::NLOCAL + j];
// }
// else
// {
// std::cout << std::setw(12) << "0";
// }
// }
// std::cout << std::endl;
// }
// return;
//}

// be called in force_lo.cpp
template <>
Expand All @@ -193,6 +194,7 @@ void Force_LCAO<double>::ftable(const bool isforce,
#endif
TGint<double>::type& gint,
const TwoCenterBundle& two_center_bundle,
const LCAO_Orbitals& orb,
const Parallel_Orbitals& pv,
const K_Vectors* kv,
Record_adj* ra)
Expand All @@ -208,7 +210,7 @@ void Force_LCAO<double>::ftable(const bool isforce,

// allocate DSloc_x, DSloc_y, DSloc_z
// allocate DHloc_fixed_x, DHloc_fixed_y, DHloc_fixed_z
this->allocate(pv, fsr, two_center_bundle);
this->allocate(pv, fsr, two_center_bundle, orb);

// calculate the force related to 'energy density matrix'.
this->cal_fedm(isforce, isstress, fsr, ucell, dm, psi, pv, pelec, foverlap, soverlap);
Expand All @@ -218,7 +220,7 @@ void Force_LCAO<double>::ftable(const bool isforce,
this->cal_fvnl_dbeta(dm,
pv,
ucell,
GlobalC::ORB,
orb,
*(two_center_bundle.overlap_orb_beta),
GlobalC::GridD,
isforce,
Expand All @@ -235,7 +237,7 @@ void Force_LCAO<double>::ftable(const bool isforce,
const std::vector<std::vector<double>>& dm_gamma = dm->get_DMK_vector();

// when deepks_scf is on, the init pdm should be same as the out pdm, so we should not recalculate the pdm
//GlobalC::ld.cal_projected_DM(dm, ucell, GlobalC::ORB, GlobalC::GridD);
//GlobalC::ld.cal_projected_DM(dm, ucell, orb, GlobalC::GridD);

GlobalC::ld.cal_descriptor(ucell.nat);

Expand All @@ -244,7 +246,7 @@ void Force_LCAO<double>::ftable(const bool isforce,
DeePKS_domain::cal_f_delta_gamma(
dm_gamma,
ucell,
GlobalC::ORB,
orb,
GlobalC::GridD,
*this->ParaV,
GlobalC::ld.lmaxd,
Expand Down
Loading

0 comments on commit 8f86cad

Please sign in to comment.