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'
166 if (me.eq.king .or. .not. out1file) &
167 write (iout,*) "Calling chainbuild"
171 end subroutine exec_MD
172 !---------------------------------------------------------------------------
173 subroutine exec_MREMD
174 use MPI_data !include 'COMMON.SETUP'
175 use control_data !include 'COMMON.CONTROL'
176 use io_units !include 'COMMON.IOUNITS'
178 use REMD_data !include 'COMMON.REMD'
179 use geometry, only:chainbuild
183 ! include 'DIMENSIONS'
190 call alloc_MREMD_arrays
192 if (me.eq.king .or. .not. out1file) &
193 write (iout,*) "Calling chainbuild"
195 if (me.eq.king .or. .not. out1file) &
196 write (iout,*) "Calling REMD"
206 end subroutine exec_MREMD
207 !-----------------------------------------------------------------------------
208 subroutine exec_eeval_or_minim
209 use MPI_data !include 'COMMON.SETUP'
210 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
211 use io_units !include 'COMMON.IOUNITS'
213 ! use energy !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
214 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
215 ! use REMD !include 'COMMON.REMD'
216 ! use MD !include 'COMMON.MD'
221 use geometry, only:chainbuild
223 use compare, only:alloc_compare_arrays,hairpin,secondary2,rms_nac_nnc
224 use minimm, only:minimize,minim_dc,sc_move
228 ! implicit real*8 (a-h,o-z)
229 ! include 'DIMENSIONS'
234 !el common /srutu/ icall
235 real(kind=8) :: energy_(0:n_ene)
236 real(kind=8) :: energy_long(0:n_ene),energy_short(0:n_ene)
237 real(kind=8) :: varia(6*nres) !(maxvar) (maxvar=6*maxres)
238 real(kind=8) :: time00, evals, etota, etot, time_ene, time1
239 integer :: nharp,nft_sc,iretcode,nfun
240 integer,dimension(4,nres/3) :: iharp !(4,nres/3)(4,maxres/3)
242 real(kind=8) :: rms,frac,frac_nn,co
245 call alloc_compare_arrays
246 if (indpdb.eq.0) then
248 write(iout,*) 'Warning: Calling chainbuild'
254 ! write(iout,*)"in exec_eeval or minimim",split_ene
256 ! write(iout,*)"cccccc",j,(c(i,j),i=1,3)
257 ! write(iout,*)"dcccccc",j,(dc(i,j),i=1,3)
260 ! write(iout,*)"in exec_eeval or minimim"
262 print *,"Processor",myrank," after chainbuild"
265 call etotal_long(energy_long)
266 write (iout,*) "Printing long range energy"
267 call enerprint(energy_long)
269 call etotal_short(energy_short)
270 write (iout,*) "Printing short range energy"
271 call enerprint(energy_short)
273 energy_(i)=energy_long(i)+energy_short(i)
274 write (iout,*) i,energy_long(i),energy_short(i),energy_(i)
276 write (iout,*) "Printing long+short range energy"
277 call enerprint(energy_)
282 time_ene=MPI_Wtime()-time00
284 write (iout,*) "Time for energy evaluation",time_ene
285 print *,"after etotal"
288 call enerprint(energy_)
289 call hairpin(.true.,nharp,iharp)
290 call secondary2(.true.)
294 print *, 'Calling OVERLAP_SC'
295 call overlap_sc(fail)
299 call sc_move(2,nres-1,10,1d10,nft_sc,etot)
300 print *,'SC_move',nft_sc,etot
301 write(iout,*) 'SC_move',nft_sc,etot
305 print *, 'Calling MINIM_DC'
309 ! call check_ecartint !el
310 call minim_dc(etot,iretcode,nfun)
311 ! call check_ecartint !el
313 if (indpdb.ne.0) then
315 write(iout,*) 'Warning: Calling chainbuild'
318 call geom_to_var(nvar,varia)
319 print *,'Calling MINIMIZE.'
324 ! call exec_checkgrad !el
325 call minimize(etot,varia,iretcode,nfun)
327 ! call exec_checkgrad !el
329 print *,'SUMSL return code is',iretcode,' eval ',nfun
331 evals=nfun/(MPI_WTIME()-time1)
333 print *,'# eval/s',evals
334 print *,'refstr=',refstr
335 call hairpin(.true.,nharp,iharp)
336 call secondary2(.true.)
339 call enerprint(energy_)
342 call briefout(0,etot)
343 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
344 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
345 write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
346 write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
348 print *,'refstr=',refstr
349 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
350 call briefout(0,etot)
352 if (outpdb) call pdbout(etot,titel(:32),ipdb)
353 if (outmol2) call mol2out(etot,titel(:32))
354 !elwrite(iout,*) "after exec_eeval_or_minim"
356 end subroutine exec_eeval_or_minim
357 !-----------------------------------------------------------------------------
358 subroutine exec_regularize
359 ! use MPI_data !include 'COMMON.SETUP'
360 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
361 use io_units !include 'COMMON.IOUNITS'
363 use energy_data !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
364 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
365 ! use REMD !include 'COMMON.REMD'
366 ! use MD !include 'COMMON.MD'
370 ! implicit real*8 (a-h,o-z)
371 ! include 'DIMENSIONS'
375 real(kind=8) :: energy_(0:n_ene)
377 real(kind=8) :: rms,frac,frac_nn,co
380 call alloc_compare_arrays
384 call regularize(nct-nnt+1,etot,rms,cref(1,nnt,1),iretcode)
386 energy_(0)=energy_(0)-energy_(14)
388 call enerprint(energy_)
390 call briefout(0,etot)
391 if (outpdb) call pdbout(etot,titel(:32),ipdb)
392 if (outmol2) call mol2out(etot,titel(:32))
393 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
394 write (iout,'(a,i3)') 'SUMSL return code:',iretcode
396 end subroutine exec_regularize
397 !-----------------------------------------------------------------------------
398 subroutine exec_thread
399 ! use MPI_data !include 'COMMON.SETUP'
402 ! include 'DIMENSIONS'
406 call alloc_compare_arrays
409 end subroutine exec_thread
410 !-----------------------------------------------------------------------------
412 ! use MPI_data !include 'COMMON.SETUP'
413 use control_data !include 'COMMON.CONTROL'
418 ! implicit real*8 (a-h,o-z)
419 ! include 'DIMENSIONS'
420 character(len=10) :: nodeinfo
421 real(kind=8) :: varia(6*nres) !(maxvar) (maxvar=6*maxres)
426 call alloc_MCM_arrays
430 if (modecalc.eq.3) then
436 if (modecalc.eq.3) then
446 end subroutine exec_MC
447 !-----------------------------------------------------------------------------
448 subroutine exec_mult_eeval_or_minim
449 use MPI_data !include 'COMMON.SETUP'
450 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
451 use io_units !include 'COMMON.IOUNITS'
453 use energy_data !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
454 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
455 ! use REMD !include 'COMMON.REMD'
456 ! use MD !include 'COMMON.MD'
458 use geometry, only:chainbuild,geom_to_var,int_from_cart1,var_to_geom
459 use energy, only:etotal,enerprint
460 use compare, only:rms_nac_nnc
461 use minimm, only:minimize!,minim_mcmf
462 ! implicit real*8 (a-h,o-z)
463 ! include 'DIMENSIONS'
465 use minimm, only:minim_mcmf
468 integer :: ierror,ierr
470 real(kind=8),dimension(mpi_status_size) :: muster
474 real(kind=8) :: varia(6*nres) !(maxvar) (maxvar=6*maxres)
475 integer,dimension(6) :: ind
476 real(kind=8) :: energy_(0:n_ene)
478 real(kind=8) :: etot,ene0
479 integer :: mm,imm,nft,n,iconf,nmin,i,iretcode,nfun,it,&
481 real(kind=8) :: rms,frac,frac_nn,co,time,ene
491 open(intin,file=intinname,status='old')
492 write (istat,'(a5,20a12)')"# ",&
493 (wname(print_order(i)),i=1,nprint_ene)
495 write (istat,'(a5,20a12)')"# ",&
496 (ename(print_order(i)),i=1,nprint_ene),&
497 "ETOT total","RMSD","nat.contact","nnt.contact"
499 write (istat,'(a5,20a12)')"# ",&
500 (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
506 read (intin,'(e15.10,e15.5)',end=1100,err=1100) time,ene
507 call read_x(intin,*11)
509 ! Broadcast the order to compute internal coordinates to the slaves.
511 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
513 call int_from_cart1(.false.)
515 read (intin,'(i5)',end=1100,err=1100) iconf
516 call read_angles(intin,*11)
517 call geom_to_var(nvar,varia)
518 write(iout,*) 'Warning: Calling chainbuild'
521 write (iout,'(a,i7)') 'Conformation #',iconf
523 call briefout(iconf,energy_(0))
524 call enerprint(energy_)
527 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
528 write (istat,'(i5,20(f12.3))') iconf,&
529 (energy_(print_order(i)),i=1,nprint_ene),etot,&
533 write (istat,'(i5,16(f12.3))') iconf,&
534 (energy_(print_order(i)),i=1,nprint_ene),etot
550 if (mm.lt.nodes) then
552 read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
553 call read_x(intin,*11)
555 ! Broadcast the order to compute internal coordinates to the slaves.
557 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
559 call int_from_cart1(.false.)
561 read (intin,'(i5)',end=11,err=11) iconf
562 call read_angles(intin,*11)
563 call geom_to_var(nvar,varia)
564 write(iout,*) 'Warning: Calling chainbuild'
567 write (iout,'(a,i7)') 'Conformation #',iconf
577 call mpi_send(ind,6,mpi_integer,mm,idint,CG_COMM,&
579 call mpi_send(varia,nvar,mpi_double_precision,mm,&
581 call mpi_send(ene0,1,mpi_double_precision,mm,&
583 ! print *,'task ',n,' sent to worker ',mm,nvar
585 call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,&
587 man=muster(mpi_source)
588 ! print *,'receiving result from worker ',man,' (',iii1,iii,')'
589 call mpi_recv(varia,nvar,mpi_double_precision,&
590 man,idreal,CG_COMM,muster,ierr)
591 call mpi_recv(ene,1,&
592 mpi_double_precision,man,idreal,&
594 call mpi_recv(ene0,1,&
595 mpi_double_precision,man,idreal,&
597 ! print *,'result received from worker ',man,' sending now'
599 call var_to_geom(nvar,varia)
600 write(iout,*) 'Warning: Calling chainbuild'
606 write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
609 call enerprint(energy_)
610 call briefout(it,etot)
611 ! if (minim) call briefout(it,etot)
613 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
614 write (istat,'(i5,19(f12.3))') iconf,&
615 (energy_(print_order(i)),i=1,nprint_ene),etot,&
618 write (istat,'(i5,15(f12.3))') iconf,&
619 (energy_(print_order(i)),i=1,nprint_ene),etot
624 read (intin,'(e15.10,e15.5)',end=1101,err=1101) time,ene
625 call read_x(intin,*11)
627 ! Broadcast the order to compute internal coordinates to the slaves.
629 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
631 call int_from_cart1(.false.)
633 read (intin,'(i5)',end=1101,err=1101) iconf
634 call read_angles(intin,*11)
635 call geom_to_var(nvar,varia)
636 write(iout,*) 'Warning: Calling chainbuild'
647 call mpi_send(ind,6,mpi_integer,man,idint,CG_COMM,&
649 call mpi_send(varia,nvar,mpi_double_precision,man,&
651 call mpi_send(ene0,1,mpi_double_precision,man,&
653 nf_mcmf=nf_mcmf+ind(4)
659 call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,&
661 man=muster(mpi_source)
662 call mpi_recv(varia,nvar,mpi_double_precision,&
663 man,idreal,CG_COMM,muster,ierr)
664 call mpi_recv(ene,1,&
665 mpi_double_precision,man,idreal,&
667 call mpi_recv(ene0,1,&
668 mpi_double_precision,man,idreal,&
671 call var_to_geom(nvar,varia)
672 write(iout,*) 'Warning: Calling chainbuild'
678 write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
681 call enerprint(energy_)
682 call briefout(it,etot)
684 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
685 write (istat,'(i5,19(f12.3))') iconf,&
686 (energy_(print_order(i)),i=1,nprint_ene),etot,&
689 write (istat,'(i5,15(f12.3))') iconf,&
690 (energy_(print_order(i)),i=1,nprint_ene),etot
702 call mpi_send(ind,6,mpi_integer,i,idint,CG_COMM,&
707 open(intin,file=intinname,status='old')
708 write (istat,'(a5,20a12)')"# ",&
709 (wname(print_order(i)),i=1,nprint_ene)
710 write (istat,'("# ",20(1pe12.4))') &
711 (weights(print_order(i)),i=1,nprint_ene)
713 write (istat,'(a5,20a12)')"# ",&
714 (ename(print_order(i)),i=1,nprint_ene),&
715 "ETOT total","RMSD","nat.contact","nnt.contact"
717 write (istat,'(a5,14a12)')"# ",&
718 (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
722 read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
723 call read_x(intin,*11)
725 ! Broadcast the order to compute internal coordinates to the slaves.
727 call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
729 call int_from_cart1(.false.)
731 read (intin,'(i5)',end=11,err=11) iconf
732 call read_angles(intin,*11)
733 call geom_to_var(nvar,varia)
734 write(iout,*) 'Warning: Calling chainbuild'
737 write (iout,'(a,i7)') 'Conformation #',iconf
738 if (minim) call minimize(etot,varia,iretcode,nfun)
742 call enerprint(energy_)
743 if (minim) call briefout(it,etot)
745 call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
746 write (istat,'(i5,18(f12.3))') iconf,&
747 (energy_(print_order(i)),i=1,nprint_ene),&
748 etot,rms,frac,frac_nn,co
751 write (istat,'(i5,14(f12.3))') iconf,&
752 (energy_(print_order(i)),i=1,nprint_ene),etot
758 end subroutine exec_mult_eeval_or_minim
759 !-----------------------------------------------------------------------------
760 subroutine exec_checkgrad
761 ! use MPI_data !include 'COMMON.SETUP'
762 use control_data !include 'COMMON.CONTROL''COMMON.TIME1''COMMON.NAMES''COMMON.HEADER'
763 use io_units !include 'COMMON.IOUNITS'
764 !el use energy_data, only:icall !include 'COMMON.INTERACT''COMMON.CONTACTS''COMMON.VAR''COMMON.FFIELD' 'COMMON.SBRIDGE'
765 use geometry_data !include 'COMMON.GEO''COMMON.CHAIN'
766 ! use REMD !include 'COMMON.REMD'
767 use MD_data !include 'COMMON.MD'
768 use io_base, only:intout
769 use io_config, only:read_fragments
774 ! implicit real*8 (a-h,o-z)
775 ! include 'DIMENSIONS'
780 !el common /srutu/ icall
781 real(kind=8) :: energy_(0:max_ene)
785 ! vbld(i)=vbld(i)+ran_number(-0.1d0,0.1d0)
786 ! if (itype(i).ne.10)
787 ! & vbld(i+nres)=vbld(i+nres)+ran_number(-0.001d0,0.001d0)
789 if (indpdb.eq.0) then
790 write(iout,*) 'Warning: Calling chainbuild'
795 ! dc(j,i)=dc(j,i)+ran_number(-0.2d0,0.2d0)
799 ! if (itype(i).ne.10) then
801 ! dc(j,i+nres)=dc(j,i+nres)+ran_number(-0.2d0,0.2d0)
806 ! dc(j,0)=ran_number(-0.2d0,0.2d0)
816 call etotal(energy_(0))
818 call enerprint(energy_(0))
819 write (iout,*) "Uconst",Uconst," Uconst_back",uconst_back
820 print *,'icheckgrad=',icheckgrad
821 goto (10,20,30) icheckgrad
822 10 call check_ecartint
824 20 call check_cartgrad
828 end subroutine exec_checkgrad
829 !-----------------------------------------------------------------------------
833 use io_config, only:map_read
836 call alloc_map_arrays
840 end subroutine exec_map
841 !-----------------------------------------------------------------------------
844 use io_units !include 'COMMON.IOUNITS'
850 ! include 'DIMENSIONS'
851 ! Conformational Space Annealling programmed by Jooyoung Lee.
852 ! This method works only with parallel machines!
854 call alloc_CSA_arrays
857 write (iout,*) "CSA works on parallel machines only"
860 end subroutine exec_CSA
861 !-----------------------------------------------------------------------------
862 subroutine exec_softreg
863 use io_units !include 'COMMON.IOUNITS'
864 use control_data !include 'COMMON.CONTROL'
866 use io_base, only:intout,briefout
867 use geometry, only:chainbuild
871 ! include 'DIMENSIONS'
872 real(kind=8) :: energy_(0:n_ene)
874 real(kind=8) :: rms,frac,frac_nn,co,etot
877 call alloc_compare_arrays
878 write(iout,*) 'Warning: Calling chainbuild'
881 call enerprint(energy_)
882 if (.not.lsecondary) then
883 write(iout,*) 'Calling secondary structure recognition'
884 call secondary2(debug)
886 write(iout,*) 'Using secondary structure supplied in pdb'
893 call enerprint(energy_)
895 call briefout(0,etot)
896 call secondary2(.true.)
897 if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
899 end subroutine exec_softreg
900 !-----------------------------------------------------------------------------
902 !-----------------------------------------------------------------------------
904 subroutine ergastulum
906 ! implicit real*8 (a-h,o-z)
907 ! include 'DIMENSIONS'
910 use MDyn, only:setup_fricmat
911 use REMD, only:fricmat_mult,ginv_mult
915 ! include 'COMMON.SETUP'
916 ! include 'COMMON.DERIV'
917 ! include 'COMMON.VAR'
918 ! include 'COMMON.IOUNITS'
919 ! include 'COMMON.FFIELD'
920 ! include 'COMMON.INTERACT'
921 ! include 'COMMON.MD'
922 ! include 'COMMON.TIME1'
923 real(kind=8),dimension(6*nres) :: z,d_a_tmp !(maxres6) maxres6=6*maxres
924 real(kind=8) :: edum(0:n_ene),time_order(0:10)
925 !el real(kind=8),dimension(2*nres,2*nres) :: Gcopy !(maxres2,maxres2) maxres2=2*maxres
926 !el common /przechowalnia/ Gcopy
930 real(kind=8) :: time00
931 integer :: iorder,i,j,nres2,ierr,ierror
934 if(.not.allocated(Gcopy)) allocate(Gcopy(nres2,nres2))
936 if(.not.allocated(Ginv)) allocate(Ginv(nres2,nres2)) !(maxres2,maxres2)
939 if(.not.allocated(ng_start)) allocate(ng_start(0:nfgtasks-1))
940 if(.not.allocated(ng_counts)) allocate(ng_counts(0:nfgtasks-1))
941 if(.not.allocated(nginv_counts)) allocate(nginv_counts(0:nfgtasks-1)) !(0:MaxProcs-1)
942 if(.not.allocated(nginv_start)) allocate(nginv_start(0:nfgtasks)) !(0:MaxProcs)
944 if(.not.allocated(fricmat)) allocate(fricmat(nres2,nres2)) !maxres2=2*maxres
946 ! Workers wait for variables and NF, and NFL from the boss
948 do while (iorder.ge.0)
949 ! write (*,*) 'Processor',fg_rank,' CG group',kolor,
950 ! & ' receives order from Master'
952 call MPI_Bcast(iorder,1,MPI_INTEGER,king,FG_COMM,IERR)
953 time_Bcast=time_Bcast+MPI_Wtime()-time00
954 if (icall.gt.4 .and. iorder.ge.0) &
955 time_order(iorder)=time_order(iorder)+MPI_Wtime()-time00
958 ! & 'Processor',fg_rank,' completed receive MPI_BCAST order',iorder
959 if (iorder.eq.0) then
962 ! write (2,*) "After etotal"
963 ! write (2,*) "dimen",dimen," dimen3",dimen3
965 else if (iorder.eq.2) then
967 call etotal_short(edum)
968 ! write (2,*) "After etotal_short"
969 ! write (2,*) "dimen",dimen," dimen3",dimen3
971 else if (iorder.eq.3) then
973 call etotal_long(edum)
974 ! write (2,*) "After etotal_long"
975 ! write (2,*) "dimen",dimen," dimen3",dimen3
977 else if (iorder.eq.1) then
979 ! write (2,*) "After sum_gradient"
980 ! write (2,*) "dimen",dimen," dimen3",dimen3
982 else if (iorder.eq.4) then
983 call ginv_mult(z,d_a_tmp)
984 else if (iorder.eq.5) then
985 ! Setup MD things for a slave
986 dimen=(nct-nnt+1)+nside
987 dimen1=(nct-nnt)+(nct-nnt+1)
989 ! write (2,*) "dimen",dimen," dimen3",dimen3
991 call int_bounds(dimen,igmult_start,igmult_end)
992 igmult_start=igmult_start-1
993 call MPI_Allgather(3*igmult_start,1,MPI_INTEGER,&
994 ng_start(0),1,MPI_INTEGER,FG_COMM,IERROR)
995 my_ng_count=igmult_end-igmult_start
996 call MPI_Allgather(3*my_ng_count,1,MPI_INTEGER,ng_counts(0),1,&
997 MPI_INTEGER,FG_COMM,IERROR)
998 write (2,*) "ng_start",(ng_start(i),i=0,nfgtasks-1) !sp
999 ! write (2,*) "ng_counts",(ng_counts(i),i=0,nfgtasks-1)
1000 myginv_ng_count=nres2*my_ng_count !el maxres2
1001 ! write (2,*) "igmult_start",igmult_start," igmult_end",
1002 ! & igmult_end," my_ng_count",my_ng_count
1004 call MPI_Allgather(nres2*igmult_start,1,MPI_INTEGER,& !el maxres2
1005 nginv_start(0),1,MPI_INTEGER,FG_COMM,IERROR)
1006 call MPI_Allgather(myginv_ng_count,1,MPI_INTEGER,&
1007 nginv_counts(0),1,MPI_INTEGER,FG_COMM,IERROR)
1008 ! write (2,*) "nginv_start",(nginv_start(i),i=0,nfgtasks-1)
1009 ! write (2,*) "nginv_counts",(nginv_counts(i),i=0,nfgtasks-1)
1011 ! call MPI_Barrier(FG_COMM,IERROR)
1013 call MPI_Scatterv(ginv(1,1),nginv_counts(0),&
1014 nginv_start(0),MPI_DOUBLE_PRECISION,gcopy(1,1),&
1015 myginv_ng_count,MPI_DOUBLE_PRECISION,king,FG_COMM,IERR)
1017 time_scatter_ginv=time_scatter_ginv+MPI_Wtime()-time00
1020 do j=1,2*my_ng_count
1021 ginv(j,i)=gcopy(i,j)
1024 ! write (2,*) "dimen",dimen," dimen3",dimen3
1025 ! write (2,*) "End MD setup"
1027 ! write (iout,*) "My chunk of ginv_block"
1028 ! call MATOUT2(my_ng_count,dimen3,maxres2,maxers2,ginv_block)
1029 else if (iorder.eq.6) then
1030 call int_from_cart1(.false.)
1031 else if (iorder.eq.7) then
1032 call chainbuild_cart
1033 else if (iorder.eq.8) then
1035 else if (iorder.eq.9) then
1036 call fricmat_mult(z,d_a_tmp)
1037 else if (iorder.eq.10) then
1041 write (*,*) 'Processor',fg_rank,' CG group',kolor,&
1042 ' absolute rank',myrank,' leves ERGASTULUM.'
1043 write(*,*)'Processor',fg_rank,' wait times for respective orders',&
1044 (' order[',i,']',time_order(i),i=0,10)
1046 end subroutine ergastulum