subroutine oligomer implicit none include "DIMENSIONS" include "COMMON.CHAIN" include "COMMON.INTERACT" include "COMMON.IOUNITS" integer i,j,k,l,ichain,jchain,ii,lshift(maxchain),lmin(maxchain), & nrestot double precision sumodl,sumodl_min,shift_coord(3,maxchain), & boxsize(3),cm(3,maxchain),ccm(3),ccm_prev(3) logical previous,change boxsize(1)=boxxsize boxsize(2)=boxysize boxsize(3)=boxzsize nrestot=0 do ichain=1,nchain c print *,"ichain",ichain nrestot=nrestot+chain_length(ichain) do j=1,3 cm(j,ichain)=0.0d0 enddo do i=chain_border(1,ichain),chain_border(2,ichain) c print *,i,c(:,i) do j=1,3 cm(j,ichain)=cm(j,ichain)+c(j,i) enddo enddo do j=1,3 cm(j,ichain)=cm(j,ichain)/chain_length(ichain) enddo enddo #ifdef DEBUG write (iout,*) "CM before shift" do i=1,nchain write (iout,'(i5,3f10.5)') i,(cm(j,i),j=1,3) enddo #endif DEBUG do j=1,3 lmin=0 sumodl_min=1.0d10 do i=0,3**nchain-1 ii=i do ichain=1,nchain lshift(ichain)=mod(ii,3)-1 ii=ii/3 enddo #ifdef DEBUG write (iout,'(10i3)') (lshift(ichain),ichain=1,nchain) #endif sumodl=0.0d0 do ichain=1,nchain do jchain=1,ichain-1 sumodl=sumodl+dabs(cm(j,ichain)+lshift(ichain)*boxsize(j) & -cm(j,jchain)-lshift(jchain)*boxsize(j)) enddo enddo #ifdef DEBUG write(iout,*) "sumodl",sumodl," sumodl_min",sumodl_min #endif if (sumodl.lt.sumodl_min) then sumodl_min=sumodl lmin=lshift endif enddo do ichain=1,nchain cm(j,ichain)=cm(j,ichain)+lmin(ichain)*boxsize(j) shift_coord(j,ichain)=lmin(ichain)*boxsize(j) enddo enddo do ichain=1,nchain do j=1,3 ccm(j)=ccm(j)+cm(j,ichain)*chain_length(ichain) enddo enddo do j=1,3 ccm(j)=ccm(j)/nrestot enddo #ifdef DEBUG write (iout,'(a,3f10.5)') "ccm ",(ccm(j),j=1,3) write (iout,'(a,3f10.5)') "ccm_prev",(ccm_prev(j),j=1,3) write (iout,'(a,3f10.5)') "diff ",(ccm(j)-ccm_prev(j),j=1,3) write (iout,*) "shift_coord" do i=1,nchain write (iout,'(i5,3f10.5)') i,(shift_coord(j,i),j=1,3) enddo #endif; #ifdef PREVIOUS c c Shift the system by box dimensions so that its CM be closest to the c CM of the pevious snapshot c if (previous) then do j=1,3 change=.true. do while (change) if (ccm(j)-ccm_prev(j).gt.boxsize(j)/2) then do ichain=1,nchain shift_coord(j,ichain)=shift_coord(j,ichain)-boxsize(j) enddo ccm(j)=ccm(j)-boxsize(j) else if (ccm(j)-ccm_prev(j).lt.-boxsize(j)/2) then do ichain=1,nchain shift_coord(j,ichain)=shift_coord(j,ichain)+boxsize(j) enddo ccm(j)=ccm(j)+boxsize(j) else change=.false. endif enddo enddo do j=1,3 ccm_prev(j)=ccm(j) enddo else previous=.true. do j=1,3 ccm_prev(j)=ccm(j) enddo endif #else c Simply shift the CM of the whole oligomer to the origin of the c coordinate system do ichain=1,nchain do j=1,3 shift_coord(j,ichain)=shift_coord(j,ichain)-ccm(j) enddo enddo #endif #ifdef DEBUG write (iout,*) "shift_coord" do i=1,nchain write (iout,'(i5,3f10.5)') i,(shift_coord(j,i),j=1,3) enddo #endif do ichain=1,nchain do i=chain_border1(1,ichain),chain_border1(2,ichain) do j=1,3 c(j,i)=c(j,i)+shift_coord(j,ichain) c(j,i+nres)=c(j,i+nres)+shift_coord(j,ichain) enddo enddo enddo return end