Skip to content

Commit

Permalink
final: upload all and terminate
Browse files Browse the repository at this point in the history
  • Loading branch information
@maki49
maki49 committed Sep 25, 2023
1 parent 684c770 commit 669873a
Show file tree
Hide file tree
Showing 32 changed files with 1,366 additions and 286 deletions.
8 changes: 3 additions & 5 deletions source/module_cell/module_symmetry/symmetry_basic.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1068,11 +1068,9 @@ void Symmetry_Basic::atom_ordering_new(double *posi, const int natom, int *subin
ModuleBase::heapsort(nxequal, weighted_func, tmpindex.data());
this->order_atoms(&posi[i * 3], nxequal, tmpindex.data());
//rearange subindex using tmpindex
for (int j = 0; j < nxequal; ++j)
{
tmpindex[j] = subindex[i + tmpindex[j]];
subindex[i + j] = tmpindex[j];
}
for (int j = 0; j < nxequal; ++j)tmpindex[j] = subindex[i + tmpindex[j]];
for (int j = 0; j < nxequal; ++j)subindex[i + j] = tmpindex[j];

}
i=ix_right;
}
Expand Down
19 changes: 12 additions & 7 deletions source/module_esolver/esolver_ks_lcao.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -485,10 +485,13 @@ void ESolver_KS_LCAO::eachiterinit(const int istep, const int iter)

#ifdef __EXX
// calculate exact-exchange
if (GlobalC::exx_info.info_global.cal_exx && !GlobalC::exx_info.info_global.separate_loop// && exx_lri_complex->two_level_step
&& !GlobalV::GAMMA_ONLY_LOCAL && ModuleSymmetry::Symmetry::symm_flag == 1)
this->exc->restore_dm(this->UHM, this->LOC, kv, GlobalC::ucell, *psi, symm, pelec->wg);
if (GlobalC::exx_info.info_ri.real_number)
this->exd->exx_eachiterinit(this->LOC, *(this->p_chgmix), iter);
this->exd->exx_eachiterinit(this->LOC, *(this->p_chgmix), this->symm, iter);
else
this->exc->exx_eachiterinit(this->LOC, *(this->p_chgmix), iter);
this->exc->exx_eachiterinit(this->LOC, *(this->p_chgmix), this->symm, iter);
#endif

if (GlobalV::dft_plus_u)
Expand Down Expand Up @@ -551,9 +554,9 @@ void ESolver_KS_LCAO::hamilt2density(int istep, int iter, double ethr)

#ifdef __EXX
if (GlobalC::exx_info.info_ri.real_number)
this->exd->exx_hamilt2density(*this->pelec, *this->LOWF.ParaV);
this->exd->exx_hamilt2density(*this->pelec, *this->LOWF.ParaV, symm);
else
this->exc->exx_hamilt2density(*this->pelec, *this->LOWF.ParaV);
this->exc->exx_hamilt2density(*this->pelec, *this->LOWF.ParaV, symm);
#endif

// if DFT+U calculation is needed, this function will calculate
Expand Down Expand Up @@ -822,7 +825,7 @@ void ESolver_KS_LCAO::afterscf(const int istep)
// rpa_interface.rpa_exx_lcao().info.files_abfs = GlobalV::rpa_orbitals;
// rpa_interface.out_for_RPA(*(this->LOWF.ParaV), *(this->psi), this->LOC, this->pelec);
RPA_LRI<double> rpa_lri_double(GlobalC::exx_info.info_ri);
rpa_lri_double.cal_postSCF_exx(this->LOC, MPI_COMM_WORLD, kv, *this->LOWF.ParaV);
rpa_lri_double.cal_postSCF_exx(this->LOC, MPI_COMM_WORLD, kv, *this->LOWF.ParaV, symm);
rpa_lri_double.init(MPI_COMM_WORLD, kv);
rpa_lri_double.out_for_RPA(*(this->LOWF.ParaV), *(this->psi), this->pelec);
}
Expand All @@ -847,10 +850,12 @@ void ESolver_KS_LCAO::afterscf(const int istep)
bool ESolver_KS_LCAO::do_after_converge(int& iter)
{
#ifdef __EXX
if (!GlobalV::GAMMA_ONLY_LOCAL && ModuleSymmetry::Symmetry::symm_flag == 1)
this->exc->restore_dm(this->UHM, this->LOC, kv, GlobalC::ucell, *psi, symm, pelec->wg);
if (GlobalC::exx_info.info_ri.real_number)
return this->exd->exx_after_converge(*this->p_hamilt, this->LM, this->LOC, kv, iter);
return this->exd->exx_after_converge(*this->p_hamilt, this->LM, this->LOC, kv, symm, iter);
else
return this->exc->exx_after_converge(*this->p_hamilt, this->LM, this->LOC, kv, iter);
return this->exc->exx_after_converge(*this->p_hamilt, this->LM, this->LOC, kv, symm, iter);
#endif // __EXX
return true;
}
Expand Down
4 changes: 2 additions & 2 deletions source/module_hamilt_lcao/hamilt_lcaodft/FORCE_k.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -132,7 +132,7 @@ void Force_LCAO_k::ftable_k(const bool isforce,
GlobalC::ld.check_v_delta_k(pv->nnr);
for (int ik = 0; ik < kv.nks; ik++)
{
uhm.LM->folding_fixedH(ik, kv.kvec_d);
uhm.LM->folding_fixedH(kv.kvec_d[ik]);
}
GlobalC::ld.cal_e_delta_band_k(loc.dm_k,
kv.nks);
Expand Down Expand Up @@ -265,7 +265,7 @@ void Force_LCAO_k::allocate_k(const Parallel_Orbitals& pv,
for (int ik = 0; ik < nks; ik++)
{
this->UHM->genH.LM->zeros_HSk('S');
this->UHM->genH.LM->folding_fixedH(ik, kvec_d, 1);
this->UHM->genH.LM->folding_fixedH(kvec_d[ik], 1);
bool bit = false; // LiuXh, 2017-03-21
ModuleIO::saving_HS(0,
this->UHM->genH.LM->Hloc2.data(),
Expand Down
7 changes: 4 additions & 3 deletions source/module_hamilt_lcao/hamilt_lcaodft/LCAO_matrix.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -24,9 +24,7 @@ class LCAO_Matrix

// folding the fixed Hamiltonian (T+Vnl) if
// k-point algorithm is used.
void folding_fixedH(const int &ik,
const std::vector<ModuleBase::Vector3<double>>& kvec_d,
bool cal_syns = false);
void folding_fixedH(const ModuleBase::Vector3<double>& kvec_d, bool cal_syns = false);

Parallel_Orbitals *ParaV;

Expand Down Expand Up @@ -124,6 +122,9 @@ class LCAO_Matrix
// Records the R direct coordinates of HR and SR output, This variable will be filled with data when HR and SR files are output.
std::set<Abfs::Vector3_Order<int>> output_R_coor;

/// for restoring DM(R) when symmetry == 1 and multi - k in EXX calculation.
/// $S^{-1}'(k)S(gk)$ for each kstar. size: [nks_ibz][nsymm][nbasis*nbasis(local)]
std::vector<std::vector<std::vector<std::complex<double>>>> invSkrot_Sgk;

//========================================
// FORCE
Expand Down
9 changes: 4 additions & 5 deletions source/module_hamilt_lcao/hamilt_lcaodft/LCAO_nnr.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -318,9 +318,8 @@ int Grid_Technique::binary_search_find_R2_offset(int val, int iat) const

// be called in LCAO_Hamilt::calculate_Hk.
void LCAO_Matrix::folding_fixedH(
const int &ik,
const std::vector<ModuleBase::Vector3<double>>& kvec_d,
bool cal_syns)
const ModuleBase::Vector3<double>& kvec_d,
bool cal_syns)
{
ModuleBase::TITLE("LCAO_nnr","folding_fixedH");
ModuleBase::timer::tick("LCAO_nnr", "folding_fixedH");
Expand Down Expand Up @@ -432,8 +431,8 @@ void LCAO_Matrix::folding_fixedH(
// dR is the index of box in Crystal coordinates
//------------------------------------------------
ModuleBase::Vector3<double> dR(adjs.box[ad].x, adjs.box[ad].y, adjs.box[ad].z);
const double arg = ( kvec_d[ik] * dR ) * ModuleBase::TWO_PI;
//const double arg = ( kvec_d[ik] * GlobalC::GridD.getBox(ad) ) * ModuleBase::TWO_PI;
const double arg = (kvec_d * dR) * ModuleBase::TWO_PI;
//const double arg = ( kvec_d * GlobalC::GridD.getBox(ad) ) * ModuleBase::TWO_PI;
double sinp, cosp;
ModuleBase::libm::sincos(arg, &sinp, &cosp);
const std::complex<double> kphase = std::complex <double> ( cosp, sinp );
Expand Down
33 changes: 31 additions & 2 deletions source/module_hamilt_lcao/hamilt_lcaodft/operator_lcao/operator_lcao.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -76,13 +76,42 @@ void OperatorLCAO<std::complex<double>>::folding_fixed(const int ik, const std::
//-----------------------------------------
this->LM->zeros_HSk('S');
this->LM->zeros_HSk('T');
this->LM->folding_fixedH(ik, kvec_d);
this->LM->folding_fixedH(kvec_d[ik]);

//------------------------------------------
// Add T(k)+Vnl(k)+Vlocal(k)
// (Hloc2 += Hloc_fixed2), (std::complex matrix)
//------------------------------------------
this->LM->update_Hloc2(ik);
this->LM->update_Hloc2(ik);
// test: output HSk
GlobalV::ofs_running << "ik=" << ik << std::endl;
// GlobalV::ofs_running << "Hloc2 " << std::endl;
// for (int ir = 0;ir < GlobalV::NLOCAL;++ir)
// {
// for (int ic = 0;ic < GlobalV::NLOCAL;++ic)
// {
// GlobalV::ofs_running << this->LM->Hloc2[ir * GlobalV::NLOCAL + ic] << " ";
// }
// GlobalV::ofs_running << std::endl;
// }
// GlobalV::ofs_running << "Hloc_fixed2 " << std::endl;
// for (int ir = 0;ir < GlobalV::NLOCAL;++ir)
// {
// for (int ic = 0;ic < GlobalV::NLOCAL;++ic)
// {
// GlobalV::ofs_running << this->LM->Hloc_fixed2[ir * GlobalV::NLOCAL + ic] << " ";
// }
// GlobalV::ofs_running << std::endl;
// }
GlobalV::ofs_running << "Sloc2 with kvec_d=" << kvec_d[ik].x << " " << kvec_d[ik].y << " " << kvec_d[ik].z << std::endl;
for (int ir = 0;ir < GlobalV::NLOCAL;++ir)
{
for (int ic = 0;ic < GlobalV::NLOCAL;++ic)
{
GlobalV::ofs_running << this->LM->Sloc2[ir * GlobalV::NLOCAL + ic] << " ";
}
GlobalV::ofs_running << std::endl;
}
ModuleBase::timer::tick("OperatorLCAO", "folding_fixed");
}

Expand Down
1 change: 1 addition & 0 deletions source/module_hamilt_lcao/module_gint/CMakeLists.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@ list(APPEND objects
gint_rho.cpp
gint_tau.cpp
gint_vl.cpp
gint_Srot.cpp
gint_k_env.cpp
gint_k_sparse.cpp
gint_k_sparse1.cpp
Expand Down
52 changes: 42 additions & 10 deletions source/module_hamilt_lcao/module_gint/gint.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,16 +23,20 @@ void Gint::cal_gint(Gint_inout *inout)
if(inout->job==Gint_Tools::job_type::rho) ModuleBase::TITLE("Gint_interface","cal_gint_rho");
if(inout->job==Gint_Tools::job_type::tau) ModuleBase::TITLE("Gint_interface","cal_gint_tau");
if(inout->job==Gint_Tools::job_type::force) ModuleBase::TITLE("Gint_interface","cal_gint_force");
if(inout->job==Gint_Tools::job_type::force_meta) ModuleBase::TITLE("Gint_interface","cal_gint_force_meta");
if (inout->job == Gint_Tools::job_type::force_meta) ModuleBase::TITLE("Gint_interface", "cal_gint_force_meta");
if (inout->job == Gint_Tools::job_type::srot) ModuleBase::TITLE("Gint_interface", "cal_gint_srot");
if (inout->job == Gint_Tools::job_type::srotk) ModuleBase::TITLE("Gint_interface", "cal_gint_srotk");

if(inout->job==Gint_Tools::job_type::vlocal) ModuleBase::timer::tick("Gint_interface", "cal_gint_vlocal");
if(inout->job==Gint_Tools::job_type::vlocal_meta) ModuleBase::timer::tick("Gint_interface","cal_gint_vlocal_meta");
if(inout->job==Gint_Tools::job_type::rho) ModuleBase::timer::tick("Gint_interface","cal_gint_rho");
if(inout->job==Gint_Tools::job_type::tau) ModuleBase::timer::tick("Gint_interface","cal_gint_tau");
if(inout->job==Gint_Tools::job_type::force) ModuleBase::timer::tick("Gint_interface","cal_gint_force");
if(inout->job==Gint_Tools::job_type::force_meta) ModuleBase::timer::tick("Gint_interface","cal_gint_force_meta");
if (inout->job == Gint_Tools::job_type::srot) ModuleBase::timer::tick("Gint_interface", "cal_gint_srot");
if (inout->job == Gint_Tools::job_type::srotk) ModuleBase::timer::tick("Gint_interface", "cal_gint_srotk");

const int max_size = this->gridt->max_atom;
const int max_size = this->gridt->max_atom;
const int LD_pool = max_size*GlobalC::ucell.nwmax;
const int lgd = this->gridt->lgd;
const int nnrg = this->gridt->nnrg;
Expand All @@ -55,7 +59,8 @@ void Gint::cal_gint(Gint_inout *inout)
// it's a uniform grid to save orbital values, so the delta_r is a constant.
const double delta_r = GlobalC::ORB.dr_uniform;

if((inout->job==Gint_Tools::job_type::vlocal || inout->job==Gint_Tools::job_type::vlocal_meta) && !GlobalV::GAMMA_ONLY_LOCAL)
if ((inout->job == Gint_Tools::job_type::vlocal || inout->job == Gint_Tools::job_type::vlocal_meta || inout->job == Gint_Tools::job_type::srot)
&& !GlobalV::GAMMA_ONLY_LOCAL)
{
if(!pvpR_alloc_flag)
{
Expand All @@ -81,8 +86,9 @@ void Gint::cal_gint(Gint_inout *inout)
//perpare auxiliary arrays to store thread-specific values
#ifdef _OPENMP
double* pvpR_thread;
if(inout->job==Gint_Tools::job_type::vlocal || inout->job==Gint_Tools::job_type::vlocal_meta)
{
if (inout->job == Gint_Tools::job_type::vlocal || inout->job == Gint_Tools::job_type::vlocal_meta
|| inout->job == Gint_Tools::job_type::srot)
{
if(!GlobalV::GAMMA_ONLY_LOCAL)
{
pvpR_thread = new double[nnrg];
Expand Down Expand Up @@ -247,8 +253,21 @@ void Gint::cal_gint(Gint_inout *inout)
#endif
delete[] vldr3;
delete[] vkdr3;
}
} // int grid_index
}
else if (inout->job == Gint_Tools::job_type::srot)
{
#ifdef _OPENMP
this->gint_kernel_Srot(inout->ginv, na_grid, grid_index, delta_r, dv, LD_pool, pvpR_thread);
#else
this->gint_kernel_Srot(inout->ginv, na_grid, grid_index, delta_r, dv, LD_pool, this->pvpR_reduced[0]);
#endif
}
else if (inout->job == Gint_Tools::job_type::srotk)
{
this->gint_kernel_Srot(na_grid, grid_index, delta_r, dv, LD_pool, inout);
}

} // int grid_index

#ifdef _OPENMP
if(inout->job==Gint_Tools::job_type::vlocal || inout->job==Gint_Tools::job_type::vlocal_meta)
Expand Down Expand Up @@ -284,7 +303,17 @@ void Gint::cal_gint(Gint_inout *inout)
{
inout->svl_dphi[0]+=svl_dphi_thread;
}
}
}

if (inout->job == Gint_Tools::job_type::srot)
{
#pragma omp critical(gint_k)
for (int innrg = 0; innrg < nnrg; innrg++)
{
pvpR_reduced[0][innrg] += pvpR_thread[innrg];
}
delete[] pvpR_thread;
}
#endif
} // end of #pragma omp parallel

Expand All @@ -300,8 +329,11 @@ void Gint::cal_gint(Gint_inout *inout)
if(inout->job==Gint_Tools::job_type::rho) ModuleBase::timer::tick("Gint_interface","cal_gint_rho");
if(inout->job==Gint_Tools::job_type::tau) ModuleBase::timer::tick("Gint_interface","cal_gint_tau");
if(inout->job==Gint_Tools::job_type::force) ModuleBase::timer::tick("Gint_interface","cal_gint_force");
if(inout->job==Gint_Tools::job_type::force_meta) ModuleBase::timer::tick("Gint_interface","cal_gint_force_meta");
return;
if (inout->job == Gint_Tools::job_type::force_meta) ModuleBase::timer::tick("Gint_interface", "cal_gint_force_meta");
if (inout->job == Gint_Tools::job_type::srot) ModuleBase::timer::tick("Gint_interface", "cal_gint_srot");
if (inout->job == Gint_Tools::job_type::srotk) ModuleBase::timer::tick("Gint_interface", "cal_gint_srotk");

return;
}

void Gint::prep_grid(
Expand Down
22 changes: 22 additions & 0 deletions source/module_hamilt_lcao/module_gint/gint.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -196,6 +196,28 @@ class Gint
double** dpsiz_dm,
double* rho);

//------------------------------------------------------
// in gint_Skrot.cpp
//------------------------------------------------------
///calculate < phi_0|g^{-1}phi_R > for exx-symmetry
void gint_kernel_Srot(
const ModuleBase::Matrix3& ginv,
const int na_grid,
const int grid_index,
const double delta_r,
const double dv,
const int LD_pool,
double* pvpR_reduced);

/// calculate $<\phi_{k_1}|\phi_{k_2}> = \int{dr} \sum_{R_1}\phi_{R_1}(r)exp(-ik_1*R_1) \sum_{R_2}\phi_{R_2}(r)exp(ik_2*R_2)$
void gint_kernel_Srot(
const int na_grid,
const int grid_index,
const double delta_r,
const double dv,
const int LD_pool,
Gint_inout* inout);

// dimension: [GlobalC::LNNR.nnrg]
// save the < phi_0i | V | phi_Rj > in sparse H matrix.
bool pvpR_alloc_flag = false;
Expand Down
Loading

0 comments on commit 669873a

Please sign in to comment.