1 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
5 C Program to carry out conformational search of proteins in an united-residue C
8 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
9 implicit real*8 (a-h,o-z)
15 include 'COMMON.SETUP'
17 include 'COMMON.TIME1'
18 include 'COMMON.INTERACT'
19 include 'COMMON.NAMES'
21 include 'COMMON.HEADER'
22 include 'COMMON.CONTROL'
23 include 'COMMON.CONTACTS'
24 include 'COMMON.CHAIN'
26 include 'COMMON.IOUNITS'
27 include 'COMMON.FFIELD'
30 include 'COMMON.SBRIDGE'
31 double precision hrtime,mintime,sectime
32 character*64 text_mode_calc(-2:14) /'test',
33 & 'SC rotamer distribution',
34 & 'Energy evaluation or minimization',
35 & 'Regularization of PDB structure',
36 & 'Threading of a sequence on PDB structures',
37 & 'Monte Carlo (with minimization) ',
38 & 'Energy minimization of multiple conformations',
39 & 'Checking energy gradient',
40 & 'Entropic sampling Monte Carlo (with minimization)',
45 & 'Soft regularization of PDB structure',
46 & 'Mesoscopic molecular dynamics (MD) ',
48 & 'Replica exchange molecular dynamics (REMD)'/
51 c call memmon_print_usage()
55 & write(iout,*)'### LAST MODIFIED 03/28/12 23:29 by czarek'
56 if (me.eq.king) call cinfo
57 C Read force field parameters and job setup data
59 if (me.eq.king .or. .not. out1file) then
61 & text_mode_calc(modecalc)(:ilen(text_mode_calc(modecalc))),
63 if (minim) write (iout,'(a)')
64 & 'Conformations will be energy-minimized.'
65 write (iout,'(80(1h*)/)')
69 if (modecalc.eq.-2) then
72 else if (modecalc.eq.-1) then
73 write(iout,*) "call check_sc_map next"
78 if (fg_rank.gt.0) then
79 C Fine-grain slaves just do energy and gradient components.
80 call ergastulum ! slave workhouse in Latin
83 if (modecalc.eq.0) then
84 call exec_eeval_or_minim
85 else if (modecalc.eq.1) then
87 else if (modecalc.eq.2) then
89 else if (modecalc.eq.3 .or. modecalc .eq.6) then
91 else if (modecalc.eq.4) then
92 call exec_mult_eeval_or_minim
93 else if (modecalc.eq.5) then
95 else if (ModeCalc.eq.7) then
97 else if (ModeCalc.eq.8) then
99 else if (modecalc.eq.11) then
101 else if (modecalc.eq.12) then
103 else if (modecalc.eq.14) then
106 write (iout,'(a)') 'This calculation type is not supported',
112 if (fg_rank.eq.0) call finish_task
113 c call memmon_print_usage()
115 call print_detailed_timing
117 call MPI_Finalize(ierr)
120 call dajczas(tcpu(),hrtime,mintime,sectime)
121 stop '********** Program terminated normally.'
124 c--------------------------------------------------------------------------
130 include 'COMMON.SETUP'
131 include 'COMMON.CONTROL'
132 include 'COMMON.IOUNITS'
133 if (me.eq.king .or. .not. out1file)
134 & write (iout,*) "Calling chainbuild"
139 c---------------------------------------------------------------------------
140 subroutine exec_MREMD
144 include 'COMMON.SETUP'
145 include 'COMMON.CONTROL'
146 include 'COMMON.IOUNITS'
147 include 'COMMON.REMD'
148 if (me.eq.king .or. .not. out1file)
149 & write (iout,*) "Calling chainbuild"
151 if (me.eq.king .or. .not. out1file)
152 & write (iout,*) "Calling REMD"
162 write (iout,*) "MREMD works on parallel machines only"
166 c---------------------------------------------------------------------------
167 subroutine exec_eeval_or_minim
168 implicit real*8 (a-h,o-z)
173 include 'COMMON.SETUP'
174 include 'COMMON.TIME1'
175 include 'COMMON.INTERACT'
176 include 'COMMON.NAMES'
178 include 'COMMON.HEADER'
179 include 'COMMON.CONTROL'
180 include 'COMMON.CONTACTS'
181 include 'COMMON.CHAIN'
183 include 'COMMON.IOUNITS'
184 include 'COMMON.FFIELD'
185 include 'COMMON.REMD'
187 include 'COMMON.SBRIDGE'
189 double precision energy(0:n_ene)
190 double precision energy_long(0:n_ene),energy_short(0:n_ene)
191 double precision varia(maxvar)
192 if (indpdb.eq.0) call chainbuild
200 print *,"Processor",myrank," after chainbuild"
202 call etotal_long(energy_long(0))
203 write (iout,*) "Printing long range energy"
204 call enerprint(energy_long(0))
205 call etotal_short(energy_short(0))
206 write (iout,*) "Printing short range energy"
207 call enerprint(energy_short(0))
209 energy(i)=energy_long(i)+energy_short(i)
210 write (iout,*) i,energy_long(i),energy_short(i),energy(i)
212 write (iout,*) "Printing long+short range energy"
213 call enerprint(energy(0))
215 call etotal(energy(0))
217 time_ene=MPI_Wtime()-time00
219 time_ene=tcpu()-time00
221 write (iout,*) "Time for energy evaluation",time_ene
222 print *,"after etotal"
225 call enerprint(energy(0))
226 call hairpin(.true.,nharp,iharp)
227 call secondary2(.true.)
231 print *, 'Calling OVERLAP_SC'
232 call overlap_sc(fail)
236 call sc_move(2,nres-1,10,1d10,nft_sc,etot)
237 print *,'SC_move',nft_sc,etot
238 write(iout,*) 'SC_move',nft_sc,etot
242 print *, 'Calling MINIM_DC'
248 call minim_dc(etot,iretcode,nfun)
250 if (indpdb.ne.0) then
254 call geom_to_var(nvar,varia)
255 print *,'Calling MINIMIZE.'
261 call minimize(etot,varia,iretcode,nfun)
263 print *,'SUMSL return code is',iretcode,' eval ',nfun
265 evals=nfun/(MPI_WTIME()-time1)
267 evals=nfun/(tcpu()-time1)
269 print *,'# eval/s',evals
270 print *,'refstr=',refstr
271 call hairpin(.true.,nharp,iharp)
272 call secondary2(.true.)
273 call etotal(energy(0))
275 call enerprint(energy(0))
278 call briefout(0,etot)
279 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
280 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
281 write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
282 write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
284 print *,'refstr=',refstr
285 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
286 call briefout(0,etot)
288 if (outpdb) call pdbout(etot,titel(:32),ipdb)
289 if (outmol2) call mol2out(etot,titel(:32))
292 c---------------------------------------------------------------------------
293 subroutine exec_regularize
294 implicit real*8 (a-h,o-z)
299 include 'COMMON.SETUP'
300 include 'COMMON.TIME1'
301 include 'COMMON.INTERACT'
302 include 'COMMON.NAMES'
304 include 'COMMON.HEADER'
305 include 'COMMON.CONTROL'
306 include 'COMMON.CONTACTS'
307 include 'COMMON.CHAIN'
309 include 'COMMON.IOUNITS'
310 include 'COMMON.FFIELD'
311 include 'COMMON.REMD'
313 include 'COMMON.SBRIDGE'
314 double precision energy(0:n_ene)
319 call regularize(nct-nnt+1,etot,rms,cref(1,nnt),iretcode)
320 call etotal(energy(0))
321 energy(0)=energy(0)-energy(14)
323 call enerprint(energy(0))
325 call briefout(0,etot)
326 if (outpdb) call pdbout(etot,titel(:32),ipdb)
327 if (outmol2) call mol2out(etot,titel(:32))
328 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
329 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
332 c---------------------------------------------------------------------------
333 subroutine exec_thread
338 include "COMMON.SETUP"
342 c---------------------------------------------------------------------------
344 implicit real*8 (a-h,o-z)
346 character*10 nodeinfo
347 double precision varia(maxvar)
351 include "COMMON.SETUP"
352 include 'COMMON.CONTROL'
356 if (modecalc.eq.3) then
362 if (modecalc.eq.3) then
373 c---------------------------------------------------------------------------
374 subroutine exec_mult_eeval_or_minim
375 implicit real*8 (a-h,o-z)
379 dimension muster(mpi_status_size)
381 include 'COMMON.SETUP'
382 include 'COMMON.TIME1'
383 include 'COMMON.INTERACT'
384 include 'COMMON.NAMES'
386 include 'COMMON.HEADER'
387 include 'COMMON.CONTROL'
388 include 'COMMON.CONTACTS'
389 include 'COMMON.CHAIN'
391 include 'COMMON.IOUNITS'
392 include 'COMMON.FFIELD'
393 include 'COMMON.REMD'
395 include 'COMMON.SBRIDGE'
396 double precision varia(maxvar)
398 double precision energy(0:n_ene)
408 open(intin,file=intinname,status='old')
409 write (istat,'(a5,30a12)')"# ",
410 & (wname(print_order(i)),i=1,nprint_ene)
412 write (istat,'(a5,30a12)')"# ",
413 & (ename(print_order(i)),i=1,nprint_ene),
414 & "ETOT total","RMSD","nat.contact","nnt.contact","cont.order"
416 write (istat,'(a5,30a12)')"# ",
417 & (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
423 read (intin,'(e15.10,e15.5)',end=1100,err=1100) time,ene
424 call read_x(intin,*11)
426 c Broadcast the order to compute internal coordinates to the slaves.
428 & call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
430 call int_from_cart1(.false.)
432 read (intin,'(i5)',end=1100,err=1100) iconf
433 call read_angles(intin,*11)
434 call geom_to_var(nvar,varia)
437 write (iout,'(a,i7)') 'Conformation #',iconf
438 call etotal(energy(0))
439 call briefout(iconf,energy(0))
440 call enerprint(energy(0))
443 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
444 write (istat,'(i5,30(f12.3))') iconf,
445 & (energy(print_order(i)),i=1,nprint_ene),etot,
446 & rms,frac,frac_nn,co
449 write (istat,'(i5,30(f12.3))') iconf,
450 & (energy(print_order(i)),i=1,nprint_ene),etot
466 if (mm.lt.nodes) then
468 read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
469 call read_x(intin,*11)
471 c Broadcast the order to compute internal coordinates to the slaves.
473 & call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
475 call int_from_cart1(.false.)
477 read (intin,'(i5)',end=11,err=11) iconf
478 call read_angles(intin,*11)
479 call geom_to_var(nvar,varia)
484 write (iout,*) 'Conformation #',iconf,' read'
493 call mpi_send(ind,6,mpi_integer,mm,idint,CG_COMM,
495 call mpi_send(varia,nvar,mpi_double_precision,mm,
496 * idreal,CG_COMM,ierr)
497 call mpi_send(ene0,1,mpi_double_precision,mm,
498 * idreal,CG_COMM,ierr)
499 c print *,'task ',n,' sent to worker ',mm,nvar
501 call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,
502 * CG_COMM,muster,ierr)
503 man=muster(mpi_source)
504 c print *,'receiving result from worker ',man,' (',iii1,iii,')'
505 call mpi_recv(varia,nvar,mpi_double_precision,
506 * man,idreal,CG_COMM,muster,ierr)
508 * mpi_double_precision,man,idreal,
509 * CG_COMM,muster,ierr)
510 call mpi_recv(ene0,1,
511 * mpi_double_precision,man,idreal,
512 * CG_COMM,muster,ierr)
513 c print *,'result received from worker ',man,' sending now'
515 call var_to_geom(nvar,varia)
517 call etotal(energy(0))
521 write (iout,*) 'Conformation #',iconf," sumsl return code ",
525 call enerprint(energy(0))
526 call briefout(it,etot)
527 c if (minim) call briefout(it,etot)
529 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
530 write (istat,'(i5,30(f12.3))') iconf,
531 & (energy(print_order(i)),i=1,nprint_ene),etot,
532 & rms,frac,frac_nn,co
534 write (istat,'(i5,30(f12.3))') iconf,
535 & (energy(print_order(i)),i=1,nprint_ene),etot
540 read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
541 call read_x(intin,*11)
543 c Broadcast the order to compute internal coordinates to the slaves.
545 & call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
547 call int_from_cart1(.false.)
549 read (intin,'(i5)',end=11,err=11) iconf
550 call read_angles(intin,*11)
551 call geom_to_var(nvar,varia)
555 write (iout,*) 'Conformation #',iconf,' read'
563 call mpi_send(ind,6,mpi_integer,man,idint,CG_COMM,
565 call mpi_send(varia,nvar,mpi_double_precision,man,
566 * idreal,CG_COMM,ierr)
567 call mpi_send(ene0,1,mpi_double_precision,man,
568 * idreal,CG_COMM,ierr)
569 nf_mcmf=nf_mcmf+ind(4)
575 call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,
576 * CG_COMM,muster,ierr)
577 man=muster(mpi_source)
578 call mpi_recv(varia,nvar,mpi_double_precision,
579 * man,idreal,CG_COMM,muster,ierr)
581 * mpi_double_precision,man,idreal,
582 * CG_COMM,muster,ierr)
583 call mpi_recv(ene0,1,
584 * mpi_double_precision,man,idreal,
585 * CG_COMM,muster,ierr)
587 call var_to_geom(nvar,varia)
589 call etotal(energy(0))
593 write (iout,*) 'Conformation #',iconf," sumsl return code ",
597 call enerprint(energy(0))
598 call briefout(it,etot)
600 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
601 write (istat,'(i5,30(f12.3))') iconf,
602 & (energy(print_order(i)),i=1,nprint_ene),etot,
603 & rms,frac,frac_nn,co
605 write (istat,'(i5,30(f12.3))') iconf,
606 & (energy(print_order(i)),i=1,nprint_ene),etot
618 call mpi_send(ind,6,mpi_integer,i,idint,CG_COMM,
623 open(intin,file=intinname,status='old')
624 write (istat,'(a5,20a12)')"# ",
625 & (wname(print_order(i)),i=1,nprint_ene)
626 write (istat,'("# ",20(1pe12.4))')
627 & (weights(print_order(i)),i=1,nprint_ene)
629 write (istat,'(a5,20a12)')"# ",
630 & (ename(print_order(i)),i=1,nprint_ene),
631 & "ETOT total","RMSD","nat.contact","nnt.contact"
633 write (istat,'(a5,14a12)')"# ",
634 & (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
638 read (intin,'(e15.10,e15.5)',end=1100,err=1100) time,ene
639 call read_x(intin,*11)
641 c Broadcast the order to compute internal coordinates to the slaves.
643 & call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
645 call int_from_cart1(.false.)
647 read (intin,'(i5)',end=1100,err=1100) iconf
648 call read_angles(intin,*11)
649 call geom_to_var(nvar,varia)
652 write (iout,'(a,i7)') 'Conformation #',iconf
653 if (minim) call minimize(etot,varia,iretcode,nfun)
654 call etotal(energy(0))
657 call enerprint(energy(0))
658 if (minim) call briefout(it,etot)
660 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
661 write (istat,'(i5,18(f12.3))') iconf,
662 & (energy(print_order(i)),i=1,nprint_ene),
663 & etot,rms,frac,frac_nn,co
666 write (istat,'(i5,14(f12.3))') iconf,
667 & (energy(print_order(i)),i=1,nprint_ene),etot
675 c---------------------------------------------------------------------------
676 subroutine exec_checkgrad
677 implicit real*8 (a-h,o-z)
682 include 'COMMON.SETUP'
683 include 'COMMON.TIME1'
684 include 'COMMON.INTERACT'
685 include 'COMMON.NAMES'
687 include 'COMMON.HEADER'
688 include 'COMMON.CONTROL'
689 include 'COMMON.CONTACTS'
690 include 'COMMON.CHAIN'
692 include 'COMMON.IOUNITS'
693 include 'COMMON.FFIELD'
694 include 'COMMON.REMD'
696 include 'COMMON.SBRIDGE'
698 double precision energy(0:max_ene)
700 c vbld(i)=vbld(i)+ran_number(-0.1d0,0.1d0)
701 c if (itype(i).ne.10)
702 c & vbld(i+nres)=vbld(i+nres)+ran_number(-0.001d0,0.001d0)
704 if (indpdb.eq.0) call chainbuild
707 c dc(j,i)=dc(j,i)+ran_number(-0.2d0,0.2d0)
711 c if (itype(i).ne.10) then
713 c dc(j,i+nres)=dc(j,i+nres)+ran_number(-0.2d0,0.2d0)
718 c dc(j,0)=ran_number(-0.2d0,0.2d0)
725 call rescale_weights(t_bath)
730 call etotal(energy(0))
732 call enerprint(energy(0))
733 write (iout,*) "Uconst",Uconst," Uconst_back",uconst_back
734 print *,'icheckgrad=',icheckgrad
735 goto (10,20,30) icheckgrad
736 10 call check_ecartint
738 20 call check_cartgrad
743 c---------------------------------------------------------------------------
750 c---------------------------------------------------------------------------
756 include 'COMMON.IOUNITS'
757 C Conformational Space Annealling programmed by Jooyoung Lee.
758 C This method works only with parallel machines!
761 write (iout,*) "CSA is not supported in this version"
763 csa write (iout,*) "CSA works on parallel machines only"
764 write (iout,*) "CSA is not supported in this version"
768 c---------------------------------------------------------------------------
769 subroutine exec_softreg
770 implicit real*8 (a-h,o-z)
772 include 'COMMON.IOUNITS'
773 include 'COMMON.CONTROL'
774 double precision energy(0:max_ene)
775 logical debug /.false./
777 call etotal(energy(0))
778 call enerprint(energy(0))
779 if (.not.lsecondary) then
780 write(iout,*) 'Calling secondary structure recognition'
781 call secondary2(debug)
783 write(iout,*) 'Using secondary structure supplied in pdb'
788 call etotal(energy(0))
790 call enerprint(energy(0))
792 call briefout(0,etot)
793 call secondary2(.true.)
794 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)