2 C Program to cluster united-residue MCM results.
9 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
12 include 'COMMON.TIME1'
13 include 'COMMON.INTERACT'
14 include 'COMMON.NAMES'
16 include 'COMMON.HEADER'
17 include 'COMMON.CONTROL'
18 include 'COMMON.CHAIN'
20 include 'COMMON.CLUSTER'
21 include 'COMMON.IOUNITS'
23 logical printang(max_cut)
24 integer printpdb(max_cut)
25 integer printmol2(max_cut)
27 REAL CRIT(maxconf),MEMBR(maxconf)
28 REAL CRITVAL(maxconf-1)
29 INTEGER IA(maxconf),IB(maxconf)
30 INTEGER ICLASS(maxconf,maxconf-1),HVALS(maxconf-1)
31 INTEGER IORDER(maxconf-1),HEIGHT(maxconf-1)
34 DIMENSION NN(maxconf),DISNN(maxconf)
36 integer i,j,k,l,m,n,len,lev,idum,ii,ind,ioffset,jj,icut,ncon,
38 double precision t1,t2,tcpu,difconf
40 double precision varia(maxvar)
41 double precision hrtime,mintime,sectime
45 call MPI_Init( IERROR )
46 call MPI_Comm_rank( MPI_COMM_WORLD, me, IERROR )
47 call MPI_Comm_size( MPI_COMM_WORLD, nprocs, IERROR )
50 write(iout,*) "SEVERE ERROR - Can't initialize MPI."
51 call mpi_finalize(ierror)
54 if (nprocs.gt.MaxProcs+1) then
55 write (2,*) "Error - too many processors",
57 write (2,*) "Increase MaxProcs and recompile"
58 call MPI_Finalize(IERROR)
68 c if (refstr) call read_ref_structure(*30)
76 print *,'MAIN: nnt=',nnt,' nct=',nct
78 write (iout,*) "1,4 SCSC repulsive interactions sacled down by 10"
84 IF (RCUTOFF(I).LT.0.0) THEN
85 RCUTOFF(I)=ABS(RCUTOFF(I))
91 write (iout,*) 'Number of cutoffs:',NCUT
92 write (iout,*) 'Cutoff values:'
94 WRITE(IOUT,'(8HRCUTOFF(,I2,2H)=,F8.1,2i2)')ICUT,RCUTOFF(ICUT),
95 & printpdb(icut),printmol2(icut)
103 call read_coords(ncon,*20)
104 write (iout,*) 'from read_coords: ncon',ncon
106 write (iout,*) "nT",nT
108 write (iout,*) "iT",iT
110 call work_partition(.true.,ncon)
113 call probabl(iT,ncon_work,ncon,*20)
115 if (ncon_work.lt.2) then
116 write (iout,*) "Too few conformations; clustering skipped"
120 ndis=ncon_work*(ncon_work-1)/2
121 call work_partition(.true.,ndis)
127 WRITE (iout,'(A80)') TITEL
130 C CALCULATE DISTANCES
132 call daread_ccoords(1,ncon_work)
135 if (mod(i,100).eq.0) print *,'Calculating RMS i=',i
138 c(l,k)=allcart(l,k,i)
147 IND=IOFFSET(NCON_work,I,J)
149 if (ind.ge.indstart(me) .and. ind.le.indend(me)) then
152 DISS(IND1)=DIFCONF(I,J)
153 c write (iout,'(2i4,i10,f10.5)') i,j,ind,DISS(IND)
160 WRITE (iout,'(/a,1pe14.5,a/)')
161 & 'Time for distance calculation:',T2-T1,' sec.'
163 PRINT '(a)','End of distance computation'
166 call MPI_Gatherv(diss(1),scount(me),MPI_REAL,diss(1),
167 & scount(0),idispl(0),MPI_REAL,Master,MPI_COMM_WORLD, IERROR)
168 if (me.eq.master) then
170 open(80,file='/tmp/distance',form='unformatted')
177 IND=IOFFSET(NCON,I,J)
178 write (jrms,'(2i5,2f10.5)') i,j,diss(IND),
179 & energy(j)-energy(i)
184 C Print out the RMS deviation matrix.
186 if (print_dist) CALL DISTOUT(NCON_work)
188 C call hierarchical clustering HC from F. Murtagh
192 write(iout,*) "-------------------------------------------"
193 write(iout,*) "HIERARCHICAL CLUSTERING using"
195 write(iout,*) "WARD'S MINIMUM VARIANCE METHOD"
196 elseif (iopt.eq.2) then
197 write(iout,*) "SINGLE LINK METHOD"
198 elseif (iopt.eq.3) then
199 write(iout,*) "COMPLETE LINK METHOD"
200 elseif (iopt.eq.4) then
201 write(iout,*) "AVERAGE LINK (OR GROUP AVERAGE) METHOD"
202 elseif (iopt.eq.5) then
203 write(iout,*) "MCQUITTY'S METHOD"
204 elseif (iopt.eq.6) then
205 write(iout,*) "MEDIAN (GOWER'S) METHOD"
206 elseif (iopt.eq.7) then
207 write(iout,*) "CENTROID METHOD"
209 write(iout,*) "IOPT=",iopt," IS INVALID, use 1-7"
210 write(*,*) "IOPT=",iopt," IS INVALID, use 1-7"
214 write(iout,*) "hc.f by F. Murtagh, ESA/ESO/STECF, Garching"
215 write(iout,*) "February 1986"
216 write(iout,*) "References:"
217 write(iout,*) "1. Multidimensional clustering algorithms"
218 write(iout,*) " Fionn Murtagh"
219 write(iout,*) " Vienna : Physica-Verlag, 1985."
220 write(iout,*) "2. Multivariate data analysis"
221 write(iout,*) " Fionn Murtagh and Andre Heck"
222 write(iout,*) " Kluwer Academic Publishers, 1987"
223 write(iout,*) "-------------------------------------------"
227 write (iout,*) "The TOTFREE array"
229 write (iout,'(2i5,f10.5)') i,list_conf(i),totfree(i)
233 CALL HC(N,M,LEN,IOPT,IA,IB,CRIT,MEMBR,NN,DISNN,FLAG,DISS)
235 write (iout,*) "n",n," ncon_work",ncon_work," lev",lev
237 write (iout,*) "Too few conformations to cluster."
240 CALL HCASS(N,IA,IB,CRIT,LEV,ICLASS,HVALS,IORDER,CRITVAL,HEIGHT)
241 c CALL HCDEN(LEV,IORDER,HEIGHT,CRITVAL)
250 licz(iclass(j,i))=licz(iclass(j,i))+1
251 nconf(iclass(j,i),licz(iclass(j,i)))=j
252 c write (iout,*) i,j,iclass(j,i),licz(iclass(j,i)),
253 c & nconf(iclass(j,i),licz(iclass(j,i)))
259 IF (HEIGHT(L).EQ.IDUM) GOTO 190
262 write(IOUT,*) "i+1",i+1," idum",idum," critval",CRITVAL(IDUM),
263 & " icut",icut," cutoff",rcutoff(icut)
264 IF (CRITVAL(IDUM).LT. RCUTOFF(ICUT)) THEN
265 WRITE (iout,'(/a,f10.5)') 'AT CUTOFF:',rcutoff(icut)
266 write (iout,'(a,f8.2)') 'Maximum distance found:',
268 CALL SRTCLUST(ICUT,ncon_work,iT)
270 CALL WRTCLUST(ncon_work,ICUT,PRINTANG,PRINTPDB,PRINTMOL2,iT)
272 if (icut.gt.ncut) goto 191
279 licz(iclass(j,i))=licz(iclass(j,i))+1
280 nconf(iclass(j,i),licz(iclass(j,i)))=j
281 c write (iout,*) i,j,iclass(j,i),licz(iclass(j,i)),
282 c & nconf(iclass(j,i),licz(iclass(j,i)))
283 cd print *,j,iclass(j,i),
284 cd & licz(iclass(j,i)),nconf(iclass(j,i),licz(iclass(j,i)))
295 WRITE (iout,'(/a,1pe14.5,a/)')
296 & 'Total time for clustering:',T2-T1,' sec.'
304 close(icbase,status="delete")
306 call MPI_Finalize(MPI_COMM_WORLD,IERROR)
308 stop '********** Program terminated normally.'
309 20 write (iout,*) "Error reading coordinates"
311 call MPI_Finalize(MPI_COMM_WORLD,IERROR)
314 30 write (iout,*) "Error reading reference structure"
316 call MPI_Finalize(MPI_COMM_WORLD,IERROR)
320 c---------------------------------------------------------------------------
321 double precision function difconf(icon,jcon)
324 include 'sizesclu.dat'
325 include 'COMMON.CONTROL'
326 include 'COMMON.CLUSTER'
327 include 'COMMON.CHAIN'
328 include 'COMMON.INTERACT'
330 include 'COMMON.IOUNITS'
332 double precision przes(3),obrot(3,3)
333 double precision xx(3,maxres2),yy(3,maxres2)
334 integer i,ii,j,icon,jcon
341 xx(j,ii)=allcart(j,i,jcon)
346 c if (itype(i).ne.10) then
349 xx(j,ii)=allcart(j,i+nres,jcon)
350 yy(j,ii)=cref(j,i+nres)
354 call fitsq(rms,xx(1,1),yy(1,1),ii,przes,obrot,non_conv)
358 c(j,i)=allcart(j,i,jcon)
361 call fitsq(rms,c(1,nstart),cref(1,nstart),nend-nstart+1,przes,
365 print *,'error, rms^2 = ',rms,icon,jcon
368 if (non_conv) print *,non_conv,icon,jcon
372 C------------------------------------------------------------------------------
373 subroutine distout(ncon)
376 include 'sizesclu.dat'
379 include 'COMMON.IOUNITS'
380 include 'COMMON.CLUSTER'
381 integer i,j,k,jlim,jlim1,nlim,ind,ioffset
384 write (iout,'(a)') 'The distance matrix'
386 nlim=min0(i+ncol-1,ncon)
387 write (iout,1000) (k,k=i,nlim)
388 write (iout,'(8h--------,10a)') ('-------',k=i,nlim)
389 1000 format (/8x,10(i4,3x))
390 1020 format (/1x,80(1h-)/)
401 IND=IOFFSET(NCON,j,k)
403 IND=IOFFSET(NCON,k,j)
406 write (iout,1010) j,(b(k),k=1,jlim-i+1)
409 1010 format (i5,3x,10(f6.2,1x))