unres Adam's changes

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

diff --git a/source/unres/src-HCD-5D/lagrangian_lesyng.F b/source/unres/src-HCD-5D/lagrangian_lesyng.F
index f57a432..6dd113b 100644 (file)
--- a/source/unres/src-HCD-5D/lagrangian_lesyng.F
+++ b/source/unres/src-HCD-5D/lagrangian_lesyng.F
@@ -8,7 +8,7 @@ c-------------------------------------------------------------------------
        include 'DIMENSIONS'
 #ifdef MPI
        include 'mpif.h'
-       integer time00
+       double precision time00
 #endif
        include 'COMMON.VAR'
        include 'COMMON.CHAIN'
@@ -55,7 +55,7 @@ c-------------------------------------------------------------------------
 #endif
 #ifdef FIVEDIAG
       call grad_transform
-      d_a=0.0d0
+      d_a(:,:2*nres)=0.0d0
       if (lprn) then
         write (iout,*) "Potential forces backbone"
         do i=1,nres
@@ -459,9 +459,9 @@ c          write (iout,*) "i",i," itype",itype(i),ntyp1
           endif
         enddo
       enddo
-      DMorig=DM
-      DU1orig=DU1
-      DU2orig=DU2
+      DMorig(:2*nres)=DM(:2*nres)
+      DU1orig(:2*nres)=DU1(:2*nres)
+      DU2orig(:2*nres)=DU2(:2*nres)
       if (gmatout) then
       write (iout,*)"The upper part of the five-diagonal inertia matrix"
       endif
@@ -886,15 +886,25 @@ c---------------------------------------------------------------------------
 c---------------------------------------------------------------------------
       subroutine fivediaginv_mult(ndim,forces,d_a_vec)
       implicit none
+#ifdef MPI
+      include 'mpif.h'
+#endif
       include 'DIMENSIONS'
       include 'COMMON.CHAIN'
       include 'COMMON.IOUNITS'
       include 'COMMON.LAGRANGE.5diag'
       include 'COMMON.INTERACT'
+      include 'COMMON.VAR'
       integer ndim
       double precision forces(3*ndim),accel(3,0:maxres2),rs(ndim),
      &  xsolv(ndim),d_a_vec(6*nres)
       integer i,j,ind,ichain,n,iposc,innt,inct,inct_prev
+#ifdef TIMING
+      include 'COMMON.TIME1'
+      double precision time01
+      time01=MPI_Wtime()
+#endif
+      accel(:,:2*nres)=0.0d0
       do j=1,3
 Compute accelerations in Calpha and SC
         do ichain=1,nchain
@@ -903,9 +913,18 @@ Compute accelerations in Calpha and SC
           innt=chain_border(1,ichain)
           inct=chain_border(2,ichain)
           do i=iposc,iposc+n-1
-            rs(i)=forces(3*(i-1)+j)
+            rs(i-iposc+1)=forces(3*(i-1)+j)
           enddo
+#ifdef DEBUG
+          write (iout,*) "j",j," chain",ichain
+          write (iout,*) "rs"
+          write (iout,'(f10.5)') (rs(i),i=1,n)
+#endif
           call FDISYS (n,DM(iposc),DU1(iposc),DU2(iposc),rs,xsolv)
+#ifdef DEBUG
+          write (iout,*) "xsolv"
+          write (iout,'(f10.5)') (xsolv(i),i=1,n)
+#endif
           ind=1
           do i=innt,inct
             if (itype(i).eq.10)then
@@ -919,7 +938,7 @@ Compute accelerations in Calpha and SC
           enddo
         enddo
       enddo
-C Conevert d_a to virtual-bon-vector basis
+C Convert d_a to virtual-bon-vector basis
 #ifdef DEBUG
       write (iout,*) "accel in CA-SC basis"
       do i=1,nres
@@ -944,6 +963,7 @@ C Conevert d_a to virtual-bon-vector basis
         end if
       enddo
       accel(:,nres)=0.0d0
+      accel(:,nct)=0.0d0
       accel(:,2*nres)=0.0d0
       if (nnt.gt.1) then
         accel(:,0)=accel(:,1)
@@ -984,9 +1004,12 @@ C Conevert d_a to virtual-bon-vector basis
           ind=ind+3
         endif
       enddo
+#ifdef TIMING
+        time_ginvmult=time_ginvmult+MPI_Wtime()-time01
+#endif
 #ifdef DEBUG
       write (iout,*) "d_a_vec"
-      write (iout,'(3f10.5)') (d_a_vec(j),j=1,dimen3)
+      write (iout,'(3f10.5)') (d_a_vec(j),j=1,3*(nct-nnt+nside))
 #endif
       return
       end