Adam 7/30/2014

[unres.git] / source / wham / src / energy_p_new.F

diff --git a/source/wham/src/energy_p_new.F b/source/wham/src/energy_p_new.F
index bf5c3e7..9b69cf7 100644 (file)
--- a/source/wham/src/energy_p_new.F
+++ b/source/wham/src/energy_p_new.F
@@ -2869,16 +2869,24 @@ C Evaluate bridge-strain energy and its gradient in virtual-bond and SC vectors.
 C
       implicit real*8 (a-h,o-z)
       include 'DIMENSIONS'
+      include 'DIMENSIONS.ZSCOPT'
       include 'COMMON.SBRIDGE'
       include 'COMMON.CHAIN'
       include 'COMMON.DERIV'
       include 'COMMON.VAR'
       include 'COMMON.INTERACT'
       include 'COMMON.IOUNITS'
+      include 'COMMON.NAMES'
       dimension ggg(3)
       ehpb=0.0D0
+#ifdef DEBUG
+      do i=1,nres
+        write (iout,'(a4,2x,i4,3f10.5,5x,3f10.5)') restyp(itype(i)),i,
+     &      (c(j,i),j=1,3),(c(j,i+nres),j=1,3)
+      enddo
 cd      write(iout,*)'edis: nhpb=',nhpb,' fbr=',fbr
 cd      write(iout,*)'link_start=',link_start,' link_end=',link_end
+#endif
       if (link_end.eq.0) return
       do i=link_start,link_end
 C If ihpb(i) and jhpb(i) > NRES, this is a SC-SC distance, otherwise a
@@ -2893,8 +2901,10 @@ C iii and jjj point to the residues for which the distance is assigned.
           iii=ii
           jjj=jj
         endif
-c        write (iout,*) "i",i," ii",ii," iii",iii," jj",jj," jjj",jjj,
-c     &    dhpb(i),dhpb1(i),forcon(i)
+#ifdef DEBUG
+        write (iout,*) "i",i," ii",ii," iii",iii," jj",jj," jjj",jjj,
+     &    dhpb(i),dhpb1(i),forcon(i)
+#endif
 C 24/11/03 AL: SS bridges handled separately because of introducing a specific
 C    distance and angle dependent SS bond potential.
         if (ii.gt.nres .and. itype(iii).eq.1 .and. itype(jjj).eq.1) then
@@ -2907,8 +2917,10 @@ c Restraints from contact prediction
           if (dhpb1(i).gt.0.0d0) then
             ehpb=ehpb+2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i))
             fac=forcon(i)*gnmr1prim(dd,dhpb(i),dhpb1(i))/dd
-c            write (iout,*) "beta nmr",
-c     &        dd,2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i))
+#ifdef DEBUG
+            write (iout,*) "beta nmr",
+     &        dd,2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i))
+#endif
           else
             dd=dist(ii,jj)
             rdis=dd-dhpb(i)
@@ -2916,7 +2928,9 @@ C Get the force constant corresponding to this distance.
             waga=forcon(i)
 C Calculate the contribution to energy.
             ehpb=ehpb+waga*rdis*rdis
-c            write (iout,*) "beta reg",dd,waga*rdis*rdis
+#ifdef DEBUG
+            write (iout,*) "beta reg",dd,waga*rdis*rdis
+#endif
 C
 C Evaluate gradient.
 C
@@ -2940,15 +2954,19 @@ C target distance.
           if (dhpb1(i).gt.0.0d0) then
             ehpb=ehpb+2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i))
             fac=forcon(i)*gnmr1prim(dd,dhpb(i),dhpb1(i))/dd
-c            write (iout,*) "alph nmr",
-c     &        dd,2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i))
+#ifdef DEBUG
+            write (iout,*) "alph nmr",
+     &        dd,2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i))
+#endif
           else
             rdis=dd-dhpb(i)
 C Get the force constant corresponding to this distance.
             waga=forcon(i)
 C Calculate the contribution to energy.
             ehpb=ehpb+waga*rdis*rdis
-c            write (iout,*) "alpha reg",dd,waga*rdis*rdis
+#ifdef DEBUG
+            write (iout,*) "alpha reg",dd,waga*rdis*rdis
+#endif
 C
 C Evaluate gradient.
 C
@@ -4463,24 +4481,16 @@ c     lprn=.true.
         itori=itortyp(itype(i-2))
         itori1=itortyp(itype(i-1))
         itori2=itortyp(itype(i))
-c        iblock=1
-c        if (iabs(itype(i+1)).eq.20) iblock=2
         phii=phi(i)
         phii1=phi(i+1)
         gloci1=0.0D0
         gloci2=0.0D0
 C Regular cosine and sine terms
-c c       do j=1,ntermd_1(itori,itori1,itori2,iblock)
-c          v1cij=v1c(1,j,itori,itori1,itori2,iblock)
-c          v1sij=v1s(1,j,itori,itori1,itori2,iblock)
-c          v2cij=v1c(2,j,itori,itori1,itori2,iblock)
-c          v2sij=v1s(2,j,itori,itori1,itori2,iblock)
-       do j=1,ntermd_1(itori,itori1,itori2)
+        do j=1,ntermd_1(itori,itori1,itori2)
           v1cij=v1c(1,j,itori,itori1,itori2)
           v1sij=v1s(1,j,itori,itori1,itori2)
           v2cij=v1c(2,j,itori,itori1,itori2)
           v2sij=v1s(2,j,itori,itori1,itori2)
-
           cosphi1=dcos(j*phii)
           sinphi1=dsin(j*phii)
           cosphi2=dcos(j*phii1)
@@ -4491,12 +4501,7 @@ c          v2sij=v1s(2,j,itori,itori1,itori2,iblock)
           gloci2=gloci2+j*(v2sij*cosphi2-v2cij*sinphi2)
         enddo
         do k=2,ntermd_2(itori,itori1,itori2)
-c        do k=2,ntermd_2(itori,itori1,itori2,iblock)
           do l=1,k-1
-c            v1cdij = v2c(k,l,itori,itori1,itori2,iblock)
-c            v2cdij = v2c(l,k,itori,itori1,itori2,iblock)
-c            v1sdij = v2s(k,l,itori,itori1,itori2,iblock)
-c            v2sdij = v2s(l,k,itori,itori1,itori2,iblock)
             v1cdij = v2c(k,l,itori,itori1,itori2)
             v2cdij = v2c(l,k,itori,itori1,itori2)
             v1sdij = v2s(k,l,itori,itori1,itori2)