Skip to content

Commit

Permalink
New printing
Browse files Browse the repository at this point in the history
  • Loading branch information
ahurta92 committed Jul 1, 2024
1 parent e59abae commit 1e9571c
Show file tree
Hide file tree
Showing 2 changed files with 35 additions and 22 deletions.
40 changes: 21 additions & 19 deletions src/apps/molresponse/FrequencyResponse.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -570,7 +570,7 @@ Tensor<double> QuadraticResponse::compute_beta_tensor(World &world, const X_spac
truncate(world, dipole_vectors, true);

std::vector<int> indices_A{0, 0, 0, 1, 1, 1, 2, 2, 2, 0};
std::vector<int> indices_BC{0, 1, 2, 0, 1, 2, 0, 1, 2, 0};
std::vector<int> indices_BC{0, 1, 2, 0, 1, 2, 0, 1, 2, 3};

Tensor<double> beta(10);
// for each vector in dipole vectors
Expand Down Expand Up @@ -867,64 +867,64 @@ std::tuple<X_space, X_space, X_space, X_space, X_space, X_space> QuadraticRespon
const X_space &zeta_bc_right, const X_space &zeta_cb_left, const X_space &zeta_cb_right,
const X_space &phi0)
{

if (r_params.print_level() >= 1)
{
molresponse::start_timer(world);
}
auto bphi0c = compute_exchange_term(world, B, phi0, C);
auto zeta_bc = compute_exchange_term(world, zeta_bc_left, zeta_bc_right, phi0);
if (r_params.print_level() >= 1)
{
molresponse::end_timer(world, "k: bphi0c");
molresponse::end_timer(world, "k: zeta_bc");
}
if (r_params.print_level() >= 1)
{
molresponse::start_timer(world);
}
auto cphi0B = compute_exchange_term(world, C, phi0, B);
auto zeta_cb = compute_exchange_term(world, zeta_cb_left, zeta_cb_right, phi0);
if (r_params.print_level() >= 1)
{
molresponse::end_timer(world, "k: cphi0B");
molresponse::end_timer(world, "k: zeta_cb");
}

if (r_params.print_level() >= 1)
{
molresponse::start_timer(world);
}
auto bphi0_phi0 = compute_exchange_term(world, B, phi0, phi0);
auto bphi0c = compute_exchange_term(world, B, phi0, C);
if (r_params.print_level() >= 1)
{
molresponse::end_timer(world, "k: bphi0_phi0");
molresponse::end_timer(world, "k: bphi0c");
}
if (r_params.print_level() >= 1)
{
molresponse::start_timer(world);
}
auto cphi0_phi0 = compute_exchange_term(world, C, phi0, phi0);
auto cphi0B = compute_exchange_term(world, C, phi0, B);
if (r_params.print_level() >= 1)
{
molresponse::end_timer(world, "k: cphi0_phi0");
molresponse::end_timer(world, "k: cphi0B");
}

if (r_params.print_level() >= 1)
{
molresponse::start_timer(world);
}
auto zeta_bc = compute_exchange_term(world, zeta_bc_left, zeta_bc_right, phi0);
auto bphi0_phi0 = compute_exchange_term(world, B, phi0, phi0);
if (r_params.print_level() >= 1)
{
molresponse::end_timer(world, "k: zeta_bc");
molresponse::end_timer(world, "k: bphi0_phi0");
}
if (r_params.print_level() >= 1)
{
molresponse::start_timer(world);
}
auto zeta_cb = compute_exchange_term(world, zeta_cb_left, zeta_cb_right, phi0);
auto cphi0_phi0 = compute_exchange_term(world, C, phi0, phi0);
if (r_params.print_level() >= 1)
{
molresponse::end_timer(world, "k: zeta_cb");
molresponse::end_timer(world, "k: cphi0_phi0");
}


world.gop.fence();

return {bphi0c, cphi0B, bphi0_phi0, cphi0_phi0, zeta_bc, zeta_cb};
Expand All @@ -944,12 +944,12 @@ X_space QuadraticResponse::compute_second_order_perturbation_terms_v2(World &wor
{
molresponse::start_timer(world);
}
auto g_bphi0c = 2.0 * j_bphi0c - k_bphi0c;
auto g_cphi0B = 2.0 * J_cphi0B - K_cphi0B;

auto g_zeta_bc = 2.0 * J_zeta_bc - K_zeta_bc;
auto g_zeta_cb = 2.0 * J_zeta_cb - K_zeta_cb;

auto g_bphi0c = 2.0 * j_bphi0c - k_bphi0c;
auto g_cphi0B = 2.0 * J_cphi0B - K_cphi0B;

auto g1b = 2.0 * J_bphi0_phi0 - K_bphi0_phi0;
auto g1c = 2.0 * J_cphi0_phi0 - K_cphi0_phi0;

Expand All @@ -962,7 +962,7 @@ X_space QuadraticResponse::compute_second_order_perturbation_terms_v2(World &wor
{
molresponse::start_timer(world);
}
auto [VB, VC] = dipole_perturbation(world, phi0, phi0);
auto [VBphi0, VCphi0] = dipole_perturbation(world, phi0, phi0);
auto [vbxc, vcxb] = dipole_perturbation(world, C, B);
if (r_params.print_level() >= 1)
{
Expand All @@ -973,7 +973,7 @@ X_space QuadraticResponse::compute_second_order_perturbation_terms_v2(World &wor
auto FB_C = vbxc + g_bphi0c;
auto FC_B = vcxb + g_cphi0B;

auto [fb_C, fc_B] = compute_first_order_fock_matrix_terms_v2(world, B, C, g1b, g1c, VB, VC, phi0);
auto [fb_C, fc_B] = compute_first_order_fock_matrix_terms_v2(world, B, C, g1b, g1c, VBphi0, VBphi0, phi0);


// now project terms
Expand All @@ -990,6 +990,8 @@ X_space QuadraticResponse::compute_second_order_perturbation_terms_v2(World &wor
g_zeta_cb = -1.0 * oop_apply(g_zeta_cb, apply_projector, false);
FB_C = -1.0 * oop_apply(FB_C, apply_projector, false);
FC_B = -1.0 * oop_apply(FC_B, apply_projector, false);

world.gop.fence();
if (r_params.print_level() >= 1)
{
molresponse::end_timer(world, "projecting terms");
Expand Down
17 changes: 14 additions & 3 deletions src/apps/molresponse/maddft/response_manager.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -1143,13 +1143,24 @@ class ResponseCalcManager

quad_calculation.set_x_data(world, omegas, restarts);
auto beta_abc = quad_calculation.compute_beta_v2(world);
// print the beta values for the frequency combination

// Make a table printing the beta value for each direction
//[XXX,XYY,XZZ,YXX,YYY,YZZ,ZXX,ZYY,ZZZ,XYZ]

vector<std::string> directions{"XXX", "XYY", "XZZ", "YXX", "YYY", "YZZ", "ZXX", "ZYY", "ZZZ", "XYZ"};

if (world.rank() == 0)
{
::print("Beta values for omega_a = ", omega_a, " omega_b = ", omega_b, " omega_c = ", omega_c);
::print(beta_abc);
for (int i = 0; i < 10; i++)
{
if (world.rank() == 0)
{
::print(directions[i], " : ", beta_abc[i]);
}
}
}


nlohmann::ordered_json beta_entry;

std::array<double, 10> beta_vector{};
Expand Down

0 comments on commit 1e9571c

Please sign in to comment.