Skip to content

Commit

Permalink
[hubbard] Use enum type to set up the method for generating the hubba…
Browse files Browse the repository at this point in the history
…rd subspace (#928)

* [hubbard] Use enum type to set up the method used to generate the hubbard subspace
  • Loading branch information
mtaillefumier authored Nov 22, 2023
1 parent 743f52d commit 8166bae
Show file tree
Hide file tree
Showing 22 changed files with 1,377 additions and 855 deletions.
38 changes: 32 additions & 6 deletions src/api/sirius_api.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -851,7 +851,7 @@ sirius_set_parameters(void* const* handler__, int const* lmax_apw__, int const*
}
if (hubbard_full_orthogonalization__ != nullptr) {
if (*hubbard_full_orthogonalization__) {
sim_ctx.cfg().hubbard().full_orthogonalization(true);
sim_ctx.cfg().hubbard().hubbard_subspace_method("full_orthogonalization");
}
}

Expand All @@ -862,12 +862,17 @@ sirius_set_parameters(void* const* handler__, int const* lmax_apw__, int const*
if (hubbard_orbitals__ != nullptr) {
std::string s(hubbard_orbitals__);
std::transform(s.begin(), s.end(), s.begin(), ::tolower);
bool jump = false;
if (s == "ortho-atomic") {
sim_ctx.cfg().hubbard().orthogonalize(true);
sim_ctx.cfg().hubbard().full_orthogonalization(true);
sim_ctx.cfg().hubbard().hubbard_subspace_method("full_orthogonalization");
jump = true;
}
if (s == "norm-atomic") {
sim_ctx.cfg().hubbard().normalize(true);
sim_ctx.cfg().hubbard().hubbard_subspace_method("normalize");
jump = true;
}
if (!jump) {
sim_ctx.cfg().hubbard().hubbard_subspace_method(s);
}
}
if (fft_grid_size__ != nullptr) {
Expand Down Expand Up @@ -2240,8 +2245,29 @@ sirius_set_atom_type_hubbard(void* const* handler__, char const* label__, int co
auto& sim_ctx = get_sim_ctx(handler__);
auto& type = sim_ctx.unit_cell().atom_type(std::string(label__));
type.hubbard_correction(true);
type.add_hubbard_orbital(*n__, *l__, *occ__, *U__, J__[1], J__, *alpha__, *beta__, *J0__,
std::vector<double>(), true);
if (type.file_name().empty()) {
type.add_hubbard_orbital(*n__, *l__, *occ__, *U__, J__[1], J__, *alpha__, *beta__, *J0__,
std::vector<double>(), true);
} else {
// we use a an external file containing the potential
// information which means that we do not have all information
// yet.
//
// let's build the local hubbard section instead
//

json elem;

elem["atom_type"] = label__;
elem["n"] = *n__;
elem["l"] = *l__;
elem["total_initial_occupancy"] = *occ__;
elem["U"] = *U__;
elem["J"] = *J__;
elem["alpha"] = *alpha__;
elem["beta"] = *beta__;
sim_ctx.cfg().hubbard().local().append(elem);
}
},
error_code__);
}
Expand Down
34 changes: 5 additions & 29 deletions src/context/config.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -1587,41 +1587,17 @@ class config_t
: dict_(dict__)
{
}
/// If true, orthogonalization is applied to Hubbard orbitals.
inline auto orthogonalize() const
/// Method to use for generating the hubbard subspace. [none] bare hubbard orbitals, [full_orthogonaliation] use all atomic wave functions to generate the hubbard subspace, [normalize] normalize the original hubbard wave functions, [orthogonalize] orthogonalize the hubbard wave functions
inline auto hubbard_subspace_method() const
{
return dict_.at("/hubbard/orthogonalize"_json_pointer).get<bool>();
return dict_.at("/hubbard/hubbard_subspace_method"_json_pointer).get<std::string>();
}
inline void orthogonalize(bool orthogonalize__)
inline void hubbard_subspace_method(std::string hubbard_subspace_method__)
{
if (dict_.contains("locked")) {
throw std::runtime_error(locked_msg);
}
dict_["/hubbard/orthogonalize"_json_pointer] = orthogonalize__;
}
/// If true, all atomic orbitals from all atoms are used to orthogonalize the hubbard subspace
inline auto full_orthogonalization() const
{
return dict_.at("/hubbard/full_orthogonalization"_json_pointer).get<bool>();
}
inline void full_orthogonalization(bool full_orthogonalization__)
{
if (dict_.contains("locked")) {
throw std::runtime_error(locked_msg);
}
dict_["/hubbard/full_orthogonalization"_json_pointer] = full_orthogonalization__;
}
/// If true, normalization is applied to Hubbard orbitals.
inline auto normalize() const
{
return dict_.at("/hubbard/normalize"_json_pointer).get<bool>();
}
inline void normalize(bool normalize__)
{
if (dict_.contains("locked")) {
throw std::runtime_error(locked_msg);
}
dict_["/hubbard/normalize"_json_pointer] = normalize__;
dict_["/hubbard/hubbard_subspace_method"_json_pointer] = hubbard_subspace_method__;
}
/// If true, simplified version of Hubbard correction is used.
inline auto simplified() const
Expand Down
24 changes: 10 additions & 14 deletions src/context/input_schema.json
Original file line number Diff line number Diff line change
Expand Up @@ -854,20 +854,16 @@
"type": "object",
"title": "Hubbard U correction",
"properties": {
"orthogonalize": {
"type": "boolean",
"default": false,
"title": "If true, orthogonalization is applied to Hubbard orbitals."
},
"full_orthogonalization": {
"type": "boolean",
"default": false,
"title": "If true, all atomic orbitals from all atoms are used to orthogonalize the hubbard subspace"
},
"normalize": {
"type": "boolean",
"default": false,
"title": "If true, normalization is applied to Hubbard orbitals."
"hubbard_subspace_method": {
"type": "string",
"default": "none",
"enum": [
"none",
"full_orthogonalization",
"normalize",
"orthogonalize"
],
"title": "Method to use for generating the hubbard subspace. [none] bare hubbard orbitals, [full_orthogonaliation] use all atomic wave functions to generate the hubbard subspace, [normalize] normalize the original hubbard wave functions, [orthogonalize] orthogonalize the hubbard wave functions"
},
"simplified": {
"type": "boolean",
Expand Down
16 changes: 8 additions & 8 deletions src/hubbard/hubbard_occupancies_derivatives.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,7 @@ build_phi_hub_s_psi_deriv(Simulation_context const& ctx__, int nbnd__, int nawf_
int offset_in_wf = atomic_wf_offset__[ia] + type.indexb_wfs().index_of(rf_index(idxr_wf));
int offset_in_hwf = hubbard_wf_offset__[ia] + type.indexb_hub().index_of(e.idxrf);

if (ctx__.cfg().hubbard().full_orthogonalization()) {
if (ctx__.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
/* compute \sum_{m} d/d r_{alpha} O^{-1/2}_{m,i} <phi_atomic_{m} | S | psi_{jk} > */
la::wrap(la::lib_t::blas)
.gemm('C', 'N', mmax, nbnd__, nawf__, &la::constant<std::complex<double>>::one(),
Expand Down Expand Up @@ -174,7 +174,7 @@ Hubbard::compute_occupancies_derivatives(K_point<double>& kp__, Q_operator<doubl
std::unique_ptr<la::dmatrix<std::complex<double>>> inv_sqrt_O;
std::unique_ptr<la::dmatrix<std::complex<double>>> evec_O;
std::vector<double> eval_O;
if (ctx_.cfg().hubbard().full_orthogonalization()) {
if (ctx_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
ovlp = la::dmatrix<std::complex<double>>(nawf, nawf);
wf::inner(ctx_.spla_context(), mt, wf::spin_range(0), phi_atomic, wf::band_range(0, nawf), phi_atomic_S,
wf::band_range(0, nawf), ovlp, 0, 0);
Expand Down Expand Up @@ -231,7 +231,7 @@ Hubbard::compute_occupancies_derivatives(K_point<double>& kp__, Q_operator<doubl

/* compute < d phi_atomic / d r_{alpha} | S | phi_atomic >
* used to compute derivative of the inverse square root of the overlap matrix */
if (ctx_.cfg().hubbard().full_orthogonalization()) {
if (ctx_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
grad_phi_atomic_s_phi_atomic[x] = la::dmatrix<std::complex<double>>(nawf, nawf);
wf::inner(ctx_.spla_context(), mt, wf::spin_range(0), *s_phi_atomic_tmp, wf::band_range(0, nawf),
phi_atomic, wf::band_range(0, nawf), grad_phi_atomic_s_phi_atomic[x], 0, 0);
Expand All @@ -249,7 +249,7 @@ Hubbard::compute_occupancies_derivatives(K_point<double>& kp__, Q_operator<doubl

/* compute <phi_atomic | S | psi_{ik} > */
std::array<la::dmatrix<std::complex<double>>, 2> phi_atomic_s_psi;
if (ctx_.cfg().hubbard().full_orthogonalization()) {
if (ctx_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
for (int ispn = 0; ispn < ctx_.num_spins(); ispn++) {
phi_atomic_s_psi[ispn] = la::dmatrix<std::complex<double>>(nawf, kp__.num_occupied_bands(ispn));
/* compute < phi_atomic | S | psi_{ik} > for all atoms */
Expand Down Expand Up @@ -312,7 +312,7 @@ Hubbard::compute_occupancies_derivatives(K_point<double>& kp__, Q_operator<doubl
/* from O = <phi | S | phi > we get
* O' = <phi' | S | phi> + <phi | S' |phi> + <phi | S | phi'> */

if (ctx_.cfg().hubbard().full_orthogonalization()) {
if (ctx_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
/* <phi | S' | phi> */
wf::inner(ctx_.spla_context(), mt, wf::spin_range(0), phi_atomic, wf::band_range(0, nawf),
*phi_atomic_tmp, wf::band_range(0, nawf), ovlp, 0, 0);
Expand Down Expand Up @@ -463,7 +463,7 @@ Hubbard::compute_occupancies_stress_derivatives(K_point<double>& kp__, Q_operato
std::unique_ptr<la::dmatrix<std::complex<double>>> inv_sqrt_O;
std::unique_ptr<la::dmatrix<std::complex<double>>> evec_O;
std::vector<double> eval_O;
if (ctx_.cfg().hubbard().full_orthogonalization()) {
if (ctx_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
ovlp = la::dmatrix<std::complex<double>>(nawf, nawf);
wf::inner(ctx_.spla_context(), mt, wf::spin_range(0), phi_atomic, wf::band_range(0, nawf), phi_atomic_S,
wf::band_range(0, nawf), ovlp, 0, 0);
Expand All @@ -477,7 +477,7 @@ Hubbard::compute_occupancies_stress_derivatives(K_point<double>& kp__, Q_operato

/* compute <phi_atomic | S | psi_{ik} > */
std::array<la::dmatrix<std::complex<double>>, 2> phi_atomic_s_psi;
if (ctx_.cfg().hubbard().full_orthogonalization()) {
if (ctx_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
for (int ispn = 0; ispn < ctx_.num_spins(); ispn++) {
phi_atomic_s_psi[ispn] = la::dmatrix<std::complex<double>>(nawf, kp__.num_occupied_bands(ispn));
/* compute < phi_atomic | S | psi_{ik} > for all atoms */
Expand Down Expand Up @@ -508,7 +508,7 @@ Hubbard::compute_occupancies_stress_derivatives(K_point<double>& kp__, Q_operato
sirius::apply_S_operator_strain_deriv(mt, 3 * nu + mu, bp_gen, bp_coeffs, bp_strain_gen, bp_strain_coeffs,
phi_atomic, q_op__, *ds_phi_atomic);

if (ctx_.cfg().hubbard().full_orthogonalization()) {
if (ctx_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
/* compute <phi_atomic | ds/d epsilon | phi_atomic> */
wf::inner(ctx_.spla_context(), mt, wf::spin_range(0), phi_atomic, wf::band_range(0, nawf),
*ds_phi_atomic, wf::band_range(0, nawf), ovlp, 0, 0);
Expand Down
4 changes: 2 additions & 2 deletions src/k_point/k_point.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -238,7 +238,7 @@ K_point<T>::generate_hubbard_orbitals()
sirius::apply_S_operator<T, std::complex<T>>(mem, wf::spin_range(0), wf::band_range(0, nwf), bp_gen, bp_coeffs,
*atomic_wave_functions_, q_op.get(), *atomic_wave_functions_S_);

if (ctx_.cfg().hubbard().full_orthogonalization()) {
if (ctx_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
/* save phi and sphi */

wf_tmp = std::make_unique<wf::Wave_functions<T>>(gkvec_, wf::num_mag_dims(0), wf::num_bands(nwf),
Expand Down Expand Up @@ -319,7 +319,7 @@ K_point<T>::generate_hubbard_orbitals()
}
}
/* restore phi and sphi */
if (ctx_.cfg().hubbard().full_orthogonalization()) {
if (ctx_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {

wf::copy(memory_t::host, *wf_tmp, wf::spin_index(0), wf::band_range(0, nwf), *atomic_wave_functions_,
wf::spin_index(0), wf::band_range(0, nwf));
Expand Down
16 changes: 9 additions & 7 deletions src/unit_cell/atom_type.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -358,13 +358,15 @@ Atom_type::print_info(std::ostream& out__) const
}
out__ << lo_descriptors_hub_[i];
}

bool orthogonalize_ = parameters_.cfg().hubbard().hubbard_subspace_method() == "orthogonalize";
bool full_orthogonalization_ =
parameters_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization";
bool normalize_ = parameters_.cfg().hubbard().hubbard_subspace_method() == "normalize";
out__ << std::endl;
out__ << " orthogonalize : " << boolstr(parameters_.cfg().hubbard().orthogonalize())
<< std::endl
<< " normalize : " << boolstr(parameters_.cfg().hubbard().normalize())
<< std::endl
<< " full_orthogonalization : "
<< boolstr(parameters_.cfg().hubbard().full_orthogonalization()) << std::endl
out__ << " orthogonalize : " << boolstr(orthogonalize_) << std::endl
<< " normalize : " << boolstr(normalize_) << std::endl
<< " full_orthogonalization : " << boolstr(full_orthogonalization_) << std::endl
<< " simplified : " << boolstr(parameters_.cfg().hubbard().simplified())
<< std::endl;
}
Expand Down Expand Up @@ -858,7 +860,7 @@ Atom_type::read_hubbard_input()
}
}

if (parameters_.cfg().hubbard().full_orthogonalization()) {
if (parameters_.cfg().hubbard().hubbard_subspace_method() == "full_orthogonalization") {
this->hubbard_correction_ = true;
if (lo_descriptors_hub_.empty()) {
for (int s = 0; s < (int)ps_atomic_wfs_.size(); s++) {
Expand Down
Loading

0 comments on commit 8166bae

Please sign in to comment.