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
61 if (me.eq.king .or. .not. out1file) then
63 & text_mode_calc(modecalc)(:ilen(text_mode_calc(modecalc))),
65 if (minim) write (iout,'(a)')
66 & 'Conformations will be energy-minimized.'
67 write (iout,'(80(1h*)/)')
71 if (modecalc.eq.-2) then
74 else if (modecalc.eq.-1) then
75 write(iout,*) "call check_sc_map next"
80 if (fg_rank.gt.0) then
81 C Fine-grain slaves just do energy and gradient components.
82 call ergastulum ! slave workhouse in Latin
85 if (modecalc.eq.0) then
86 call exec_eeval_or_minim
87 else if (modecalc.eq.1) then
89 else if (modecalc.eq.2) then
91 else if (modecalc.eq.3 .or. modecalc .eq.6) then
93 else if (modecalc.eq.4) then
94 call exec_mult_eeval_or_minim
95 else if (modecalc.eq.5) then
97 else if (ModeCalc.eq.7) then
99 else if (ModeCalc.eq.8) then
101 else if (modecalc.eq.11) then
103 else if (modecalc.eq.12) then
105 else if (modecalc.eq.14) then
108 write (iout,'(a)') 'This calculation type is not supported',
114 if (fg_rank.eq.0) call finish_task
115 c call memmon_print_usage()
117 call print_detailed_timing
119 call MPI_Finalize(ierr)
122 call dajczas(tcpu(),hrtime,mintime,sectime)
123 stop '********** Program terminated normally.'
126 c--------------------------------------------------------------------------
132 include 'COMMON.SETUP'
133 include 'COMMON.CONTROL'
134 include 'COMMON.IOUNITS'
135 if (me.eq.king .or. .not. out1file)
136 & write (iout,*) "Calling chainbuild"
141 c---------------------------------------------------------------------------
142 subroutine exec_MREMD
146 include 'COMMON.SETUP'
147 include 'COMMON.CONTROL'
148 include 'COMMON.IOUNITS'
149 include 'COMMON.REMD'
150 if (me.eq.king .or. .not. out1file)
151 & write (iout,*) "Calling chainbuild"
153 if (me.eq.king .or. .not. out1file)
154 & write (iout,*) "Calling REMD"
164 write (iout,*) "MREMD works on parallel machines only"
168 c---------------------------------------------------------------------------
169 subroutine exec_eeval_or_minim
170 implicit real*8 (a-h,o-z)
175 include 'COMMON.SETUP'
176 include 'COMMON.TIME1'
177 include 'COMMON.INTERACT'
178 include 'COMMON.NAMES'
180 include 'COMMON.HEADER'
181 include 'COMMON.CONTROL'
182 include 'COMMON.CONTACTS'
183 include 'COMMON.CHAIN'
185 include 'COMMON.IOUNITS'
186 include 'COMMON.FFIELD'
187 include 'COMMON.REMD'
189 include 'COMMON.SBRIDGE'
191 double precision energy(0:n_ene)
192 double precision energy_long(0:n_ene),energy_short(0:n_ene)
193 double precision varia(maxvar)
194 if (indpdb.eq.0) call chainbuild
202 print *,"Processor",myrank," after chainbuild"
204 call etotal_long(energy_long(0))
205 write (iout,*) "Printing long range energy"
206 call enerprint(energy_long(0))
207 call etotal_short(energy_short(0))
208 write (iout,*) "Printing short range energy"
209 call enerprint(energy_short(0))
211 energy(i)=energy_long(i)+energy_short(i)
212 write (iout,*) i,energy_long(i),energy_short(i),energy(i)
214 write (iout,*) "Printing long+short range energy"
215 call enerprint(energy(0))
217 call etotal(energy(0))
219 time_ene=MPI_Wtime()-time00
221 time_ene=tcpu()-time00
223 write (iout,*) "Time for energy evaluation",time_ene
224 print *,"after etotal"
227 call enerprint(energy(0))
228 call hairpin(.true.,nharp,iharp)
229 call secondary2(.true.)
233 print *, 'Calling OVERLAP_SC'
234 call overlap_sc(fail)
238 call sc_move(2,nres-1,10,1d10,nft_sc,etot)
239 print *,'SC_move',nft_sc,etot
240 write(iout,*) 'SC_move',nft_sc,etot
244 print *, 'Calling MINIM_DC'
250 call minim_dc(etot,iretcode,nfun)
252 if (indpdb.ne.0) then
256 call geom_to_var(nvar,varia)
257 print *,'Calling MINIMIZE.'
263 call minimize(etot,varia,iretcode,nfun)
265 print *,'SUMSL return code is',iretcode,' eval ',nfun
267 evals=nfun/(MPI_WTIME()-time1)
269 evals=nfun/(tcpu()-time1)
271 print *,'# eval/s',evals
272 print *,'refstr=',refstr
273 call hairpin(.true.,nharp,iharp)
274 call secondary2(.true.)
275 call etotal(energy(0))
277 call enerprint(energy(0))
280 call briefout(0,etot)
281 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
282 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
283 write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
284 write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
286 print *,'refstr=',refstr
287 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
288 call briefout(0,etot)
290 if (outpdb) call pdbout(etot,titel(:32),ipdb)
291 if (outmol2) call mol2out(etot,titel(:32))
294 c---------------------------------------------------------------------------
295 subroutine exec_regularize
296 implicit real*8 (a-h,o-z)
301 include 'COMMON.SETUP'
302 include 'COMMON.TIME1'
303 include 'COMMON.INTERACT'
304 include 'COMMON.NAMES'
306 include 'COMMON.HEADER'
307 include 'COMMON.CONTROL'
308 include 'COMMON.CONTACTS'
309 include 'COMMON.CHAIN'
311 include 'COMMON.IOUNITS'
312 include 'COMMON.FFIELD'
313 include 'COMMON.REMD'
315 include 'COMMON.SBRIDGE'
316 double precision energy(0:n_ene)
321 call regularize(nct-nnt+1,etot,rms,cref(1,nnt),iretcode)
322 call etotal(energy(0))
323 energy(0)=energy(0)-energy(14)
325 call enerprint(energy(0))
327 call briefout(0,etot)
328 if (outpdb) call pdbout(etot,titel(:32),ipdb)
329 if (outmol2) call mol2out(etot,titel(:32))
330 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
331 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
334 c---------------------------------------------------------------------------
335 subroutine exec_thread
340 include "COMMON.SETUP"
344 c---------------------------------------------------------------------------
346 implicit real*8 (a-h,o-z)
348 character*10 nodeinfo
349 double precision varia(maxvar)
353 include "COMMON.SETUP"
354 include 'COMMON.CONTROL'
358 if (modecalc.eq.3) then
364 if (modecalc.eq.3) then
375 c---------------------------------------------------------------------------
376 subroutine exec_mult_eeval_or_minim
377 implicit real*8 (a-h,o-z)
381 dimension muster(mpi_status_size)
383 include 'COMMON.SETUP'
384 include 'COMMON.TIME1'
385 include 'COMMON.INTERACT'
386 include 'COMMON.NAMES'
388 include 'COMMON.HEADER'
389 include 'COMMON.CONTROL'
390 include 'COMMON.CONTACTS'
391 include 'COMMON.CHAIN'
393 include 'COMMON.IOUNITS'
394 include 'COMMON.FFIELD'
395 include 'COMMON.REMD'
397 include 'COMMON.SBRIDGE'
398 double precision varia(maxvar)
400 double precision energy(0:n_ene)
413 call xdrfopen_(ixdrf,intinname, "r", iret)
415 call xdrfopen(ixdrf,intinname, "r", iret)
418 open(intin,file=intinname,status='old')
420 write (istat,'(a5,30a12)')"# ",
421 & (wname(print_order(i)),i=1,nprint_ene)
423 write (istat,'(a5,30a12)')"# ",
424 & (ename(print_order(i)),i=1,nprint_ene),
425 & "ETOT total","RMSD","nat.contact","nnt.contact","cont.order"
427 write (istat,'(a5,30a12)')"# ",
428 & (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
434 call read_cx(ixdrf,*1100)
437 read (intin,'(i5)',end=1100,err=1100) iconf
438 call read_angles(intin,*11)
439 call geom_to_var(nvar,varia)
442 write (iout,'(/a,i7)') 'Conformation #',iconf
443 call etotal(energy(0))
444 call briefout(iconf,energy(0))
445 call enerprint(energy(0))
448 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
449 write (istat,'(i5,30(f12.3))') iconf,
450 & (energy(print_order(i)),i=1,nprint_ene),etot,
451 & rms,frac,frac_nn,co
454 write (istat,'(i5,30(f12.3))') iconf,
455 & (energy(print_order(i)),i=1,nprint_ene),etot
471 if (mm.lt.nodes) then
473 call read_cx(ixdrf,*11)
475 call geom_to_var(nvar,varia)
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 call read_cx(ixdrf,*11)
542 call geom_to_var(nvar,varia)
544 read (intin,'(i5)',end=11,err=11) iconf
545 call read_angles(intin,*11)
546 call geom_to_var(nvar,varia)
550 write (iout,*) 'Conformation #',iconf,' read'
558 call mpi_send(ind,6,mpi_integer,man,idint,CG_COMM,
560 call mpi_send(varia,nvar,mpi_double_precision,man,
561 * idreal,CG_COMM,ierr)
562 call mpi_send(ene0,1,mpi_double_precision,man,
563 * idreal,CG_COMM,ierr)
564 nf_mcmf=nf_mcmf+ind(4)
570 call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,
571 * CG_COMM,muster,ierr)
572 man=muster(mpi_source)
573 call mpi_recv(varia,nvar,mpi_double_precision,
574 * man,idreal,CG_COMM,muster,ierr)
576 * mpi_double_precision,man,idreal,
577 * CG_COMM,muster,ierr)
578 call mpi_recv(ene0,1,
579 * mpi_double_precision,man,idreal,
580 * CG_COMM,muster,ierr)
582 call var_to_geom(nvar,varia)
584 call etotal(energy(0))
588 write (iout,*) 'Conformation #',iconf," sumsl return code ",
592 call enerprint(energy(0))
593 call briefout(it,etot)
595 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
596 write (istat,'(i5,30(f12.3))') iconf,
597 & (energy(print_order(i)),i=1,nprint_ene),etot,
598 & rms,frac,frac_nn,co
600 write (istat,'(i5,30(f12.3))') iconf,
601 & (energy(print_order(i)),i=1,nprint_ene),etot
613 call mpi_send(ind,6,mpi_integer,i,idint,CG_COMM,
618 open(intin,file=intinname,status='old')
619 write (istat,'(a5,20a12)')"# ",
620 & (wname(print_order(i)),i=1,nprint_ene)
621 write (istat,'("# ",20(1pe12.4))')
622 & (weights(print_order(i)),i=1,nprint_ene)
624 write (istat,'(a5,20a12)')"# ",
625 & (ename(print_order(i)),i=1,nprint_ene),
626 & "ETOT total","RMSD","nat.contact","nnt.contact"
628 write (istat,'(a5,14a12)')"# ",
629 & (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
633 read (intin,'(e15.10,e15.5)',end=1100,err=1100) time,ene
634 call read_x(intin,*11)
635 call int_from_cart1(.false.)
637 read (intin,'(i5)',end=1100,err=1100) iconf
638 call read_angles(intin,*11)
639 call geom_to_var(nvar,varia)
642 write (iout,'(a,i7)') 'Conformation #',iconf
643 if (minim) call minimize(etot,varia,iretcode,nfun)
644 call etotal(energy(0))
647 call enerprint(energy(0))
648 if (minim) call briefout(it,etot)
650 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
651 write (istat,'(i5,18(f12.3))') iconf,
652 & (energy(print_order(i)),i=1,nprint_ene),
653 & etot,rms,frac,frac_nn,co
656 write (istat,'(i5,14(f12.3))') iconf,
657 & (energy(print_order(i)),i=1,nprint_ene),etot
665 c---------------------------------------------------------------------------
666 subroutine exec_checkgrad
667 implicit real*8 (a-h,o-z)
672 include 'COMMON.SETUP'
673 include 'COMMON.TIME1'
674 include 'COMMON.INTERACT'
675 include 'COMMON.NAMES'
677 include 'COMMON.HEADER'
678 include 'COMMON.CONTROL'
679 include 'COMMON.CONTACTS'
680 include 'COMMON.CHAIN'
682 include 'COMMON.IOUNITS'
683 include 'COMMON.FFIELD'
684 include 'COMMON.REMD'
686 include 'COMMON.SBRIDGE'
688 double precision energy(0:max_ene)
690 c vbld(i)=vbld(i)+ran_number(-0.1d0,0.1d0)
691 c if (itype(i).ne.10)
692 c & vbld(i+nres)=vbld(i+nres)+ran_number(-0.001d0,0.001d0)
694 if (indpdb.eq.0) call chainbuild
697 c dc(j,i)=dc(j,i)+ran_number(-0.2d0,0.2d0)
701 c if (itype(i).ne.10) then
703 c dc(j,i+nres)=dc(j,i+nres)+ran_number(-0.2d0,0.2d0)
708 c dc(j,0)=ran_number(-0.2d0,0.2d0)
715 call rescale_weights(t_bath)
720 call etotal(energy(0))
722 call enerprint(energy(0))
723 write (iout,*) "Uconst",Uconst," Uconst_back",uconst_back
724 print *,'icheckgrad=',icheckgrad
725 goto (10,20,30) icheckgrad
726 10 call check_ecartint
728 20 call check_cartgrad
733 c---------------------------------------------------------------------------
740 c---------------------------------------------------------------------------
746 include 'COMMON.IOUNITS'
747 C Conformational Space Annealling programmed by Jooyoung Lee.
748 C This method works only with parallel machines!
751 write (iout,*) "CSA is not supported in this version"
753 csa write (iout,*) "CSA works on parallel machines only"
754 write (iout,*) "CSA is not supported in this version"
758 c---------------------------------------------------------------------------
759 subroutine exec_softreg
760 implicit real*8 (a-h,o-z)
762 include 'COMMON.IOUNITS'
763 include 'COMMON.CONTROL'
764 double precision energy(0:max_ene)
765 logical debug /.false./
767 call etotal(energy(0))
768 call enerprint(energy(0))
769 if (.not.lsecondary) then
770 write(iout,*) 'Calling secondary structure recognition'
771 call secondary2(debug)
773 write(iout,*) 'Using secondary structure supplied in pdb'
778 call etotal(energy(0))
780 call enerprint(energy(0))
782 call briefout(0,etot)
783 call secondary2(.true.)
784 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)