box shift in wham from Adam
authorCezary Czaplewski <czarek@chem.univ.gda.pl>
Wed, 25 Jul 2018 09:13:36 +0000 (11:13 +0200)
committerCezary Czaplewski <czarek@chem.univ.gda.pl>
Wed, 25 Jul 2018 09:13:36 +0000 (11:13 +0200)
source/wham/src-M/CMakeLists.txt
source/wham/src-M/cxread.F
source/wham/src-M/oligomer.F [new file with mode: 0644]

index 8e71493..56defe7 100644 (file)
@@ -62,6 +62,7 @@ set(UNRES_WHAM_M_SRC0
        define_pairs.f
        mysort.f
         ssMD.F
+        oligomer.F
 )
 
 set(UNRES_WHAM_M_PP_SRC
@@ -90,6 +91,7 @@ set(UNRES_WHAM_M_PP_SRC
        wham_multparm.F
        xread.F
        proc_proc.c
+        oligomer.F
 ) 
 
 
index c0319f6..15777b1 100644 (file)
@@ -168,6 +168,20 @@ c      call flush(iout)
         enddo
       enddo
 
+c Box shift  
+      call oligomer
+      do i=1,nres
+        do j=1,3
+          xoord(j,i)=c(j,i)
+        enddo
+      enddo  
+      do i=1,nct-nnt+1
+        do j=1,3
+          xoord(j,i+nres)=c(j,i+nres+nnt-1)
+        enddo
+      enddo  
+c end change 
+
       if (islice.gt.0 .and. islice.le.nslice .and. (.not.separate_parset
      &    .or. iset.eq.myparm)) then
         ii=ii+1
diff --git a/source/wham/src-M/oligomer.F b/source/wham/src-M/oligomer.F
new file mode 100644 (file)
index 0000000..017bafc
--- /dev/null
@@ -0,0 +1,94 @@
+      subroutine oligomer
+      implicit none
+      include "DIMENSIONS"
+      include "COMMON.CHAIN"
+      include "COMMON.INTERACT"
+      include "COMMON.IOUNITS"
+      integer nchain,i,ii,ipi,ipj,ipmin,j,jmin,k,ix,iy,iz,
+     &  ixmin,iymin,izmin
+      logical newchain 
+      integer ichain(2,20),iper(20),iaux
+      double precision dchain,dchainmin,cmchain(3,20)
+      nchain=1
+      newchain=.false.
+      ichain(1,nchain)=1
+      do i=2,nres
+        if (itype(i).eq.ntyp1) then
+          if (.not.newchain) then
+            ichain(2,nchain)=i
+            if (i.lt.nres) nchain=nchain+1
+            newchain=.true.
+          else 
+            newchain=.false.
+            ichain(1,nchain)=i
+          endif
+        endif   
+      enddo
+#ifdef DEBUG
+      write (iout,*) "Chains"
+      do i=1,nchain
+        write (iout,*) i,ichain(1,i),ichain(2,i)
+      enddo
+#endif
+      cmchain=0.0d0
+      do i=1,nchain
+        ii=0
+        do j=ichain(1,i),ichain(2,i)
+          if (itype(j).eq.ntyp1) cycle
+          ii=ii+1
+          do k=1,3
+            cmchain(k,i)=cmchain(k,i)+c(k,j)
+          enddo
+        enddo
+        do k=1,3
+          cmchain(k,i)=cmchain(k,i)/ii
+        enddo
+      enddo
+      do i=1,nchain
+        iper(i)=i
+      enddo
+      do i=1,nchain-1
+        dchainmin=1.0d10
+        do j=i+1,nchain
+          ipi=iper(i)
+          ipj=iper(j)
+          do ix=-1,1
+            do iy=-1,1
+              do iz=-1,1
+                dchain=(cmchain(1,ipj)-cmchain(1,ipi)+ix*boxxsize)**2+
+     &                 (cmchain(2,ipj)-cmchain(2,ipi)+iy*boxysize)**2+
+     &                 (cmchain(3,ipj)-cmchain(3,ipi)+iz*boxzsize)**2
+c                write (iout,*) "i",i," ipi",ipi," j",j," ipj",ipj," d",
+c     &               dsqrt(dchain)," dmin",dsqrt(dchainmin)," jmin",jmin
+                if (dchain.lt.dchainmin) then
+                  dchainmin=dchain
+                  ixmin=ix
+                  iymin=iy
+                  izmin=iz
+                  jmin=j
+                endif
+              enddo
+            enddo
+          enddo
+        enddo
+        if (ixmin.eq.0 .and. iymin.eq.0 .and. izmin.eq.0) cycle
+        ipj=iper(jmin)
+        cmchain(1,ipj)=cmchain(1,ipj)+ixmin*boxxsize
+        cmchain(2,ipj)=cmchain(2,ipj)+iymin*boxysize
+        cmchain(3,ipj)=cmchain(3,ipj)+izmin*boxzsize
+        do k=ichain(1,ipj),ichain(2,ipj)
+          c(1,k)=c(1,k)+ixmin*boxxsize
+          c(2,k)=c(2,k)+iymin*boxysize
+          c(3,k)=c(3,k)+izmin*boxzsize
+          c(1,k+nres)=c(1,k+nres)+ixmin*boxxsize
+          c(2,k+nres)=c(2,k+nres)+iymin*boxysize
+          c(3,k+nres)=c(3,k+nres)+izmin*boxzsize
+        enddo
+c        write (iout,*) "jmin",jmin," ipj",ipj,
+c     &   " ixmin",ixmin," iymin",iymin," izmin",izmin
+        iaux=iper(i+1)
+        iper(i+1)=iper(jmin)
+        iper(jmin)=iaux
+      enddo
+      return
+      end