Merge branch 'prerelease-3.2.1' of mmka.chem.univ.gda.pl:unres into prerelease-3.2.1
[unres.git] / source / unres / src_MD / energy_p_new_barrier.F
index 5ae5a43..47583a4 100644 (file)
@@ -3023,6 +3023,9 @@ C
 C Loop over i,i+2 and i,i+3 pairs of the peptide groups
 C
       do i=iturn3_start,iturn3_end
 C Loop over i,i+2 and i,i+3 pairs of the peptide groups
 C
       do i=iturn3_start,iturn3_end
+C        if (itype(i).eq.21 .or. itype(i+1).eq.21
+C     &  .or. itype(i+2).eq.21 .or. itype(i+3).eq.21.or.itype(i+4).eq.21)
+C     &  cycle
         dxi=dc(1,i)
         dyi=dc(2,i)
         dzi=dc(3,i)
         dxi=dc(1,i)
         dyi=dc(2,i)
         dzi=dc(3,i)
@@ -3038,6 +3041,10 @@ C
         num_cont_hb(i)=num_conti
       enddo
       do i=iturn4_start,iturn4_end
         num_cont_hb(i)=num_conti
       enddo
       do i=iturn4_start,iturn4_end
+C        if (itype(i).eq.21 .or. itype(i+1).eq.21
+C     &  .or. itype(i+2).eq.21 .or. itype(i+3).eq.21.or.itype(i+4).eq.21
+C     &  .or. itype(i+5).eq.21)
+C     & cycle
         dxi=dc(1,i)
         dyi=dc(2,i)
         dzi=dc(3,i)
         dxi=dc(1,i)
         dyi=dc(2,i)
         dzi=dc(3,i)
@@ -3056,6 +3063,8 @@ c
 c Loop over all pairs of interacting peptide groups except i,i+2 and i,i+3
 c
       do i=iatel_s,iatel_e
 c Loop over all pairs of interacting peptide groups except i,i+2 and i,i+3
 c
       do i=iatel_s,iatel_e
+C          if (itype(i).eq.21 .or. itype(i+1).eq.21
+C     &.or.itype(i+2)) cycle
         dxi=dc(1,i)
         dyi=dc(2,i)
         dzi=dc(3,i)
         dxi=dc(1,i)
         dyi=dc(2,i)
         dzi=dc(3,i)
@@ -3068,6 +3077,8 @@ c
 c        write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i)
         num_conti=num_cont_hb(i)
         do j=ielstart(i),ielend(i)
 c        write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i)
         num_conti=num_cont_hb(i)
         do j=ielstart(i),ielend(i)
+C          if (itype(j).eq.21 .or. itype(j+1).eq.21
+C     &.or.itype(j+2)) cycle
           call eelecij(i,j,ees,evdw1,eel_loc)
         enddo ! j
         num_cont_hb(i)=num_conti
           call eelecij(i,j,ees,evdw1,eel_loc)
         enddo ! j
         num_cont_hb(i)=num_conti
@@ -4775,7 +4786,7 @@ C
      & cosph1(maxsingle),sinph1(maxsingle),cosph2(maxsingle),
      & sinph2(maxsingle),cosph1ph2(maxdouble,maxdouble),
      & sinph1ph2(maxdouble,maxdouble)
      & cosph1(maxsingle),sinph1(maxsingle),cosph2(maxsingle),
      & sinph2(maxsingle),cosph1ph2(maxdouble,maxdouble),
      & sinph1ph2(maxdouble,maxdouble)
-      logical lprn /.false./, lprn1 /.true./
+      logical lprn /.false./, lprn1 /.false./
       etheta=0.0D0
       do i=ithet_start,ithet_end
         if ((itype(i-1).eq.ntyp1).or.(itype(i-2).eq.ntyp1).or.
       etheta=0.0D0
       do i=ithet_start,ithet_end
         if ((itype(i-1).eq.ntyp1).or.(itype(i-2).eq.ntyp1).or.
@@ -4789,7 +4800,8 @@ C
           coskt(k)=dcos(k*theti2)
           sinkt(k)=dsin(k*theti2)
         enddo
           coskt(k)=dcos(k*theti2)
           sinkt(k)=dsin(k*theti2)
         enddo
-        if (i.gt.3) then
+C        if (i.gt.3) then
+         if (i.gt.3 .and. itype(i-3).ne.ntyp1) then
 #ifdef OSF
           phii=phi(i)
           if (phii.ne.phii) phii=150.0
 #ifdef OSF
           phii=phi(i)
           if (phii.ne.phii) phii=150.0
           enddo
         else
           phii=0.0d0
           enddo
         else
           phii=0.0d0
-          ityp1=nthetyp+1
+          ityp1=ithetyp(itype(i-2))
           do k=1,nsingle
             cosph1(k)=0.0d0
             sinph1(k)=0.0d0
           enddo 
         endif
           do k=1,nsingle
             cosph1(k)=0.0d0
             sinph1(k)=0.0d0
           enddo 
         endif
-        if (i.lt.nres) then
+        if ((i.lt.nres).and. itype(i+1).ne.ntyp1) then
 #ifdef OSF
           phii1=phi(i+1)
           if (phii1.ne.phii1) phii1=150.0
 #ifdef OSF
           phii1=phi(i+1)
           if (phii1.ne.phii1) phii1=150.0
@@ -4824,7 +4836,7 @@ C
           enddo
         else
           phii1=0.0d0
           enddo
         else
           phii1=0.0d0
-          ityp3=nthetyp+1
+          ityp3=ithetyp(itype(i))
           do k=1,nsingle
             cosph2(k)=0.0d0
             sinph2(k)=0.0d0
           do k=1,nsingle
             cosph2(k)=0.0d0
             sinph2(k)=0.0d0
@@ -4930,11 +4942,9 @@ C
           enddo
         enddo
 10      continue
           enddo
         enddo
 10      continue
-c        lprn1=.true.
         if (lprn1) write (iout,'(a4,i2,3f8.1,9h ethetai ,f10.5)') 
      &  'ebe', i,theta(i)*rad2deg,phii*rad2deg,
      &   phii1*rad2deg,ethetai
         if (lprn1) write (iout,'(a4,i2,3f8.1,9h ethetai ,f10.5)') 
      &  'ebe', i,theta(i)*rad2deg,phii*rad2deg,
      &   phii1*rad2deg,ethetai
-c        lprn1=.false.
         etheta=etheta+ethetai
         if (i.gt.3) gloc(i-3,icg)=gloc(i-3,icg)+wang*dephii
         if (i.lt.nres) gloc(i-2,icg)=gloc(i-2,icg)+wang*dephii1
         etheta=etheta+ethetai
         if (i.gt.3) gloc(i-3,icg)=gloc(i-3,icg)+wang*dephii
         if (i.lt.nres) gloc(i-2,icg)=gloc(i-2,icg)+wang*dephii1
@@ -5958,17 +5968,19 @@ c      lprn=.true.
 c      write (iout,*) "EBACK_SC_COR",iphi_start,iphi_end,nterm_sccor
       esccor=0.0D0
       do i=itau_start,itau_end
 c      write (iout,*) "EBACK_SC_COR",iphi_start,iphi_end,nterm_sccor
       esccor=0.0D0
       do i=itau_start,itau_end
+C        do i=42,42
         esccor_ii=0.0D0
         if ((itype(i-2).eq.ntyp1).or.(itype(i-1).eq.ntyp1)) cycle
         isccori=isccortyp(itype(i-2))
         isccori1=isccortyp(itype(i-1))
         phii=phi(i)
         esccor_ii=0.0D0
         if ((itype(i-2).eq.ntyp1).or.(itype(i-1).eq.ntyp1)) cycle
         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)
 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)
-
+C        print *,i,tauangle(1,i)
         
 c        do intertyp=1,3 !intertyp
         do intertyp=2,2 !intertyp
         
 c        do intertyp=1,3 !intertyp
         do intertyp=2,2 !intertyp
@@ -5996,6 +6008,7 @@ c   3 = SC...Ca...Ca...SCi
           esccor=esccor+v1ij*cosphi+v2ij*sinphi
           gloci=gloci+j*(v2ij*cosphi-v1ij*sinphi)
         enddo
           esccor=esccor+v1ij*cosphi+v2ij*sinphi
           gloci=gloci+j*(v2ij*cosphi-v1ij*sinphi)
         enddo
+C        print *,i,tauangle(1,i),gloci
         gloc_sc(intertyp,i-3,icg)=gloc_sc(intertyp,i-3,icg)+wsccor*gloci
 c        write (iout,*) "WTF",intertyp,i,itype(i),v1ij*cosphi+v2ij*sinphi
 c     &gloc_sc(intertyp,i-3,icg)
         gloc_sc(intertyp,i-3,icg)=gloc_sc(intertyp,i-3,icg)+wsccor*gloci
 c        write (iout,*) "WTF",intertyp,i,itype(i),v1ij*cosphi+v2ij*sinphi
 c     &gloc_sc(intertyp,i-3,icg)