ctest src_MD-M ref ene and missing files

[unres.git] / source / unres / src_MD-M / intcartderiv.F

diff --git a/source/unres/src_MD-M/intcartderiv.F b/source/unres/src_MD-M/intcartderiv.F
index 94b3989..a890c0a 100644 (file)
--- a/source/unres/src_MD-M/intcartderiv.F
+++ b/source/unres/src_MD-M/intcartderiv.F
@@ -46,14 +46,14 @@ c We need dtheta(:,:,i-1) to compute dphi(:,:,i)
         cost=dcos(theta(i))
        sint=sqrt(1-cost*cost)
         do j=1,3
-C         if (itype(i-1).ne.21) then
           dcostheta(j,1,i)=-(dc_norm(j,i-1)+cost*dc_norm(j,i-2))/
      &   vbld(i-1)
-          if (itype(i-1).ne.21)  dtheta(j,1,i)=-dcostheta(j,1,i)/sint
+c          if (itype(i-1).ne.ntyp1)
+          dtheta(j,1,i)=-dcostheta(j,1,i)/sint
           dcostheta(j,2,i)=-(dc_norm(j,i-2)+cost*dc_norm(j,i-1))/
      &   vbld(i)
-          if (itype(i-1).ne.21)  dtheta(j,2,i)=-dcostheta(j,2,i)/sint
-C          endif
+c          if (itype(i-1).ne.ntyp1)
+          dtheta(j,2,i)=-dcostheta(j,2,i)/sint
         enddo
       enddo
 #if defined(MPI) && defined(PARINTDER)
@@ -101,7 +101,8 @@ c conventional formulas around 0 and 180.
 #else
       do i=4,nres      
 #endif
-c        if (itype(i-1).eq.21 .or. itype(i-2).eq.21 ) cycle
+c        if (itype(i-2).eq.ntyp1.or. itype(i-1).eq.ntyp1
+c     &      .or. itype(i).eq.ntyp1 .or. itype(i-3).eq.ntyp1) cycle
 c the conventional case
         sint=dsin(theta(i))
        sint1=dsin(theta(i-1))
@@ -126,8 +127,8 @@ c    Obtaining the gamma derivatives from sine derivative
             ctgt=cost/sint
             ctgt1=cost1/sint1
             cosg_inv=1.0d0/cosg
-            if (itype(i-1).ne.21 .and. itype(i-2).ne.21) then
-           dsinphi(j,1,i)=-sing*ctgt1*dtheta(j,1,i-1)
+c            if (itype(i-1).ne.ntyp1 .and. itype(i-2).ne.ntyp1) then
+      dsinphi(j,1,i)=-sing*ctgt1*dtheta(j,1,i-1)
      &        -(fac0*vp1(j)+sing*dc_norm(j,i-3))*vbld_inv(i-2)
             dphi(j,1,i)=cosg_inv*dsinphi(j,1,i)
             dsinphi(j,2,i)=
@@ -138,13 +139,13 @@ c    Obtaining the gamma derivatives from sine derivative
      &        +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i)
 c     &        +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
             dphi(j,3,i)=cosg_inv*dsinphi(j,3,i)
-            endif
+c            endif
 c Bug fixed 3/24/05 (AL)
         enddo                                              
 c   Obtaining the gamma derivatives from cosine derivative
         else
            do j=1,3
-           if (itype(i-1).ne.21 .and. itype(i-2).ne.21) then
+c           if (itype(i-1).ne.ntyp1 .and. itype(i-2).ne.ntyp1) then
            dcosphi(j,1,i)=fac1*dcostheta(j,1,i-1)+fac3*
      &    dcostheta(j,1,i-1)-fac0*(dc_norm(j,i-1)-scalp*
      &     dc_norm(j,i-3))/vbld(i-2)
@@ -157,10 +158,16 @@ c   Obtaining the gamma derivatives from cosine derivative
      &    dcostheta(j,2,i)-fac0*(dc_norm(j,i-3)-scalp*
      &     dc_norm(j,i-1))/vbld(i)
            dphi(j,3,i)=-1/sing*dcosphi(j,3,i)       
-           endif
+c           endif
          enddo
         endif                                                                                           
       enddo
+      do i=1,nres-1
+        do j=1,3
+        dc_norm2(j,i+nres)=-dc_norm(j,i+nres)
+cc       write(iout,*) dc_norm2(j,i-2+nres),"dcnorm"
+        enddo
+      enddo
 Calculate derivative of Tauangle
 #ifdef PARINTDER
       do i=itau_start,itau_end
@@ -179,10 +186,6 @@ c the conventional case
         cost=dcos(theta(i))
         cost1=dcos(omicron(2,i-1))
         cosg=dcos(tauangle(1,i))
-        do j=1,3
-        dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres)
-cc       write(iout,*) dc_norm2(j,i-2+nres),"dcnorm"
-        enddo
         scalp=scalar(dc_norm2(1,i-2+nres),dc_norm(1,i-1))
         fac0=1.0d0/(sint1*sint)
         fac1=cost*fac0
@@ -208,7 +211,7 @@ c    Obtaining the gamma derivatives from sine derivative
             dsintau(j,1,2,i)=
      &        -sing*(ctgt1*domicron(j,2,1,i-1)+ctgt*dtheta(j,1,i))
      &        -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
-c            write(iout,*) "dsintau", dsintau(j,1,2,i)
+c            write(iout,*) "dsintau", dsintau(j,1,1,i),dsintau(j,1,2,i)
             dtauangle(j,1,2,i)=cosg_inv*dsintau(j,1,2,i)
 c Bug fixed 3/24/05 (AL)
             dsintau(j,1,3,i)=-sing*ctgt*dtheta(j,2,i)
@@ -235,7 +238,7 @@ c         write (iout,*) "else",i
          enddo
         endif
 c        do k=1,3                 
-c        write(iout,*) "tu",i,k,(dtauangle(j,1,k,i),j=1,3)        
+c        write(iout,*) "tu",1,i,k,(dtauangle(j,1,k,i),j=1,3)        
 c        enddo                
       enddo
 CC Second case Ca...Ca...Ca...SC
@@ -390,7 +393,7 @@ c   Derivatives of side-chain angles alpha and omega
 #else
         do i=2,nres-1          
 #endif
-          if(itype(i).ne.10 .and. itype(i).ne.21) then   
+          if(itype(i).ne.10 .and. itype(i).ne.ntyp1) then        
              fac5=1.0d0/dsqrt(2*(1+dcos(theta(i+1))))
              fac6=fac5/vbld(i)
              fac7=fac5*fac5
@@ -624,7 +627,7 @@ c Check alpha gradient
       write (iout,*) 
      & "Analytical (upper) and numerical (lower) gradient of alpha"
       do i=2,nres-1
-       if((itype(i).ne.10).and.(itype(i).ne.ntyp1)) then
+       if(itype(i).ne.10) then
                    do j=1,3
              dcji=dc(j,i-1)
                      dc(j,i-1)=dcji+aincr
@@ -660,7 +663,7 @@ c     Check omega gradient
       write (iout,*) 
      & "Analytical (upper) and numerical (lower) gradient of omega"
       do i=2,nres-1
-       if((itype(i).ne.10).and.(itype(i).ne.ntyp1)) then
+       if(itype(i).ne.10) then
                    do j=1,3
              dcji=dc(j,i-1)
                      dc(j,i-1)=dcji+aincr
@@ -694,7 +697,7 @@ c     Check omega gradient
       enddo
       return
       end
-
+c------------------------------------------------------------
       subroutine chainbuild_cart
       implicit real*8 (a-h,o-z)
       include 'DIMENSIONS'
@@ -738,6 +741,7 @@ c        call flush(iout)
 #endif
       do j=1,3
         c(j,1)=dc(j,0)
+c        c(j,1)=c(j,1)
       enddo
       do i=2,nres
         do j=1,3
@@ -749,6 +753,7 @@ c        call flush(iout)
           c(j,i+nres)=c(j,i)+dc(j,i+nres)
         enddo
       enddo
+C      print *,'tutu'
 c      write (iout,*) "CHAINBUILD_CART"
 c      call cartprint
       call int_from_cart1(.false.)