dihedral class 2 bug fix

Hello.

Here is the bug fix for dihedral_style class2 which seems to fix the
sign problem (and another small problem).

(Test results after the bug fix
The test-system consists of 4 consecutively bonded atoms, initially at
a +45 degree phi angle. The in.test input script, I use a phi1
parameter of -60 degrees. All other terms and coefficients were set
to zero. After equilibration, the molecule adopts a -60 degree angle
as it should. The energy values range from 0 to 20.0 as they should
(K1=10.0), with the 0 minima at -60.)

Cheers
Andrew

dihedral_class2.cpp (29.7 KB)

test_system_class2.tar.gz (1.47 KB)

Does it also work w/out the minimum_image() calls
you added? Those should no longer be necessary
in any bonded interaction.

Steve

Does it also work w/out the minimum_image() calls
you added? Those should no longer be necessary
in any bonded interaction.

Hi
It -seems- to work fine without them. I'm attaching this code without
those lines, and with all of the annoying flags I added removed. I
also made this change to "dihedral_table.h" as well. (See attached.
"dihedral_table.cpp" does not need to be modified.)

After removing these lines, I tested the system moving part of the
molecule to the other side of the boundary, (in a way that should
cause the phi angle to change sign). The molecule still behaves
exactly as before, at least as far as I can tell visually. Test
details are at the end of this email. (I attached the messy testing
files I used).

Although not a big deal, some source files do continue to invoke calls
to domain->minimum_image(). They are:

pair_hbond_dreiding_lj.cpp
pair_hbond_dreiding_morse.cpp
pair_lj_cut_coul_long_tip4p.cpp
pppm_tip4p.cpp

(I didn't check the other optional packages)

------Here are the ugly testing details. Feel free to skip them:---------

For reference, in all tests, I am using a cubic boundary box of size 40

  0.000000 40.0000 xlo xhi
  0.000000 40.0000 ylo yhi
  0.000000 40.0000 zlo zhi

I don't know if I tested this correctly, but I tried two different
versions of of the 4-atom test system, each with one of the atoms on
the other side of the boundary. (The 4th atom. See data file
excerpts below.) (If the dihedral code had been looking at the wrong
image, then this should have had a huge effect on the phi angle,
changing it from a positive to a negative angle, causing the restoring
force to pull it in the opposite direction.) However in all cases,
the initial frames of animation pull on this atom in the correct
direction, eventually coming to rest at the -60 degree minima, as it
did in the earlier tests.

-- excerpt from data file version #1 --

Atoms

    1 1 2 0.0 11.00 11.00 0.50
    2 1 1 0.0 10.00 12.00 0.50
    3 1 1 0.0 12.00 12.00 0.50
    4 1 2 0.0 13.00 11.00 39.5 # <-- this atom crosses the boundary

-- excerpt from data file version #2 --

Atoms

    1 1 2 0.0 11.00 11.00 0.50
    2 1 1 0.0 10.00 12.00 0.50
    3 1 1 0.0 12.00 12.00 0.50
    4 1 2 0.0 13.00 11.00 -0.5 # <-- this atom crosses the boundary

-- original version --

Atoms

    1 1 2 0.0 11.00 11.00 10.0
    2 1 1 0.0 10.00 12.00 10.0
    3 1 1 0.0 12.00 12.00 10.0
    4 1 2 0.0 13.00 11.00 9.0

dihedral_class2.cpp (29.1 KB)

dihedral_table.h (10.7 KB)

testclass2_at_boundary.tar.gz (2.19 KB)

Although not a big deal, some source files do continue to invoke calls
to domain->minimum_image(). They are:

pair_hbond_dreiding_lj.cpp
pair_hbond_dreiding_morse.cpp
pair_lj_cut_coul_long_tip4p.cpp
pppm_tip4p.cpp

That's b/c they do something different. They look
up a bonded partner by its atom ID. That may
return a far-away image of the close atom,
so they need to use minimum_image().
The bond/angle/dihedral files no longer
do that.

Steve