CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C U N R E S C C C C Program to carry out conformational search of proteins in an united-residue C C approximation. C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC implicit real*8 (a-h,o-z) include 'DIMENSIONS' #ifdef MPI include 'mpif.h' include 'COMMON.SETUP' #endif include 'COMMON.TIME1' include 'COMMON.INTERACT' include 'COMMON.NAMES' include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' include 'COMMON.FFIELD' include 'COMMON.REMD' include 'COMMON.MD' include 'COMMON.SBRIDGE' double precision hrtime,mintime,sectime character*64 text_mode_calc(-2:14) /'test', & 'SC rotamer distribution', & 'Energy evaluation or minimization', & 'Regularization of PDB structure', & 'Threading of a sequence on PDB structures', & 'Monte Carlo (with minimization) ', & 'Energy minimization of multiple conformations', & 'Checking energy gradient', & 'Entropic sampling Monte Carlo (with minimization)', & 'Energy map', & 'CSA calculations', & 'Not used 9', & 'Not used 10', & 'Soft regularization of PDB structure', & 'Mesoscopic molecular dynamics (MD) ', & 'Not used 13', & 'Replica exchange molecular dynamics (REMD)'/ external ilen c call memmon_print_usage() call init_task if (me.eq.king) & write(iout,*)'### LAST MODIFIED 4/25/08 7:29PM by adam' if (me.eq.king) call cinfo C Read force field parameters and job setup data call readrtns C if (me.eq.king .or. .not. out1file) then write (iout,'(2a/)') & text_mode_calc(modecalc)(:ilen(text_mode_calc(modecalc))), & ' calculation.' if (minim) write (iout,'(a)') & 'Conformations will be energy-minimized.' write (iout,'(80(1h*)/)') endif call flush(iout) C if (modecalc.eq.-2) then call test stop else if (modecalc.eq.-1) then write(iout,*) "call check_sc_map next" call check_bond stop endif #ifdef MPI if (fg_rank.gt.0) then C Fine-grain slaves just do energy and gradient components. call ergastulum ! slave workhouse in Latin else #endif if (modecalc.eq.0) then call exec_eeval_or_minim else if (modecalc.eq.1) then call exec_regularize else if (modecalc.eq.2) then call exec_thread else if (modecalc.eq.3 .or. modecalc .eq.6) then call exec_MC else if (modecalc.eq.4) then call exec_mult_eeval_or_minim else if (modecalc.eq.5) then call exec_checkgrad else if (ModeCalc.eq.7) then call exec_map else if (ModeCalc.eq.8) then call exec_CSA else if (modecalc.eq.11) then call exec_softreg else if (modecalc.eq.12) then call exec_MD else if (modecalc.eq.14) then #ifdef MPI call exec_MREMD #else write (iout,*) "Need a parallel version to run MREMD." stop #endif else write (iout,'(a)') 'This calculation type is not supported', & ModeCalc endif #ifdef MPI endif C Finish task. if (fg_rank.eq.0) call finish_task c call memmon_print_usage() #ifdef TIMING call print_detailed_timing #endif call MPI_Finalize(ierr) stop 'Bye Bye...' #else call dajczas(tcpu(),hrtime,mintime,sectime) stop '********** Program terminated normally.' #endif end c-------------------------------------------------------------------------- subroutine exec_MD include 'DIMENSIONS' #ifdef MPI include "mpif.h" #endif include 'COMMON.SETUP' include 'COMMON.CONTROL' include 'COMMON.IOUNITS' c if (me.eq.king .or. .not. out1file) then c write (iout,*) "Calling chainbuild" c call flush(iout) c endif call chainbuild c if (me.eq.king .or. .not. out1file) then c write (iout,*) "Calling MD" c call flush(iout) c endif call MD return end c--------------------------------------------------------------------------- #ifdef MPI subroutine exec_MREMD include 'DIMENSIONS' #ifdef MPI include "mpif.h" #endif include 'COMMON.SETUP' include 'COMMON.CONTROL' include 'COMMON.IOUNITS' include 'COMMON.REMD' if (me.eq.king .or. .not. out1file) & write (iout,*) "Calling chainbuild" call chainbuild if (me.eq.king .or. .not. out1file) & write (iout,*) "Calling REMD" if (remd_mlist) then call MREMD else do i=1,nrep remd_m(i)=1 enddo call MREMD endif return end #endif c--------------------------------------------------------------------------- subroutine exec_eeval_or_minim implicit real*8 (a-h,o-z) include 'DIMENSIONS' #ifdef MPI include 'mpif.h' #endif include 'COMMON.SETUP' include 'COMMON.TIME1' include 'COMMON.INTERACT' include 'COMMON.NAMES' include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' include 'COMMON.FFIELD' include 'COMMON.REMD' include 'COMMON.MD' include 'COMMON.SBRIDGE' common /srutu/ icall double precision energy(0:n_ene) double precision energy_long(0:n_ene),energy_short(0:n_ene) double precision varia(maxvar) if (indpdb.eq.0) call chainbuild #ifdef MPI time00=MPI_Wtime() #else time00=tcpu() #endif call chainbuild_cart if (split_ene) then print *,"Processor",myrank," after chainbuild" icall=1 call etotal_long(energy_long(0)) write (iout,*) "Printing long range energy" call enerprint(energy_long(0)) call etotal_short(energy_short(0)) write (iout,*) "Printing short range energy" call enerprint(energy_short(0)) do i=0,n_ene energy(i)=energy_long(i)+energy_short(i) write (iout,*) i,energy_long(i),energy_short(i),energy(i) enddo write (iout,*) "Printing long+short range energy" call enerprint(energy(0)) endif call etotal(energy(0)) #ifdef MPI time_ene=MPI_Wtime()-time00 #else time_ene=tcpu()-time00 #endif write (iout,*) "Time for energy evaluation",time_ene print *,"after etotal" etota = energy(0) etot =etota call enerprint(energy(0)) call hairpin(.true.,nharp,iharp) call secondary2(.true.) if (minim) then crc overlap test if (overlapsc) then print *, 'Calling OVERLAP_SC' call overlap_sc(fail) endif if (searchsc) then call sc_move(2,nres-1,10,1d10,nft_sc,etot) print *,'SC_move',nft_sc,etot write(iout,*) 'SC_move',nft_sc,etot endif if (dccart) then print *, 'Calling MINIM_DC' #ifdef MPI time1=MPI_WTIME() #else time1=tcpu() #endif call minim_dc(etot,iretcode,nfun) if(iretcode.eq.8) call check_ecartint else if (indpdb.ne.0) then call bond_regular call chainbuild endif call geom_to_var(nvar,varia) print *,'Calling MINIMIZE.' #ifdef MPI time1=MPI_WTIME() #else time1=tcpu() #endif call minimize(etot,varia,iretcode,nfun) endif print *,'SUMSL return code is',iretcode,' eval ',nfun #ifdef MPI evals=nfun/(MPI_WTIME()-time1) #else evals=nfun/(tcpu()-time1) #endif print *,'# eval/s',evals print *,'refstr=',refstr call hairpin(.true.,nharp,iharp) call secondary2(.true.) call etotal(energy(0)) etot = energy(0) call enerprint(energy(0)) call intout call briefout(0,etot) if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.) write (iout,'(a,i3)') 'SUMSL return code:',iretcode write (iout,'(a,i20)') '# of energy evaluations:',nfun+1 write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals else print *,'refstr=',refstr if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.) call briefout(0,etot) endif if (outpdb) call pdbout(etot,titel(:32),ipdb) if (outmol2) call mol2out(etot,titel(:32)) return end c--------------------------------------------------------------------------- subroutine exec_regularize implicit real*8 (a-h,o-z) include 'DIMENSIONS' #ifdef MPI include 'mpif.h' #endif include 'COMMON.SETUP' include 'COMMON.TIME1' include 'COMMON.INTERACT' include 'COMMON.NAMES' include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' include 'COMMON.FFIELD' include 'COMMON.REMD' include 'COMMON.MD' include 'COMMON.SBRIDGE' double precision energy(0:n_ene) call gen_dist_constr call sc_conf call intout call regularize(nct-nnt+1,etot,rms,cref(1,nnt,1),iretcode) call etotal(energy(0)) energy(0)=energy(0)-energy(14) etot=energy(0) call enerprint(energy(0)) call intout call briefout(0,etot) if (outpdb) call pdbout(etot,titel(:32),ipdb) if (outmol2) call mol2out(etot,titel(:32)) if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.) write (iout,'(a,i3)') 'SUMSL return code:',iretcode return end c--------------------------------------------------------------------------- subroutine exec_thread include 'DIMENSIONS' #ifdef MP include "mpif.h" #endif include "COMMON.SETUP" call thread_seq return end c--------------------------------------------------------------------------- subroutine exec_MC implicit real*8 (a-h,o-z) include 'DIMENSIONS' character*10 nodeinfo double precision varia(maxvar) #ifdef MPI include "mpif.h" #endif include "COMMON.SETUP" include 'COMMON.CONTROL' call mcm_setup if (minim) then #ifdef MPI if (modecalc.eq.3) then call do_mcm(ipar) else call entmcm endif #else if (modecalc.eq.3) then call do_mcm(ipar) else call entmcm endif #endif else call monte_carlo endif return end c--------------------------------------------------------------------------- subroutine exec_mult_eeval_or_minim implicit real*8 (a-h,o-z) include 'DIMENSIONS' #ifdef MPI include 'mpif.h' dimension muster(mpi_status_size) #endif include 'COMMON.SETUP' include 'COMMON.TIME1' include 'COMMON.INTERACT' include 'COMMON.NAMES' include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' include 'COMMON.FFIELD' include 'COMMON.REMD' include 'COMMON.MD' include 'COMMON.SBRIDGE' double precision varia(maxvar) dimension ind(6) double precision energy(0:n_ene) logical eof eof=.false. #ifdef MPI if(me.ne.king) then call minim_mcmf return endif close (intin) open(intin,file=intinname,status='old') write (istat,'(a5,20a12)')"# ", & (wname(print_order(i)),i=1,nprint_ene) if (refstr) then write (istat,'(a5,20a12)')"# ", & (ename(print_order(i)),i=1,nprint_ene), & "ETOT total","RMSD","nat.contact","nnt.contact" else write (istat,'(a5,20a12)')"# ", & (ename(print_order(i)),i=1,nprint_ene),"ETOT total" endif if (.not.minim) then do while (.not. eof) if (read_cart) then read (intin,'(e15.10,e15.5)',end=1100,err=1100) time,ene call read_x(intin,*11) #ifdef MPI c Broadcast the order to compute internal coordinates to the slaves. if (nfgtasks.gt.1) & call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR) #endif call int_from_cart1(.false.) else read (intin,'(i5)',end=1100,err=1100) iconf call read_angles(intin,*11) call geom_to_var(nvar,varia) call chainbuild endif write (iout,'(a,i7)') 'Conformation #',iconf call etotal(energy(0)) call briefout(iconf,energy(0)) call enerprint(energy(0)) etot=energy(0) if (refstr) then call rms_nac_nnc(rms,frac,frac_nn,co,.true.) write (istat,'(i5,20(f12.3))') iconf, & (energy(print_order(i)),i=1,nprint_ene),etot, & rms,frac,frac_nn,co cjlee end else write (istat,'(i5,16(f12.3))') iconf, & (energy(print_order(i)),i=1,nprint_ene),etot endif enddo 1100 continue goto 1101 endif mm=0 imm=0 nft=0 ene0=0.0d0 n=0 iconf=0 c do n=1,nzsc do while (.not. eof) mm=mm+1 if (mm.lt.nodes) then if (read_cart) then read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene call read_x(intin,*11) #ifdef MPI c Broadcast the order to compute internal coordinates to the slaves. if (nfgtasks.gt.1) & call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR) #endif call int_from_cart1(.false.) else read (intin,'(i5)',end=11,err=11) iconf call read_angles(intin,*11) call geom_to_var(nvar,varia) call chainbuild endif write (iout,'(a,i7)') 'Conformation #',iconf n=n+1 imm=imm+1 ind(1)=1 ind(2)=n ind(3)=0 ind(4)=0 ind(5)=0 ind(6)=0 ene0=0.0d0 call mpi_send(ind,6,mpi_integer,mm,idint,CG_COMM, * ierr) call mpi_send(varia,nvar,mpi_double_precision,mm, * idreal,CG_COMM,ierr) call mpi_send(ene0,1,mpi_double_precision,mm, * idreal,CG_COMM,ierr) c print *,'task ',n,' sent to worker ',mm,nvar else call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint, * CG_COMM,muster,ierr) man=muster(mpi_source) c print *,'receiving result from worker ',man,' (',iii1,iii,')' call mpi_recv(varia,nvar,mpi_double_precision, * man,idreal,CG_COMM,muster,ierr) call mpi_recv(ene,1, * mpi_double_precision,man,idreal, * CG_COMM,muster,ierr) call mpi_recv(ene0,1, * mpi_double_precision,man,idreal, * CG_COMM,muster,ierr) c print *,'result received from worker ',man,' sending now' call var_to_geom(nvar,varia) call chainbuild call etotal(energy(0)) iconf=ind(2) write (iout,*) write (iout,*) write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5) etot=energy(0) call enerprint(energy(0)) call briefout(it,etot) c if (minim) call briefout(it,etot) if (refstr) then call rms_nac_nnc(rms,frac,frac_nn,co,.true.) write (istat,'(i5,19(f12.3))') iconf, & (energy(print_order(i)),i=1,nprint_ene),etot, & rms,frac,frac_nn,co else write (istat,'(i5,15(f12.3))') iconf, & (energy(print_order(i)),i=1,nprint_ene),etot endif imm=imm-1 if (read_cart) then read (intin,'(e15.10,e15.5)',end=1101,err=1101) time,ene call read_x(intin,*11) #ifdef MPI c Broadcast the order to compute internal coordinates to the slaves. if (nfgtasks.gt.1) & call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR) #endif call int_from_cart1(.false.) else read (intin,'(i5)',end=1101,err=1101) iconf call read_angles(intin,*11) call geom_to_var(nvar,varia) call chainbuild endif n=n+1 imm=imm+1 ind(1)=1 ind(2)=n ind(3)=0 ind(4)=0 ind(5)=0 ind(6)=0 call mpi_send(ind,6,mpi_integer,man,idint,CG_COMM, * ierr) call mpi_send(varia,nvar,mpi_double_precision,man, * idreal,CG_COMM,ierr) call mpi_send(ene0,1,mpi_double_precision,man, * idreal,CG_COMM,ierr) nf_mcmf=nf_mcmf+ind(4) nmin=nmin+1 endif enddo 11 continue do j=1,imm call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint, * CG_COMM,muster,ierr) man=muster(mpi_source) call mpi_recv(varia,nvar,mpi_double_precision, * man,idreal,CG_COMM,muster,ierr) call mpi_recv(ene,1, * mpi_double_precision,man,idreal, * CG_COMM,muster,ierr) call mpi_recv(ene0,1, * mpi_double_precision,man,idreal, * CG_COMM,muster,ierr) call var_to_geom(nvar,varia) call chainbuild call etotal(energy(0)) iconf=ind(2) write (iout,*) write (iout,*) write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5) etot=energy(0) call enerprint(energy(0)) call briefout(it,etot) if (refstr) then call rms_nac_nnc(rms,frac,frac_nn,co,.true.) write (istat,'(i5,19(f12.3))') iconf, & (energy(print_order(i)),i=1,nprint_ene),etot, & rms,frac,frac_nn,co else write (istat,'(i5,15(f12.3))') iconf, & (energy(print_order(i)),i=1,nprint_ene),etot endif nmin=nmin+1 enddo 1101 continue do i=1, nodes-1 ind(1)=0 ind(2)=0 ind(3)=0 ind(4)=0 ind(5)=0 ind(6)=0 call mpi_send(ind,6,mpi_integer,i,idint,CG_COMM, * ierr) enddo #else close (intin) open(intin,file=intinname,status='old') write (istat,'(a5,20a12)')"# ", & (wname(print_order(i)),i=1,nprint_ene) write (istat,'("# ",20(1pe12.4))') & (weights(print_order(i)),i=1,nprint_ene) if (refstr) then write (istat,'(a5,20a12)')"# ", & (ename(print_order(i)),i=1,nprint_ene), & "ETOT total","RMSD","nat.contact","nnt.contact" else write (istat,'(a5,14a12)')"# ", & (ename(print_order(i)),i=1,nprint_ene),"ETOT total" endif do while (.not. eof) if (read_cart) then read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene call read_x(intin,*11) #ifdef MPI c Broadcast the order to compute internal coordinates to the slaves. if (nfgtasks.gt.1) & call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR) #endif call int_from_cart1(.false.) else read (intin,'(i5)',end=11,err=11) iconf call read_angles(intin,*11) call geom_to_var(nvar,varia) call chainbuild endif write (iout,'(a,i7)') 'Conformation #',iconf if (minim) call minimize(etot,varia,iretcode,nfun) call etotal(energy(0)) etot=energy(0) call enerprint(energy(0)) if (minim) call briefout(it,etot) if (refstr) then call rms_nac_nnc(rms,frac,frac_nn,co,.true.) write (istat,'(i5,18(f12.3))') iconf, & (energy(print_order(i)),i=1,nprint_ene), & etot,rms,frac,frac_nn,co cjlee end else write (istat,'(i5,14(f12.3))') iconf, & (energy(print_order(i)),i=1,nprint_ene),etot endif enddo 11 continue #endif return end c--------------------------------------------------------------------------- subroutine exec_checkgrad implicit real*8 (a-h,o-z) include 'DIMENSIONS' #ifdef MPI include 'mpif.h' #endif include 'COMMON.SETUP' include 'COMMON.TIME1' include 'COMMON.INTERACT' include 'COMMON.NAMES' include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' include 'COMMON.FFIELD' include 'COMMON.REMD' include 'COMMON.MD' include 'COMMON.SBRIDGE' common /srutu/ icall double precision energy(0:n_ene) c do i=2,nres c vbld(i)=vbld(i)+ran_number(-0.1d0,0.1d0) c if (itype(i).ne.10) c & vbld(i+nres)=vbld(i+nres)+ran_number(-0.001d0,0.001d0) c enddo if (indpdb.eq.0) call chainbuild c do i=0,nres c do j=1,3 c dc(j,i)=dc(j,i)+ran_number(-0.2d0,0.2d0) c enddo c enddo c do i=1,nres-1 c if (itype(i).ne.10) then c do j=1,3 c dc(j,i+nres)=dc(j,i+nres)+ran_number(-0.2d0,0.2d0) c enddo c endif c enddo c do j=1,3 c dc(j,0)=ran_number(-0.2d0,0.2d0) c enddo usampl=.true. totT=1.d0 eq_time=0.0d0 call read_fragments call chainbuild_cart call cartprint call intout icall=1 call etotal(energy(0)) etot = energy(0) call enerprint(energy(0)) write (iout,*) "Uconst",Uconst," Uconst_back",uconst_back print *,'icheckgrad=',icheckgrad goto (10,20,30) icheckgrad 10 call check_ecartint write(iout,*) "kupadupa" call check_ecartint return 20 call check_cartgrad return 30 call check_eint return end c--------------------------------------------------------------------------- subroutine exec_map C Energy maps call map_read call map return end c--------------------------------------------------------------------------- subroutine exec_CSA #ifdef MPI include "mpif.h" #endif include 'DIMENSIONS' include 'COMMON.IOUNITS' C Conformational Space Annealling programmed by Jooyoung Lee. C This method works only with parallel machines! #ifdef MPI call together #else write (iout,*) "CSA works on parallel machines only" #endif return end c--------------------------------------------------------------------------- subroutine exec_softreg include 'DIMENSIONS' include 'COMMON.IOUNITS' include 'COMMON.CONTROL' double precision energy(0:n_ene) call chainbuild call etotal(energy(0)) call enerprint(energy(0)) if (.not.lsecondary) then write(iout,*) 'Calling secondary structure recognition' call secondary2(debug) else write(iout,*) 'Using secondary structure supplied in pdb' endif call softreg call etotal(energy(0)) etot=energy(0) call enerprint(energy(0)) call intout call briefout(0,etot) call secondary2(.true.) if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.) return end