1 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
5 ! Program to carry out conformational search of proteins in an united-residue !
8 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
13 use control, only:tcpu
14 use io_base, only:ilen
15 use geometry, only:chainbuild
16 use control, only:dajczas
19 use compare, only: test
25 ! implicit real*8 (a-h,o-z)
26 ! include 'DIMENSIONS'
28 use MPI_data ! include 'COMMON.SETUP'
33 use MPI_data, only: me,king
36 ! include 'COMMON.TIME1'
37 ! include 'COMMON.INTERACT'
38 ! include 'COMMON.NAMES'
39 ! include 'COMMON.GEO'
40 ! include 'COMMON.HEADER'
41 ! include 'COMMON.CONTROL'
42 ! include 'COMMON.CONTACTS'
43 ! include 'COMMON.CHAIN'
44 ! include 'COMMON.VAR'
45 ! include 'COMMON.IOUNITS'
46 ! include 'COMMON.FFIELD'
47 ! include 'COMMON.REMD'
49 ! include 'COMMON.SBRIDGE'
51 real(kind=8) :: hrtime,mintime,sectime
52 character(len=64) :: text_mode_calc(-2:14)
53 text_mode_calc(-2) = 'test'
54 text_mode_calc(-1) = 'cos'
55 text_mode_calc(0) = 'Energy evaluation or minimization'
56 text_mode_calc(1) = 'Regularization of PDB structure'
57 text_mode_calc(2) = 'Threading of a sequence on PDB structures'
58 text_mode_calc(3) = 'Monte Carlo (with minimization) '
59 text_mode_calc(4) = 'Energy minimization of multiple conformations'
60 text_mode_calc(5) = 'Checking energy gradient'
61 text_mode_calc(6) = 'Entropic sampling Monte Carlo (with minimization)'
62 text_mode_calc(7) = 'Energy map'
63 text_mode_calc(8) = 'CSA calculations'
64 text_mode_calc(9) = 'Not used 9'
65 text_mode_calc(10) = 'Not used 10'
66 text_mode_calc(11) = 'Soft regularization of PDB structure'
67 text_mode_calc(12) = 'Mesoscopic molecular dynamics (MD) '
68 text_mode_calc(13) = 'Not used 13'
69 text_mode_calc(14) = 'Replica exchange molecular dynamics (REMD)'
71 ! call memmon_print_usage()
75 write(iout,*)'### LAST MODIFIED 09/03/15 15:32PM by EL'
76 if (me.eq.king) call cinfo
77 ! Read force field parameters and job setup data
80 write (iout,*) "After readrtns"
83 if (me.eq.king .or. .not. out1file) then
84 write (iout,'(2a/)') &
85 text_mode_calc(modecalc)(:ilen(text_mode_calc(modecalc))), &
87 if (minim) write (iout,'(a)') &
88 'Conformations will be energy-minimized.'
89 write (iout,'(80(1h*)/)')
93 if (modecalc.eq.-2) then
96 else if (modecalc.eq.-1) then
97 write(iout,*) "call check_sc_map next"
101 !elwrite(iout,*)"!!!!!!!!!!!!!!!!! in unres"
104 if (fg_rank.gt.0) then
105 ! Fine-grain slaves just do energy and gradient components.
106 call ergastulum ! slave workhouse in Latin
109 if (modecalc.eq.0) then
110 call exec_eeval_or_minim
111 else if (modecalc.eq.1) then
113 else if (modecalc.eq.2) then
115 else if (modecalc.eq.3 .or. modecalc .eq.6) then
117 else if (modecalc.eq.4) then
118 call exec_mult_eeval_or_minim
119 else if (modecalc.eq.5) then
121 else if (ModeCalc.eq.7) then
123 else if (ModeCalc.eq.8) then
125 else if (modecalc.eq.11) then
127 else if (modecalc.eq.12) then
130 else if (modecalc.eq.14) then
133 write (iout,'(a)') 'This calculation type is not supported',&
140 if (fg_rank.eq.0) call finish_task
141 ! call memmon_print_usage()
143 call print_detailed_timing
145 call MPI_Finalize(ierr)
148 call dajczas(tcpu(),hrtime,mintime,sectime)
149 stop '********** Program terminated normally.'
153 !-----------------------------------------------------------------------------
155 !-----------------------------------------------------------------------------
157 use MPI_data !include 'COMMON.SETUP'
158 use control_data !include 'COMMON.CONTROL'
159 use geometry, only:chainbuild,chainbuild_cart
161 use io_units !include 'COMMON.IOUNITS'
162 use compare, only:alloc_compare_arrays
165 ! include 'DIMENSIONS'
171 print *,'After MD alloc'
172 call alloc_compare_arrays
173 print *,'After compare alloc'
174 if (me.eq.king .or. .not. out1file) &
175 write (iout,*) "Calling chainbuild"
183 end subroutine exec_MD
184 !---------------------------------------------------------------------------
185 subroutine exec_MREMD
186 use MPI_data !include 'COMMON.SETUP'
187 use control_data !include 'COMMON.CONTROL'
188 use io_units !include 'COMMON.IOUNITS'
190 use REMD_data !include 'COMMON.REMD'
191 use geometry, only:chainbuild
193 use compare, only:alloc_compare_arrays
196 ! include 'DIMENSIONS'
203 call alloc_MREMD_arrays
204 call alloc_compare_arrays
205 ! if (me.eq.king .or. .not. out1file) &
206 ! write (iout,*) "Calling chainbuild"
208 if (me.eq.king .or. .not. out1file) &
209 write (iout,*) "Calling REMD",remd_mlist,nrep
219 end subroutine exec_MREMD
220 !-----------------------------------------------------------------------------
221 subroutine exec_eeval_or_minim
222 use MPI_data !include 'COMMON.SETUP'
223 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
224 use io_units !include 'COMMON.IOUNITS'
226 ! use energy !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
227 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
228 ! use REMD !include 'COMMON.REMD'
229 ! use MD !include 'COMMON.MD'
232 use MD_data, only: iset
234 use geometry, only:chainbuild
236 use compare, only:alloc_compare_arrays,hairpin,secondary2,rms_nac_nnc
237 use minimm, only:minimize,minim_dc,sc_move
241 ! implicit real*8 (a-h,o-z)
242 ! include 'DIMENSIONS'
247 !el common /srutu/ icall
248 real(kind=8) :: energy_(0:n_ene)
249 real(kind=8) :: energy_long(0:n_ene),energy_short(0:n_ene)
250 real(kind=8) :: varia(6*nres) !(maxvar) (maxvar=6*maxres)
251 real(kind=8) :: time00, evals, etota, etot, time_ene, time1
252 integer :: nharp,nft_sc,iretcode,nfun
253 integer,dimension(4,nres) :: iharp !(4,nres/3)(4,maxres/3)
255 real(kind=8) :: rms,frac,frac_nn,co
258 if (iset.eq.0) iset=1
259 call alloc_compare_arrays
260 if ((indpdb.eq.0).and.(.not.read_cart)) then
262 write(iout,*) 'Warning: Calling chainbuild'
268 ! write(iout,*)"in exec_eeval or minimim",split_ene
270 ! write(iout,*)"cccccc",j,(c(i,j),i=1,3)
271 ! write(iout,*)"dcccccc",j,(dc(i,j),i=1,3)
274 ! write(iout,*)"in exec_eeval or minimim"
276 print *,"Processor",myrank," after chainbuild"
279 call etotal_long(energy_long)
280 write (iout,*) "Printing long range energy"
281 call enerprint(energy_long)
283 call etotal_short(energy_short)
284 write (iout,*) "Printing short range energy"
285 call enerprint(energy_short)
287 energy_(i)=energy_long(i)+energy_short(i)
288 write (iout,*) i,energy_long(i),energy_short(i),energy_(i)
290 write (iout,*) "Printing long+short range energy"
291 call enerprint(energy_)
296 time_ene=MPI_Wtime()-time00
298 write (iout,*) "Time for energy evaluation",time_ene
299 print *,"after etotal"
302 call enerprint(energy_)
304 call hairpin(.true.,nharp,iharp)
305 call secondary2(.true.)
309 print *,"overlap",searchsc,overlapsc
311 print *, 'Calling OVERLAP_SC'
312 call overlap_sc(fail)
316 call sc_move(2,nres-1,10,1d10,nft_sc,etot)
317 print *,'SC_move',nft_sc,etot
318 write(iout,*) 'SC_move',nft_sc,etot
322 print *, 'Calling MINIM_DC'
326 ! call check_ecartint !el
327 call minim_dc(etot,iretcode,nfun)
328 ! call check_ecartint !el
330 if (indpdb.ne.0) then
332 write(iout,*) 'Warning: Calling chainbuild'
335 call geom_to_var(nvar,varia)
336 print *,'Calling MINIMIZE.'
341 ! call exec_checkgrad !el
342 call minimize(etot,varia,iretcode,nfun)
344 ! call exec_checkgrad !el
346 print *,'SUMSL return code is',iretcode,' eval ',nfun
348 evals=nfun/(MPI_WTIME()-time1)
350 print *,'# eval/s',evals
351 print *,'refstr=',refstr
352 ! call hairpin(.true.,nharp,iharp)
353 ! call secondary2(.true.)
356 call enerprint(energy_)
359 call briefout(0,etot)
360 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
361 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
362 write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
363 write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
365 print *,'refstr=',refstr,frac,frac_nn,co
366 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
367 print *,"after rms_nac_ncc"
368 call briefout(0,etot)
370 if (outpdb) call pdbout(etot,titel(:32),ipdb)
371 if (outmol2) call mol2out(etot,titel(:32))
372 write(iout,*) "after exec_eeval_or_minim"
374 end subroutine exec_eeval_or_minim
375 !-----------------------------------------------------------------------------
376 subroutine exec_regularize
377 ! use MPI_data !include 'COMMON.SETUP'
378 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
379 use io_units !include 'COMMON.IOUNITS'
381 use energy_data !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
382 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
383 ! use REMD !include 'COMMON.REMD'
384 ! use MD !include 'COMMON.MD'
388 ! implicit real*8 (a-h,o-z)
389 ! include 'DIMENSIONS'
393 real(kind=8) :: energy_(0:n_ene)
395 real(kind=8) :: rms,frac,frac_nn,co
398 call alloc_compare_arrays
402 call regularize(nct-nnt+1,etot,rms,cref(1,nnt,1),iretcode)
404 energy_(0)=energy_(0)-energy_(14)
406 call enerprint(energy_)
408 call briefout(0,etot)
409 if (outpdb) call pdbout(etot,titel(:32),ipdb)
410 if (outmol2) call mol2out(etot,titel(:32))
411 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
412 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
414 end subroutine exec_regularize
415 !-----------------------------------------------------------------------------
416 subroutine exec_thread
417 ! use MPI_data !include 'COMMON.SETUP'
420 ! include 'DIMENSIONS'
424 call alloc_compare_arrays
427 end subroutine exec_thread
428 !-----------------------------------------------------------------------------
430 ! use MPI_data !include 'COMMON.SETUP'
431 use control_data !include 'COMMON.CONTROL'
436 ! implicit real*8 (a-h,o-z)
437 ! include 'DIMENSIONS'
438 character(len=10) :: nodeinfo
439 real(kind=8) :: varia(6*nres) !(maxvar) (maxvar=6*maxres)
444 call alloc_MCM_arrays
448 if (modecalc.eq.3) then
454 if (modecalc.eq.3) then
464 end subroutine exec_MC
465 !-----------------------------------------------------------------------------
466 subroutine exec_mult_eeval_or_minim
467 use MPI_data !include 'COMMON.SETUP'
468 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
469 use io_units !include 'COMMON.IOUNITS'
471 use energy_data !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
472 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
473 ! use REMD !include 'COMMON.REMD'
474 ! use MD !include 'COMMON.MD'
476 use geometry, only:chainbuild,geom_to_var,int_from_cart1,var_to_geom
477 use energy, only:etotal,enerprint
478 use compare, only:rms_nac_nnc
479 use minimm, only:minimize!,minim_mcmf
480 ! implicit real*8 (a-h,o-z)
481 ! include 'DIMENSIONS'
483 use minimm, only:minim_mcmf
486 integer :: ierror,ierr
488 real(kind=8),dimension(mpi_status_size) :: muster
492 real(kind=8) :: varia(6*nres) !(maxvar) (maxvar=6*maxres)
493 integer,dimension(6) :: ind
494 real(kind=8) :: energy_(0:n_ene)
496 real(kind=8) :: etot,ene0
497 integer :: mm,imm,nft,n,iconf,nmin,i,iretcode,nfun,it,&
499 real(kind=8) :: rms,frac,frac_nn,co,time,ene
509 open(intin,file=intinname,status='old')
510 write (istat,'(a5,20a12)')"# ",&
511 (wname(print_order(i)),i=1,nprint_ene)
513 write (istat,'(a5,20a12)')"# ",&
514 (ename(print_order(i)),i=1,nprint_ene),&
515 "ETOT total","RMSD","nat.contact","nnt.contact"
517 write (istat,'(a5,20a12)')"# ",&
518 (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
524 read (intin,'(e15.10,e15.5)',end=1100,err=1100) time,ene
525 call read_x(intin,*11)
527 ! Broadcast the order to compute internal coordinates to the slaves.
529 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
531 call int_from_cart1(.false.)
533 read (intin,'(i5)',end=1100,err=1100) iconf
534 call read_angles(intin,*11)
535 call geom_to_var(nvar,varia)
536 write(iout,*) 'Warning: Calling chainbuild1'
539 write (iout,'(a,i7)') 'Conformation #',iconf
541 call briefout(iconf,energy_(0))
542 call enerprint(energy_)
545 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
546 write (istat,'(i5,20(f12.3))') iconf,&
547 (energy_(print_order(i)),i=1,nprint_ene),etot,&
551 write (istat,'(i5,16(f12.3))') iconf,&
552 (energy_(print_order(i)),i=1,nprint_ene),etot
568 if (mm.lt.nodes) then
570 read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
571 call read_x(intin,*11)
573 ! Broadcast the order to compute internal coordinates to the slaves.
575 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
577 call int_from_cart1(.false.)
579 read (intin,'(i5)',end=11,err=11) iconf
580 call read_angles(intin,*11)
581 call geom_to_var(nvar,varia)
582 write(iout,*) 'Warning: Calling chainbuild2'
585 write (iout,'(a,i7)') 'Conformation #',iconf
595 call mpi_send(ind,6,mpi_integer,mm,idint,CG_COMM,&
597 call mpi_send(varia,nvar,mpi_double_precision,mm,&
599 call mpi_send(ene0,1,mpi_double_precision,mm,&
601 ! print *,'task ',n,' sent to worker ',mm,nvar
603 call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,&
605 man=muster(mpi_source)
606 ! print *,'receiving result from worker ',man,' (',iii1,iii,')'
607 call mpi_recv(varia,nvar,mpi_double_precision,&
608 man,idreal,CG_COMM,muster,ierr)
609 call mpi_recv(ene,1,&
610 mpi_double_precision,man,idreal,&
612 call mpi_recv(ene0,1,&
613 mpi_double_precision,man,idreal,&
615 ! print *,'result received from worker ',man,' sending now'
617 call var_to_geom(nvar,varia)
618 write(iout,*) 'Warning: Calling chainbuild3'
624 write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
627 call enerprint(energy_)
628 call briefout(it,etot)
629 ! if (minim) call briefout(it,etot)
631 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
632 write (istat,'(i5,19(f12.3))') iconf,&
633 (energy_(print_order(i)),i=1,nprint_ene),etot,&
636 write (istat,'(i5,15(f12.3))') iconf,&
637 (energy_(print_order(i)),i=1,nprint_ene),etot
642 read (intin,'(e15.10,e15.5)',end=1101,err=1101) time,ene
643 call read_x(intin,*11)
645 ! Broadcast the order to compute internal coordinates to the slaves.
647 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
649 call int_from_cart1(.false.)
651 read (intin,'(i5)',end=1101,err=1101) iconf
652 call read_angles(intin,*11)
653 call geom_to_var(nvar,varia)
654 write(iout,*) 'Warning: Calling chainbuild4'
665 call mpi_send(ind,6,mpi_integer,man,idint,CG_COMM,&
667 call mpi_send(varia,nvar,mpi_double_precision,man,&
669 call mpi_send(ene0,1,mpi_double_precision,man,&
671 nf_mcmf=nf_mcmf+ind(4)
677 call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,&
679 man=muster(mpi_source)
680 call mpi_recv(varia,nvar,mpi_double_precision,&
681 man,idreal,CG_COMM,muster,ierr)
682 call mpi_recv(ene,1,&
683 mpi_double_precision,man,idreal,&
685 call mpi_recv(ene0,1,&
686 mpi_double_precision,man,idreal,&
689 call var_to_geom(nvar,varia)
690 write(iout,*) 'Warning: Calling chainbuild5'
696 write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
699 call enerprint(energy_)
700 call briefout(it,etot)
702 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
703 write (istat,'(i5,19(f12.3))') iconf,&
704 (energy_(print_order(i)),i=1,nprint_ene),etot,&
707 write (istat,'(i5,15(f12.3))') iconf,&
708 (energy_(print_order(i)),i=1,nprint_ene),etot
720 call mpi_send(ind,6,mpi_integer,i,idint,CG_COMM,&
725 open(intin,file=intinname,status='old')
726 write (istat,'(a5,20a12)')"# ",&
727 (wname(print_order(i)),i=1,nprint_ene)
728 write (istat,'("# ",20(1pe12.4))') &
729 (weights(print_order(i)),i=1,nprint_ene)
731 write (istat,'(a5,20a12)')"# ",&
732 (ename(print_order(i)),i=1,nprint_ene),&
733 "ETOT total","RMSD","nat.contact","nnt.contact"
735 write (istat,'(a5,14a12)')"# ",&
736 (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
740 read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
741 call read_x(intin,*11)
743 ! Broadcast the order to compute internal coordinates to the slaves.
745 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
747 call int_from_cart1(.false.)
749 read (intin,'(i5)',end=11,err=11) iconf
750 call read_angles(intin,*11)
751 call geom_to_var(nvar,varia)
752 write(iout,*) 'Warning: Calling chainbuild5'
755 write (iout,'(a,i7)') 'Conformation #',iconf
756 if (minim) call minimize(etot,varia,iretcode,nfun)
760 call enerprint(energy_)
761 if (minim) call briefout(it,etot)
763 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
764 write (istat,'(i5,18(f12.3))') iconf,&
765 (energy_(print_order(i)),i=1,nprint_ene),&
766 etot,rms,frac,frac_nn,co
769 write (istat,'(i5,14(f12.3))') iconf,&
770 (energy_(print_order(i)),i=1,nprint_ene),etot
776 end subroutine exec_mult_eeval_or_minim
777 !-----------------------------------------------------------------------------
778 subroutine exec_checkgrad
779 ! use MPI_data !include 'COMMON.SETUP'
780 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
781 use io_units !include 'COMMON.IOUNITS'
782 !el use energy_data, only:icall !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
783 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
784 ! use REMD !include 'COMMON.REMD'
785 use MD_data !include 'COMMON.MD'
786 use io_base, only:intout
787 use io_config, only:read_fragments
792 ! implicit real*8 (a-h,o-z)
793 ! include 'DIMENSIONS'
798 !el common /srutu/ icall
799 real(kind=8) :: energy_(0:max_ene)
803 ! vbld(i)=vbld(i)+ran_number(-0.1d0,0.1d0)
804 ! if (itype(i).ne.10)
805 ! & vbld(i+nres)=vbld(i+nres)+ran_number(-0.001d0,0.001d0)
807 if (indpdb.eq.0) then
808 write(iout,*) 'Warning: Calling chainbuild'
813 ! dc(j,i)=dc(j,i)+ran_number(-0.2d0,0.2d0)
817 ! if (itype(i).ne.10) then
819 ! dc(j,i+nres)=dc(j,i+nres)+ran_number(-0.2d0,0.2d0)
824 ! dc(j,0)=ran_number(-0.2d0,0.2d0)
836 write (iout,*) "before etotal"
837 call etotal(energy_(0))
839 call enerprint(energy_(0))
840 write (iout,*) "Uconst",Uconst," Uconst_back",uconst_back
841 print *,'icheckgrad=',icheckgrad
842 goto (10,20,30) icheckgrad
843 10 call check_ecartint
846 20 call check_cartgrad
851 end subroutine exec_checkgrad
852 !-----------------------------------------------------------------------------
856 use io_config, only:map_read
859 call alloc_map_arrays
863 end subroutine exec_map
864 !-----------------------------------------------------------------------------
867 use io_units !include 'COMMON.IOUNITS'
873 ! include 'DIMENSIONS'
874 ! Conformational Space Annealling programmed by Jooyoung Lee.
875 ! This method works only with parallel machines!
877 call alloc_CSA_arrays
880 write (iout,*) "CSA works on parallel machines only"
883 end subroutine exec_CSA
884 !-----------------------------------------------------------------------------
885 subroutine exec_softreg
886 use io_units !include 'COMMON.IOUNITS'
887 use control_data !include 'COMMON.CONTROL'
889 use io_base, only:intout,briefout
890 use geometry, only:chainbuild
894 ! include 'DIMENSIONS'
895 real(kind=8) :: energy_(0:n_ene)
897 real(kind=8) :: rms,frac,frac_nn,co,etot
900 call alloc_compare_arrays
901 write(iout,*) 'Warning: Calling chainbuild'
904 call enerprint(energy_)
905 if (.not.lsecondary) then
906 write(iout,*) 'Calling secondary structure recognition'
907 call secondary2(debug)
909 write(iout,*) 'Using secondary structure supplied in pdb'
916 call enerprint(energy_)
918 call briefout(0,etot)
919 call secondary2(.true.)
920 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
922 end subroutine exec_softreg
923 !-----------------------------------------------------------------------------
925 !-----------------------------------------------------------------------------
927 subroutine ergastulum
929 ! implicit real*8 (a-h,o-z)
930 ! include 'DIMENSIONS'
933 use MDyn, only:setup_fricmat
935 use REMD, only:fricmat_mult,ginv_mult
940 ! include 'COMMON.SETUP'
941 ! include 'COMMON.DERIV'
942 ! include 'COMMON.VAR'
943 ! include 'COMMON.IOUNITS'
944 ! include 'COMMON.FFIELD'
945 ! include 'COMMON.INTERACT'
946 ! include 'COMMON.MD'
947 ! include 'COMMON.TIME1'
948 real(kind=8),dimension(6*nres) :: z,d_a_tmp !(maxres6) maxres6=6*maxres
949 real(kind=8) :: edum(0:n_ene),time_order(0:10)
950 !el real(kind=8),dimension(2*nres,2*nres) :: Gcopy !(maxres2,maxres2) maxres2=2*maxres
951 !el common /przechowalnia/ Gcopy
955 real(kind=8) :: time00
956 integer :: iorder,i,j,nres2,ierr,ierror
959 if(.not.allocated(Gcopy)) allocate(Gcopy(nres2,nres2))
961 if(.not.allocated(Ginv)) allocate(Ginv(nres2,nres2)) !(maxres2,maxres2)
964 if(.not.allocated(ng_start)) allocate(ng_start(0:nfgtasks-1))
965 if(.not.allocated(ng_counts)) allocate(ng_counts(0:nfgtasks-1))
966 if(.not.allocated(nginv_counts)) allocate(nginv_counts(0:nfgtasks-1)) !(0:MaxProcs-1)
967 if(.not.allocated(nginv_start)) allocate(nginv_start(0:nfgtasks)) !(0:MaxProcs)
969 if(.not.allocated(fricmat)) allocate(fricmat(nres2,nres2)) !maxres2=2*maxres
971 ! Workers wait for variables and NF, and NFL from the boss
973 do while (iorder.ge.0)
974 ! write (*,*) 'Processor',fg_rank,' CG group',kolor,
975 ! & ' receives order from Master'
977 call MPI_Bcast(iorder,1,MPI_INTEGER,king,FG_COMM,IERR)
978 time_Bcast=time_Bcast+MPI_Wtime()-time00
979 if (icall.gt.4 .and. iorder.ge.0) &
980 time_order(iorder)=time_order(iorder)+MPI_Wtime()-time00
983 ! & 'Processor',fg_rank,' completed receive MPI_BCAST order',iorder
984 if (iorder.eq.0) then
987 ! write (2,*) "After etotal"
988 ! write (2,*) "dimen",dimen," dimen3",dimen3
990 else if (iorder.eq.2) then
992 call etotal_short(edum)
993 ! write (2,*) "After etotal_short"
994 ! write (2,*) "dimen",dimen," dimen3",dimen3
996 else if (iorder.eq.3) then
998 call etotal_long(edum)
999 ! write (2,*) "After etotal_long"
1000 ! write (2,*) "dimen",dimen," dimen3",dimen3
1002 else if (iorder.eq.1) then
1004 ! write (2,*) "After sum_gradient"
1005 ! write (2,*) "dimen",dimen," dimen3",dimen3
1008 else if (iorder.eq.4) then
1009 call ginv_mult(z,d_a_tmp)
1010 else if (iorder.eq.5) then
1011 ! Setup MD things for a slave
1012 dimen=(nct-nnt+1)+nside
1013 dimen1=(nct-nnt)+(nct-nnt+1)
1015 ! write (2,*) "dimen",dimen," dimen3",dimen3
1017 call int_bounds(dimen,igmult_start,igmult_end)
1018 igmult_start=igmult_start-1
1019 call MPI_Allgather(3*igmult_start,1,MPI_INTEGER,&
1020 ng_start(0),1,MPI_INTEGER,FG_COMM,IERROR)
1021 my_ng_count=igmult_end-igmult_start
1022 call MPI_Allgather(3*my_ng_count,1,MPI_INTEGER,ng_counts(0),1,&
1023 MPI_INTEGER,FG_COMM,IERROR)
1024 write (2,*) "ng_start",(ng_start(i),i=0,nfgtasks-1) !sp
1025 ! write (2,*) "ng_counts",(ng_counts(i),i=0,nfgtasks-1)
1026 myginv_ng_count=nres2*my_ng_count !el maxres2
1027 ! write (2,*) "igmult_start",igmult_start," igmult_end",
1028 ! & igmult_end," my_ng_count",my_ng_count
1030 call MPI_Allgather(nres2*igmult_start,1,MPI_INTEGER,& !el maxres2
1031 nginv_start(0),1,MPI_INTEGER,FG_COMM,IERROR)
1032 call MPI_Allgather(myginv_ng_count,1,MPI_INTEGER,&
1033 nginv_counts(0),1,MPI_INTEGER,FG_COMM,IERROR)
1034 ! write (2,*) "nginv_start",(nginv_start(i),i=0,nfgtasks-1)
1035 ! write (2,*) "nginv_counts",(nginv_counts(i),i=0,nfgtasks-1)
1037 ! call MPI_Barrier(FG_COMM,IERROR)
1039 call MPI_Scatterv(ginv(1,1),nginv_counts(0),&
1040 nginv_start(0),MPI_DOUBLE_PRECISION,gcopy(1,1),&
1041 myginv_ng_count,MPI_DOUBLE_PRECISION,king,FG_COMM,IERR)
1043 time_scatter_ginv=time_scatter_ginv+MPI_Wtime()-time00
1046 do j=1,2*my_ng_count
1047 ginv(j,i)=gcopy(i,j)
1050 ! write (2,*) "dimen",dimen," dimen3",dimen3
1051 ! write (2,*) "End MD setup"
1053 ! write (iout,*) "My chunk of ginv_block"
1054 ! call MATOUT2(my_ng_count,dimen3,maxres2,maxers2,ginv_block)
1056 else if (iorder.eq.6) then
1057 call int_from_cart1(.false.)
1058 else if (iorder.eq.7) then
1059 call chainbuild_cart
1060 else if (iorder.eq.8) then
1063 else if (iorder.eq.9) then
1064 call fricmat_mult(z,d_a_tmp)
1066 else if (iorder.eq.10) then
1070 write (*,*) 'Processor',fg_rank,' CG group',kolor,&
1071 ' absolute rank',myrank,' leves ERGASTULUM.'
1072 write(*,*)'Processor',fg_rank,' wait times for respective orders',&
1073 (' order[',i,']',time_order(i),i=0,10)
1075 end subroutine ergastulum