diff --git a/doc/tests/scripts/loop_mm_sys_creation.py b/doc/tests/scripts/loop_mm_sys_creation.py
index a3f29d39973a814c4b52b708f6dd15d8701da9b9..642c4dac3b7403b5eb7a389608c6b28eeaa5336c 100644
--- a/doc/tests/scripts/loop_mm_sys_creation.py
+++ b/doc/tests/scripts/loop_mm_sys_creation.py
@@ -12,24 +12,20 @@ num_residues = len(res_list)
# get all atom positions for full protein
all_atoms = loop.AllAtomPositions(res_list)
-# get res. indices for loop followed by rest
-# -> here we add everything to rest but in practice this could be
-# just a subset of residues relevant to the loop of interest
-loop_start = 10
-loop_length = 6
-res_indices = list(range(loop_start, loop_start + loop_length)) \
- + list(range(0, loop_start)) \
- + list(range(loop_start + loop_length, num_residues))
+# here full structure in res_indices but in practice this could
+# be just a subset of residues relevant to the loop of interest
+res_indices = range(num_residues)
+# define two loops (indices into res_indices list)
+loop_start_indices = [10, 20]
+loop_lengths = [6, 4]
# define which of the res_indices is terminal
-is_n_ter = [False] * num_residues
-is_c_ter = [False] * num_residues
-is_n_ter[res_indices.index(0)] = True
-is_c_ter[res_indices.index(num_residues-1)] = True
+is_n_ter = [True] + [False] * (num_residues - 1)
+is_c_ter = [False] * (num_residues - 1) + [True]
# get disulfid bridges
disulfid_bridges = mm_sys.GetDisulfidBridges(all_atoms, res_indices)
# setup MM system
-mm_sys.SetupSystem(all_atoms, res_indices, loop_length,
- is_n_ter, is_c_ter, disulfid_bridges)
+mm_sys.SetupSystem(all_atoms, res_indices, loop_start_indices,
+ loop_lengths, is_n_ter, is_c_ter, disulfid_bridges)
# run simulation
sim = mm_sys.GetSimulation()
diff --git a/loop/doc/mm_system_creation.rst b/loop/doc/mm_system_creation.rst
index fbfd17382a7e7aa116915ec33ae878d0059eeb55..4cea3afc29cd59e1239391d46f13a14161cdc944 100644
--- a/loop/doc/mm_system_creation.rst
+++ b/loop/doc/mm_system_creation.rst
@@ -11,6 +11,15 @@ The example below showcases the creation and use of an MM system:
.. literalinclude:: ../../../tests/doc/scripts/loop_mm_sys_creation.py
+.. note::
+
+ To simplify use when sidechains may be missing and the region of interest for
+ the MM system has to be determined, it might be better to use the simpler
+ :class:`promod3.modelling.SidechainReconstructor` and
+ :class:`promod3.modelling.AllAtomRelaxer` classes. Even if all the sidechains
+ are available, those classes will be helpful since the overhead to check
+ sidechains without reconstructing them is minimal.
+
Create MM systems for loops
--------------------------------------------------------------------------------
@@ -61,6 +70,8 @@ Create MM systems for loops
.. method:: SetupSystem(all_pos, res_indices, loop_length, is_n_ter, \
is_c_ter, disulfid_bridges)
+ SetupSystem(all_pos, res_indices, loop_start_indices, \
+ loop_lengths, is_n_ter, is_c_ter, disulfid_bridges)
Setup a new system for the loop / surrounding defined by *all_pos* and
*res_indices*. Positions are taken from *all_pos[res_indices[i]]* and each
@@ -72,11 +83,16 @@ Create MM systems for loops
approximation. Similarly, we do not generate OXT atoms for residues followed
by a gap unless they are specifically labelled as C-terminal.
- Residues *res_indices[i]* for *i* < *loop_length* are loop, rest is
- surrounding. Residues *res_indices[0]* / *res_indices[loop_length-1]* are
- assumed to be N- / C-stems of loop. The stems of the loop are kept fixed (N,
+ Each loop is defined by an entry in *loop_start_indices* and *loop_lengths*.
+ For loop *i_loop*, *res_indices[loop_start_indices[i_loop]]* is the N-stem
+ and *res_indices[loop_start_indices[i_loop] + loop_lengths[i_loop] - 1]* is
+ the C-stem of the loop. The loop indices are expected to be contiguous in
+ *res_indices* between the stems. The stems of the loop are kept fixed (N,
CA, CB for N-stem / CA, CB, C, O for C-stem) unless they are terminal
- residues.
+ residues. Overlapping loops are merged and 0-length loops are removed.
+
+ If *loop_length* is given, a single loop with start index 0 and the given
+ loop length is defined.
If a system was setup previously, it is overwritten completely.
@@ -88,8 +104,12 @@ Create MM systems for loops
:type all_pos: :class:`AllAtomPositions`
:param res_indices: Residue indices into *all_pos*.
:type res_indices: :class:`list` of :class:`int`
- :param loop_length: Length of loop.
+ :param loop_length: Length of single loop (incl. stems).
:type loop_length: :class:`int`
+ :param loop_start_indices: Start indices of loops.
+ :type loop_start_indices: :class:`list` of :class:`int`
+ :param loop_lengths: Lengths of loops (incl. stems).
+ :type loop_lengths: :class:`list` of :class:`int`
:param is_n_ter: For each *res_indices[i]*, *is_n_ter[i]* defines whether
that residue is N-terminal.
:type is_n_ter: :class:`list` of :class:`bool`
@@ -102,10 +122,11 @@ Create MM systems for loops
:type disulfid_bridges: :class:`list` of :class:`tuple` with two
:class:`int`
- :raises: :exc:`~exceptions.RuntimeError` if *loop_length* is longer than
- *res_indices*, *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.
+ :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.
.. method:: UpdatePositions(all_pos, res_indices)
@@ -130,14 +151,21 @@ Create MM systems for loops
:meth:`SetupSystem`. This function then extracts these positions back into a
compatible output storage.
- :param loop_pos: Storage for loop positions. Residues *loop_pos[i]* and
- *all_pos[res_indices[i]]* (of :meth:`SetupSystem`) must
- have matching amino acid types for *i* < *loop_length*.
+ If *loop_pos* is passed, only the loops are stored. The loop residues are
+ stored contiguously and the loops are stored in the order given by
+ :meth:`GetLoopStartIndices` / :meth:`GetLoopLengths`. The *loop_pos* storage
+ must have at least :meth:`GetNumLoopResidues` residues and must have
+ matching amino acid types with respect to the loop residues passed as
+ *all_pos* of :meth:`SetupSystem`.
+
+ If *out_pos* and *res_indices* are passed, residues must have matching amino
+ acid types for *out_pos[res_indices[i]]* and *all_pos[res_indices[i]]* of
+ :meth:`SetupSystem`. Only loop residues with *i* in the ranges defined by
+ :meth:`GetLoopStartIndices` / :meth:`GetLoopLengths` are touched.
+
+ :param loop_pos: Reduced storage only covering loop positions.
:type loop_pos: :class:`AllAtomPositions`
- :param out_pos: Storage for loop positions. Residues
- *out_pos[res_indices[i]]* and *all_pos[res_indices[i]]* (of
- :meth:`SetupSystem`) must have matching amino acid types for
- *i* < *loop_length*.
+ :param out_pos: Storage for loop positions linked to *res_indices*.
:type out_pos: :class:`AllAtomPositions`
:param res_indices: Residue indices into *out_pos*.
:type res_indices: :class:`list` of :class:`int`
@@ -151,6 +179,26 @@ Create MM systems for loops
MM simulations.
:rtype: :class:`ost.mol.mm.Simulation`
+ .. method:: GetNumResidues
+
+ :return: Number of residues of current simulation setup.
+ :rtype: :class:`int`
+
+ .. method:: GetNumLoopResidues
+
+ :return: Number of loop residues (incl. stems) of current simulation setup.
+ :rtype: :class:`int`
+
+ .. method:: GetLoopStartIndices
+
+ :return: Start indices of loops (see :meth:`SetupSystem`).
+ :rtype: :class:`list` of :class:`int`
+
+ .. method:: GetLoopLengths
+
+ :return: Lengths of loops (see :meth:`SetupSystem`).
+ :rtype: :class:`list` of :class:`int`
+
.. method:: GetForcefieldAminoAcids
:return: Forcefield-specific amino acid type for each residue of last
diff --git a/loop/src/mm_system_creator.cc b/loop/src/mm_system_creator.cc
index c071411790bb073461f89b5ea3a92f260ab145b2..5ffe7a97b1b6bcfbfc2f190eb31d0bbe9170822f 100644
--- a/loop/src/mm_system_creator.cc
+++ b/loop/src/mm_system_creator.cc
@@ -197,7 +197,8 @@ MmSystemCreator::MmSystemCreator(ForcefieldLookupPtr ff_lookup,
bool kill_electrostatics,
Real nonbonded_cutoff,
bool inaccurate_pot_energy)
- : aa_lookup_(AminoAcidLookup::GetInstance()) {
+ : num_loop_residues_(0)
+ , aa_lookup_(AminoAcidLookup::GetInstance()) {
// store settings and check for CPU support
ff_lookup_ = ff_lookup;
fix_surrounding_hydrogens_ = fix_surrounding_hydrogens;