bugfix in shield FGPROC>1

[unres.git] / source / unres / src_MD-M / energy_p_new-sep_barrier.F

diff --git a/source/unres/src_MD-M/energy_p_new-sep_barrier.F b/source/unres/src_MD-M/energy_p_new-sep_barrier.F
index 0a638ef..44384a6 100644 (file)
--- a/source/unres/src_MD-M/energy_p_new-sep_barrier.F
+++ b/source/unres/src_MD-M/energy_p_new-sep_barrier.F
@@ -6,7 +6,7 @@ C      if(r.lt.r_cut-rlamb) then
 C        sscale=1.0d0
 C      else if(r.le.r_cut.and.r.ge.r_cut-rlamb) then
 C        gamm=(r-(r_cut-rlamb))/rlamb
-        sscalelip=1.0d0+r*r*(2*r-3.0d0)
+        sscalelip=1.0d0+r*r*(2.0d0*r-3.0d0)
 C      else
 C        sscale=0d0
 C      endif
@@ -20,7 +20,7 @@ C     if(r.lt.r_cut-rlamb) then
 C        sscagrad=0.0d0
 C      else if(r.le.r_cut.and.r.ge.r_cut-rlamb) then
 C        gamm=(r-(r_cut-rlamb))/rlamb
-        sscagradlip=r*(6*r-6.0d0)
+        sscagradlip=r*(6.0d0*r-6.0d0)
 C      else
 C        sscagrad=0.0d0
 C      endif
@@ -2683,7 +2683,6 @@ C Uncomment following three lines for Ca-p interactions
        endif
 
           rrij=1.0D0/(xj*xj+yj*yj+zj*zj)
-
           sss=sscale(1.0d0/(dsqrt(rrij)))
           sssgrad=sscagrad(1.0d0/(dsqrt(rrij)))
           if (sss.lt.1.0d0) then
@@ -2704,7 +2703,7 @@ C Calculate contributions to the gradient in the virtual-bond and SC vectors.
 C
              
             fac=-(evdwij+e1)*rrij*(1.0d0-sss)
-            fac=fac-(evdwij)*sssgrad*dsqrt(rrij)
+            fac=fac-(evdwij)*sssgrad*dsqrt(rrij)/expon
             ggg(1)=xj*fac
             ggg(2)=yj*fac
             ggg(3)=zj*fac
@@ -2795,7 +2794,7 @@ C Uncomment following three lines for Ca-p interactions
           if (xj.lt.0) xj=xj+boxxsize
           yj=mod(yj,boxysize)
           if (yj.lt.0) yj=yj+boxysize
-          zj=mod(zi,boxzsize)
+          zj=mod(zj,boxzsize)
           if (zj.lt.0) zj=zj+boxzsize
       dist_init=(xj-xi)**2+(yj-yi)**2+(zj-zi)**2
       xj_safe=xj
@@ -2831,7 +2830,7 @@ C Uncomment following three lines for Ca-p interactions
           rrij=1.0D0/(xj*xj+yj*yj+zj*zj)
           sss=sscale(1.0d0/(dsqrt(rrij)))
           sssgrad=sscagrad(1.0d0/(dsqrt(rrij)))
-          if (sss.gt.0.0d0) then
+            if (sss.gt.0.0d0) then
 
             fac=rrij**expon2
             e1=fac*fac*aad(itypj,iteli)
@@ -2845,11 +2844,12 @@ C Uncomment following three lines for Ca-p interactions
             evdw2=evdw2+evdwij*sss
             if (energy_dec) write (iout,'(a6,2i5,2(0pf7.3))')
      &          'evdw2',i,j,sss,evdwij
+
 C
 C Calculate contributions to the gradient in the virtual-bond and SC vectors.
 C
             fac=-(evdwij+e1)*rrij*sss
-            fac=fac+(evdwij)*sssgrad*dsqrt(rrij)
+            fac=fac+(evdwij)*sssgrad*dsqrt(rrij)/expon
             ggg(1)=xj*fac
             ggg(2)=yj*fac
             ggg(3)=zj*fac