Skip to content

Commit

Permalink
Roll back changes made to semicanonicalize
Browse files Browse the repository at this point in the history
  • Loading branch information
@brianz98
brianz98 committed Nov 15, 2024
1 parent bb89c72 commit d878f0d
Show file tree
Hide file tree
Showing 5 changed files with 11 additions and 119 deletions.
4 changes: 0 additions & 4 deletions forte/api/orbital_api.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -88,12 +88,8 @@ void export_SemiCanonical(py::module& m) {
"Semicanonicalize the orbitals and transform the integrals and reference")
.def("Ua", &SemiCanonical::Ua, "Return the alpha rotation matrix")
.def("Ub", &SemiCanonical::Ub, "Return the alpha rotation matrix")
.def("Ua_c", &SemiCanonical::Ua_c, "Return the alpha rotation matrix in the core space")
.def("Ub_c", &SemiCanonical::Ub_c, "Return the beta rotation matrix in the core space")
.def("Ua_t", &SemiCanonical::Ua_t, "Return the alpha rotation matrix in the active space")
.def("Ub_t", &SemiCanonical::Ub_t, "Return the beta rotation matrix in the active space")
.def("Ua_v", &SemiCanonical::Ua_v, "Return the alpha rotation matrix in the virtual space")
.def("Ub_v", &SemiCanonical::Ub_v, "Return the beta rotation matrix in the virtual space")
.def("fix_orbital_success", &SemiCanonical::fix_orbital_success,
"Return if the orbital ordering and phases are fixed successfully");
}
Expand Down
1 change: 0 additions & 1 deletion forte/mrdsrg-spin-integrated/master_mrdsrg.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -707,7 +707,6 @@ void MASTER_DSRG::deGNO_ints_full(const std::string& name, double& H0, BlockedTe
L1h.block("CC").iterate(
[&](const std::vector<size_t>& i, double& value) { value = i[0] == i[1] ? 1.0 : 0.0; });


// scalar from H1
double scalar1 = 0.0;
scalar1 -= H1["ji"] * L1h["ij"];
Expand Down
74 changes: 0 additions & 74 deletions forte/orbital-helpers/semi_canonicalize.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -64,26 +64,16 @@ SemiCanonical::SemiCanonical(std::shared_ptr<MOSpaceInfo> mo_space_info,
void SemiCanonical::startup() {
nirrep_ = mo_space_info_->nirrep();
nmopi_ = mo_space_info_->dimension("ALL");
ncore_ = mo_space_info_->size("CORE");
nact_ = mo_space_info_->size("ACTIVE");
nvirt_ = mo_space_info_->size("VIRTUAL");

// Prepare orbital rotation matrix, which transforms all MOs
Ua_ = std::make_shared<psi::Matrix>("Ua", nmopi_, nmopi_);
Ub_ = std::make_shared<psi::Matrix>("Ub", nmopi_, nmopi_);

// Prepare orbital rotation matrix, which transforms only core MOs
Ua_c_ = ambit::Tensor::build(ambit::CoreTensor, "Ua", {ncore_, ncore_});
Ub_c_ = ambit::Tensor::build(ambit::CoreTensor, "Ub", {ncore_, ncore_});

// Prepare orbital rotation matrix, which transforms only active MOs
Ua_t_ = ambit::Tensor::build(ambit::CoreTensor, "Ua", {nact_, nact_});
Ub_t_ = ambit::Tensor::build(ambit::CoreTensor, "Ub", {nact_, nact_});

// Prepare orbital rotation matrix, which transforms only virtual MOs
Ua_v_ = ambit::Tensor::build(ambit::CoreTensor, "Ua", {nvirt_, nvirt_});
Ub_v_ = ambit::Tensor::build(ambit::CoreTensor, "Ub", {nvirt_, nvirt_});

// Initialize U to identity
set_U_to_identity();

Expand Down Expand Up @@ -122,23 +112,11 @@ void SemiCanonical::set_U_to_identity() {
Ua_->identity();
Ub_->identity();

Ua_c_.iterate(
[&](const std::vector<size_t>& i, double& value) { value = (i[0] == i[1]) ? 1.0 : 0.0; });

Ub_c_.iterate(
[&](const std::vector<size_t>& i, double& value) { value = (i[0] == i[1]) ? 1.0 : 0.0; });

Ua_t_.iterate(
[&](const std::vector<size_t>& i, double& value) { value = (i[0] == i[1]) ? 1.0 : 0.0; });

Ub_t_.iterate(
[&](const std::vector<size_t>& i, double& value) { value = (i[0] == i[1]) ? 1.0 : 0.0; });

Ua_v_.iterate(
[&](const std::vector<size_t>& i, double& value) { value = (i[0] == i[1]) ? 1.0 : 0.0; });

Ub_v_.iterate(
[&](const std::vector<size_t>& i, double& value) { value = (i[0] == i[1]) ? 1.0 : 0.0; });
}

void SemiCanonical::semicanonicalize(std::shared_ptr<RDMs> rdms, bool build_fock,
Expand Down Expand Up @@ -299,64 +277,12 @@ void SemiCanonical::build_transformation_matrices(const bool& semi) {
if (!inactive_mix_)
fix_orbital_success_ = ints_->fix_orbital_phases(Ua_, true);

// fill in Ua_c_
fill_Ucore(Ua_, Ua_c_);
// fill in Ua_t_
fill_Uactv(Ua_, Ua_t_);
// fill in Ua_v_
fill_Uvirt(Ua_, Ua_v_);

// pass to Ub
Ub_->copy(Ua_);
Ub_c_.copy(Ua_c_);
Ub_t_.copy(Ua_t_);
Ub_v_.copy(Ua_v_);
}

void SemiCanonical::fill_Ucore(const std::shared_ptr<psi::Matrix>& U, ambit::Tensor& Ut) {
auto core_names = mo_space_info_->composite_spaces_def().at("CORE");
auto& Ut_data = Ut.data();

for (const std::string& name : core_names) {
auto size = mo_space_info_->size(name);
if (size == 0)
continue;

auto pos = mo_space_info_->pos_in_space(name, "CORE");
auto relative_mos = mo_space_info_->relative_mo(name);
for (size_t p = 0; p < size; ++p) {
const auto& [hp, np] = relative_mos[p];
for (size_t q = 0; q < size; ++q) {
const auto& [hq, nq] = relative_mos[q];
if (hp != hq)
continue;
Ut_data[pos[p] * ncore_ + pos[q]] = U->get(hp, np, nq);
}
}
}
}

void SemiCanonical::fill_Uvirt(const std::shared_ptr<psi::Matrix>& U, ambit::Tensor& Ut) {
auto virtual_names = mo_space_info_->composite_spaces_def().at("VIRTUAL");
auto& Ut_data = Ut.data();

for (const std::string& name : virtual_names) {
auto size = mo_space_info_->size(name);
if (size == 0)
continue;

auto pos = mo_space_info_->pos_in_space(name, "VIRTUAL");
auto relative_mos = mo_space_info_->relative_mo(name);
for (size_t p = 0; p < size; ++p) {
const auto& [hp, np] = relative_mos[p];
for (size_t q = 0; q < size; ++q) {
const auto& [hq, nq] = relative_mos[q];
if (hp != hq)
continue;
Ut_data[pos[p] * nvirt_ + pos[q]] = U->get(hp, np, nq);
}
}
}
}

void SemiCanonical::fill_Uactv(const std::shared_ptr<psi::Matrix>& U, ambit::Tensor& Ut) {
Expand Down
32 changes: 0 additions & 32 deletions forte/orbital-helpers/semi_canonicalize.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -109,22 +109,12 @@ class SemiCanonical {
/// @return the beta rotation matrix
std::shared_ptr<psi::Matrix> Ub() { return Ub_; }

/// @return the alpha rotation matrix in the core space
ambit::Tensor Ua_c() const { return Ua_c_.clone(); }
/// @return the beta rotation matrix in the core space
ambit::Tensor Ub_c() const { return Ub_c_.clone(); }

/// @return the alpha rotation matrix in the active space
ambit::Tensor Ua_t() const { return Ua_t_.clone(); }

/// @return the beta rotation matrix in the active space
ambit::Tensor Ub_t() const { return Ub_t_.clone(); }

/// @return the alpha rotation matrix in the virtual space
ambit::Tensor Ua_v() const { return Ua_v_.clone(); }
/// @return the beta rotation matrix in the virtual space
ambit::Tensor Ub_v() const { return Ub_v_.clone(); }

/// @return if the orbital ordering and phases are fixed successfully
bool fix_orbital_success() const { return fix_orbital_success_; }

Expand Down Expand Up @@ -158,15 +148,9 @@ class SemiCanonical {
/// Offset of GAS orbitals within ACTIVE
std::map<std::string, psi::Dimension> actv_offsets_;

/// Number of core MOs
size_t ncore_;

/// Number of active MOs
size_t nact_;

/// Number of virtual MOs
size_t nvirt_;

/// Number of irreps
size_t nirrep_;

Expand All @@ -175,21 +159,11 @@ class SemiCanonical {
/// Unitary matrix for beta orbital rotation
std::shared_ptr<psi::Matrix> Ub_;

/// Unitary matrix for alpha orbital rotation in the core space
ambit::Tensor Ua_c_;
/// Unitary matrix for beta orbital rotation in the core space
ambit::Tensor Ub_c_;

/// Unitary matrix for alpha orbital rotation in the active space
ambit::Tensor Ua_t_;
/// Unitary matrix for beta orbital rotation in the active space
ambit::Tensor Ub_t_;

/// Unitary matrix for alpha orbital rotation in the virtual space
ambit::Tensor Ua_v_;
/// Unitary matrix for beta orbital rotation in the virtual space
ambit::Tensor Ub_v_;

/// Set Ua_, Ub_, Ua_t_, and Ub_t_ to identity
void set_U_to_identity();

Expand All @@ -211,15 +185,9 @@ class SemiCanonical {
/// Builds unitary matrices used to diagonalize diagonal blocks of Fock
void build_transformation_matrices(const bool& semi);

/// Fill ambit::Tensor Ua_c_ (Ub_c_) using std::shared_ptr<psi::Matrix> Ua_ (Ub_)
void fill_Ucore(const std::shared_ptr<psi::Matrix>& U, ambit::Tensor& Ut);

/// Fill ambit::Tensor Ua_t_ (Ub_t_) using std::shared_ptr<psi::Matrix> Ua_ (Ub_)
void fill_Uactv(const std::shared_ptr<psi::Matrix>& U, ambit::Tensor& Ut);

/// Fill ambit::Tensor Ua_v_ (Ub_v_) using std::shared_ptr<psi::Matrix> Ua_ (Ub_)
void fill_Uvirt(const std::shared_ptr<psi::Matrix>& U, ambit::Tensor& Ut);

/// Successfully fix the orbital ordering and phases
bool fix_orbital_success_;
};
Expand Down
19 changes: 11 additions & 8 deletions forte/proc/dsrg.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -249,7 +249,7 @@ def compute_energy(self):
if not self.Heff_implemented:
self.relax_maxiter = 0

if self.options.get_bool("FULL_HBAR") and self.solver_type in ["MRDSRG","SA-MRDSRG","SA_MRDSRG"]:
if self.options.get_bool("FULL_HBAR") and self.solver_type in ["MRDSRG","SA-MRDSRG","SA_MRDSRG"] and self.relax_maxiter == 0:
psi4.core.print_out("\n =>** Saving Full Hbar (unrelaxed) **<=\n")
Hbar1, Hbar2 = self.dsrg_solver.compute_Heff_full()
np.savez("save_Hbar", **Hbar1, **Hbar2)
Expand All @@ -267,7 +267,7 @@ def compute_energy(self):

del Mbar0, Mbar1, Mbar2

if self.options.get_bool("FULL_HBAR_DEGNO") and self.solver_type in ["SA-MRDSRG","SA_MRDSRG","MRDSRG"]:
if self.options.get_bool("FULL_HBAR_DEGNO") and self.solver_type in ["MRDSRG","SA-MRDSRG","SA_MRDSRG"] and self.relax_maxiter == 0:
psi4.core.print_out("\n =>** Saving Full Hbar in de-normal-ordered basis (unrelaxed) **<=\n")
Hbar1, Hbar2 = self.dsrg_solver.compute_Heff_full_degno()
np.savez("save_Hbar_degno", **Hbar1, **Hbar2)
Expand All @@ -283,16 +283,14 @@ def compute_energy(self):
# However, the ForteIntegrals object and the dipole integrals always refer to the current semi-canonical basis.
# so to compute the dipole moment correctly, we need to make the RDMs and orbital basis consistent

if self.options.get_bool("FULL_HBAR") and n == self.relax_maxiter - 1:
psi4.core.print_out("\n =>** Saving Full Hbar (relax) **<=\n")
Heff = self.dsrg_solver.compute_Heff_full()
Heff_dict = forte.Heff_dict(Heff)
np.savez("save_Hbar", **Heff_dict)
if self.options.get_bool("FULL_HBAR") and self.solver_type in ["MRDSRG","SA-MRDSRG","SA_MRDSRG"] and n == self.relax_maxiter - 1:
psi4.core.print_out("\n =>** Saving Full Hbar (relaxed) **<=\n")
Hbar1, Hbar2 = self.dsrg_solver.compute_Heff_full()
np.savez("save_Hbar", **Hbar1, **Hbar2)

if self.solver_type not in ["MRDSRG_SO", "MRDSRG-SO"]:
psi4.core.print_out("\n =>** Getting dipole integral **<=\n")
Mbar0 = self.dsrg_solver.compute_Mbar0_full()
print(Mbar0)
np.save("Mbar0", Mbar0)
Mbar1 = self.dsrg_solver.compute_Mbar1_full()
Mbar2 = self.dsrg_solver.compute_Mbar2_full()
Expand All @@ -302,6 +300,11 @@ def compute_energy(self):
np.savez(f"Mbar2_{i}", **Mbar2[i])

del Mbar0, Mbar1, Mbar2

if self.options.get_bool("FULL_HBAR_DEGNO") and self.solver_type in ["MRDSRG","SA-MRDSRG","SA_MRDSRG"] and n == self.relax_maxiter - 1:
psi4.core.print_out("\n =>** Saving Full Hbar in de-normal-ordered basis (relaxed) **<=\n")
Hbar1, Hbar2 = self.dsrg_solver.compute_Heff_full_degno()
np.savez("save_Hbar_degno", **Hbar1, **Hbar2)

ints_dressed = self.dsrg_solver.compute_Heff_actv()
if self.fno_pt2_Heff_shift is not None:
Expand Down

0 comments on commit d878f0d

Please sign in to comment.