rozgrzebany shield

diff --git a/source/unres/src_MD-M/energy_p_new_barrier.F b/source/unres/src_MD-M/energy_p_new_barrier.F
index 82d6a3f..5cb9a2a 100644 (file)
--- a/source/unres/src_MD-M/energy_p_new_barrier.F
+++ b/source/unres/src_MD-M/energy_p_new_barrier.F
@@ -3708,8 +3708,8 @@ C 12/26/95 - for the evaluation of multi-body H-bonding interactions
           if (shield_mode.gt.0) then
 C          fac_shield(i)=0.4
 C          fac_shield(j)=0.6
-          el1=el1*fac_shield(i)*fac_shield(j)
-          el2=el2*fac_shield(i)*fac_shield(j)
+          el1=el1*fac_shield(i)**2*fac_shield(j)**2
+          el2=el2*fac_shield(i)**2*fac_shield(j)**2
           eesij=(el1+el2)
           ees=ees+eesij
           else
@@ -3756,9 +3756,10 @@ C          print *,i,j
            iresshield=shield_list(ilist,i)
            do k=1,3
            rlocshield=grad_shield_side(k,ilist,i)*eesij/fac_shield(i)
+     &      *2.0
            gshieldx(k,iresshield)=gshieldx(k,iresshield)+
      &              rlocshield
-     & +grad_shield_loc(k,ilist,i)*eesij/fac_shield(i)
+     & +grad_shield_loc(k,ilist,i)*eesij/fac_shield(i)*2.0
             gshieldc(k,iresshield-1)=gshieldc(k,iresshield-1)+rlocshield
 C           gshieldc_loc(k,iresshield)=gshieldc_loc(k,iresshield)
 C     & +grad_shield_loc(k,ilist,i)*eesij/fac_shield(i)
@@ -3781,9 +3782,10 @@ C              endif
            iresshield=shield_list(ilist,j)
            do k=1,3
            rlocshield=grad_shield_side(k,ilist,j)*eesij/fac_shield(j)
+     &     *2.0
            gshieldx(k,iresshield)=gshieldx(k,iresshield)+
      &              rlocshield
-     & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j)
+     & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j)*2.0
            gshieldc(k,iresshield-1)=gshieldc(k,iresshield-1)+rlocshield
 
 C     & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j)
@@ -3806,13 +3808,13 @@ C              endif
 
           do k=1,3
             gshieldc(k,i)=gshieldc(k,i)+
-     &              grad_shield(k,i)*eesij/fac_shield(i)
+     &              grad_shield(k,i)*eesij/fac_shield(i)*2.0
             gshieldc(k,j)=gshieldc(k,j)+
-     &              grad_shield(k,j)*eesij/fac_shield(j)
+     &              grad_shield(k,j)*eesij/fac_shield(j)*2.0
             gshieldc(k,i-1)=gshieldc(k,i-1)+
-     &              grad_shield(k,i)*eesij/fac_shield(i)
+     &              grad_shield(k,i)*eesij/fac_shield(i)*2.0
             gshieldc(k,j-1)=gshieldc(k,j-1)+
-     &              grad_shield(k,j)*eesij/fac_shield(j)
+     &              grad_shield(k,j)*eesij/fac_shield(j)*2.0
 
            enddo
            endif
@@ -3931,7 +3933,7 @@ cd        print '(2i3,2(3(1pd14.5),3x))',i,j,(dcosb(k),k=1,3),
 cd   &          (dcosg(k),k=1,3)
           do k=1,3
             ggg(k)=(ecosb*dcosb(k)+ecosg*dcosg(k))*
-     &      fac_shield(i)*fac_shield(j)
+     &      fac_shield(i)**2*fac_shield(j)**2
           enddo
 c          do k=1,3
 c            ghalf=0.5D0*ggg(k)
@@ -3952,11 +3954,11 @@ C                     print *,"before22", gelc_long(1,i), gelc_long(1,j)
             gelc(k,i)=gelc(k,i)
      &           +((ecosa*(dc_norm(k,j)-cosa*dc_norm(k,i))
      &           + ecosb*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1))
-     &           *fac_shield(i)*fac_shield(j)   
+     &           *fac_shield(i)**2*fac_shield(j)**2   
             gelc(k,j)=gelc(k,j)
      &           +((ecosa*(dc_norm(k,i)-cosa*dc_norm(k,j))
      &           + ecosg*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1))
-     &           *fac_shield(i)*fac_shield(j)
+     &           *fac_shield(i)**2*fac_shield(j)**2
             gelc_long(k,j)=gelc_long(k,j)+ggg(k)
             gelc_long(k,i)=gelc_long(k,i)-ggg(k)
           enddo
@@ -7302,6 +7304,7 @@ c------------------------------------------------------------------------------
       include 'COMMON.DERIV'
       include 'COMMON.INTERACT'
       include 'COMMON.CONTACTS'
+      include 'COMMON.SHIELD'
       double precision gx(3),gx1(3)
       logical lprn
       lprn=.false.
@@ -7720,6 +7723,7 @@ C This subroutine calculates multi-body contributions to hydrogen-bonding
       include 'COMMON.CONTACTS'
       include 'COMMON.CHAIN'
       include 'COMMON.CONTROL'
+      include 'COMMON.SHIELD'
       double precision gx(3),gx1(3)
       integer num_cont_hb_old(maxres)
       logical lprn,ldone
@@ -8022,6 +8026,7 @@ cd               write (iout,*) "grij_hb_cont i1",grij_hb_cont(:,jj,i1)
                 call calc_eello(i,jp,i+1,jp1,jj,kk)
                 if (wcorr4.gt.0.0d0) 
      &            ecorr=ecorr+eello4(i,jp,i+1,jp1,jj,kk)
+CC     &            *fac_shield(i)**2*fac_shield(j)**2
                   if (energy_dec.and.wcorr4.gt.0.0d0) 
      1                 write (iout,'(a6,4i5,0pf7.3)')
      2                'ecorr4',i,j,i+1,j1,eello4(i,jp,i+1,jp1,jj,kk)
@@ -8141,6 +8146,7 @@ c------------------------------------------------------------------------------
       include 'COMMON.DERIV'
       include 'COMMON.INTERACT'
       include 'COMMON.CONTACTS'
+      include 'COMMON.SHIELD'
       double precision gx(3),gx1(3)
       logical lprn
       lprn=.false.
@@ -8151,7 +8157,8 @@ c------------------------------------------------------------------------------
       ees0mij=ees0m(jj,i)
       ees0mkl=ees0m(kk,k)
       ekont=eij*ekl
-      ees=-(coeffp*ees0pij*ees0pkl+coeffm*ees0mij*ees0mkl)
+      ees=-(coeffp*ees0pij*ees0pkl+coeffm*ees0mij*ees0mkl)*
+     & fac_shield(i)**2*fac_shield(j)**2
 cd    ees=-(coeffp*ees0pkl+coeffm*ees0mkl)
 C Following 4 lines for diagnostics.
 cd    ees0pkl=0.0D0