- subroutine thread_seq
-C Thread the sequence through a database of known structures
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
-#ifdef MPI
- include 'mpif.h'
-#endif
- include 'COMMON.CONTROL'
- include 'COMMON.CHAIN'
- include 'COMMON.DBASE'
- include 'COMMON.INTERACT'
- include 'COMMON.VAR'
- include 'COMMON.THREAD'
- include 'COMMON.FFIELD'
- include 'COMMON.SBRIDGE'
- include 'COMMON.HEADER'
- include 'COMMON.IOUNITS'
- include 'COMMON.TIME1'
- include 'COMMON.CONTACTS'
- include 'COMMON.MCM'
- include 'COMMON.NAMES'
-#ifdef MPI
- include 'COMMON.INFO'
- integer ThreadId,ThreadType,Kwita
-#endif
- double precision varia(maxvar)
- double precision przes(3),obr(3,3)
- double precision time_for_thread
- logical found_pattern,non_conv
- character*32 head_pdb
- double precision energia(0:n_ene)
- n_ene_comp=nprint_ene
-C
-C Body
-C
-#ifdef MPI
- if (me.eq.king) then
- do i=1,nctasks
- nsave_part(i)=0
- enddo
- endif
- nacc_tot=0
-#endif
- Kwita=0
- close(igeom)
- close(ipdb)
- close(istat)
- do i=1,maxthread
- do j=1,14
- ener0(j,i)=0.0D0
- ener(j,i)=0.0D0
- enddo
- enddo
- nres0=nct-nnt+1
- ave_time_for_thread=0.0D0
- max_time_for_thread=0.0D0
-cd print *,'nthread=',nthread,' nseq=',nseq,' nres0=',nres0
- nthread=nexcl+nthread
- do ithread=1,nthread
- found_pattern=.false.
- itrial=0
- do while (.not.found_pattern)
- itrial=itrial+1
- if (itrial.gt.1000) then
- write (iout,'(/a/)') 'Too many attempts to find pattern.'
- nthread=ithread-1
-#ifdef MPI
- call recv_stop_sig(Kwita)
- call send_stop_sig(-3)
-#endif
- goto 777
- endif
-C Find long enough chain in the database
- ii=iran_num(1,nseq)
- nres_t=nres_base(1,ii)
-C Select the starting position to thread.
- print *,'nseq',nseq,' ii=',ii,' nres_t=',
- & nres_t,' nres0=',nres0
- if (nres_t.ge.nres0) then
- ist=iran_num(0,nres_t-nres0)
-#ifdef MPI
- if (Kwita.eq.0) call recv_stop_sig(Kwita)
- if (Kwita.lt.0) then
- write (iout,*) 'Stop signal received. Terminating.'
- write (*,*) 'Stop signal received. Terminating.'
- nthread=ithread-1
- write (*,*) 'ithread=',ithread,' nthread=',nthread
- goto 777
- endif
- call pattern_receive
-#endif
- do i=1,nexcl
- if (iexam(1,i).eq.ii .and. iexam(2,i).eq.ist) goto 10
- enddo
- found_pattern=.true.
- endif
-C If this point is reached, the pattern has not yet been examined.
- 10 continue
-c print *,'found_pattern:',found_pattern
- enddo
- nexcl=nexcl+1
- iexam(1,nexcl)=ii
- iexam(2,nexcl)=ist
-#ifdef MPI
- if (Kwita.eq.0) call recv_stop_sig(Kwita)
- if (Kwita.lt.0) then
- write (iout,*) 'Stop signal received. Terminating.'
- nthread=ithread-1
- write (*,*) 'ithread=',ithread,' nthread=',nthread
- goto 777
- endif
- call pattern_send
-#endif
- ipatt(1,ithread)=ii
- ipatt(2,ithread)=ist
-#ifdef MPI
- write (iout,'(/80(1h*)/a,i4,a,i5,2a,i3,a,i3,a,i3/)')
- & 'Processor:',me,' Attempt:',ithread,
- & ' pattern: ',str_nam(ii),nres_base(2,ii),':',nres_base(3,ii),
- & ' start at res.',ist+1
- write (*,'(a,i4,a,i5,2a,i3,a,i3,a,i3)') 'Processor:',me,
- & ' Attempt:',ithread,
- & ' pattern: ',str_nam(ii),nres_base(2,ii),':',nres_base(3,ii),
- & ' start at res.',ist+1
-#else
- write (iout,'(/80(1h*)/a,i5,2a,i3,a,i3,a,i3/)')
- & 'Attempt:',ithread,
- & ' pattern: ',str_nam(ii),nres_base(2,ii),':',nres_base(3,ii),
- & ' start at res.',ist+1
- write (*,'(a,i5,2a,i3,a,i3,a,i3)')
- & 'Attempt:',ithread,
- & ' pattern: ',str_nam(ii),nres_base(2,ii),':',nres_base(3,ii),
- & ' start at res.',ist+1
-#endif
- ipattern=ii
-C Copy coordinates from the database.
- ist=ist-(nnt-1)
- do i=nnt,nct
- do j=1,3
- c(j,i)=cart_base(j,i+ist,ii)
-c cref(j,i)=c(j,i)
- enddo
-cd write (iout,'(a,i4,3f10.5)') restyp(itype(i)),i,(c(j,i),j=1,3)
- enddo
-cd call fitsq(rms,c(1,nnt),cref(1,nnt),nct-nnt+1,przes,obr,
-cd non_conv)
-cd write (iout,'(a,f10.5)')
-cd & 'Initial RMS deviation from reference structure:',rms
- if (itype(nres).eq.21) then
- do j=1,3
- dcj=c(j,nres-2)-c(j,nres-3)
- c(j,nres)=c(j,nres-1)+dcj
- c(j,2*nres)=c(j,nres)
- enddo
- endif
- if (itype(1).eq.21) then
- do j=1,3
- dcj=c(j,4)-c(j,3)
- c(j,1)=c(j,2)-dcj
- c(j,nres+1)=c(j,1)
- enddo
- endif
- call int_from_cart(.false.,.false.)
-cd print *,'Exit INT_FROM_CART.'
-cd print *,'nhpb=',nhpb
- do i=nss+1,nhpb
- ii=ihpb(i)
- jj=jhpb(i)
- dhpb(i)=dist(ii,jj)
-c write (iout,'(2i5,2f10.5)') ihpb(i),jhpb(i),dhpb(i),forcon(i)
- enddo
-c stop 'End generate'
-C Generate SC conformations.
- call sc_conf
-c call intout
-#ifdef MPI
-cd print *,'Processor:',me,': exit GEN_SIDE.'
-#else
-cd print *,'Exit GEN_SIDE.'
-#endif
-C Calculate initial energy.
- call chainbuild
- call etotal(energia(0))
- etot=energia(0)
- do i=1,n_ene_comp
- ener0(i,ithread)=energia(i)
- enddo
- ener0(n_ene_comp+1,ithread)=energia(0)
- if (refstr) then
- call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
- ener0(n_ene_comp+3,ithread)=contact_fract(ncont,ncont_ref,
- & icont,icont_ref)
- ener0(n_ene_comp+2,ithread)=rms
- ener0(n_ene_comp+4,ithread)=frac
- ener0(n_ene_comp+5,ithread)=frac_nn
- endif
- ener0(n_ene_comp+3,ithread)=0.0d0
-C Minimize energy.
-#ifdef MPI
- print*,'Processor:',me,' ithread=',ithread,' Start REGULARIZE.'
-#else
- print*,'ithread=',ithread,' Start REGULARIZE.'
-#endif
- curr_tim=tcpu()
- call regularize(nct-nnt+1,etot,rms,
- & cart_base(1,ist+nnt,ipattern),iretcode)
- curr_tim1=tcpu()
- time_for_thread=curr_tim1-curr_tim
- ave_time_for_thread=
- & ((ithread-1)*ave_time_for_thread+time_for_thread)/ithread
- if (time_for_thread.gt.max_time_for_thread)
- & max_time_for_thread=time_for_thread
-#ifdef MPI
- print *,'Processor',me,': Exit REGULARIZE.'
- if (WhatsUp.eq.2) then
- write (iout,*)
- & 'Sufficient number of confs. collected. Terminating.'
- nthread=ithread-1
- goto 777
- else if (WhatsUp.eq.-1) then
- nthread=ithread-1
- write (iout,*) 'Time up in REGULARIZE. Call SEND_STOP_SIG.'
- if (Kwita.eq.0) call recv_stop_sig(Kwita)
- call send_stop_sig(-2)
- goto 777
- else if (WhatsUp.eq.-2) then
- nthread=ithread-1
- write (iout,*) 'Timeup signal received. Terminating.'
- goto 777
- else if (WhatsUp.eq.-3) then
- nthread=ithread-1
- write (iout,*) 'Error stop signal received. Terminating.'
- goto 777
- endif
-#else
- print *,'Exit REGULARIZE.'
- if (iretcode.eq.11) then
- write (iout,'(/a/)')
- &'******* Allocated time exceeded in SUMSL. The program will stop.'
- nthread=ithread-1
- goto 777
- endif
-#endif
- head_pdb=titel(:24)//':'//str_nam(ipattern)
- if (outpdb) call pdbout(etot,head_pdb,ipdb)
- if (outmol2) call mol2out(etot,head_pdb)
-c call intout
- call briefout(ithread,etot)
- link_end0=link_end
- link_end=min0(link_end,nss)
- write (iout,*) 'link_end=',link_end,' link_end0=',link_end0,
- & ' nss=',nss
- call etotal(energia(0))
-c call enerprint(energia(0))
- link_end=link_end0
-cd call chainbuild
-cd call fitsq(rms,c(1,nnt),cref(1,nnt),nct-nnt+1,przes,obr,non_conv)
-cd write (iout,'(a,f10.5)')
-cd & 'RMS deviation from reference structure:',dsqrt(rms)
- do i=1,n_ene_comp
- ener(i,ithread)=energia(i)
- enddo
- ener(n_ene_comp+1,ithread)=energia(0)
- ener(n_ene_comp+3,ithread)=rms
- if (refstr) then
- call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
- ener(n_ene_comp+2,ithread)=rms
- ener(n_ene_comp+4,ithread)=frac
- ener(n_ene_comp+5,ithread)=frac_nn
- endif
- call write_stat_thread(ithread,ipattern,ist)
-c write (istat,'(i4,2x,a8,i4,11(1pe14.5),2(0pf8.3),f8.5)')
-c & ithread,str_nam(ipattern),ist+1,(ener(k,ithread),k=1,11),
-c & (ener(k,ithread),k=12,14)
-#ifdef MPI
- if (me.eq.king) then
- nacc_tot=nacc_tot+1
- call pattern_receive
- call receive_MCM_info
- if (nacc_tot.ge.nthread) then
- write (iout,*)
- & 'Sufficient number of conformations collected nacc_tot=',
- & nacc_tot,'. Stopping other processors and terminating.'
- write (*,*)
- & 'Sufficient number of conformations collected nacc_tot=',
- & nacc_tot,'. Stopping other processors and terminating.'
- call recv_stop_sig(Kwita)
- if (Kwita.eq.0) call send_stop_sig(-1)
- nthread=ithread
- goto 777
- endif
- else
- call send_MCM_info(2)
- endif
-#endif
- if (timlim-curr_tim1-safety .lt. max_time_for_thread) then
- write (iout,'(/2a)')
- & '********** There would be not enough time for another thread. ',
- & 'The program will stop.'
- write (*,'(/2a)')
- & '********** There would be not enough time for another thread. ',
- & 'The program will stop.'
- write (iout,'(a,1pe14.4/)')
- & 'Elapsed time for last threading step: ',time_for_thread
- nthread=ithread
-#ifdef MPI
- call recv_stop_sig(Kwita)
- call send_stop_sig(-2)
-#endif
- goto 777
- else
- curr_tim=curr_tim1
- write (iout,'(a,1pe14.4)')
- & 'Elapsed time for this threading step: ',time_for_thread
- endif
-#ifdef MPI
- if (Kwita.eq.0) call recv_stop_sig(Kwita)
- if (Kwita.lt.0) then
- write (iout,*) 'Stop signal received. Terminating.'
- write (*,*) 'Stop signal received. Terminating.'
- nthread=ithread
- write (*,*) 'nthread=',nthread,' ithread=',ithread
- goto 777
- endif
-#endif
- enddo
-#ifdef MPI
- call send_stop_sig(-1)
-#endif
- 777 continue
-#ifdef MPI
-C Any messages left for me?
- call pattern_receive
- if (Kwita.eq.0) call recv_stop_sig(Kwita)
-#endif
- call write_thread_summary
-#ifdef MPI
- if (king.eq.king) then
- Kwita=1
- do while (Kwita.ne.0 .or. nacc_tot.ne.0)
- Kwita=0
- nacc_tot=0
- call recv_stop_sig(Kwita)
- call receive_MCM_info
- enddo
- do iproc=1,nprocs-1
- call receive_thread_results(iproc)
- enddo
- call write_thread_summary
- else
- call send_thread_results
- endif
-#endif
- return
- end
-c--------------------------------------------------------------------------
- subroutine write_thread_summary
-C Thread the sequence through a database of known structures
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
-#ifdef MPI
- include 'mpif.h'
-#endif
- include 'COMMON.CONTROL'
- include 'COMMON.CHAIN'
- include 'COMMON.DBASE'
- include 'COMMON.INTERACT'
- include 'COMMON.VAR'
- include 'COMMON.THREAD'
- include 'COMMON.FFIELD'
- include 'COMMON.SBRIDGE'
- include 'COMMON.HEADER'
- include 'COMMON.NAMES'
- include 'COMMON.IOUNITS'
- include 'COMMON.TIME1'
-#ifdef MPI
- include 'COMMON.INFO'
-#endif
- dimension ip(maxthread)
- double precision energia(0:n_ene)
- write (iout,'(30x,a/)')
- & ' *********** Summary threading statistics ************'
- write (iout,'(a)') 'Initial energies:'
- write (iout,'(a4,2x,a12,14a14,3a8)')
- & 'No','seq',(ename(print_order(i)),i=1,nprint_ene),'ETOT',
- & 'RMSnat','NatCONT','NNCONT','RMS'
-C Energy sort patterns
- do i=1,nthread
- ip(i)=i
- enddo
- do i=1,nthread-1
- enet=ener(n_ene-1,ip(i))
- jj=i
- do j=i+1,nthread
- if (ener(n_ene-1,ip(j)).lt.enet) then
- jj=j
- enet=ener(n_ene-1,ip(j))
- endif
- enddo
- if (jj.ne.i) then
- ipj=ip(jj)
- ip(jj)=ip(i)
- ip(i)=ipj
- endif
- enddo
- do ik=1,nthread
- i=ip(ik)
- ii=ipatt(1,i)
- ist=nres_base(2,ii)+ipatt(2,i)
- do kk=1,n_ene_comp
- energia(i)=ener0(kk,i)
- enddo
- etot=ener0(n_ene_comp+1,i)
- rmsnat=ener0(n_ene_comp+2,i)
- rms=ener0(n_ene_comp+3,i)
- frac=ener0(n_ene_comp+4,i)
- frac_nn=ener0(n_ene_comp+5,i)
-
- if (refstr) then
- write (iout,'(i4,2x,a8,i4,14(1pe14.5),0pf8.3,f8.5,f8.5,f8.3)')
- & i,str_nam(ii),ist+1,
- & (energia(print_order(kk)),kk=1,nprint_ene),
- & etot,rmsnat,frac,frac_nn,rms
- else
- write (iout,'(i4,2x,a8,i4,14(1pe14.5),0pf8.3)')
- & i,str_nam(ii),ist+1,
- & (energia(print_order(kk)),kk=1,nprint_ene),etot
- endif
- enddo
- write (iout,'(//a)') 'Final energies:'
- write (iout,'(a4,2x,a12,17a14,3a8)')
- & 'No','seq',(ename(print_order(kk)),kk=1,nprint_ene),'ETOT',
- & 'RMSnat','NatCONT','NNCONT','RMS'
- do ik=1,nthread
- i=ip(ik)
- ii=ipatt(1,i)
- ist=nres_base(2,ii)+ipatt(2,i)
- do kk=1,n_ene_comp
- energia(kk)=ener(kk,ik)
- enddo
- etot=ener(n_ene_comp+1,i)
- rmsnat=ener(n_ene_comp+2,i)
- rms=ener(n_ene_comp+3,i)
- frac=ener(n_ene_comp+4,i)
- frac_nn=ener(n_ene_comp+5,i)
- write (iout,'(i4,2x,a8,i4,14(1pe14.5),0pf8.3,f8.5,f8.5,f8.3)')
- & i,str_nam(ii),ist+1,
- & (energia(print_order(kk)),kk=1,nprint_ene),
- & etot,rmsnat,frac,frac_nn,rms
- enddo
- write (iout,'(/a/)') 'IEXAM array:'
- write (iout,'(i5)') nexcl
- do i=1,nexcl
- write (iout,'(2i5)') iexam(1,i),iexam(2,i)
- enddo
- write (iout,'(/a,1pe14.4/a,1pe14.4/)')
- & 'Max. time for threading step ',max_time_for_thread,
- & 'Average time for threading step: ',ave_time_for_thread
- return
- end
-c----------------------------------------------------------------------------
- subroutine sc_conf
-C Sample (hopefully) optimal SC orientations given backcone conformation.
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'COMMON.CHAIN'
- include 'COMMON.DBASE'
- include 'COMMON.INTERACT'
- include 'COMMON.VAR'
- include 'COMMON.THREAD'
- include 'COMMON.FFIELD'
- include 'COMMON.SBRIDGE'
- include 'COMMON.HEADER'
- include 'COMMON.GEO'
- include 'COMMON.IOUNITS'
- double precision varia(maxvar)
- common /srutu/ icall
- double precision energia(0:n_ene)
- logical glycine,fail
- maxsample=10
- link_end0=link_end
- link_end=min0(link_end,nss)
- do i=nnt,nct
- if (itype(i).ne.10) then
-cd print *,'i=',i,' itype=',itype(i),' theta=',theta(i+1)
- call gen_side(itype(i),theta(i+1),alph(i),omeg(i),fail)
- endif
- enddo
- call chainbuild
- call etotal(energia(0))
- do isample=1,maxsample
-C Choose a non-glycine side chain.
- glycine=.true.
- do while(glycine)
- ind_sc=iran_num(nnt,nct)
- glycine=(itype(ind_sc).eq.10)
- enddo
- alph0=alph(ind_sc)
- omeg0=omeg(ind_sc)
- call gen_side(itype(ind_sc),theta(ind_sc+1),alph(ind_sc),
- & omeg(ind_sc),fail)
- call chainbuild
- call etotal(energia(0))
-cd write (iout,'(a,i5,a,i4,2(a,f8.3),2(a,1pe14.5))')
-cd & 'Step:',isample,' SC',ind_sc,' alpha',alph(ind_sc)*rad2deg,
-cd & ' omega',omeg(ind_sc)*rad2deg,' old energy',e0,' new energy',e1
- e1=0.0d0
- if (e0.le.e1) then
- alph(ind_sc)=alph0
- omeg(ind_sc)=omeg0
- else
- e0=e1
- endif
- enddo
- link_end=link_end0
- return
- end
-c---------------------------------------------------------------------------
- subroutine write_stat_thread(ithread,ipattern,ist)
- implicit real*8 (a-h,o-z)
- include "DIMENSIONS"
- include "COMMON.CONTROL"
- include "COMMON.IOUNITS"
- include "COMMON.THREAD"
- include "COMMON.FFIELD"
- include "COMMON.DBASE"
- include "COMMON.NAMES"
- double precision energia(0:n_ene)
-
-#if defined(AIX) || defined(PGI)
- open(istat,file=statname,position='append')
-#else
- open(istat,file=statname,access='append')
-#endif
- do i=1,n_ene_comp
- energia(i)=ener(i,ithread)
- enddo
- etot=ener(n_ene_comp+1,ithread)
- rmsnat=ener(n_ene_comp+2,ithread)
- rms=ener(n_ene_comp+3,ithread)
- frac=ener(n_ene_comp+4,ithread)
- frac_nn=ener(n_ene_comp+5,ithread)
- write (istat,'(i4,2x,a8,i4,14(1pe14.5),0pf8.3,f8.5,f8.5,f8.3)')
- & ithread,str_nam(ipattern),ist+1,
- & (energia(print_order(i)),i=1,nprint_ene),
- & etot,rmsnat,frac,frac_nn,rms
- close (istat)
- return
- end