From ba8f40a2336651a42eb42e9d3379d62aa11bb5d7 Mon Sep 17 00:00:00 2001 From: Cezary Czaplewski Date: Sat, 12 Aug 2017 04:42:03 +0200 Subject: [PATCH] multichain cluster can select number of clusters instead of cutoff --- source/cluster/wham/src-M/COMMON.CLUSTER | 4 +-- source/cluster/wham/src-M/main_clust.F | 54 +++++++++++++++++++++++------- source/cluster/wham/src-M/readrtns.F | 6 ++-- source/cluster/wham/src-M/wrtclust.f | 17 ++++++---- 4 files changed, 58 insertions(+), 23 deletions(-) diff --git a/source/cluster/wham/src-M/COMMON.CLUSTER b/source/cluster/wham/src-M/COMMON.CLUSTER index 4fedc52..a82948d 100644 --- a/source/cluster/wham/src-M/COMMON.CLUSTER +++ b/source/cluster/wham/src-M/COMMON.CLUSTER @@ -4,11 +4,11 @@ real*4 diss,allcart double precision enetb,entfac,totfree,energy,rmstb integer ncut,ngr,licz,nconf,iass,icc,mult,list_conf, - & nss_all,ihpb_all,jhpb_all,iass_tot,iscore,nprop + & nss_all,ihpb_all,jhpb_all,iass_tot,iscore,nprop,nclust common /clu/ diss(maxdist),energy(0:maxconf), & enetb(max_ene,maxconf),ecut, & entfac(maxconf),totfree(0:maxconf),totfree_gr(maxgr), - & rcutoff(max_cut+1),ncut,min_var,tree,plot_tree,lgrp + & rcutoff(max_cut+1),ncut,nclust,min_var,tree,plot_tree,lgrp common /clu1/ ngr,licz(maxgr),nconf(maxgr,maxingr),iass(maxgr), & iass_tot(maxgr,max_cut),list_conf(maxconf) common /alles/ allcart(3,maxres2,maxconf),rmstb(maxconf), diff --git a/source/cluster/wham/src-M/main_clust.F b/source/cluster/wham/src-M/main_clust.F index a3f8094..9835c01 100644 --- a/source/cluster/wham/src-M/main_clust.F +++ b/source/cluster/wham/src-M/main_clust.F @@ -23,7 +23,7 @@ C logical printang(max_cut) integer printpdb(max_cut) integer printmol2(max_cut) - character*240 lineh + character*240 lineh,scrachdir2d REAL CRIT(maxconf),MEMBR(maxconf) REAL CRITVAL(maxconf-1) INTEGER IA(maxconf),IB(maxconf) @@ -34,12 +34,13 @@ C DIMENSION NN(maxconf),DISNN(maxconf) LOGICAL FLAG(maxconf) integer i,j,k,l,m,n,len,lev,idum,ii,ind,ioffset,jj,icut,ncon, - & it,ncon_work,ind1,kkk, ijk + & it,ncon_work,ind1,kkk, ijk,is,ie,ilen double precision t1,t2,tcpu,difconf double precision varia(maxvar) double precision hrtime,mintime,sectime logical eof + external ilen #ifdef MPI call MPI_Init( IERROR ) call MPI_Comm_rank( MPI_COMM_WORLD, me, IERROR ) @@ -79,6 +80,12 @@ c write (iout,*) "after permut" c call flush(iout) print *,'MAIN: nnt=',nnt,' nct=',nct + if (nclust.gt.0) then + PRINTANG(1)=.TRUE. + PRINTPDB(1)=outpdb + printmol2(1)=outmol2 + ncut=0 + else DO I=1,NCUT PRINTANG(I)=.FALSE. PRINTPDB(I)=0 @@ -90,12 +97,21 @@ c call flush(iout) printmol2(i)=outmol2 ENDIF ENDDO + endif + if (ncut.gt.0) then write (iout,*) 'Number of cutoffs:',NCUT write (iout,*) 'Cutoff values:' DO ICUT=1,NCUT WRITE(IOUT,'(8HRCUTOFF(,I2,2H)=,F8.1,2i2)')ICUT,RCUTOFF(ICUT), & printpdb(icut),printmol2(icut) ENDDO + else if (nclust.gt.0) then + write (iout,'("Number of clusters requested",i5)') nclust + else + if (me.eq.Master) + & write (iout,*) "ERROR: Either nclust or ncut must be >0" + stop + endif DO I=1,NRES-3 MULT(I)=1 ENDDO @@ -180,7 +196,8 @@ c write (iout,'(2i4,i10,f10.5)') i,j,ind,DISS(IND) & scount(0),idispl(0),MPI_REAL,Master,MPI_COMM_WORLD, IERROR) if (me.eq.master) then #endif - open(80,file='/tmp/distance',form='unformatted') + scrachdir2d=scratchdir(:ilen(scratchdir))//'distance' + open(80,file=scrachdir2d,form='unformatted') do i=1,ndis write(80) diss(i) enddo @@ -253,29 +270,39 @@ C CALL HCASS(N,IA,IB,CRIT,LEV,ICLASS,HVALS,IORDER,CRITVAL,HEIGHT) c CALL HCDEN(LEV,IORDER,HEIGHT,CRITVAL) +c 3/3/16 AL: added explicit number of cluters + if (nclust.gt.0) then + is=nclust-1 + ie=nclust-1 + icut=1 + else + is=1 + ie=lev-1 + endif do i=1,maxgr licz(i)=0 enddo icut=1 - i=1 - NGR=i+1 + i=is + NGR=is+1 do j=1,n licz(iclass(j,i))=licz(iclass(j,i))+1 nconf(iclass(j,i),licz(iclass(j,i)))=j c write (iout,*) i,j,iclass(j,i),licz(iclass(j,i)), c & nconf(iclass(j,i),licz(iclass(j,i))) enddo - do i=1,lev-1 - +c do i=1,lev-1 + do i=is,ie idum=lev-i DO L=1,LEV IF (HEIGHT(L).EQ.IDUM) GOTO 190 ENDDO 190 IDUM=L - write(IOUT,*) "i+1",i+1," idum",idum," critval",CRITVAL(IDUM), - & " icut",icut," cutoff",rcutoff(icut) - IF (CRITVAL(IDUM).LT. RCUTOFF(ICUT)) THEN - WRITE (iout,'(/a,f10.5)') 'AT CUTOFF:',rcutoff(icut) +c write(IOUT,*) "i+1",i+1," idum",idum," critval",CRITVAL(IDUM), +c & " icut",icut," cutoff",rcutoff(icut) + IF (nclust.gt.0.or.CRITVAL(IDUM).LT. RCUTOFF(ICUT)) THEN + if (nclust.le.0) + & WRITE (iout,'(/a,f10.5)') 'AT CUTOFF:',rcutoff(icut) write (iout,'(a,f8.2)') 'Maximum distance found:', & CRITVAL(IDUM) CALL SRTCLUST(ICUT,ncon_work,iT) @@ -288,9 +315,10 @@ c & nconf(iclass(j,i),licz(iclass(j,i))) do l=1,maxgr licz(l)=0 enddo + ii=i-is+1 do j=1,n - licz(iclass(j,i))=licz(iclass(j,i))+1 - nconf(iclass(j,i),licz(iclass(j,i)))=j + licz(iclass(j,ii))=licz(iclass(j,ii))+1 + nconf(iclass(j,ii),licz(iclass(j,ii)))=j c write (iout,*) i,j,iclass(j,i),licz(iclass(j,i)), c & nconf(iclass(j,i),licz(iclass(j,i))) cd print *,j,iclass(j,i), diff --git a/source/cluster/wham/src-M/readrtns.F b/source/cluster/wham/src-M/readrtns.F index ae8f8a8..a3feaa5 100644 --- a/source/cluster/wham/src-M/readrtns.F +++ b/source/cluster/wham/src-M/readrtns.F @@ -91,7 +91,8 @@ C long axis of side chain min_var=(index(controlcard,'MINVAR').gt.0) plot_tree=(index(controlcard,'PLOT_TREE').gt.0) punch_dist=(index(controlcard,'PUNCH_DIST').gt.0) - call readi(controlcard,'NCUT',ncut,1) + call readi(controlcard,'NCUT',ncut,0) + call readi(controlcard,'NCLUST',nclust,5) call readi(controlcard,'SYM',symetr,1) write (iout,*) 'sym', symetr call readi(controlcard,'NSTART',nstart,0) @@ -102,7 +103,8 @@ C long axis of side chain lgrp=(index(controlcard,'LGRP').gt.0) caonly=(index(controlcard,'CA_ONLY').gt.0) print_dist=(index(controlcard,'PRINT_DIST').gt.0) - call multreada(controlcard,'CUTOFF',rcutoff,ncut,-1.0d0) + if (ncut.gt.0) + & call multreada(controlcard,'CUTOFF',rcutoff,ncut,-1.0d0) call readi(controlcard,'IOPT',iopt,2) lside = index(controlcard,"SIDE").gt.0 efree = index(controlcard,"EFREE").gt.0 diff --git a/source/cluster/wham/src-M/wrtclust.f b/source/cluster/wham/src-M/wrtclust.f index f2f3eb7..917d4f3 100644 --- a/source/cluster/wham/src-M/wrtclust.f +++ b/source/cluster/wham/src-M/wrtclust.f @@ -84,12 +84,14 @@ C 12/8/93 Estimation of "diameters" of the subsequent families. ave_dim=0.0 amax_dim=0.0 c write (iout,*) "ecut",ecut + emin=totfree(nconf(igr,1)) +c write (2,*) "emin",emin," ecut",ecut do i=2,licz(igr) ii=nconf(igr,i) if (totfree(ii)-emin .gt. ecut) goto 10 do j=1,i-1 jj=nconf(igr,j) - if (jj.eq.1) exit +c if (jj.eq.1) exit if (ii.lt.jj) then ind=ioffset(ncon,ii,jj) else @@ -112,9 +114,12 @@ c & list_conf(jj),curr_dist & '; average distance in the family:',ave_dim rmsave(igr)=0.0d0 qpart=0.0d0 + emin=totfree(nconf(igr,1)) do i=1,licz(igr) icon=nconf(igr,i) - boltz=dexp(-totfree(icon)) + boltz=dexp(-totfree(icon)+emin) +c write (2,*) "igr",igr," i",i," icon",icon," totfree", +c & totfree(icon)," emin",emin," boltz",boltz," rms",rmstb(icon) rmsave(igr)=rmsave(igr)+boltz*rmstb(icon) qpart=qpart+boltz enddo @@ -191,8 +196,8 @@ c Write conformations of the family i to PDB files c write (iout,*) i,ncon_out,nconf(i,ncon_out), c & totfree(nconf(i,ncon_out)),emin1,ecut enddo - write (iout,*) "ncon_out",ncon_out - call flush(iout) +c write (iout,*) "ncon_out",ncon_out +c call flush(iout) do j=1,nres tempfac(1,j)=5.0d0 tempfac(2,j)=5.0d0 @@ -263,9 +268,9 @@ c create InsightII command file for their displaying in different colors 200 FORMAT (/'FAMILY ',I4,' WITH TOTAL FREE ENERGY',1pE15.5, & ' CONTAINS ',I4,' CONFORMATION(S): ') c 300 FORMAT ( 8(I4,F6.1)) - 300 FORMAT (5(I4,1pe12.3)) + 300 FORMAT (5(I6,1pe12.3)) 400 FORMAT (//'ASSIGNMENT OF CONSECUTIVE CONFORMATIONS TO FAMILIES:') - 500 FORMAT (8(2I4,2X)) + 500 FORMAT (8(I6,I4,2X)) 600 FORMAT ('REMARK FAMILY',I4,' CONFORMATION',I4,' ENERGY ',E15.6) RETURN END -- 1.7.9.5