From c7719f4843c775b8663c1cc7ad7997ed7d45e9f8 Mon Sep 17 00:00:00 2001
From: Gerardo Tauriello <gerardo.tauriello@unibas.ch>
Date: Tue, 19 Sep 2017 18:38:32 +0200
Subject: [PATCH] SCHWED-1740: handle loops of length 0 in MM system creator

---
 loop/doc/mm_system_creation.rst      | 10 +++----
 loop/src/mm_system_creator.cc        |  4 +--
 loop/tests/test_mm_system_creator.cc | 42 ++++++++++++++++++++++------
 3 files changed, 39 insertions(+), 17 deletions(-)

diff --git a/loop/doc/mm_system_creation.rst b/loop/doc/mm_system_creation.rst
index 4cea3afc..7ad0c03a 100644
--- a/loop/doc/mm_system_creation.rst
+++ b/loop/doc/mm_system_creation.rst
@@ -38,7 +38,7 @@ Create MM systems for loops
   :type ff_lookup:  :class:`ForcefieldLookup`
   :param fix_surrounding_hydrogens: If False, the hydrogens of the surrounding
                                     residues can move to improve H-bond building
-                                    (turned off by default as it only has a minor
+                                    (True by default as it only has a minor
                                     impact on the result and a big one on
                                     performance).
   :type fix_surrounding_hydrogens:  :class:`bool`
@@ -123,10 +123,10 @@ Create MM systems for loops
                              :class:`int`
 
     :raises: :exc:`~exceptions.RuntimeError` if loops out of bounds with respect
-             to *res_indices*, all loops have length 0, *loop_start_indices* /
-             *loop_lengths* or *res_indices* / *is_n_ter* / *is_c_ter* have
-             inconsistent lengths, or if any *all_pos[res_indices[i]]* is
-             invalid or has unset heavy atoms.
+             to *res_indices*, *loop_start_indices* / *loop_lengths* or
+             *res_indices* / *is_n_ter* / *is_c_ter* have inconsistent lengths,
+             or if any *all_pos[res_indices[i]]* is invalid or has unset heavy
+             atoms.
 
   .. method:: UpdatePositions(all_pos, res_indices)
 
diff --git a/loop/src/mm_system_creator.cc b/loop/src/mm_system_creator.cc
index 2395048d..fdda4908 100644
--- a/loop/src/mm_system_creator.cc
+++ b/loop/src/mm_system_creator.cc
@@ -319,9 +319,7 @@ void MmSystemCreator::SetupSystem(const AllAtomPositions& all_pos,
   loop_start_indices_ = loop_start_indices;
   loop_lengths_ = loop_lengths;
   IdxHandler::ResolveOverlaps(loop_start_indices_, loop_lengths_, true);
-  if (loop_lengths_.empty()) {
-    throw promod3::Error("No valid loops in MmSystemCreator::SetupSystem!");
-  }
+  // note: ok to not have any loops here (should not bother anyone really...)
   num_loop_residues_ = 0;
   for (uint i_loop = 0; i_loop < loop_lengths_.size(); ++i_loop) {
     num_loop_residues_ += loop_lengths_[i_loop];
diff --git a/loop/tests/test_mm_system_creator.cc b/loop/tests/test_mm_system_creator.cc
index 743cac93..90713693 100644
--- a/loop/tests/test_mm_system_creator.cc
+++ b/loop/tests/test_mm_system_creator.cc
@@ -586,6 +586,32 @@ void CheckLoopMinimize(const AllAtomPositions& all_pos,
       }
     }
   }
+
+  ////////////////////////////////////////////////////////////////////////////
+  // try loop of length 0 (nothing should move!)
+  ////////////////////////////////////////////////////////////////////////////
+  // build system
+  loop_lengths[0] = 0;
+  mm_sim.SetupSystem(all_pos, res_indices, loop_start_indices, loop_lengths,
+                     is_n_ter, is_c_ter, disulfid_bridges);
+  // run it for a bit
+  sim = mm_sim.GetSimulation();
+  pot_ref_loop = sim->GetPotentialEnergy();
+  sim->ApplySD(0.01, 25);
+  if (fix_surrounding_hydrogens) {
+    // nothing should move...
+    BOOST_CHECK_CLOSE(sim->GetPotentialEnergy(), pot_ref_loop, 1);
+  } else {
+    // hydrogens should move and get better...
+    BOOST_CHECK(sim->GetPotentialEnergy() < pot_ref_loop);
+  }
+  // check heavy atoms (none should move)
+  new_pos = all_pos.Copy();
+  mm_sim.ExtractLoopPositions(*new_pos, res_indices);
+  for (uint i = 0; i < all_pos.GetNumAtoms(); ++i) {
+    const Real dist = geom::Distance(all_pos.GetPos(i), new_pos->GetPos(i));
+    BOOST_CHECK_EQUAL(dist, Real(0));
+  }
 }
 
 } // anon ns
@@ -633,7 +659,7 @@ BOOST_AUTO_TEST_CASE(test_mm_system_creator_loop_minimize) {
   ost::mol::EntityHandle test_ent = LoadTestStructure("data/1CRN.pdb");
   AllAtomPositions all_pos(test_ent.GetResidueList());
   // check loop minimization with and w/o inaccurate_pot_e
-  CheckLoopMinimize(all_pos, true, false, 8);
+  CheckLoopMinimize(all_pos, false, false, 8);
   CheckLoopMinimize(all_pos, false, true, 5);
   CheckLoopMinimize(all_pos, true, false, 8);
   CheckLoopMinimize(all_pos, true, true, 5);
@@ -749,10 +775,9 @@ BOOST_AUTO_TEST_CASE(test_mm_system_creator_loop_exceptions) {
   // check/break multi-loop
   std::vector<uint> loop_start_indices;
   std::vector<uint> loop_lengths;
-  BOOST_CHECK_THROW(mm_sim.SetupSystem(all_pos, res_indices, loop_start_indices,
-                                       loop_lengths, is_n_ter, is_c_ter,
-                                       dis_bridges),
-                    promod3::Error);
+  BOOST_CHECK_NO_THROW(mm_sim.SetupSystem(all_pos, res_indices,
+                                          loop_start_indices, loop_lengths,
+                                          is_n_ter, is_c_ter, dis_bridges));
   loop_start_indices.push_back(1);
   loop_lengths.push_back(5);
   BOOST_CHECK_THROW(mm_sim.SetupSystem(all_pos, res_indices, loop_start_indices,
@@ -760,10 +785,9 @@ BOOST_AUTO_TEST_CASE(test_mm_system_creator_loop_exceptions) {
                                        dis_bridges),
                     promod3::Error);
   loop_lengths[0] = 0;
-  BOOST_CHECK_THROW(mm_sim.SetupSystem(all_pos, res_indices, loop_start_indices,
-                                       loop_lengths, is_n_ter, is_c_ter,
-                                       dis_bridges),
-                    promod3::Error);
+  BOOST_CHECK_NO_THROW(mm_sim.SetupSystem(all_pos, res_indices,
+                                          loop_start_indices, loop_lengths,
+                                          is_n_ter, is_c_ter, dis_bridges));
   loop_start_indices.push_back(1);
   BOOST_CHECK_THROW(mm_sim.SetupSystem(all_pos, res_indices, loop_start_indices,
                                        loop_lengths, is_n_ter, is_c_ter,
-- 
GitLab