<html>
  <head>
    <meta content="text/html; charset=ISO-8859-1"
      http-equiv="Content-Type">
  </head>
  <body text="#000000" bgcolor="#FFFFFF">
    <div class="moz-cite-prefix">Hi,<br>
      <br>
      I fully understand your sentiments.<br>
      But I think this solution, although fine for personal use, it to
      ad-hoc to include in Gromacs. Currently the couple-moltype
      decouples a whole molecule, which is a very clear operation
      (although already somewhat complicated to understand due to all
      the different coupling mechanisms). I certainly don't want to
      replace this by a index group option. Adding an index group option
      is also messy.<br>
      <br>
      I looked at the domain decomposition code and there is nothing
      molecule specific to how the bonded interactions are assigned. So
      on the mdrun level there is nothing that prevents the working of
      inter-molecular interactions. The only thing that needs to be done
      is adding something like an [ intermolecular-interactions ]
      section to the topology file processing and the topology data
      structure. Having the current processing code go through this is a
      matter of adding a few lines. I estimate that implementing this
      would take me an hour or so, if we stick to simple, absolute
      system (not molecule) atom indices in this section. I can do this
      within a week, promise :) The only disadvantage is that the ligand
      will not automatically end up in the same periodic image as the
      protein.<br>
      <br>
      I will discuss this with our local developers to get some
      feedback.<br>
      <br>
      Cheers,<br>
      <br>
      Berk<br>
      <br>
      On 08/19/2013 08:31 PM, David Mobley wrote:<br>
    </div>
    <blockquote
cite="mid:CAPCBi1_16k692Wc7FR=M5QshX2iFcYFrv8UC8XN_fvJtzm0emA@mail.gmail.com"
      type="cite">
      <meta http-equiv="Content-Type" content="text/html;
        charset=ISO-8859-1">
      <div dir="ltr">Hi, Berk,<br>
        <div class="gmail_quote">
          <div dir="ltr">
            <div>
              <div>
                <div><br>
                </div>
                I agree that this solution would be the "proper" one and
                would be nice. However, it's been a long time coming
                (this is an issue I've been raising for approximately 5
                years, and there's no progress; normally no one even
                answers). Currently, for the restraints I need to use, I
                *have* to merge the ligand into my protein topology,
                inconvenient or not. Doing this also means there are
                some features I simply can't use (couple-moltype for
                example). <br>
                <br>
                It strikes me that Floris's solution is a good interim
                one, in that it at least solves most of the typical use
                cases we're facing, and has the advantage of being
                *already implemented*. <br>
                <br>
              </div>
              Is there any way we could get this into the main GROMACS,
              at least until someone gets time to implementing the
              "proper" solution (which could be years from now)?<br>
              <br>
            </div>
            Thanks!<span class="HOEnZb"><font color="#888888"><br>
                David<br>
                <br>
              </font></span></div>
          <div class="gmail_extra">
            <div>
              <div class="h5"><br>
                <br>
                <div class="gmail_quote">
                  On Mon, Aug 19, 2013 at 1:10 AM, Berk Hess <span
                    dir="ltr">&lt;<a moz-do-not-send="true"
                      href="mailto:hess@kth.se" target="_blank">hess@kth.se</a>&gt;</span>
                  wrote:<br>
                  <blockquote class="gmail_quote" style="margin:0 0 0
                    .8ex;border-left:1px #ccc solid;padding-left:1ex">Hi,<br>
                    <br>
                    I fully understood this point, maybe my reply wasn't
                    well structured.<br>
                    I was trying to argue that the proper solution would
                    be to implement inter-molecular restraints, instead
                    of introducing an option which could me misused.
                    Also merging a ligand topology into your protein
                    topology is very inconvenient. This problem with the
                    proper solution is, of course, that someone will
                    need to implement inter-molecular
                    restraints/potentials.<br>
                    <br>
                    One issue used to be that we corrected molecules for
                    PBC before calculating bonded interactions. But in
                    4.6 this is usually not faster and not used any
                    more. That makes it easier to treat intra- and
                    inter-molecular interactions the same way at the
                    mdrun level.<br>
                    <br>
                    Cheers,<br>
                    <br>
                    Berk
                    <div>
                      <div><br>
                        <br>
                        On 08/19/2013 09:55 AM, Floris Buelens wrote:<br>
                        <blockquote class="gmail_quote" style="margin:0
                          0 0 .8ex;border-left:1px #ccc
                          solid;padding-left:1ex">
                          Hi Berk,<br>
                          <br>
                          I think some clarification is needed. What
                          we're talking about is only relevant in the
                          context of non-bonded interactions, typically
                          between a small molecule and a protein. To
                          access the binding affinity of a ligand
                          through alchemical methods, it's useful to
                          switch off only the interactions of the ligand
                          with its environment while maintaining
                          intra-ligand potentials - this is what's
                          provided by the 'couple-moltype' code path.<br>
                          <br>
                          The limitation that we're trying to circumvent
                          comes from the fact that as you scale down
                          interactions with the environment, the ligand
                          is no longer held in the binding site. To
                          counteract this, it's necessary to perturb in
                          restraining potentials as the intramolecular
                          nonbonded interactions are perturbed out. The
                          most practical method makes use of a single
                          distance restraint, two angle and three
                          dihedral restraints, whose effect can later be
                          accounted for analytically.<br>
                          <br>
                          The current code only allows decoupling of a
                          whole logical 'molecule' as specified in the
                          topology file. This precludes the use of
                          regular (fully perturbation-aware) bonded
                          interactions. As far as I'm aware, the pull
                          code doesn't provide what's required.
                          Inter-molecular bonded interactions would be a
                          great general solution but presumably won't
                          show up any time soon.<br>
                          <br>
                          The workaround here is instead to represent
                          two physical molecules (e.g. protein and
                          ligand) as a single logical molecule (in the
                          topology file), so regular bonded potentials
                          can be applied. Current Gromacs doesn't allow
                          the decoupling ('couple-moltype') code to be
                          used in this scenario. My modification allows
                          you to target the ligand by its residue name
                          and get this functionality back.<br>
                          <br>
                          Decoupling a single residue of a multi-residue
                          chain is indeed probably not correct in this
                          framework (I did call it 'weird' in my last
                          message :-) ). An extra check would block
                          users from trying this.<br>
                          <br>
                          I hope that clarifies the problem we're trying
                          to solve. I agree this will be useful a
                          relatively small number of users, but on the
                          other hand it's a very unobtrusive,
                          self-contained modification which can maintain
                          full backwards compatibility.<br>
                          <br>
                          thanks,<br>
                          <br>
                          Floris<br>
                          <br>
                          <br>
                          ----- Original Message -----<br>
                          From: Berk Hess &lt;<a moz-do-not-send="true"
                            href="mailto:hess@kth.se" target="_blank">hess@kth.se</a>&gt;<br>
                          To: Floris Buelens &lt;<a
                            moz-do-not-send="true"
                            href="mailto:floris_buelens@yahoo.com"
                            target="_blank">floris_buelens@yahoo.com</a>&gt;;
                          Discussion list for GROMACS development &lt;<a
                            moz-do-not-send="true"
                            href="mailto:gmx-developers@gromacs.org"
                            target="_blank">gmx-developers@gromacs.org</a>&gt;<br>
                          Cc:<br>
                          Sent: Monday, 19 August 2013, 9:03<br>
                          Subject: Re: [gmx-developers] free energy
                          calculations with restraints<br>
                          <br>
                          Hi,<br>
                          <br>
                          There is the rotational pull code, but that
                          might not provide the exact<br>
                          functionality you want.<br>
                          Already for some time we have been discussing
                          inter-molecular<br>
                          interactions, by specifying molecule type,
                          molecule index and atom<br>
                          index, but there are no concrete plans for
                          implementing this yet.<br>
                          <br>
                          An option for decoupling part of a molecule
                          indeed sound useful. But in<br>
                          practice you always need to replace that part
                          by something else, at<br>
                          least a hydrogen, and modify some potentials
                          of connecting atoms, so I<br>
                          don't know how generally useful such an option
                          is.<br>
                          <br>
                          Cheers,<br>
                          <br>
                          Berk<br>
                          <br>
                          On 08/19/2013 07:56 AM, Floris Buelens wrote:<br>
                          <blockquote class="gmail_quote"
                            style="margin:0 0 0 .8ex;border-left:1px
                            #ccc solid;padding-left:1ex">
                            I suppose that would make sense. The changes
                            are fairly minor (more or less just the
                            function convert_moltype_couple and the four
                            functions it calls) and mainly consist of a
                            bit of extra juggling with nonbonded
                            exclusions during preprocessing. The current
                            functionality (decouple a whole molecule)
                            could be handled as a special case of the
                            new code. A check for bonds to the rest of
                            the structure could stop people from trying
                            weird stuff like decoupling residues from a
                            chain.<br>
                            <br>
                            <br>
                            <br>
                            <br>
                            ________________________________<br>
                            From: David Mobley &lt;<a
                              moz-do-not-send="true"
                              href="mailto:dmobley@gmail.com"
                              target="_blank">dmobley@gmail.com</a>&gt;<br>
                            To: Floris Buelens &lt;<a
                              moz-do-not-send="true"
                              href="mailto:floris_buelens@yahoo.com"
                              target="_blank">floris_buelens@yahoo.com</a>&gt;<br>
                            Cc: Discussion list for GROMACS development
                            &lt;<a moz-do-not-send="true"
                              href="mailto:gmx-developers@gromacs.org"
                              target="_blank">gmx-developers@gromacs.org</a>&gt;<br>
                            Sent: Friday, 16 August 2013, 18:30<br>
                            Subject: Re: [gmx-developers] free energy
                            calculations with restraints<br>
                            <br>
                            <br>
                            <br>
                            This does sound useful, though it would be
                            more useful if this could be implemented
                            into the main GROMACS rather than a separate
                            code (since otherwise it will go away as
                            GROMACS is further developed).<br>
                            <br>
                            Would this be a possibility going forward,
                            devels?<br>
                            <br>
                            Thanks!<br>
                            <br>
                            <br>
                            <br>
                            <br>
                            On Fri, Aug 16, 2013 at 3:07 AM, Floris
                            Buelens &lt;<a moz-do-not-send="true"
                              href="mailto:floris_buelens@yahoo.com"
                              target="_blank">floris_buelens@yahoo.com</a>&gt;
                            wrote:<br>
                            <br>
                            Hi David,<br>
                            <blockquote class="gmail_quote"
                              style="margin:0 0 0 .8ex;border-left:1px
                              #ccc solid;padding-left:1ex">
                              I have the same requirement as you, but
                              I've gone about it slightly differently.
                              I've hacked the couple-moltype code path
                              to allow decoupling of a specific residue
                              (identified by name) instead of a
                              molecule. This allows the residue of
                              interest to be part of another molecule
                              block so you can set up distance, angle
                              and dihedral restraints using regular
                              bonded interactions with full perturbation
                              support.<br>
                              <br>
                              If this is useful to you or to anyone
                              else, let me know and I'll be happy to
                              share.<br>
                              <br>
                              best,<br>
                              <br>
                              Floris<br>
                              <br>
                              ________________________________<br>
                              From: David Mobley &lt;<a
                                moz-do-not-send="true"
                                href="mailto:dmobley@gmail.com"
                                target="_blank">dmobley@gmail.com</a>&gt;<br>
                              <br>
                              To: Discussion list for GROMACS
                              development &lt;<a moz-do-not-send="true"
                                href="mailto:gmx-developers@gromacs.org"
                                target="_blank">gmx-developers@gromacs.org</a>&gt;<br>
                              Sent: Friday, 14 June 2013, 21:47<br>
                              Subject: [gmx-developers] free energy
                              calculations with restraints<br>
                              <br>
                              <br>
                              <br>
                              <br>
                              Hi, devs,<br>
                              <br>
                              <br>
                              I'm writing with an issue relating to the
                              interplay of new free energy features with
                              restraints.<br>
                              <br>
                              I'm very much appreciating some of the new
                              free energy features in gromacs, such as
                              the 'couple-moltype' option which provides
                              a way to set up decoupling or annihilation
                              of a specific molecule via free energy
                              calculations without having to edit the
                              topology file directly. This is especially
                              great when it comes to decoupling --
                              charge decoupling was not previously
                              possible via topology file editing, and
                              vdW decoupling took substantial
                              manipulation of the topology file.<br>
                              <br>
                              However, for binding free energies, my
                              work employs orientational restraints
                              between the protein and ligand. I need to
                              be able to impose both dihedral and angle
                              restraints on angles between the protein
                              and ligand. Currently, I do this using
                              angle-restraints and dihedral-restraints.
                              This requires that both the protein and
                              ligand be within the same logical
                              'molecule', which (unfortunately) means
                              that I can't make use of the new free
                              energy features above, since
                              couple-moltype has to apply to a whole
                              molecule, not just part of a molecule.<br>
                              <br>
                              So, my I see two possible solutions to the
                              problem, and hence have these questions:<br>
                              1) Can dihedral and angle restraints be
                              applied via the pull code? If not, are
                              there any plans to add that?<br>
                              2) Alternatively, what about modifying the
                              restraints code so it uses (or at least
                              optionally allows) absolute atom
                              numbering, rather than numbering within a
                              specific molecule, thus allowing
                              restraints (angle-restraints and
                              dihedral-restraints) to be applied between
                              'molecules'?<br>
                              <br>
                              Thanks!<br>
                              David<br>
                              <br>
                              <br>
                              <br>
                              --<br>
                              David Mobley<br>
                              <a moz-do-not-send="true"
                                href="mailto:dmobley@gmail.com"
                                target="_blank">dmobley@gmail.com</a><br>
                              <a moz-do-not-send="true"
                                href="tel:949-385-2436"
                                value="+19493852436" target="_blank">949-385-2436</a><br>
                              <br>
                              <br>
                              --<br>
                              gmx-developers mailing list<br>
                              <a moz-do-not-send="true"
                                href="mailto:gmx-developers@gromacs.org"
                                target="_blank">gmx-developers@gromacs.org</a><br>
                              <a moz-do-not-send="true"
                                href="http://lists.gromacs.org/mailman/listinfo/gmx-developers"
                                target="_blank">http://lists.gromacs.org/mailman/listinfo/gmx-developers</a><br>
                              Please don't post (un)subscribe requests
                              to the list. Use the<br>
                              www interface or send it to <a
                                moz-do-not-send="true"
                                href="mailto:gmx-developers-request@gromacs.org"
                                target="_blank">gmx-developers-request@gromacs.org</a>.<br>
                              <br>
                            </blockquote>
                          </blockquote>
                        </blockquote>
                        <br>
                        -- <br>
                        gmx-developers mailing list<br>
                        <a moz-do-not-send="true"
                          href="mailto:gmx-developers@gromacs.org"
                          target="_blank">gmx-developers@gromacs.org</a><br>
                        <a moz-do-not-send="true"
                          href="http://lists.gromacs.org/mailman/listinfo/gmx-developers"
                          target="_blank">http://lists.gromacs.org/mailman/listinfo/gmx-developers</a><br>
                        Please don't post (un)subscribe requests to the
                        list. Use the www interface or send it to <a
                          moz-do-not-send="true"
                          href="mailto:gmx-developers-request@gromacs.org"
                          target="_blank">gmx-developers-request@gromacs.org</a>.<br>
                      </div>
                    </div>
                  </blockquote>
                </div>
                <br>
                <br clear="all">
                <br>
              </div>
            </div>
            <div class="im">-- <br>
              <div dir="ltr">David Mobley
                <div>Assistant Professor</div>
                <div>Department of Pharmaceutical Sciences<br>
                </div>
                <div>Department of Chemistry<br>
                </div>
                <div>3134B Natural Sciences I</div>
                <div>University of California, Irvine</div>
                <div>Irvine, CA 92697</div>
                <div><a moz-do-not-send="true"
                    href="mailto:dmobley@uci.edu" target="_blank">dmobley@uci.edu</a></div>
                <div>work <a moz-do-not-send="true"
                    href="tel:%28949%29%20824-6383" value="+19498246383"
                    target="_blank">(949) 824-6383</a></div>
                <div>cell <a moz-do-not-send="true"
                    href="tel:%28949%29%20385-2436" value="+19493852436"
                    target="_blank">(949) 385-2436</a></div>
                <div style="text-align:left">
                  <font color="#666666" face="Arial, sans-serif"><br>
                  </font></div>
              </div>
            </div>
          </div>
        </div>
        <br>
        <br clear="all">
        <br>
        -- <br>
        David Mobley<br>
        <a moz-do-not-send="true" href="mailto:dmobley@gmail.com"
          target="_blank">dmobley@gmail.com</a><br>
        949-385-2436<br>
      </div>
      <br>
      <fieldset class="mimeAttachmentHeader"></fieldset>
      <br>
    </blockquote>
    <br>
  </body>
</html>