Fix dihedral definition, add pbc_shortest_vectors method

- Dihedral definition was incorrect and giving negative results for pbc.
  Correct it. pbc_shortest_vectors indices need to be inverted.
- Update test correcting this. Now pbc dihedral gives correct with
  respect to ASE, i.e, ASE gives 300, here we get -60, probably ASE range
  is 0 to 360. Here it is (-180, 180).
- Add method for pbc_shortes_vectors with system and index of atoms.

src/getters.jl

@@ -180,13 +180,13 @@ function dihedral(
         frpos2 = pos2' * icell
         frpos3 = pos3' * icell
         frpos4 = pos4' * icell
-        vec1 = pbc_shortest_vectors(cell, frpos3, frpos4)[1,1,:]  # Axis
+        vec1 = pbc_shortest_vectors(cell, frpos1, frpos2)[1,1,:]  # Axis
         vec2 = pbc_shortest_vectors(cell, frpos2, frpos3)[1,1,:]
-        vec3 = pbc_shortest_vectors(cell, frpos1, frpos2)[1,1,:]
+        vec3 = pbc_shortest_vectors(cell, frpos3, frpos4)[1,1,:]
     else
-        vec1 = pos3 - pos4
-        vec2 = pos2 - pos3
-        vec3 = pos1 - pos2
+        vec1 = pos2 - pos1
+        vec2 = pos3 - pos2
+        vec3 = pos4 - pos3
     end
     p23 = vec2 × vec3
     p12 = vec1 × vec2

src/pbc_utils.jl

@@ -24,6 +24,10 @@ function dot_2d_mod(a,b)
 end
 
 """
+    pbc_shortest_vectors(
+        system::AbstractSystem, at1::Int, at2::Int,
+        return_dists::Val{RD}=Val(false), return_vects::Val{RV}=Val(true)
+    ) where {RD, RV}
     pbc_shortest_vectors(
         system::AbstractSystem, pos1::T, pos2::T, 
         return_dists::Val{RD}=Val(false), return_vects::Val{RV}=Val(true)
@@ -32,6 +36,15 @@ end
 Obtain the shortest vectors between position vectors, `pos1` and `pos2`,
 taking into account PBC of the given `system`.
 """
+function pbc_shortest_vectors(
+        system::AbstractSystem, at1::Int, at2::Int,
+        return_dists::Val{RD}=Val(false), return_vects::Val{RV}=Val(true)
+    ) where {RD, RV}
+    pos1 = position(system, at1)
+    pos2 = position(system, at2)
+
+    return pbc_shortest_vectors(system, pos1, pos2, return_dists, return_vects)[1,1,:]
+end
 function pbc_shortest_vectors(
         system::AbstractSystem, pos1::T, pos2::U,
         return_dists::Val{RD}=Val(false), return_vects::Val{RV}=Val(true)

test/runtests.jl

@@ -33,10 +33,9 @@ using AtomsIO
             @test angle(si, 3, 2, 5; pbc=false) ≈ 58.51784589u"°"
         end
 
-        # This needs to be verified for pbc=true. Seems inconsistent with ASE
         @testset "dihedral" begin
             @test dihedral(si, 1,5,2,6; pbc=false) ≈ 0.0u"°" atol=1e-8u"°"
-            @test dihedral(si, 5,6,1,4) ≈ 60.0u"°" atol=1e-8u"°"
+            @test dihedral(si, 5,6,1,4) ≈ -60.0u"°" atol=1e-8u"°"
             @test dihedral(si, 5,6,1,4; pbc=false) ≈ 30.0u"°" atol=1e-8u"°"
         end