Adam's unres update

[unres.git] / source / unres / src-HCD-5D / energy_p_new-sep_barrier.F

diff --git a/source/unres/src-HCD-5D/energy_p_new-sep_barrier.F b/source/unres/src-HCD-5D/energy_p_new-sep_barrier.F
index c4e54bc..f92aebb 100644 (file)
--- a/source/unres/src-HCD-5D/energy_p_new-sep_barrier.F
+++ b/source/unres/src-HCD-5D/energy_p_new-sep_barrier.F
@@ -94,9 +94,9 @@ c      write(iout,*)'Entering ELJ nnt=',nnt,' nct=',nct,' expon=',expon
       gg_lipi=0.0d0
       gg_lipj=0.0d0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      do ikont=g_listscsc_start_long,g_listscsc_end_long
+        i=newcontlisti_long(ikont)
+        j=newcontlistj_long(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -219,9 +219,9 @@ c      write(iout,*)'Entering ELJ nnt=',nnt,' nct=',nct,' expon=',expon
       gg_lipi=0.0d0
       gg_lipj=0.0d0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      do ikont=g_listscsc_start_short,g_listscsc_end_short
+        i=newcontlisti_short(ikont)
+        j=newcontlistj_short(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -340,9 +340,9 @@ c     print *,'Entering ELJK nnt=',nnt,' nct=',nct,' expon=',expon
       gg_lipi=0.0d0
       gg_lipj=0.0d0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      do ikont=g_listscsc_start_long,g_listscsc_end_long
+        i=newcontlisti_long(ikont)
+        j=newcontlistj_long(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -462,9 +462,9 @@ c     print *,'Entering ELJK nnt=',nnt,' nct=',nct,' expon=',expon
       gg_lipi=0.0d0
       gg_lipj=0.0d0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      do ikont=g_listscsc_start_short,g_listscsc_end_short
+        i=newcontlisti_short(ikont)
+        j=newcontlistj_short(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -586,9 +586,9 @@ c     else
 c     endif
       ind=0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      do ikont=g_listscsc_start_long,g_listscsc_end_long
+        i=newcontlisti_long(ikont)
+        j=newcontlistj_long(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -734,9 +734,9 @@ c     else
 c     endif
       ind=0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      do ikont=g_listscsc_start_short,g_listscsc_end_short
+        i=newcontlisti_short(ikont)
+        j=newcontlistj_short(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -872,7 +872,7 @@ C
      & xj_temp,yj_temp,zj_temp,dist_temp,sig,rij_shift,faclip
       double precision dist,sscale,sscagrad,sscagradlip,sscalelip
       double precision subchap,sss1,sssgrad1
-      double precision boxshift
+      double precision boxshift,rij1
       evdw=0.0D0
 ccccc      energy_dec=.false.
 c     print *,'Entering EGB nnt=',nnt,' nct=',nct,' expon=',expon
@@ -882,9 +882,12 @@ c     print *,'Entering EGB nnt=',nnt,' nct=',nct,' expon=',expon
 c     if (icall.eq.0) lprn=.false.
       ind=0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      if (energy_dec)
+     & write(2,*) "g_listscsc_start_long,g_listscsc_end_long",
+     & g_listscsc_start_long,g_listscsc_end_long
+      do ikont=g_listscsc_start_long,g_listscsc_end_long
+        i=newcontlisti_long(ikont)
+        j=newcontlistj_long(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -940,9 +943,13 @@ c            write (iout,*) "i",i," j", j," itype",itype(i),itype(j)
             dzj=dc_norm(3,nres+j)
             rrij=1.0D0/(xj*xj+yj*yj+zj*zj)
             rij=dsqrt(rrij)
-            sss1=sscale(1.0d0/rij,r_cut_int)
+            rij1=1.0d0/rij
+c            sss1=sscale(1.0d0/rij,r_cut_int)
+            sss1=sscale(rij1,r_cut_int)
             if (sss1.eq.0.0d0) cycle
-            sss=sscale(1.0d0/(rij*sigmaii(itypi,itypj)),r_cut_respa)
+            rij1=rij1/sigmaii(itypi,itypj)
+            sss=sscale(rij1,r_cut_respa)
+c            sss=sscale(1.0d0/(rij*sigmaii(itypi,itypj)),r_cut_respa)
             if (sss.lt.1.0d0) then
 C Calculate angle-dependent terms of energy and contributions to their
 C derivatives.
@@ -1026,6 +1033,7 @@ C This subroutine calculates the interaction energy of nonbonded side chains
 C assuming the Gay-Berne potential of interaction.
 C
       implicit none
+      include 'mpif.h'
       include 'DIMENSIONS'
       include 'COMMON.GEO'
       include 'COMMON.VAR'
@@ -1038,6 +1046,7 @@ C
       include 'COMMON.CALC'
       include 'COMMON.CONTROL'
       include "COMMON.SPLITELE"
+      include 'COMMON.TIME1'
       logical lprn
       double precision evdw
       integer itypi,itypj,itypi1,iint,ind,ikont
@@ -1046,6 +1055,8 @@ C
      & sslipj,ssgradlipj,ssgradlipi,sig,rij_shift,faclip
       double precision dist,sscale,sscagrad,sscagradlip,sscalelip
       double precision boxshift
+      double precision time01
+c      time01=MPI_Wtime()
       evdw=0.0D0
 ccccc      energy_dec=.false.
 c     print *,'Entering EGB nnt=',nnt,' nct=',nct,' expon=',expon
@@ -1055,9 +1066,12 @@ c     print *,'Entering EGB nnt=',nnt,' nct=',nct,' expon=',expon
 c     if (icall.eq.0) lprn=.false.
       ind=0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      if (energy_dec)
+     & write(2,*) "g_listscsc_start_short,g_listscsc_end_short",
+     & g_listscsc_start_short,g_listscsc_end_short
+      do ikont=g_listscsc_start_short,g_listscsc_end_short
+        i=newcontlisti_short(ikont)
+        j=newcontlistj_short(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -1119,8 +1133,8 @@ c     &        (2.0d0-sslipi-sslipj)/2.0d0
             rrij=1.0D0/(xj*xj+yj*yj+zj*zj)
             rij=dsqrt(rrij)
             sss=sscale(1.0d0/(rij*sigmaii(itypi,itypj)),r_cut_respa)
-          sssgrad=sscagrad((1.0d0/rij)/sigmaii(itypi,itypj),r_cut_respa)
             if (sss.gt.0.0d0) then
+          sssgrad=sscagrad((1.0d0/rij)/sigmaii(itypi,itypj),r_cut_respa)
 
 C Calculate angle-dependent terms of energy and contributions to their
 C derivatives.
@@ -1189,6 +1203,7 @@ C Calculate angular part of the gradient.
 c          enddo      ! j
 c        enddo        ! iint
       enddo          ! i
+c      time_evdw_short=time_evdw_short+MPI_Wtime()-time01
 c      write (iout,*) "Number of loop steps in EGB:",ind
 cccc      energy_dec=.false.
       return
@@ -1230,9 +1245,9 @@ c     print *,'Entering EGB nnt=',nnt,' nct=',nct,' expon=',expon
 c     if (icall.eq.0) lprn=.true.
       ind=0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      do ikont=g_listscsc_start_long,g_listscsc_end_long
+        i=newcontlisti_long(ikont)
+        j=newcontlistj_long(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -1392,9 +1407,9 @@ c     print *,'Entering EGB nnt=',nnt,' nct=',nct,' expon=',expon
 c     if (icall.eq.0) lprn=.true.
       ind=0
 c      do i=iatsc_s,iatsc_e
-      do ikont=g_listscsc_start,g_listscsc_end
-        i=newcontlisti(ikont)
-        j=newcontlistj(ikont)
+      do ikont=g_listscsc_start_short,g_listscsc_end_short
+        i=newcontlisti_short(ikont)
+        j=newcontlistj_short(ikont)
         itypi=iabs(itype(i))
         if (itypi.eq.ntyp1) cycle
         itypi1=iabs(itype(i+1))
@@ -1742,6 +1757,9 @@ c
 c Loop over all pairs of interacting peptide groups except i,i+2 and i,i+3
 c
 c      do i=iatel_s,iatel_e
+      if (energy_dec)
+     & write(iout,*) "g_listpp_start,g_listpp_end",
+     & g_listpp_start,g_listpp_end
       do ikont=g_listpp_start,g_listpp_end
         i=newcontlistppi(ikont)
         j=newcontlistppj(ikont)
@@ -2812,9 +2830,12 @@ c      write (iout,*) "iatel_s_vdw",iatel_s_vdw,
 c     & " iatel_e_vdw",iatel_e_vdw
 c      call flush(iout)
 c      do i=iatel_s_vdw,iatel_e_vdw
-      do ikont=g_listpp_vdw_start,g_listpp_vdw_end
-        i=newcontlistpp_vdwi(ikont)
-        j=newcontlistpp_vdwj(ikont)
+      if (energy_dec)
+     & write(iout,*) "g_listpp_vdw_start_short,g_listpp_vdw_end_short",
+     & g_listpp_vdw_start_short,g_listpp_vdw_end_short
+      do ikont=g_listpp_vdw_start_short,g_listpp_vdw_end_short
+        i=newcontlistpp_vdwi_short(ikont)
+        j=newcontlistpp_vdwj_short(ikont)
         if (itype(i).eq.ntyp1.or. itype(i+1).eq.ntyp1) cycle
         dxi=dc(1,i)
         dyi=dc(2,i)
@@ -2939,9 +2960,12 @@ c      if (lprint_short)
 c     &  write (iout,*) 'ESCP_LONG iatscp_s=',iatscp_s,
 c     & ' iatscp_e=',iatscp_e
 c      do i=iatscp_s,iatscp_e
-      do ikont=g_listscp_start,g_listscp_end
-        i=newcontlistscpi(ikont)
-        j=newcontlistscpj(ikont)
+      if (energy_dec)
+     & write(iout,*)"g_listscp_start_long,g_listscp_end_long",
+     & g_listscp_start_long,g_listscp_end_long
+      do ikont=g_listscp_start_long,g_listscp_end_long
+        i=newcontlistscpi_long(ikont)
+        j=newcontlistscpj_long(ikont)
         if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) cycle
         iteli=itel(i)
         xi=0.5D0*(c(1,i)+c(1,i+1))
@@ -3065,14 +3089,20 @@ C
       double precision ggg(3)
       double precision sscale,sscagrad
       double precision boxshift
+      integer ikont
       evdw2=0.0D0
       evdw2_14=0.0d0
 cd    print '(a)','Enter ESCP'
 c      if (lprint_short) 
 c     &  write (iout,*) 'ESCP_SHORT iatscp_s=',iatscp_s,
 c     & ' iatscp_e=',iatscp_e
-      if (energy_dec) write (iout,*) "escp_short:",r_cut_int,rlamb
-      do i=iatscp_s,iatscp_e
+c      if (energy_dec) write (iout,*) "escp_short:",r_cut_int,rlamb
+      if (energy_dec)
+     & write(iout,*) "g_listscp_start_short,g_listscp_end_short",
+     & g_listscp_start_short,g_listscp_end_short
+      do ikont=g_listscp_start_short,g_listscp_end_short
+        i=newcontlistscpi_short(ikont)
+        j=newcontlistscpj_short(ikont)
         if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) cycle
         iteli=itel(i)
         xi=0.5D0*(c(1,i)+c(1,i+1))
@@ -3083,9 +3113,9 @@ c     & ' iatscp_e=',iatscp_e
 c        if (lprint_short) 
 c     &    write (iout,*) "i",i," itype",itype(i),itype(i+1),
 c     &     " nscp_gr",nscp_gr(i)   
-        do iint=1,nscp_gr(i)
-
-        do j=iscpstart(i,iint),iscpend(i,iint)
+c        do iint=1,nscp_gr(i)
+c
+c        do j=iscpstart(i,iint),iscpend(i,iint)
           itypj=iabs(itype(j))
 c        if (lprint_short)
 c     &    write (iout,*) "j",j," itypj",itypj
@@ -3149,9 +3179,9 @@ c             gradx_scp(k,j)=gradx_scp(k,j)+ggg(k)
               gvdwc_scp(k,j)=gvdwc_scp(k,j)+ggg(k)
             enddo
           endif
-        enddo
+c        enddo
 
-        enddo ! iint
+c        enddo ! iint
       enddo ! i
       do i=1,nct
         do j=1,3