From f33e2b993602aee12deae872b234fd1a061200d8 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Conrad=20H=C3=BCbler?= Date: Sun, 26 May 2024 15:18:02 +0200 Subject: [PATCH] stuff MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: Conrad Hübler --- src/capabilities/confscan.cpp | 20 +++++-------- src/capabilities/confscan.h | 6 ++-- src/capabilities/rmsd.cpp | 56 ++++++++++++++++++++--------------- src/capabilities/rmsd.h | 6 ++-- 4 files changed, 44 insertions(+), 44 deletions(-) diff --git a/src/capabilities/confscan.cpp b/src/capabilities/confscan.cpp index 0a9f311..4eea280 100644 --- a/src/capabilities/confscan.cpp +++ b/src/capabilities/confscan.cpp @@ -110,6 +110,7 @@ int ConfScanThread::execute() m_reorder_rule = m_driver->ReorderRules(); } + std::cout << m_rmsd << " "; m_driver->clear(); return 0; } @@ -225,7 +226,9 @@ void ConfScan::LoadControlJson() m_useorders = Json2KeyWord(m_defaults, "UseOrders"); m_MaxHTopoDiff = Json2KeyWord(m_defaults, "MaxHTopoDiff"); m_threads = m_defaults["threads"].get(); - m_RMSDmethod = Json2KeyWord(m_defaults, "RMSDMethod"); + m_RMSDmethod = Json2KeyWord(m_defaults, "method"); + fmt::print(fg(fmt::color::green) | fmt::emphasis::bold, "\nPermutation of atomic indices performed according to {0} \n\n", m_RMSDmethod); + if (m_RMSDmethod == "molalign") { fmt::print(fg(fmt::color::green) | fmt::emphasis::bold, "\nPlease cite the follow research report!\nJ. Chem. Inf. Model. 2023, 63, 4, 1157–1165 - DOI: 10.1021/acs.jcim.2c01187\n\n"); m_domolalign = -1; @@ -853,7 +856,7 @@ void ConfScan::CheckOnly(double sLE, double sLI, double sLH) std::string laststring; m_maxmol = m_ordered_list.size(); - json rmsd = RMSDJson; + json rmsd = m_controller; rmsd["silent"] = true; rmsd["check"] = CheckConnections(); rmsd["heavy"] = m_heavy; @@ -1031,19 +1034,10 @@ void ConfScan::Reorder(double dLE, double dLI, double dLH, bool reuse_only, bool m_rejected = 0, m_accepted = 0, m_reordered = 0, m_reordered_worked = 0, m_reordered_reused = 0; - json rmsd = RMSDJson; + json rmsd = m_controller; rmsd["silent"] = true; rmsd["reorder"] = true; - rmsd["check"] = CheckConnections(); - rmsd["heavy"] = m_heavy; - rmsd["method"] = m_RMSDmethod; - rmsd["element"] = m_rmsd_element_templates; - rmsd["update-rotation"] = m_update_rotation; - rmsd["damping"] = m_damping; - rmsd["nomunkres"] = m_nomunkres; - rmsd["molalignbin"] = m_molalign; - rmsd["molaligntol"] = m_molaligntol; - rmsd["cycles"] = m_cycles; + rmsd["threads"] = 1; std::vector cached; if (reset) cached = m_all_structures; diff --git a/src/capabilities/confscan.h b/src/capabilities/confscan.h index abc46ec..f221366 100644 --- a/src/capabilities/confscan.h +++ b/src/capabilities/confscan.h @@ -66,12 +66,12 @@ static const json ConfScanJson = { { "skip_orders", false }, { "MaxParam", -1 }, { "UseOrders", -1 }, - { "RMSDMethod", "hybrid" }, + { "method", "hybrid" }, { "MaxHTopoDiff", -1 }, { "threads", 1 }, { "RMSDElement", 7 }, { "accepted", "" }, - { "method", "" }, + { "method", "free" }, { "lastdE", -1 }, { "fewerFile", false }, { "skipinit", false }, @@ -108,7 +108,7 @@ class ConfScanThread : public CxxThread { public: ConfScanThread(const std::vector>& reorder_rules, double rmsd_threshold, int MaxHTopoDiff, bool reuse_only, const json& config) { - m_driver = new RMSDDriver(config, true); + m_driver = new RMSDDriver(config, false); m_config = config; m_reuse_only = reuse_only; m_reorder_rules = reorder_rules; diff --git a/src/capabilities/rmsd.cpp b/src/capabilities/rmsd.cpp index 9fff47c..4ade6fa 100644 --- a/src/capabilities/rmsd.cpp +++ b/src/capabilities/rmsd.cpp @@ -173,11 +173,9 @@ void RMSDDriver::LoadControlJson() m_topo = Json2KeyWord(m_defaults, "topo"); m_write = Json2KeyWord(m_defaults, "write"); m_noreorder = Json2KeyWord(m_defaults, "noreorder"); - m_moi = Json2KeyWord(m_defaults, "moi"); m_update_rotation = Json2KeyWord(m_defaults, "update-rotation"); m_split = Json2KeyWord(m_defaults, "split"); m_nofree = Json2KeyWord(m_defaults, "nofree"); - m_dmix = Json2KeyWord(m_defaults, "dmix"); m_maxtrial = Json2KeyWord(m_defaults, "maxtrial"); #pragma message("these hacks to overcome the json stuff are not nice, TODO!") try { @@ -222,7 +220,9 @@ void RMSDDriver::LoadControlJson() m_limit = Json2KeyWord(m_defaults, "limit"); } else m_method = 1; - + if (!m_silent) { + fmt::print(fg(fmt::color::green) | fmt::emphasis::bold, "\nPermutation of atomic indices performed according to {0} \n\n", m_method); + } m_costmatrix = Json2KeyWord(m_defaults, "costmatrix"); std::string order = Json2KeyWord(m_defaults, "order"); int cycles = Json2KeyWord(m_defaults, "cycles"); @@ -277,10 +277,6 @@ void RMSDDriver::start() if (m_reference.Atoms() != m_target.Atoms() || m_force_reorder) { if (!m_noreorder) { ReorderMolecule(); - m_reorder_rules = m_results.begin()->second; - m_target_reordered = ApplyOrder(m_reorder_rules, m_target); - m_target = m_target_reordered; - m_rmsd = m_results.begin()->first; // std::cout << m_rmsd << std::endl; rmsd_calculated = true; } @@ -315,6 +311,8 @@ void RMSDDriver::start() std::cout << std::endl << "RMSD calculation took " << timer.Elapsed() << " msecs." << std::endl; std::cout << "Difference in Topological Hydrogen Bond Matrix is " << m_htopo_diff << std::endl; + std::cout << std::endl + << std::endl; CheckTopology(); } if (m_swap) { @@ -504,7 +502,8 @@ void RMSDDriver::ReorderIncremental() } // m_reorder_rules = m_results.begin()->second; if (m_stored_rules.size() == 0) { - std::cout << "No new solution found, sorry" << std::endl; + if (!m_silent) + std::cout << "No new solution found, sorry" << std::endl; for (int i = 0; i < m_reference.AtomCount(); ++i) m_reorder_rules.push_back(i); } else { @@ -778,7 +777,8 @@ void RMSDDriver::ReorderMolecule() if (m_method == 4) blob.first = 1; - InsertRotation(blob); + if (m_method != 6) + InsertRotation(blob); if (m_method == 1) ReorderIncremental(); else if (m_method == 2) @@ -791,15 +791,29 @@ void RMSDDriver::ReorderMolecule() TemplateFree(); else if (m_method == 6) { if (!MolAlignLib()) { + InsertRotation(blob); + FinaliseTemplate(); + + m_target_reordered = ApplyOrder(m_reorder_rules, m_target); + m_target_aligned = m_target; + m_reorder_rules = m_results.begin()->second; + m_target_reordered = ApplyOrder(m_reorder_rules, m_target); + m_target = m_target_reordered; + m_rmsd = m_results.begin()->first; return; } } else if (m_method == 7) DistanceTemplate(); + if (m_method != 6) { + FinaliseTemplate(); - FinaliseTemplate(); - - m_target_reordered = ApplyOrder(m_reorder_rules, m_target); - m_target_aligned = m_target; + m_target_reordered = ApplyOrder(m_reorder_rules, m_target); + m_target_aligned = m_target; + m_reorder_rules = m_results.begin()->second; + m_target_reordered = ApplyOrder(m_reorder_rules, m_target); + m_target = m_target_reordered; + m_rmsd = m_results.begin()->first; + } } void RMSDDriver::FinaliseTemplate() @@ -816,7 +830,6 @@ void RMSDDriver::FinaliseTemplate() if (eq_counter > m_maxtrial || incr_counter > m_maxtrial || iter > m_maxtrial) break; if (prev_cost != -1 && prev_cost < permutation.first / 10) { - std::cout << permutation.first << std::endl; break; } prev_cost = permutation.first; @@ -1409,18 +1422,12 @@ std::vector RMSDDriver::SolveCostMatrix(Matrix& distance) } int* order = new int[dim]; assign(dim, table, order); - // new_order.clear(); - // new_order.resize(dim); for (int i = 0; i < dim; ++i) new_order[i] = order[i]; delete (table); delete (order); } else { auto result = MunkressAssign(distance); - - // new_order.clear(); - // new_order.resize(dim); - for (int i = 0; i < result.cols(); ++i) { for (int j = 0; j < result.rows(); ++j) { if (result(i, j) == 1) { @@ -1510,7 +1517,7 @@ bool RMSDDriver::MolAlignLib() m_target.writeXYZFile("molalign_tar.xyz"); FILE* FileOpen; - std::string command = m_molalign + " molaign_ref.xyz " + " molalign_tar.xyz -sort -fast -tol " + std::to_string(m_molaligntol) + " 2>&1"; + std::string command = m_molalign + " molaign_ref.xyz " + " molalign_tar.xyz -remap -fast -tol " + std::to_string(m_molaligntol) + " 2>&1"; if (!m_silent) std::cout << command << std::endl; FileOpen = popen(command.c_str(), "r"); @@ -1523,7 +1530,7 @@ bool RMSDDriver::MolAlignLib() // ok = std::string(line).find("RMSD") != std::string::npos; rndm = std::string(line).find("random") != std::string::npos; error = std::string(line).find("Error") != std::string::npos; - // printf("%s", line); + printf("%s", line); } pclose(FileOpen); @@ -1532,11 +1539,12 @@ bool RMSDDriver::MolAlignLib() fmt::print(fg(fmt::color::green) | fmt::emphasis::bold, "\nPlease cite the follow research report!\nJ. Chem. Inf. Model. 2023, 63, 4, 1157–1165 - DOI: 10.1021/acs.jcim.2c01187\n\n"); } FileIterator file("aligned.xyz", true); - file.Next(); + m_reference_centered = file.Next(); m_target_reordered = file.Next(); m_target_aligned = m_target_reordered; m_target = m_target_reordered; - // std::filesystem::remove("aligned.xyz"); + m_rmsd = CalculateRMSD(); + std::filesystem::remove("aligned.xyz"); } else { if (!rndm && !error) { fmt::print(fg(fmt::color::salmon) | fmt::emphasis::bold, "Molalign was not found. Consider getting it from\nhttps://github.com/qcuaeh/molalignlib\nEither adding the location of the binary to the path for executables or append\n-molalignbin /yourpath/molalign\nto your argument list!\n"); diff --git a/src/capabilities/rmsd.h b/src/capabilities/rmsd.h index 54d5e6b..b4d2389 100644 --- a/src/capabilities/rmsd.h +++ b/src/capabilities/rmsd.h @@ -90,12 +90,10 @@ static const json RMSDJson = { { "check", false }, { "topo", 0 }, { "write", 0 }, - { "moi", false }, { "update-rotation", false }, { "damping", 0.8 }, { "split", false }, { "nomunkres", false }, - { "dmix", -1 }, { "molalignbin", "molalign" }, { "molaligntol", 10 }, { "cycles", -1 }, @@ -304,14 +302,14 @@ class RMSDDriver : public CurcumaMethod { Molecule m_reference, m_target, m_target_original, m_reference_aligned, m_target_aligned, m_target_reordered, m_reorder_reference, m_reorder_target, m_reference_centered, m_target_centered; Geometry m_reorder_reference_geometry; - bool m_force_reorder = false, m_protons = true, m_print_intermediate = false, m_silent = false, m_moi = false; + bool m_force_reorder = false, m_protons = true, m_print_intermediate = false, m_silent = false; std::vector> m_intermediate_results; std::map> m_results, m_intermediate_cost_matrices; std::vector m_last_rmsd; std::vector m_reorder_rules; std::vector> m_stored_rules, m_intermedia_rules; std::vector m_tmp_rmsd; - double m_rmsd = 0, m_rmsd_raw = 0, m_scaling = 1.5, m_intermedia_storage = 1, m_threshold = 99, m_damping = 0.8, m_dmix = -1; + double m_rmsd = 0, m_rmsd_raw = 0, m_scaling = 1.5, m_intermedia_storage = 1, m_threshold = 99, m_damping = 0.8; bool m_check = false; bool m_check_connections = false, m_postprocess = true, m_noreorder = false, m_swap = false, m_dynamic_center = false; bool m_update_rotation = false, m_split = false, m_nofree = false;