34b7be0bbcaabeb3b03739d01928bac6bd123519
[unres.git] / source / wham / src-HCD / oligomer.F
1       subroutine oligomer
2       implicit none
3       include "DIMENSIONS"
4       include "COMMON.CHAIN"
5       include "COMMON.INTERACT"
6       include "COMMON.IOUNITS"
7       integer i,ii,ipi,ipj,ipmin,j,jmin,k,ix,iy,iz,
8      &  ixmin,iymin,izmin,ir_start,ir_end
9       integer iper(maxchain),iaux
10       double precision dchain,dchainmin,cmchain(3,20)
11       cmchain=0.0d0
12       do i=1,nchain
13         ii=0
14         do j=chain_border(1,i),chain_border(2,i)
15           if (itype(j).eq.ntyp1) cycle
16           ii=ii+1
17           do k=1,3
18             cmchain(k,i)=cmchain(k,i)+c(k,j)
19           enddo
20         enddo
21         do k=1,3
22           cmchain(k,i)=cmchain(k,i)/ii
23         enddo
24       enddo
25       do i=1,nchain
26         iper(i)=i
27       enddo
28       do i=1,nchain-1
29         dchainmin=1.0d10
30         do j=i+1,nchain
31           ipi=iper(i)
32           ipj=iper(j)
33           do ix=-1,1
34             do iy=-1,1
35               do iz=-1,1
36                 dchain=(cmchain(1,ipj)-cmchain(1,ipi)+ix*boxxsize)**2+
37      &                 (cmchain(2,ipj)-cmchain(2,ipi)+iy*boxysize)**2+
38      &                 (cmchain(3,ipj)-cmchain(3,ipi)+iz*boxzsize)**2
39 c                write (iout,*) "i",i," ipi",ipi," j",j," ipj",ipj," d",
40 c     &               dsqrt(dchain)," dmin",dsqrt(dchainmin)," jmin",jmin
41                 if (dchain.lt.dchainmin) then
42                   dchainmin=dchain
43                   ixmin=ix
44                   iymin=iy
45                   izmin=iz
46                   jmin=j
47                 endif
48               enddo
49             enddo
50           enddo
51         enddo
52         if (ixmin.eq.0 .and. iymin.eq.0 .and. izmin.eq.0) cycle
53         ipj=iper(jmin)
54         cmchain(1,ipj)=cmchain(1,ipj)+ixmin*boxxsize
55         cmchain(2,ipj)=cmchain(2,ipj)+iymin*boxysize
56         cmchain(3,ipj)=cmchain(3,ipj)+izmin*boxzsize
57         ir_start=chain_border(1,ipj)
58         if (ir_start.gt.1) ir_start=ir_start-1 
59         ir_end=chain_border(2,ipj)
60         if (ir_end.lt.nres) ir_end=ir_end+1
61         do k=ir_start,ir_end
62           c(1,k)=c(1,k)+ixmin*boxxsize
63           c(2,k)=c(2,k)+iymin*boxysize
64           c(3,k)=c(3,k)+izmin*boxzsize
65           c(1,k+nres)=c(1,k+nres)+ixmin*boxxsize
66           c(2,k+nres)=c(2,k+nres)+iymin*boxysize
67           c(3,k+nres)=c(3,k+nres)+izmin*boxzsize
68         enddo
69 c        write (iout,*) "jmin",jmin," ipj",ipj,
70 c     &   " ixmin",ixmin," iymin",iymin," izmin",izmin
71         iaux=iper(i+1)
72         iper(i+1)=iper(jmin)
73         iper(jmin)=iaux
74       enddo
75       return
76       end