Dzialajacy INTERTYP=1 i INTERTYP=2
[unres.git] / source / unres / src_MD / energy_p_new_barrier.F
index 51fbd11..3979f64 100644 (file)
@@ -483,6 +483,7 @@ cMS$ATTRIBUTES C ::  proc_proc
       include 'COMMON.CONTROL'
       include 'COMMON.TIME1'
       include 'COMMON.MAXGRAD'
+      include 'COMMON.SCCOR'
 #ifdef TIMING
 #ifdef MPI
       time01=MPI_Wtime()
@@ -493,8 +494,9 @@ cMS$ATTRIBUTES C ::  proc_proc
 #ifdef DEBUG
       write (iout,*) "sum_gradient gvdwc, gvdwx"
       do i=1,nres
-        write (iout,'(i3,3f10.5,5x,3f10.5,5x,f10.5)') 
-     &   i,(gvdwx(j,i),j=1,3),(gvdwc(j,i),j=1,3)
+        write (iout,'(i3,3f10.5,5x,3f10.5,5x,3f10.5,5x,3f10.5)') 
+     &   i,(gvdwx(j,i),j=1,3),(gvdwcT(j,i),j=1,3),(gvdwc(j,i),j=1,3),
+     &   (gvdwcT(j,i),j=1,3)
       enddo
       call flush(iout)
 #endif
@@ -527,6 +529,21 @@ c      enddo
       call flush(iout)
 #endif
 #ifdef SPLITELE
+#ifdef TSCSC
+      do i=1,nct
+        do j=1,3
+          gradbufc(j,i)=wsc*gvdwc(j,i)+wsc*wscT*gvdwcT(j,i)+
+     &                wscp*(gvdwc_scp(j,i)+gvdwc_scpp(j,i))+
+     &                welec*gelc_long(j,i)+wvdwpp*gvdwpp(j,i)+
+     &                wel_loc*gel_loc_long(j,i)+
+     &                wcorr*gradcorr_long(j,i)+
+     &                wcorr5*gradcorr5_long(j,i)+
+     &                wcorr6*gradcorr6_long(j,i)+
+     &                wturn6*gcorr6_turn_long(j,i)+
+     &                wstrain*ghpbc(j,i)
+        enddo
+      enddo 
+#else
       do i=1,nct
         do j=1,3
           gradbufc(j,i)=wsc*gvdwc(j,i)+
@@ -540,6 +557,7 @@ c      enddo
      &                wstrain*ghpbc(j,i)
         enddo
       enddo 
+#endif
 #else
       do i=1,nct
         do j=1,3
@@ -571,6 +589,16 @@ c      enddo
           gradbufc_sum(j,i)=gradbufc(j,i)
         enddo
       enddo
+c      call MPI_AllReduce(gradbufc(1,1),gradbufc_sum(1,1),3*nres,
+c     &    MPI_DOUBLE_PRECISION,MPI_SUM,FG_COMM,IERR)
+c      time_reduce=time_reduce+MPI_Wtime()-time00
+#ifdef DEBUG
+c      write (iout,*) "gradbufc_sum after allreduce"
+c      do i=1,nres
+c        write (iout,'(i3,3f10.5)') i,(gradbufc_sum(j,i),j=1,3)
+c      enddo
+c      call flush(iout)
+#endif
 #ifdef TIMING
 c      time_allreduce=time_allreduce+MPI_Wtime()-time00
 #endif
@@ -585,6 +613,17 @@ c      time_allreduce=time_allreduce+MPI_Wtime()-time00
      &                  " jgrad_end  ",jgrad_end(i),
      &                  i=igrad_start,igrad_end)
 #endif
+c
+c Obsolete and inefficient code; we can make the effort O(n) and, therefore,
+c do not parallelize this part.
+c
+c      do i=igrad_start,igrad_end
+c        do j=jgrad_start(i),jgrad_end(i)
+c          do k=1,3
+c            gradbufc(k,i)=gradbufc(k,i)+gradbufc_sum(k,j)
+c          enddo
+c        enddo
+c      enddo
       do j=1,3
         gradbufc(j,nres-1)=gradbufc_sum(j,nres)
       enddo
@@ -623,6 +662,16 @@ c      time_allreduce=time_allreduce+MPI_Wtime()-time00
           gradbufc(j,i)=gradbufc(j,i+1)+gradbufc_sum(j,i+1)
         enddo
       enddo
+c      do i=nnt,nres-1
+c        do k=1,3
+c          gradbufc(k,i)=0.0d0
+c        enddo
+c        do j=i+1,nres
+c          do k=1,3
+c            gradbufc(k,i)=gradbufc(k,i)+gradbufc(k,j)
+c          enddo
+c        enddo
+c      enddo
 #ifdef DEBUG
       write (iout,*) "gradbufc after summing"
       do i=1,nres
@@ -661,7 +710,7 @@ c      time_allreduce=time_allreduce+MPI_Wtime()-time00
           gradc(j,i,icg)=gradbufc(j,i)+welec*gelc(j,i)+
      &                wel_loc*gel_loc(j,i)+
      &                0.5d0*(wscp*gvdwc_scpp(j,i)+
-     &                welec*gelc_long(j,i)
+     &                welec*gelc_long(j,i)+
      &                wel_loc*gel_loc_long(j,i)+
      &                wcorr*gcorr_long(j,i)+
      &                wcorr5*gradcorr5_long(j,i)+
@@ -677,11 +726,20 @@ c      time_allreduce=time_allreduce+MPI_Wtime()-time00
      &                wsccor*gsccorc(j,i)
      &               +wscloc*gscloc(j,i)
 #endif
+#ifdef TSCSC
+          gradx(j,i,icg)=wsc*gvdwx(j,i)+wsc*wscT*gvdwxT(j,i)+
+     &                  wscp*gradx_scp(j,i)+
+     &                  wbond*gradbx(j,i)+
+     &                  wstrain*ghpbx(j,i)+wcorr*gradxorr(j,i)+
+     &                  wsccor*gsccorx(j,i)
+     &                 +wscloc*gsclocx(j,i)
+#else
           gradx(j,i,icg)=wsc*gvdwx(j,i)+wscp*gradx_scp(j,i)+
      &                  wbond*gradbx(j,i)+
      &                  wstrain*ghpbx(j,i)+wcorr*gradxorr(j,i)+
      &                  wsccor*gsccorx(j,i)
      &                 +wscloc*gsclocx(j,i)
+#endif
         enddo
       enddo 
 #ifdef DEBUG
@@ -698,7 +756,6 @@ c      time_allreduce=time_allreduce+MPI_Wtime()-time00
      &   +wturn3*gel_loc_turn3(i)
      &   +wturn6*gel_loc_turn6(i)
      &   +wel_loc*gel_loc_loc(i)
-     &   +wsccor*gsccor_loc(i)
       enddo
 #ifdef DEBUG
       write (iout,*) "gloc after adding corr"
@@ -764,6 +821,10 @@ c
       do i=1,nct
         gvdwc_norm=dsqrt(scalar(gvdwc(1,i),gvdwc(1,i)))
         if (gvdwc_norm.gt.gvdwc_max) gvdwc_max=gvdwc_norm
+#ifdef TSCSC
+        gvdwc_norm=dsqrt(scalar(gvdwcT(1,i),gvdwcT(1,i)))
+        if (gvdwc_norm.gt.gvdwc_max) gvdwc_max=gvdwc_norm          
+#endif
         gvdwc_scp_norm=dsqrt(scalar(gvdwc_scp(1,i),gvdwc_scp(1,i)))
         if (gvdwc_scp_norm.gt.gvdwc_scp_max) 
      &   gvdwc_scp_max=gvdwc_scp_norm
@@ -802,6 +863,10 @@ c
         if (gscloc_norm.gt.gscloc_max) gscloc_max=gscloc_norm
         gvdwx_norm=dsqrt(scalar(gvdwx(1,i),gvdwx(1,i)))
         if (gvdwx_norm.gt.gvdwx_max) gvdwx_max=gvdwx_norm
+#ifdef TSCSC
+        gvdwx_norm=dsqrt(scalar(gvdwxT(1,i),gvdwxT(1,i)))
+        if (gvdwx_norm.gt.gvdwx_max) gvdwx_max=gvdwx_norm
+#endif
         gradx_scp_norm=dsqrt(scalar(gradx_scp(1,i),gradx_scp(1,i)))
         if (gradx_scp_norm.gt.gradx_scp_max) 
      &    gradx_scp_max=gradx_scp_norm
@@ -5585,7 +5650,7 @@ C Proline-Proline pair is a special case...
      &  restyp(itype(i-2)),i-2,restyp(itype(i-1)),i-1,itori,itori1,
      &  (v1(j,itori,itori1),j=1,6),(v2(j,itori,itori1),j=1,6)
         gloc(i-3,icg)=gloc(i-3,icg)+wtor*gloci
-c       write (iout,*) 'i=',i,' gloc=',gloc(i-3,icg)
+        write (iout,*) 'i=',i,' gloc=',gloc(i-3,icg)
       enddo
 ! 6/20/98 - dihedral angle constraints
       edihcnstr=0.0d0
@@ -5700,6 +5765,7 @@ c      do i=1,ndih_constr
         else
           difi=0.0
         endif
+c        write (iout,*) "gloci", gloc(i-3,icg)
 cd        write (iout,'(2i5,4f8.3,2e14.5)') i,itori,rad2deg*phii,
 cd     &    rad2deg*phi0(i),  rad2deg*drange(i),
 cd     &    rad2deg*difi,0.25d0*ftors*difi**4,gloc(itori-3,icg)
@@ -5736,7 +5802,6 @@ c     lprn=.true.
         phii1=phi(i+1)
         gloci1=0.0D0
         gloci2=0.0D0
-C Regular cosine and sine terms
         do j=1,ntermd_1(itori,itori1,itori2)
           v1cij=v1c(1,j,itori,itori1,itori2)
           v1sij=v1s(1,j,itori,itori1,itori2)
@@ -5771,6 +5836,7 @@ C Regular cosine and sine terms
         enddo
         gloc(i-3,icg)=gloc(i-3,icg)+wtor_d*gloci1
         gloc(i-2,icg)=gloc(i-2,icg)+wtor_d*gloci2
+c        write (iout,*) "gloci", gloc(i-3,icg)
       enddo
       return
       end
@@ -5803,26 +5869,55 @@ C Set lprn=.true. for debugging
 c      lprn=.true.
 c      write (iout,*) "EBACK_SC_COR",iphi_start,iphi_end,nterm_sccor
       esccor=0.0D0
-      do i=iphi_start,iphi_end
+      do i=iphi_start-1,iphi_end+1
         esccor_ii=0.0D0
-        itori=itype(i-2)
-        itori1=itype(i-1)
+        isccori=isccortyp(itype(i-2))
+        isccori1=isccortyp(itype(i-1))
         phii=phi(i)
+cccc  Added 9 May 2012
+cc Tauangle is torsional engle depending on the value of first digit 
+c(see comment below)
+cc Omicron is flat angle depending on the value of first digit 
+c(see comment below)
+
         gloci=0.0D0
-        do j=1,nterm_sccor
-          v1ij=v1sccor(j,itori,itori1)
-          v2ij=v2sccor(j,itori,itori1)
-          cosphi=dcos(j*phii)
-          sinphi=dsin(j*phii)
+        do intertyp=1,1 !intertyp
+cc Added 09 May 2012 (Adasko)
+cc  Intertyp means interaction type of backbone mainchain correlation: 
+c   1 = SC...Ca...Ca...Ca
+c   2 = Ca...Ca...Ca...SC
+c   3 = SC...Ca...Ca...SC
+        if (((intertyp.eq.3).and.((itype(i-2).eq.10).or.
+     &      (itype(i-1).eq.10).or.(itype(i-2).eq.21).or.
+     &      (itype(i-1).eq.21)))
+     &    .or. ((intertyp.eq.1).and.((itype(i-2).eq.10)
+     &     .or.(itype(i-2).eq.21)))
+     &    .or.((intertyp.eq.2).and.((itype(i-1).eq.10).or.
+     &      (itype(i-1).eq.21)))) cycle  
+        if ((intertyp.eq.2).and.(i.le.iphi_start-1)) cycle
+        if ((intertyp.eq.1).and.(i.ge.iphi_end+1)) cycle
+        do j=1,nterm_sccor(isccori,isccori1)
+          v1ij=v1sccor(j,intertyp,isccori,isccori1)
+          v2ij=v2sccor(j,intertyp,isccori,isccori1)
+          cosphi=dcos(j*tauangle(intertyp,i))
+          sinphi=dsin(j*tauangle(intertyp,i))
           esccor=esccor+v1ij*cosphi+v2ij*sinphi
           gloci=gloci+j*(v2ij*cosphi-v1ij*sinphi)
         enddo
+        gloc_sc(intertyp,i-3,icg)=gloc_sc(intertyp,i-3,icg)+wsccor*gloci
+c        write (iout,*) "WTF",intertyp,i,itype(i),
+c     &gloc_sc(intertyp,i-3,icg)
         if (lprn)
      &  write (iout,'(2(a3,2x,i3,2x),2i3,6f8.3/26x,6f8.3/)')
      &  restyp(itype(i-2)),i-2,restyp(itype(i-1)),i-1,itori,itori1,
-     &  (v1sccor(j,itori,itori1),j=1,6),(v2sccor(j,itori,itori1),j=1,6)
+     &  (v1sccor(j,intertyp,itori,itori1),j=1,6)
+     & ,(v2sccor(j,intertyp,itori,itori1),j=1,6)
         gsccor_loc(i-3)=gsccor_loc(i-3)+gloci
+       enddo !intertyp
       enddo
+c        do i=1,nres
+c        write (iout,*) "W@T@F",  gloc_sc(1,i,icg),gloc(i,icg)
+c        enddo
       return
       end
 c----------------------------------------------------------------------------