changes in wham and unres
[unres4.git] / source / unres / geometry.f90
index fafa747..2287f0a 100644 (file)
       end subroutine alloc_geo_arrays
 !-----------------------------------------------------------------------------
 !-----------------------------------------------------------------------------
+      subroutine returnbox
+      integer :: allareout,i,j,k,nojumpval,chain_beg,mnum
+      integer :: chain_end,ireturnval
+      real*8 :: difference
+!C change suggested by Ana - end
+        j=1
+        chain_beg=1
+!C        do i=1,nres
+!C       write(*,*) 'initial', i,j,c(j,i)
+!C        enddo
+!C change suggested by Ana - begin
+        allareout=1
+!C change suggested by Ana -end
+        do i=1,nres-1
+           mnum=molnum(i)
+         if ((itype(i,mnum).eq.ntyp1_molec(mnum))&
+            .and.(itype(i+1,mnum).eq.ntyp1_molec(mnum))) then
+          chain_end=i
+          if (allareout.eq.1) then
+            ireturnval=int(c(j,i)/boxxsize)
+            if (c(j,i).le.0) ireturnval=ireturnval-1
+            do k=chain_beg,chain_end
+              c(j,k)=c(j,k)-ireturnval*boxxsize
+              c(j,k+nres)=c(j,k+nres)-ireturnval*boxxsize
+            enddo
+!C Suggested by Ana
+            if (chain_beg.eq.1) &
+            dc_old(1,0)=dc_old(1,0)-ireturnval*boxxsize
+!C Suggested by Ana -end
+           endif
+           chain_beg=i+1
+           allareout=1
+         else
+          if (int(c(j,i)/boxxsize).eq.0) allareout=0
+         endif
+        enddo
+         if (allareout.eq.1) then
+            ireturnval=int(c(j,i)/boxxsize)
+            if (c(j,i).le.0) ireturnval=ireturnval-1
+            do k=chain_beg,nres
+              c(j,k)=c(j,k)-ireturnval*boxxsize
+              c(j,k+nres)=c(j,k+nres)-ireturnval*boxxsize
+            enddo
+          endif
+!C NO JUMP 
+!C        do i=1,nres
+!C        write(*,*) 'befor no jump', i,j,c(j,i)
+!C        enddo
+        nojumpval=0
+        do i=2,nres
+           mnum=molnum(i)
+           if (itype(i,mnum).eq.ntyp1_molec(mnum)&
+              .and. itype(i-1,mnum).eq.ntyp1_molec(mnum)) then
+             difference=abs(c(j,i-1)-c(j,i))
+!C             print *,'diff', difference
+             if (difference.gt.boxxsize/2.0) then
+                if (c(j,i-1).gt.c(j,i)) then
+                  nojumpval=1
+                 else
+                   nojumpval=-1
+                 endif
+              else
+              nojumpval=0
+              endif
+              endif
+              c(j,i)=c(j,i)+nojumpval*boxxsize
+              c(j,i+nres)=c(j,i+nres)+nojumpval*boxxsize
+         enddo
+       nojumpval=0
+        do i=2,nres
+           mnum=molnum(i)
+           if (itype(i,mnum).eq.ntyp1_molec(mnum) .and. itype(i-1,mnum).eq.ntyp1_molec(mnum)) then
+             difference=abs(c(j,i-1)-c(j,i))
+             if (difference.gt.boxxsize/2.0) then
+                if (c(j,i-1).gt.c(j,i)) then
+                  nojumpval=1
+                 else
+                   nojumpval=-1
+                 endif
+              else
+              nojumpval=0
+              endif
+             endif
+              c(j,i)=c(j,i)+nojumpval*boxxsize
+              c(j,i+nres)=c(j,i+nres)+nojumpval*boxxsize
+         enddo
+
+!C        do i=1,nres
+!C        write(*,*) 'after no jump', i,j,c(j,i)
+!C        enddo
+
+!C NOW Y dimension
+!C suggesed by Ana begins
+        allareout=1
+        j=2
+        chain_beg=1
+        do i=1,nres-1
+           mnum=molnum(i)
+         if ((itype(i,mnum).eq.ntyp1_molec(mnum))&
+           .and.(itype(i+1,mnum).eq.ntyp1_molec(mnum))) then
+          chain_end=i
+          if (allareout.eq.1) then
+            ireturnval=int(c(j,i)/boxysize)
+            if (c(j,i).le.0) ireturnval=ireturnval-1
+            do k=chain_beg,chain_end
+              c(j,k)=c(j,k)-ireturnval*boxysize
+             c(j,k+nres)=c(j,k+nres)-ireturnval*boxysize
+            enddo
+!C Suggested by Ana
+            if (chain_beg.eq.1) &
+            dc_old(1,0)=dc_old(1,0)-ireturnval*boxxsize
+!C Suggested by Ana -end
+           endif
+           chain_beg=i+1
+           allareout=1
+         else
+          if (int(c(j,i)/boxysize).eq.0) allareout=0
+         endif
+        enddo
+         if (allareout.eq.1) then
+            ireturnval=int(c(j,i)/boxysize)
+            if (c(j,i).le.0) ireturnval=ireturnval-1
+            do k=chain_beg,nres
+              c(j,k)=c(j,k)-ireturnval*boxysize
+              c(j,k+nres)=c(j,k+nres)-ireturnval*boxysize
+            enddo
+          endif
+        nojumpval=0
+        do i=2,nres
+           mnum=molnum(i)
+           if (itype(i,mnum).eq.ntyp1_molec(mnum)&
+              .and. itype(i-1,mnum).eq.ntyp1_molec(mnum)) then
+             difference=abs(c(j,i-1)-c(j,i))
+             if (difference.gt.boxysize/2.0) then
+                if (c(j,i-1).gt.c(j,i)) then
+                  nojumpval=1
+                 else
+                   nojumpval=-1
+                 endif
+             else
+              nojumpval=0
+              endif
+           endif
+              c(j,i)=c(j,i)+nojumpval*boxysize
+              c(j,i+nres)=c(j,i+nres)+nojumpval*boxysize
+         enddo
+      nojumpval=0
+        do i=2,nres
+           mnum=molnum(i)
+           if (itype(i,mnum).eq.ntyp1_molec(mnum)&
+             .and. itype(i-1,mnum).eq.ntyp1) then
+             difference=abs(c(j,i-1)-c(j,i))
+             if (difference.gt.boxysize/2.0) then
+                if (c(j,i-1).gt.c(j,i)) then
+                  nojumpval=1
+                 else
+                   nojumpval=-1
+                 endif
+              else
+              nojumpval=0
+              endif
+            endif
+              c(j,i)=c(j,i)+nojumpval*boxysize
+              c(j,i+nres)=c(j,i+nres)+nojumpval*boxysize
+         enddo
+!C Now Z dimension
+!C Suggested by Ana -begins
+        allareout=1
+!C Suggested by Ana -ends
+       j=3
+        chain_beg=1
+        do i=1,nres-1
+           mnum=molnum(i)
+         if ((itype(i,mnum).eq.ntyp1_molec(mnum))&
+           .and.(itype(i+1,mnum).eq.ntyp1_molec(mnum))) then
+          chain_end=i
+          if (allareout.eq.1) then
+            ireturnval=int(c(j,i)/boxysize)
+            if (c(j,i).le.0) ireturnval=ireturnval-1
+            do k=chain_beg,chain_end
+              c(j,k)=c(j,k)-ireturnval*boxzsize
+              c(j,k+nres)=c(j,k+nres)-ireturnval*boxzsize
+            enddo
+!C Suggested by Ana
+            if (chain_beg.eq.1) dc_old(1,0)=dc_old(1,0)-ireturnval*boxxsize
+!C Suggested by Ana -end
+           endif
+           chain_beg=i+1
+           allareout=1
+         else
+          if (int(c(j,i)/boxzsize).eq.0) allareout=0
+         endif
+        enddo
+         if (allareout.eq.1) then
+            ireturnval=int(c(j,i)/boxzsize)
+            if (c(j,i).le.0) ireturnval=ireturnval-1
+            do k=chain_beg,nres
+              c(j,k)=c(j,k)-ireturnval*boxzsize
+              c(j,k+nres)=c(j,k+nres)-ireturnval*boxzsize
+            enddo
+          endif
+        nojumpval=0
+        do i=2,nres
+           mnum=molnum(i)
+           if (itype(i,mnum).eq.ntyp1_molec(mnum) .and. itype(i-1,mnum).eq.ntyp1_molec(mnum)) then
+             difference=abs(c(j,i-1)-c(j,i))
+             if (difference.gt.(boxzsize/2.0)) then
+                if (c(j,i-1).gt.c(j,i)) then
+                  nojumpval=1
+                 else
+                   nojumpval=-1
+                 endif
+              else
+              nojumpval=0
+              endif
+            endif
+              c(j,i)=c(j,i)+nojumpval*boxzsize
+              c(j,i+nres)=c(j,i+nres)+nojumpval*boxzsize
+         enddo
+       nojumpval=0
+        do i=2,nres
+           mnum=molnum(i)
+           if (itype(i,mnum).eq.ntyp1_molec(mnum) &
+            .and. itype(i-1,mnum).eq.ntyp1_molec(mnum)) then
+             difference=abs(c(j,i-1)-c(j,i))
+             if (difference.gt.boxzsize/2.0) then
+                if (c(j,i-1).gt.c(j,i)) then
+                  nojumpval=1
+                 else
+                   nojumpval=-1
+                 endif
+              else
+              nojumpval=0
+              endif
+            endif
+             c(j,i)=c(j,i)+nojumpval*boxzsize
+              c(j,i+nres)=c(j,i+nres)+nojumpval*boxzsize
+         enddo
+        do i=1,nres
+         if (molnum(i).eq.5) then
+          c(1,i)=dmod(c(1,i),boxxsize)
+          c(2,i)=dmod(c(2,i),boxysize)
+          c(3,i)=dmod(c(3,i),boxzsize)
+          c(1,i+nres)=dmod(c(1,i+nres),boxxsize)
+          c(2,i+nres)=dmod(c(2,i+nres),boxysize)
+          c(3,i+nres)=dmod(c(3,i+nres),boxzsize)
+         endif
+        enddo
+        return
+        end       subroutine returnbox
+!-------------------------------------------------------------------------------------------------------
       end module geometry