<html>
  <head>
    <meta http-equiv="Content-Type" content="text/html; charset=utf-8">
  </head>
  <body text="#000000" bgcolor="#FFFFFF">
    <div class="moz-cite-prefix">Hi,<br>
      <br>
      I'll just chime in and say that if there are triple bonds four
      atoms will, with a rather high possibility, be on a straight line
      in atomic force fields as well.<br>
      <br>
      Cheers,<br>
      <br>
      Magnus<br>
      <br>
      On 2017-10-29 09:18, Erik Lindahl wrote:<br>
    </div>
    <blockquote type="cite"
cite="mid:CAEJJM8Gwcz6-n+W54DFH6tGp5GZDoHco-Up_3f3LL3icxS2egw@mail.gmail.com">
      <div dir="ltr">Hi Pavel,
        <div><br>
        </div>
        <div>Molecular Dynamics is a chaotic process, so no matter what
          precision you use the results will diverge - it will just take
          a bit longer in double precision.</div>
        <div><br>
        </div>
        <div>The reason why the method works is that we rely on sampling
          phase space (with realistic dynamics) rather than making an
          infinite-precision prediction of how a particle/atom will move
          - remember that we anyway won't know the initial conditions
          (in particular atom velocities) perfectly. For this reason,
          single precision is virtually always accurate enough, apart
          from a few places where we do summation of forces. One good
          way to test this is to check the extent to which the total
          energy is conserved.</div>
        <div><br>
        </div>
        <div>For the dihedral formula, IIRC the reason we changed it is
          that we wanted to be able to use so-called "coarse-grained"
          force fields. For normal atomic force fields the four atoms
          can *never* be on a straight line (since the angle potentials
          are too large and will prevent that), but for coarse-grained
          force fields the dihedral angles won't correspond directly to
          four bonded atoms, and then they can all be very close to
          linear once in a blue moon.</div>
        <div><br>
        </div>
        <div>In that case, the scalar product between the two vectors
          will be *very* close to 1 or -1. Both the way you calculate
          the scalar product (where you might add a very small term to a
          large one) and the behavior of the acos() function (whose
          derivative approaches infinity at -1 and +1) can lead to very
          large numerical errors.</div>
        <div><br>
        </div>
        <div>Cheers,</div>
        <div><br>
        </div>
        <div>Erik</div>
        <div><br>
        </div>
        <div><br>
        </div>
        <div><br>
        </div>
      </div>
      <div class="gmail_extra"><br>
        <div class="gmail_quote">On Sun, Oct 29, 2017 at 12:00 AM, Pavel
          Panchekha <span dir="ltr">&lt;<a
              href="mailto:me@pavpanchekha.com" target="_blank"
              moz-do-not-send="true">me@pavpanchekha.com</a>&gt;</span>
          wrote:<br>
          <blockquote class="gmail_quote" style="margin:0 0 0
            .8ex;border-left:1px #ccc solid;padding-left:1ex">
            <div dir="ltr">
              <div class="gmail_default" style="font-size:small">Berk,
                thank you for the reply. The quality of the GROMACS
                source code is clear from reading it.</div>
              <div class="gmail_default" style="font-size:small"><br>
              </div>
              <div class="gmail_default" style="font-size:small">Thank
                you for explaining that the sign is not used except in
                computing the angle. Do you know how much error in the
                angle is acceptable? I am trying some variations on the
                dihedral angle computation and seeing answers that
                differ by 1e-8. Could that small a difference possibly
                matter?<br>
              </div>
              <div class="gmail_default" style="font-size:small"><br>
              </div>
              <div class="gmail_default" style="font-size:small">I've
                noticed that earlier versions of GROMACS used a slightly
                different formula to compute <u>dih_angle</u> (using <u>acos</u>
                instead of <u>atan2</u>). Did the switch cause any
                known changes in behavior?<br>
              </div>
            </div>
            <div class="gmail_extra"><br clear="all">
              <div>
                <div class="m_-5706029853271738848gmail_signature"
                  data-smartmail="gmail_signature">
                  <div dir="ltr">—Pavel Panchekha</div>
                </div>
              </div>
              <br>
              <div class="gmail_quote">
                <div>
                  <div class="h5">On Sat, Oct 28, 2017 at 1:06 AM, Berk
                    Hess <span dir="ltr">&lt;<a
                        href="mailto:hess@kth.se" target="_blank"
                        moz-do-not-send="true">hess@kth.se</a>&gt;</span>
                    wrote:<br>
                  </div>
                </div>
                <blockquote class="gmail_quote" style="margin:0 0 0
                  .8ex;border-left:1px #ccc solid;padding-left:1ex">
                  <div>
                    <div class="h5">
                      <div text="#000000" bgcolor="#FFFFFF">
                        <div
                          class="m_-5706029853271738848m_3495107603991124129moz-cite-prefix">Hi,<br>
                          <br>
                          In molecular dynamics it's rather likely that
                          you hit corner cases, so we pay a lot of
                          attention to numerical issues when writing and
                          reviewing code. Therefore it is unlikely that
                          there are cases that can be improved, but you
                          are welcome to look.<br>
                          <br>
                          It is good that you mention the sign in
                          dihedral angle. The return value is actually
                          never used. The sign is only used to set the
                          sign of the angle, which switches at angle=0,
                          so there is no stability issue.<br>
                          <br>
                          Cheers,<br>
                          <br>
                          Berk
                          <div>
                            <div class="m_-5706029853271738848h5"><br>
                              <br>
                              On 10/28/2017 07:16 AM, Pavel Panchekha
                              wrote:<br>
                            </div>
                          </div>
                        </div>
                        <div>
                          <div class="m_-5706029853271738848h5">
                            <blockquote type="cite">
                              <div dir="ltr">
                                <div class="gmail_default"
                                  style="font-size:small">Hello,</div>
                                <div class="gmail_default"
                                  style="font-size:small"><br>
                                </div>
                                <div class="gmail_default"
                                  style="font-size:small">I'm a student
                                  at the University of Washington
                                  working on numerical accuracy, and
                                  I've been looking into GROMACS as an
                                  example of a beloved and important
                                  large numerical program. Browsing
                                  through the code base, I was struck by
                                  the computation of the dihedral angle
                                  in <u>bonded.c</u>, at line 1569 (in
                                  GROMACS 5.1.4).</div>
                                <div class="gmail_default"
                                  style="font-size:small"><br>
                                </div>
                                <div class="gmail_default"
                                  style="font-size:small">I wonder if it
                                  is possible for the <u>sign</u>
                                  variable to get the wrong value, if
                                  the input points are close to one
                                  another or have some other reason to
                                  experience cancellation when doing the
                                  cross product. I think this is fixable
                                  using an adaptive algorithm. Could
                                  changing the sign of the angle
                                  computed by <u>dih_angle</u> affect
                                  GROMACS's results?<br>
                                </div>
                                <div class="gmail_default"
                                  style="font-size:small"><br>
                                </div>
                                <div class="gmail_default"
                                  style="font-size:small">More broadly,
                                  has anyone ever had numerical accuracy
                                  issues with GROMACS? Would there be
                                  interest in reviewing a patch that
                                  improved the accuracy (for edge case
                                  inputs) for the dihedral angle
                                  computation?<br>
                                </div>
                                <div class="gmail_default"
                                  style="font-size:small"><br
                                    clear="all">
                                </div>
                                <div>
                                  <div
class="m_-5706029853271738848m_3495107603991124129m_3853295731210088818gmail_signature"
                                    data-smartmail="gmail_signature">
                                    <div dir="ltr">—Pavel Panchekha</div>
                                  </div>
                                </div>
                              </div>
                              <br>
                              <fieldset
                                class="m_-5706029853271738848m_3495107603991124129mimeAttachmentHeader"></fieldset>
                              <br>
                            </blockquote>
                            <br>
                          </div>
                        </div>
                      </div>
                      <br>
                    </div>
                  </div>
                  <span class="HOEnZb"><font color="#888888">--<br>
                      Gromacs Developers mailing list<br>
                      <br>
                      * Please search the archive at <a
                        href="http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_List"
                        rel="noreferrer" target="_blank"
                        moz-do-not-send="true">http://www.gromacs.org/Support<wbr>/Mailing_Lists/GMX-developers_<wbr>List</a>
                      before posting!<br>
                      <br>
                      * Can't post? Read <a
                        href="http://www.gromacs.org/Support/Mailing_Lists"
                        rel="noreferrer" target="_blank"
                        moz-do-not-send="true">http://www.gromacs.org/Support<wbr>/Mailing_Lists</a><br>
                      <br>
                      * For (un)subscribe requests visit<br>
                      <a
href="https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-developers"
                        rel="noreferrer" target="_blank"
                        moz-do-not-send="true">https://maillist.sys.kth.se/ma<wbr>ilman/listinfo/gromacs.org_gmx<wbr>-developers</a>
                      or send a mail to <a
                        href="mailto:gmx-developers-request@gromacs.org"
                        target="_blank" moz-do-not-send="true">gmx-developers-request@gromacs<wbr>.org</a>.<br>
                    </font></span></blockquote>
              </div>
              <br>
            </div>
            <br>
            --<br>
            Gromacs Developers mailing list<br>
            <br>
            * Please search the archive at <a
              href="http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_List"
              rel="noreferrer" target="_blank" moz-do-not-send="true">http://www.gromacs.org/<wbr>Support/Mailing_Lists/GMX-<wbr>developers_List</a>
            before posting!<br>
            <br>
            * Can't post? Read <a
              href="http://www.gromacs.org/Support/Mailing_Lists"
              rel="noreferrer" target="_blank" moz-do-not-send="true">http://www.gromacs.org/<wbr>Support/Mailing_Lists</a><br>
            <br>
            * For (un)subscribe requests visit<br>
            <a
href="https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-developers"
              rel="noreferrer" target="_blank" moz-do-not-send="true">https://maillist.sys.kth.se/<wbr>mailman/listinfo/gromacs.org_<wbr>gmx-developers</a>
            or send a mail to <a
              href="mailto:gmx-developers-request@gromacs.org"
              moz-do-not-send="true">gmx-developers-request@<wbr>gromacs.org</a>.<br>
          </blockquote>
        </div>
        <br>
        <br clear="all">
        <div><br>
        </div>
        -- <br>
        <div class="gmail_signature" data-smartmail="gmail_signature">
          <div dir="ltr">
            <div>
              <div dir="ltr">
                <div>
                  <div dir="ltr">
                    <div>
                      <div dir="ltr">
                        <div>Erik Lindahl &lt;<a
                            href="mailto:erik.lindahl@dbb.su.se"
                            target="_blank" moz-do-not-send="true">erik.lindahl@dbb.su.se</a>&gt;</div>
                        <div>Professor of Biophysics, Dept. Biochemistry
                          &amp; Biophysics, Stockholm University</div>
                        <div>Science for Life Laboratory, Box 1031,
                          17121 Solna, Sweden</div>
                      </div>
                    </div>
                  </div>
                </div>
              </div>
            </div>
          </div>
        </div>
      </div>
      <br>
      <fieldset class="mimeAttachmentHeader"></fieldset>
      <br>
    </blockquote>
    <p><br>
    </p>
  </body>
</html>