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
81 if (me.eq.king .or. .not. out1file) then
82 write (iout,'(2a/)') &
83 text_mode_calc(modecalc)(:ilen(text_mode_calc(modecalc))), &
85 if (minim) write (iout,'(a)') &
86 'Conformations will be energy-minimized.'
87 write (iout,'(80(1h*)/)')
91 if (modecalc.eq.-2) then
94 else if (modecalc.eq.-1) then
95 write(iout,*) "call check_sc_map next"
99 !elwrite(iout,*)"!!!!!!!!!!!!!!!!! in unres"
102 if (fg_rank.gt.0) then
103 ! Fine-grain slaves just do energy and gradient components.
104 call ergastulum ! slave workhouse in Latin
107 if (modecalc.eq.0) then
108 call exec_eeval_or_minim
109 else if (modecalc.eq.1) then
111 else if (modecalc.eq.2) then
113 else if (modecalc.eq.3 .or. modecalc .eq.6) then
115 else if (modecalc.eq.4) then
116 call exec_mult_eeval_or_minim
117 else if (modecalc.eq.5) then
119 else if (ModeCalc.eq.7) then
121 else if (ModeCalc.eq.8) then
123 else if (modecalc.eq.11) then
125 else if (modecalc.eq.12) then
127 else if (modecalc.eq.14) then
130 write (iout,'(a)') 'This calculation type is not supported',&
137 if (fg_rank.eq.0) call finish_task
138 ! call memmon_print_usage()
140 call print_detailed_timing
142 call MPI_Finalize(ierr)
145 call dajczas(tcpu(),hrtime,mintime,sectime)
146 stop '********** Program terminated normally.'
150 !-----------------------------------------------------------------------------
152 !-----------------------------------------------------------------------------
154 use MPI_data !include 'COMMON.SETUP'
155 use control_data !include 'COMMON.CONTROL'
156 use geometry, only:chainbuild
158 use io_units !include 'COMMON.IOUNITS'
161 ! include 'DIMENSIONS'
167 ! if (me.eq.king .or. .not. out1file) &
168 ! write (iout,*) "Calling chainbuild"
172 end subroutine exec_MD
173 !---------------------------------------------------------------------------
174 subroutine exec_MREMD
175 use MPI_data !include 'COMMON.SETUP'
176 use control_data !include 'COMMON.CONTROL'
177 use io_units !include 'COMMON.IOUNITS'
179 use REMD_data !include 'COMMON.REMD'
180 use geometry, only:chainbuild
184 ! include 'DIMENSIONS'
191 call alloc_MREMD_arrays
193 ! if (me.eq.king .or. .not. out1file) &
194 ! write (iout,*) "Calling chainbuild"
196 if (me.eq.king .or. .not. out1file) &
197 write (iout,*) "Calling REMD"
207 end subroutine exec_MREMD
208 !-----------------------------------------------------------------------------
209 subroutine exec_eeval_or_minim
210 use MPI_data !include 'COMMON.SETUP'
211 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
212 use io_units !include 'COMMON.IOUNITS'
214 ! use energy !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
215 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
216 ! use REMD !include 'COMMON.REMD'
217 ! use MD !include 'COMMON.MD'
222 use geometry, only:chainbuild
224 use compare, only:alloc_compare_arrays,hairpin,secondary2,rms_nac_nnc
225 use minimm, only:minimize,minim_dc,sc_move
229 ! implicit real*8 (a-h,o-z)
230 ! include 'DIMENSIONS'
235 !el common /srutu/ icall
236 real(kind=8) :: energy_(0:n_ene)
237 real(kind=8) :: energy_long(0:n_ene),energy_short(0:n_ene)
238 real(kind=8) :: varia(6*nres) !(maxvar) (maxvar=6*maxres)
239 real(kind=8) :: time00, evals, etota, etot, time_ene, time1
240 integer :: nharp,nft_sc,iretcode,nfun
241 integer,dimension(4,nres/3) :: iharp !(4,nres/3)(4,maxres/3)
243 real(kind=8) :: rms,frac,frac_nn,co
246 call alloc_compare_arrays
247 if (indpdb.eq.0) then
249 write(iout,*) 'Warning: Calling chainbuild'
255 ! write(iout,*)"in exec_eeval or minimim",split_ene
257 ! write(iout,*)"cccccc",j,(c(i,j),i=1,3)
258 ! write(iout,*)"dcccccc",j,(dc(i,j),i=1,3)
261 ! write(iout,*)"in exec_eeval or minimim"
263 print *,"Processor",myrank," after chainbuild"
266 call etotal_long(energy_long)
267 write (iout,*) "Printing long range energy"
268 call enerprint(energy_long)
270 call etotal_short(energy_short)
271 write (iout,*) "Printing short range energy"
272 call enerprint(energy_short)
274 energy_(i)=energy_long(i)+energy_short(i)
275 write (iout,*) i,energy_long(i),energy_short(i),energy_(i)
277 write (iout,*) "Printing long+short range energy"
278 call enerprint(energy_)
283 time_ene=MPI_Wtime()-time00
285 write (iout,*) "Time for energy evaluation",time_ene
286 print *,"after etotal"
289 call enerprint(energy_)
290 call hairpin(.true.,nharp,iharp)
291 call secondary2(.true.)
295 print *, 'Calling OVERLAP_SC'
296 call overlap_sc(fail)
300 call sc_move(2,nres-1,10,1d10,nft_sc,etot)
301 print *,'SC_move',nft_sc,etot
302 write(iout,*) 'SC_move',nft_sc,etot
306 print *, 'Calling MINIM_DC'
310 ! call check_ecartint !el
311 call minim_dc(etot,iretcode,nfun)
312 ! call check_ecartint !el
314 if (indpdb.ne.0) then
316 write(iout,*) 'Warning: Calling chainbuild'
319 call geom_to_var(nvar,varia)
320 print *,'Calling MINIMIZE.'
325 ! call exec_checkgrad !el
326 call minimize(etot,varia,iretcode,nfun)
328 ! call exec_checkgrad !el
330 print *,'SUMSL return code is',iretcode,' eval ',nfun
332 evals=nfun/(MPI_WTIME()-time1)
334 print *,'# eval/s',evals
335 print *,'refstr=',refstr
336 call hairpin(.true.,nharp,iharp)
337 call secondary2(.true.)
340 call enerprint(energy_)
343 call briefout(0,etot)
344 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
345 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
346 write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
347 write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
349 print *,'refstr=',refstr
350 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
351 call briefout(0,etot)
353 if (outpdb) call pdbout(etot,titel(:32),ipdb)
354 if (outmol2) call mol2out(etot,titel(:32))
355 !elwrite(iout,*) "after exec_eeval_or_minim"
357 end subroutine exec_eeval_or_minim
358 !-----------------------------------------------------------------------------
359 subroutine exec_regularize
360 ! use MPI_data !include 'COMMON.SETUP'
361 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
362 use io_units !include 'COMMON.IOUNITS'
364 use energy_data !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
365 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
366 ! use REMD !include 'COMMON.REMD'
367 ! use MD !include 'COMMON.MD'
371 ! implicit real*8 (a-h,o-z)
372 ! include 'DIMENSIONS'
376 real(kind=8) :: energy_(0:n_ene)
378 real(kind=8) :: rms,frac,frac_nn,co
381 call alloc_compare_arrays
385 call regularize(nct-nnt+1,etot,rms,cref(1,nnt,1),iretcode)
387 energy_(0)=energy_(0)-energy_(14)
389 call enerprint(energy_)
391 call briefout(0,etot)
392 if (outpdb) call pdbout(etot,titel(:32),ipdb)
393 if (outmol2) call mol2out(etot,titel(:32))
394 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
395 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
397 end subroutine exec_regularize
398 !-----------------------------------------------------------------------------
399 subroutine exec_thread
400 ! use MPI_data !include 'COMMON.SETUP'
403 ! include 'DIMENSIONS'
407 call alloc_compare_arrays
410 end subroutine exec_thread
411 !-----------------------------------------------------------------------------
413 ! use MPI_data !include 'COMMON.SETUP'
414 use control_data !include 'COMMON.CONTROL'
419 ! implicit real*8 (a-h,o-z)
420 ! include 'DIMENSIONS'
421 character(len=10) :: nodeinfo
422 real(kind=8) :: varia(6*nres) !(maxvar) (maxvar=6*maxres)
427 call alloc_MCM_arrays
431 if (modecalc.eq.3) then
437 if (modecalc.eq.3) then
447 end subroutine exec_MC
448 !-----------------------------------------------------------------------------
449 subroutine exec_mult_eeval_or_minim
450 use MPI_data !include 'COMMON.SETUP'
451 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
452 use io_units !include 'COMMON.IOUNITS'
454 use energy_data !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
455 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
456 ! use REMD !include 'COMMON.REMD'
457 ! use MD !include 'COMMON.MD'
459 use geometry, only:chainbuild,geom_to_var,int_from_cart1,var_to_geom
460 use energy, only:etotal,enerprint
461 use compare, only:rms_nac_nnc
462 use minimm, only:minimize!,minim_mcmf
463 ! implicit real*8 (a-h,o-z)
464 ! include 'DIMENSIONS'
466 use minimm, only:minim_mcmf
469 integer :: ierror,ierr
471 real(kind=8),dimension(mpi_status_size) :: muster
475 real(kind=8) :: varia(6*nres) !(maxvar) (maxvar=6*maxres)
476 integer,dimension(6) :: ind
477 real(kind=8) :: energy_(0:n_ene)
479 real(kind=8) :: etot,ene0
480 integer :: mm,imm,nft,n,iconf,nmin,i,iretcode,nfun,it,&
482 real(kind=8) :: rms,frac,frac_nn,co,time,ene
492 open(intin,file=intinname,status='old')
493 write (istat,'(a5,20a12)')"# ",&
494 (wname(print_order(i)),i=1,nprint_ene)
496 write (istat,'(a5,20a12)')"# ",&
497 (ename(print_order(i)),i=1,nprint_ene),&
498 "ETOT total","RMSD","nat.contact","nnt.contact"
500 write (istat,'(a5,20a12)')"# ",&
501 (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
507 read (intin,'(e15.10,e15.5)',end=1100,err=1100) time,ene
508 call read_x(intin,*11)
510 ! Broadcast the order to compute internal coordinates to the slaves.
512 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
514 call int_from_cart1(.false.)
516 read (intin,'(i5)',end=1100,err=1100) iconf
517 call read_angles(intin,*11)
518 call geom_to_var(nvar,varia)
519 write(iout,*) 'Warning: Calling chainbuild'
522 write (iout,'(a,i7)') 'Conformation #',iconf
524 call briefout(iconf,energy_(0))
525 call enerprint(energy_)
528 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
529 write (istat,'(i5,20(f12.3))') iconf,&
530 (energy_(print_order(i)),i=1,nprint_ene),etot,&
534 write (istat,'(i5,16(f12.3))') iconf,&
535 (energy_(print_order(i)),i=1,nprint_ene),etot
551 if (mm.lt.nodes) then
553 read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
554 call read_x(intin,*11)
556 ! Broadcast the order to compute internal coordinates to the slaves.
558 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
560 call int_from_cart1(.false.)
562 read (intin,'(i5)',end=11,err=11) iconf
563 call read_angles(intin,*11)
564 call geom_to_var(nvar,varia)
565 write(iout,*) 'Warning: Calling chainbuild'
568 write (iout,'(a,i7)') 'Conformation #',iconf
578 call mpi_send(ind,6,mpi_integer,mm,idint,CG_COMM,&
580 call mpi_send(varia,nvar,mpi_double_precision,mm,&
582 call mpi_send(ene0,1,mpi_double_precision,mm,&
584 ! print *,'task ',n,' sent to worker ',mm,nvar
586 call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,&
588 man=muster(mpi_source)
589 ! print *,'receiving result from worker ',man,' (',iii1,iii,')'
590 call mpi_recv(varia,nvar,mpi_double_precision,&
591 man,idreal,CG_COMM,muster,ierr)
592 call mpi_recv(ene,1,&
593 mpi_double_precision,man,idreal,&
595 call mpi_recv(ene0,1,&
596 mpi_double_precision,man,idreal,&
598 ! print *,'result received from worker ',man,' sending now'
600 call var_to_geom(nvar,varia)
601 write(iout,*) 'Warning: Calling chainbuild'
607 write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
610 call enerprint(energy_)
611 call briefout(it,etot)
612 ! if (minim) call briefout(it,etot)
614 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
615 write (istat,'(i5,19(f12.3))') iconf,&
616 (energy_(print_order(i)),i=1,nprint_ene),etot,&
619 write (istat,'(i5,15(f12.3))') iconf,&
620 (energy_(print_order(i)),i=1,nprint_ene),etot
625 read (intin,'(e15.10,e15.5)',end=1101,err=1101) time,ene
626 call read_x(intin,*11)
628 ! Broadcast the order to compute internal coordinates to the slaves.
630 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
632 call int_from_cart1(.false.)
634 read (intin,'(i5)',end=1101,err=1101) iconf
635 call read_angles(intin,*11)
636 call geom_to_var(nvar,varia)
637 write(iout,*) 'Warning: Calling chainbuild'
648 call mpi_send(ind,6,mpi_integer,man,idint,CG_COMM,&
650 call mpi_send(varia,nvar,mpi_double_precision,man,&
652 call mpi_send(ene0,1,mpi_double_precision,man,&
654 nf_mcmf=nf_mcmf+ind(4)
660 call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,&
662 man=muster(mpi_source)
663 call mpi_recv(varia,nvar,mpi_double_precision,&
664 man,idreal,CG_COMM,muster,ierr)
665 call mpi_recv(ene,1,&
666 mpi_double_precision,man,idreal,&
668 call mpi_recv(ene0,1,&
669 mpi_double_precision,man,idreal,&
672 call var_to_geom(nvar,varia)
673 write(iout,*) 'Warning: Calling chainbuild'
679 write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
682 call enerprint(energy_)
683 call briefout(it,etot)
685 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
686 write (istat,'(i5,19(f12.3))') iconf,&
687 (energy_(print_order(i)),i=1,nprint_ene),etot,&
690 write (istat,'(i5,15(f12.3))') iconf,&
691 (energy_(print_order(i)),i=1,nprint_ene),etot
703 call mpi_send(ind,6,mpi_integer,i,idint,CG_COMM,&
708 open(intin,file=intinname,status='old')
709 write (istat,'(a5,20a12)')"# ",&
710 (wname(print_order(i)),i=1,nprint_ene)
711 write (istat,'("# ",20(1pe12.4))') &
712 (weights(print_order(i)),i=1,nprint_ene)
714 write (istat,'(a5,20a12)')"# ",&
715 (ename(print_order(i)),i=1,nprint_ene),&
716 "ETOT total","RMSD","nat.contact","nnt.contact"
718 write (istat,'(a5,14a12)')"# ",&
719 (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
723 read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
724 call read_x(intin,*11)
726 ! Broadcast the order to compute internal coordinates to the slaves.
728 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
730 call int_from_cart1(.false.)
732 read (intin,'(i5)',end=11,err=11) iconf
733 call read_angles(intin,*11)
734 call geom_to_var(nvar,varia)
735 write(iout,*) 'Warning: Calling chainbuild'
738 write (iout,'(a,i7)') 'Conformation #',iconf
739 if (minim) call minimize(etot,varia,iretcode,nfun)
743 call enerprint(energy_)
744 if (minim) call briefout(it,etot)
746 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
747 write (istat,'(i5,18(f12.3))') iconf,&
748 (energy_(print_order(i)),i=1,nprint_ene),&
749 etot,rms,frac,frac_nn,co
752 write (istat,'(i5,14(f12.3))') iconf,&
753 (energy_(print_order(i)),i=1,nprint_ene),etot
759 end subroutine exec_mult_eeval_or_minim
760 !-----------------------------------------------------------------------------
761 subroutine exec_checkgrad
762 ! use MPI_data !include 'COMMON.SETUP'
763 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
764 use io_units !include 'COMMON.IOUNITS'
765 !el use energy_data, only:icall !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
766 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
767 ! use REMD !include 'COMMON.REMD'
768 use MD_data !include 'COMMON.MD'
769 use io_base, only:intout
770 use io_config, only:read_fragments
775 ! implicit real*8 (a-h,o-z)
776 ! include 'DIMENSIONS'
781 !el common /srutu/ icall
782 real(kind=8) :: energy_(0:max_ene)
786 ! vbld(i)=vbld(i)+ran_number(-0.1d0,0.1d0)
787 ! if (itype(i).ne.10)
788 ! & vbld(i+nres)=vbld(i+nres)+ran_number(-0.001d0,0.001d0)
790 if (indpdb.eq.0) then
791 write(iout,*) 'Warning: Calling chainbuild'
796 ! dc(j,i)=dc(j,i)+ran_number(-0.2d0,0.2d0)
800 ! if (itype(i).ne.10) then
802 ! dc(j,i+nres)=dc(j,i+nres)+ran_number(-0.2d0,0.2d0)
807 ! dc(j,0)=ran_number(-0.2d0,0.2d0)
817 call etotal(energy_(0))
819 call enerprint(energy_(0))
820 write (iout,*) "Uconst",Uconst," Uconst_back",uconst_back
821 print *,'icheckgrad=',icheckgrad
822 goto (10,20,30) icheckgrad
823 10 call check_ecartint
825 20 call check_cartgrad
829 end subroutine exec_checkgrad
830 !-----------------------------------------------------------------------------
834 use io_config, only:map_read
837 call alloc_map_arrays
841 end subroutine exec_map
842 !-----------------------------------------------------------------------------
845 use io_units !include 'COMMON.IOUNITS'
851 ! include 'DIMENSIONS'
852 ! Conformational Space Annealling programmed by Jooyoung Lee.
853 ! This method works only with parallel machines!
855 call alloc_CSA_arrays
858 write (iout,*) "CSA works on parallel machines only"
861 end subroutine exec_CSA
862 !-----------------------------------------------------------------------------
863 subroutine exec_softreg
864 use io_units !include 'COMMON.IOUNITS'
865 use control_data !include 'COMMON.CONTROL'
867 use io_base, only:intout,briefout
868 use geometry, only:chainbuild
872 ! include 'DIMENSIONS'
873 real(kind=8) :: energy_(0:n_ene)
875 real(kind=8) :: rms,frac,frac_nn,co,etot
878 call alloc_compare_arrays
879 write(iout,*) 'Warning: Calling chainbuild'
882 call enerprint(energy_)
883 if (.not.lsecondary) then
884 write(iout,*) 'Calling secondary structure recognition'
885 call secondary2(debug)
887 write(iout,*) 'Using secondary structure supplied in pdb'
894 call enerprint(energy_)
896 call briefout(0,etot)
897 call secondary2(.true.)
898 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
900 end subroutine exec_softreg
901 !-----------------------------------------------------------------------------
903 !-----------------------------------------------------------------------------
905 subroutine ergastulum
907 ! implicit real*8 (a-h,o-z)
908 ! include 'DIMENSIONS'
911 use MDyn, only:setup_fricmat
912 use REMD, only:fricmat_mult,ginv_mult
916 ! include 'COMMON.SETUP'
917 ! include 'COMMON.DERIV'
918 ! include 'COMMON.VAR'
919 ! include 'COMMON.IOUNITS'
920 ! include 'COMMON.FFIELD'
921 ! include 'COMMON.INTERACT'
922 ! include 'COMMON.MD'
923 ! include 'COMMON.TIME1'
924 real(kind=8),dimension(6*nres) :: z,d_a_tmp !(maxres6) maxres6=6*maxres
925 real(kind=8) :: edum(0:n_ene),time_order(0:10)
926 !el real(kind=8),dimension(2*nres,2*nres) :: Gcopy !(maxres2,maxres2) maxres2=2*maxres
927 !el common /przechowalnia/ Gcopy
931 real(kind=8) :: time00
932 integer :: iorder,i,j,nres2,ierr,ierror
935 if(.not.allocated(Gcopy)) allocate(Gcopy(nres2,nres2))
937 if(.not.allocated(Ginv)) allocate(Ginv(nres2,nres2)) !(maxres2,maxres2)
940 if(.not.allocated(ng_start)) allocate(ng_start(0:nfgtasks-1))
941 if(.not.allocated(ng_counts)) allocate(ng_counts(0:nfgtasks-1))
942 if(.not.allocated(nginv_counts)) allocate(nginv_counts(0:nfgtasks-1)) !(0:MaxProcs-1)
943 if(.not.allocated(nginv_start)) allocate(nginv_start(0:nfgtasks)) !(0:MaxProcs)
945 if(.not.allocated(fricmat)) allocate(fricmat(nres2,nres2)) !maxres2=2*maxres
947 ! Workers wait for variables and NF, and NFL from the boss
949 do while (iorder.ge.0)
950 ! write (*,*) 'Processor',fg_rank,' CG group',kolor,
951 ! & ' receives order from Master'
953 call MPI_Bcast(iorder,1,MPI_INTEGER,king,FG_COMM,IERR)
954 time_Bcast=time_Bcast+MPI_Wtime()-time00
955 if (icall.gt.4 .and. iorder.ge.0) &
956 time_order(iorder)=time_order(iorder)+MPI_Wtime()-time00
959 ! & 'Processor',fg_rank,' completed receive MPI_BCAST order',iorder
960 if (iorder.eq.0) then
963 ! write (2,*) "After etotal"
964 ! write (2,*) "dimen",dimen," dimen3",dimen3
966 else if (iorder.eq.2) then
968 call etotal_short(edum)
969 ! write (2,*) "After etotal_short"
970 ! write (2,*) "dimen",dimen," dimen3",dimen3
972 else if (iorder.eq.3) then
974 call etotal_long(edum)
975 ! write (2,*) "After etotal_long"
976 ! write (2,*) "dimen",dimen," dimen3",dimen3
978 else if (iorder.eq.1) then
980 ! write (2,*) "After sum_gradient"
981 ! write (2,*) "dimen",dimen," dimen3",dimen3
983 else if (iorder.eq.4) then
984 call ginv_mult(z,d_a_tmp)
985 else if (iorder.eq.5) then
986 ! Setup MD things for a slave
987 dimen=(nct-nnt+1)+nside
988 dimen1=(nct-nnt)+(nct-nnt+1)
990 ! write (2,*) "dimen",dimen," dimen3",dimen3
992 call int_bounds(dimen,igmult_start,igmult_end)
993 igmult_start=igmult_start-1
994 call MPI_Allgather(3*igmult_start,1,MPI_INTEGER,&
995 ng_start(0),1,MPI_INTEGER,FG_COMM,IERROR)
996 my_ng_count=igmult_end-igmult_start
997 call MPI_Allgather(3*my_ng_count,1,MPI_INTEGER,ng_counts(0),1,&
998 MPI_INTEGER,FG_COMM,IERROR)
999 write (2,*) "ng_start",(ng_start(i),i=0,nfgtasks-1) !sp
1000 ! write (2,*) "ng_counts",(ng_counts(i),i=0,nfgtasks-1)
1001 myginv_ng_count=nres2*my_ng_count !el maxres2
1002 ! write (2,*) "igmult_start",igmult_start," igmult_end",
1003 ! & igmult_end," my_ng_count",my_ng_count
1005 call MPI_Allgather(nres2*igmult_start,1,MPI_INTEGER,& !el maxres2
1006 nginv_start(0),1,MPI_INTEGER,FG_COMM,IERROR)
1007 call MPI_Allgather(myginv_ng_count,1,MPI_INTEGER,&
1008 nginv_counts(0),1,MPI_INTEGER,FG_COMM,IERROR)
1009 ! write (2,*) "nginv_start",(nginv_start(i),i=0,nfgtasks-1)
1010 ! write (2,*) "nginv_counts",(nginv_counts(i),i=0,nfgtasks-1)
1012 ! call MPI_Barrier(FG_COMM,IERROR)
1014 call MPI_Scatterv(ginv(1,1),nginv_counts(0),&
1015 nginv_start(0),MPI_DOUBLE_PRECISION,gcopy(1,1),&
1016 myginv_ng_count,MPI_DOUBLE_PRECISION,king,FG_COMM,IERR)
1018 time_scatter_ginv=time_scatter_ginv+MPI_Wtime()-time00
1021 do j=1,2*my_ng_count
1022 ginv(j,i)=gcopy(i,j)
1025 ! write (2,*) "dimen",dimen," dimen3",dimen3
1026 ! write (2,*) "End MD setup"
1028 ! write (iout,*) "My chunk of ginv_block"
1029 ! call MATOUT2(my_ng_count,dimen3,maxres2,maxers2,ginv_block)
1030 else if (iorder.eq.6) then
1031 call int_from_cart1(.false.)
1032 else if (iorder.eq.7) then
1033 call chainbuild_cart
1034 else if (iorder.eq.8) then
1036 else if (iorder.eq.9) then
1037 call fricmat_mult(z,d_a_tmp)
1038 else if (iorder.eq.10) then
1042 write (*,*) 'Processor',fg_rank,' CG group',kolor,&
1043 ' absolute rank',myrank,' leves ERGASTULUM.'
1044 write(*,*)'Processor',fg_rank,' wait times for respective orders',&
1045 (' order[',i,']',time_order(i),i=0,10)
1047 end subroutine ergastulum