subroutine make_list(lprn,first_call,nvarr,x) implicit none include "DIMENSIONS" include "DIMENSIONS.ZSCOPT" #ifdef MPI include "mpif.h" integer IERROR,ErrCode,Status(MPI_STATUS_SIZE,10) integer req(10),msg_in(5),msg_out(5),address,size character*1 buffer(8*(2*maxstr_proc*nntyp+8000)) include "COMMON.MPI" #endif include "COMMON.WEIGHTS" include "COMMON.WEIGHTDER" include "COMMON.ENERGIES" include "COMMON.IOUNITS" include "COMMON.VMCPAR" include "COMMON.NAMES" include "COMMON.INTERACT" include "COMMON.TIME1" include "COMMON.CHAIN" include "COMMON.PROTFILES" include "COMMON.VAR" include "COMMON.GEO" include "COMMON.OPTIM" include "COMMON.CLASSES" include "COMMON.COMPAR" include "COMMON.TORSION" C Define local variables integer nvarr integer i,ii,iii,ibatch,kkk,jj,j,k,kk,l,iprot,ib,nn integer ipass_conf,istart_conf,iend_conf,Previous,Next double precision energia(0:max_ene) double precision etoti,elowesti,dene double precision tcpu_ini,tcpu_fin,tcpu double precision elowest_t(2,maxT),etot_aux,enepsjk, & emini,elowest_aux(2,maxT) integer iroof,icant external iroof,icant logical lprn,first_call logical*1 lflag(maxstr) integer ipermut(maxstr) integer list_conf_(maxstr,maxprot) double precision x(nvarr) double precision ftune_eps external ftune_eps #ifdef MPI double precision e_total_(maxstr_proc) #endif write (iout,*) "Making the worklist of conformations." write (iout,*) "enecut",(enecut(i),i=1,nprot) write (iout,*) "first_call ",first_call," nvarr",nvarr tcpu_ini = tcpu() do iprot=1,nprot do i=1,ntot(iprot) list_conf(i,iprot)=i enddo ntot_work(iprot)=ntot(iprot) enddo #ifdef MPI c c Divide the whole database between processors c Previous = me1-1 Next = me1+1 if (Previous.lt.0) Previous = Nprocs-1 if (Next.ge.Nprocs) Next = 0 call work_partition(lprn) #endif c c Loop over proteins c DO iprot=1,nprot call restore_molinfo(iprot) c write (iout,*) "Processor",me," iprot",iprot, c & " indstart",indstart(me1,iprot)," indend",indend(me1,iprot), c & " init_ene",init_ene," mod_fourier",mod_fourier(0), c & " mod_elec",mod_elec," mod_scp",mod_scp call restore_molinfo(iprot) c c Loop over the conformations of protein IPROT assigned to the current processor. c The conformations are read off a DA scratchfile and processed in passes, each c of which requires no more than MAXSTR_PROC conformations to be stored in memory c simultaneously. c if (.not.init_ene .and. mod_fourier(nloctyp).eq.0 & .and. mod_fouriertor(nloctyp).eq.0 & .and. .not. mod_elec .and. .not. mod_scp) then open (ientin,file=benefiles(iprot),status="old", & form="unformatted",access="direct",recl=lenrec_ene(iprot)) else open (icbase,file=bprotfiles(iprot),status="old", & form="unformatted",access="direct",recl=lenrec(iprot)) c Change AL 12/30/2017 if (.not. mod_other_params) & open (ientout,file=benefiles(iprot),status="old", & form="unformatted",access="direct",recl=lenrec_ene(iprot)) endif #ifdef MPI nchunk_ene(iprot) = iroof(scount(me1,iprot),maxstr) nchunk_conf(iprot) = iroof(scount(me1,iprot),maxstr_proc) write (iout,*)"Protein",iprot," energy evaluation in", & nchunk_conf(iprot)," passes." ipass_conf=0 jj=0 do i=indstart(me1,iprot),indend(me1,iprot),maxstr_proc ipass_conf=ipass_conf+1 write (iout,*) "MAKE_LIST: Pass",ipass_conf istart_conf=i iend_conf=min0(i+maxstr_proc-1,indend(me1,iprot)) #else nchunk_ene(iprot) = iroof(ntot(iprot),maxstr) nchunk_conf(iprot) = iroof(ntot(iprot),maxstr_proc) write (iout,*)"Protein",iprot," energy evaluation in", & nchunk_conf," passes." ipass_conf=0 do i=1,ntot(iprot),maxstr_proc ipass_conf=ipass_conf+1 write (iout,*) "MAKE_LIST: Pass",ipass_conf istart_conf=i iend_conf=min0(i+maxstr_proc-1,ntot(iprot)) #endif ii=0 c c Read the chunk of conformations off a DA scratchfile. c if (.not.init_ene .and. mod_fourier(nloctyp).eq.0 & .and. mod_fouriertor(nloctyp).eq.0 .and. .not. mod_elec & .and. .not. mod_scp) then c c If energy components have been pre-computed read them off a DA file. c call daread_ene(iprot,istart_conf,iend_conf) do iii=istart_conf,iend_conf ii=ii+1 jj=jj+1 if (mod_side) then enetb(ii,1,iprot)=0.0d0 do j=1,ntyp do k=1,j enetb(ii,1,iprot)=enetb(ii,1,iprot)+ftune_eps(eps(j,k))* & eneps(1,icant(j,k),ii,iprot)+ & eps(j,k)*eneps(2,icant(j,k),ii,iprot) enddo enddo endif etot_aux=0.0d0 do j=1,n_ene etot_aux=etot_aux+ww(j)*enetb(ii,j,iprot) enddo e_total(iii,iprot)=etot_aux #ifdef DEBUG write (iout,'(i5,16(1pe12.4))') iii, & (enetb(ii,j,iprot),j=1,n_ene),e_total(iii,iprot) call flush(iout) #endif enddo if (first_call .and. mod_side) then write (iout,*) "Callig x2w" call flush(iout) call x2w(nvarr,x_orig) write (iout,*) "After x2w" call flush(iout) ii=0 do iii=istart_conf,iend_conf ii=ii+1 enetb_oorig(ii,1,iprot)=0.0d0 do j=1,ntyp do k=1,j enetb_oorig(ii,1,iprot)=enetb_oorig(ii,1,iprot)+ & ftune_eps(eps(j,k))*eneps(1,icant(j,k),ii,iprot)+ & eps(j,k)*eneps(2,icant(j,k),ii,iprot) enddo enddo do j=2,n_ene enetb_oorig(ii,j,iprot)=enetb(ii,j,iprot) enddo enddo call x2w(nvarr,x) endif else call daread_ccoords(iprot,istart_conf,iend_conf) c c Compute the energies of the conformations currently in memory and compute c the lowest energies. c do iii=istart_conf,iend_conf ii=ii+1 jj=jj+1 call restore_ccoords(iprot,ii) call int_from_cart1(.false.) #ifdef DEBUG write (iout,*) "Before etotal",iii,i call flush(iout) #endif call etotal(energia(0)) #ifdef DEBUG write (iout,*) "After etotal",i call flush(iout) call enerprint(energia(0)) c write (iout,'(i5,16(1pe12.4))') i, c & (energia(j),j=1,n_ene),energia(0) call flush(iout) #endif #ifdef DEBUG write (iout,*) "Conformation:",i write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres) write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct) write (iout,'(8f10.4)') (vbld(k),k=2,nres) write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct) write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres) write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres) write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1) write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1) call enerprint(energia(0)) #endif e_total(iii,iprot)=energia(0) do j=1,n_ene enetb(ii,j,iprot)=energia(j) enddo do j=1,ntyp do k=1,j eneps(1,icant(j,k),ii,iprot)=eneps_temp(1,icant(j,k)) eneps(2,icant(j,k),ii,iprot)=eneps_temp(2,icant(j,k)) enddo enddo #ifdef DEBUG write (iout,'(i5,20(1pe12.4))') iii, & (energia(j),j=1,n_ene),energia(0),eini(iii,iprot), & entfac(iii,iprot) call flush(iout) #endif if (energia(0).ge.1.0d99) then write (iout,*) & "MAKE_LIST:CHUJ NASTAPIL in energy evaluation for", & " point",i,". Probably NaNs in some of the energy components." write (iout,*) "The components of the energy are:" call enerprint(energia(0)) write (iout,*) "Conformation:",i write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres) write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct) write (iout,'(8f10.4)') (vbld(k),k=2,nres) write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct) write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres) write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres) write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1) write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1) write (iout,*) "Calculation terminated at this point.", & " Check the database of conformations" #ifdef MPI call MPI_Abort(MPI_COMM_WORLD,ErrCode,IERROR) #endif stop "SEVERE error in energy calculation" endif enddo ! iii if (first_call) then do ii=1,iend_conf-istart_conf+1 do j=1,n_ene enetb_orig(ii,j,iprot)=enetb(ii,j,iprot) enddo enddo ii=0 call x2w(nvarr,x_orig) #ifdef DEBUG write (iout,*) "x,xorig" do k=1,nvarr write (iout,'(i5,2f10.5)') k,x(k),x_orig(k) enddo #endif do iii=istart_conf,iend_conf ii=ii+1 call restore_ccoords(iprot,ii) call int_from_cart1(.false.) call etotal(energia(0)) #ifdef DEBUG write (iout,*) "Conformation:",iii,ii write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres) write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct) write (iout,'(8f10.4)') (vbld(k),k=2,nres) write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct) write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres) write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres) write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1) write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1) call enerprint(energia(0)) #endif do j=1,n_ene enetb_oorig(ii,j,iprot)=energia(j) enddo #ifdef DEBUG write (iout,'(2i5,20(1pe12.4))') iii,ii, & (energia(j),j=1,n_ene),energia(0) write (iout,'(2i5,20(1pe12.4))') iii,ii, & (enetb(ii,j,iprot),j=1,n_ene) call flush(iout) #endif if (energia(0).ge.1.0d99) then write (iout,*) "CHUJ NASTAPIL in energy evaluation for", & " point",i,". Probably NaNs in some of the energy components." write (iout,*) "The components of the energy are:" call enerprint(energia(0)) write (iout,*) "Calculation terminated at this point.", & " Check the database of conformations" #ifdef MPI call MPI_Abort(MPI_COMM_WORLD,ErrCode,IERROR) #endif stop "SEVERE error in energy calculation" endif enddo ! iii write (iout,*) "MAKE_LIST Callig x2w" call flush(iout) call x2w(nvarr,x) write (iout,*) "After x2w" call flush(iout) endif write (iout,*) "make_list: calling dawrite_ene" write (iout,*) "istart_conf",istart_conf, & " iend_conf",iend_conf C Change AL 12/30/2017 if (.not. mod_other_params) & call dawrite_ene(iprot,istart_conf,iend_conf,ientout) #ifdef MPI c Distribute energy components through ring call MPI_Barrier(WHAM_COMM,IERROR) write (iout,*) "Processes synchronized in make_list" call flush(iout) msg_out(1)= 5*me1+1000*ipass_conf+1 msg_out(2)= 5*me1+1000*ipass_conf+2 msg_out(3)= 5*me1+1000*ipass_conf+3 msg_out(4)= 5*me1+1000*ipass_conf+4 msg_out(5)= 5*me1+1000*ipass_conf+5 do iii=1,Nprocs-1 c Send the current energy tables to the right neighbor c Receive the energy tables produced by processor kkk from the left neighbor kkk = mod(me1-iii+NProcs,Nprocs) msg_in(1)= 5*kkk+1000*ipass_conf+1 msg_in(2)= 5*kkk+1000*ipass_conf+2 msg_in(3)= 5*kkk+1000*ipass_conf+3 msg_in(4)= 5*kkk+1000*ipass_conf+4 msg_in(5)= 5*kkk+1000*ipass_conf+5 write (iout,*) "me1",me1," iii",iii," Previous",Previous, & " Next",Next," kkk",kkk write (iout,*) "msg_in",msg_in write (iout,*) "msg_out",msg_out call flush(iout) write (iout,*) "Processor",me1," Start Send and receive" call flush(iout) call MPI_Send(istart_conf,1,MPI_INTEGER,Next,msg_out(1), & WHAM_COMM,IERROR) write (iout,*) "Send",msg_out(1)," complete" call flush(iout) call MPI_Recv(istart_conf,1,MPI_INTEGER,Previous, & msg_in(1),WHAM_COMM,STATUS(1,6),IERROR) write (iout,*) "Recv",msg_in(1)," complete" call flush(iout) call MPI_Send(iend_conf,1,MPI_INTEGER,Next,msg_out(2), & WHAM_COMM,IERROR) write (iout,*) "Send",msg_out(2)," complete" call flush(iout) call MPI_Recv(iend_conf,1,MPI_INTEGER,Previous, & msg_in(2),WHAM_COMM,STATUS(1,7),IERROR) write (iout,*) "Recv",msg_in(2)," complete" call flush(iout) call MPI_Buffer_Attach(buffer(1),8*(2*maxstr_proc*n_ene+800), & IERROR) call MPI_BSend(enetb(1,1,iprot),maxstr_proc*n_ene, & MPI_DOUBLE_PRECISION,Next,msg_out(3), & WHAM_COMM,IERROR) write (iout,*) "Send",msg_out(3)," complete (enetb)" call flush(iout) call MPI_Recv(enetb(1,1,iprot),maxstr_proc*n_ene, & MPI_DOUBLE_PRECISION,Previous,msg_in(3),WHAM_COMM, & STATUS(1,8),IERROR) write (iout,*) "Recv",msg_in(3)," complete (enetb)" call flush(iout) call MPI_Buffer_Detach(address,size,IERROR) c write (iout,*) "MPI_Buffer_Detach complete (enetb)" c call flush(iout) call MPI_Buffer_Attach(buffer(1),8*(2*maxstr_proc*nntyp+800), & IERROR) c write (iout,*) "MPI_Buffer_Attach complete (eneps)" c call flush(iout) call MPI_BSend(eneps(1,1,1,iprot),2*maxstr_proc*nntyp, & MPI_DOUBLE_PRECISION,Next,msg_out(4), & WHAM_COMM,IERROR) write (iout,*) "Send",msg_out(4)," complete (eneps)" call flush(iout) call MPI_Recv(eneps(1,1,1,iprot),2*maxstr_proc*nntyp, & MPI_DOUBLE_PRECISION,Previous,msg_in(4),WHAM_COMM, & STATUS(1,9),IERROR) write (iout,*) "Recv",msg_in(4)," complete (eneps)" call flush(iout) call MPI_Buffer_Detach(address,size,IERROR) call MPI_Buffer_Attach(buffer(1),8*(2*maxstr_proc*maxnatlike & +800),IERROR) call MPI_BSend(nu(1,1,1,iprot), & maxstr_proc*maxnatlike*maxdimnat, & MPI_DOUBLE_PRECISION,Next,msg_out(5), & WHAM_COMM,IERROR) write (iout,*) "Send",msg_out(5)," complete (nu)" call flush(iout) call MPI_Recv(nu(1,1,1,iprot), & maxstr_proc*maxnatlike*maxdimnat, & MPI_DOUBLE_PRECISION,Previous,msg_in(5),WHAM_COMM, & STATUS(1,10),IERROR) write (iout,*) "Recv",msg_in(5)," complete (nu)" call flush(iout) call MPI_Buffer_Detach(address,size,IERROR) write (iout,*) "Send and receive complete" call flush(iout) if (.not. mod_other_params) then write (iout,*) "Processor",me1," calling dawrite_ene", & " istart_conf",istart_conf," iend_conf",iend_conf call flush(iout) endif if (first_call) then c write (iout,*) "assignment of enetb_orig" do ii=1,iend_conf-istart_conf+1 do j=1,n_ene enetb_orig(ii,j,iprot)=enetb(ii,j,iprot) enddo c write (iout,'(i5,20f8.2)') c & ii,(enetb_orig(ii,j,iprot),j=1,n_ene) enddo endif if (.not. mod_other_params) & call dawrite_ene(iprot,istart_conf,iend_conf,ientout) do k=1,5 msg_out(k)=msg_in(k) enddo enddo #endif endif enddo ! i if (.not.init_ene .and. mod_fourier(nloctyp).eq.0 & .and. mod_fouriertor(nloctyp).eq.0 & .and. .not. mod_elec .and. .not. mod_scp) then close (ientin) else close (icbase) c Change AL 12/30/2017 if (.not. mod_other_params) close (ientout) endif ENDDO ! iprot init_ene=.false. c Lowest free energies of structural classes DO IPROT=1,NPROT IF (NCHUNK_CONF(IPROT).EQ.1) THEN #ifdef MPI do i=1,ntot_work(iprot) i2ii(i,iprot)=0 enddo ii=0 do i=indstart(me1,iprot),indend(me1,iprot) ii=ii+1 i2ii(i,iprot)=ii enddo istart_conf=indstart(me1,iprot) iend_conf=indend(me1,iprot) #else do i=1,ntot_work(iprot) i2ii(i,iprot)=i endif istart_conf=1 iend_conf=ntot_work(iprot) #endif #ifdef DEBUG write (iout,*) "i2ii at make_list" do i=1,ntot_work(iprot) write (iout,*) "i",i," i2ii",i2ii(i,iprot) enddo call flush(iout) #endif if (.not. mod_other_params) then open (ientin,file=benefiles(iprot),status="old", & form="unformatted",access="direct",recl=lenrec_ene(iprot)) c Change AL 12/30/2017 call daread_ene(iprot,istart_conf,iend_conf) endif call emin_search(iprot) ELSE if (.not. mod_other_params) &open (ientin,file=benefiles(iprot),status="old", & form="unformatted",access="direct",recl=lenrec_ene(iprot)) ipass_conf=0 #ifdef MPI do istart_conf=indstart(me1,iprot),indend(me1,iprot),maxstr_proc iend_conf=min0(istart_conf+maxstr_proc-1,indend(me1,iprot)) #else do istart_conf=1,ntot_work(iprot),maxstr_proc iend_conf=min0(istart_conf+maxstr_proc-1,ntot_work(iprot)) #endif c c Read the chunk of energies and derivatives off a DA scratchfile. c ipass_conf=ipass_conf+1 do i=1,ntot_work(iprot) i2ii(i,iprot)=0 enddo ii=0 do i=istart_conf,iend_conf ii=ii+1 i2ii(i,iprot)=ii enddo #ifdef DEBUG write (iout,*) "ipass_conf",ipass_conf, & " istart_conf",istart_conf," iend_conf",iend_conf do i=1,ntot_work(iprot) write (iout,*) "i",i," i2ii",i2ii(i,iprot) enddo call flush(iout) #endif c Change AL 12/30/2017 if (.not. mod_other_params) & call daread_ene(iprot,istart_conf,iend_conf) call emin_search(iprot) enddo close(ientin) ENDIF ENDDO #ifdef MPI c Complete the calculation of the lowest energies over all classes and c distribute the values to all procs do iprot=1,nprot do ibatch=1,natlike(iprot)+2 #ifdef DEBUG do ib=1,nbeta(ibatch,iprot) write (iout,'(7hELOWEST,3i3,f15.3,i12)') iprot,ibatch,ib, & elowest(ib,ibatch,iprot),ind_lowest(ib,ibatch,iprot) enddo #endif do ib=1,nbeta(ibatch,iprot) elowest_aux(1,ib)=elowest(ib,ibatch,iprot) elowest_aux(2,ib)=ind_lowest(ib,ibatch,iprot) enddo call MPI_Allreduce(elowest_aux(1,1),elowest_t(1,1), & nbeta(ibatch,iprot), & MPI_2DOUBLE_PRECISION, MPI_MINLOC, Comm1, IERROR) #ifdef DEBUG do ib=1,nbeta(ibatch,iprot) write (iout,*) "beta=",betaT(ib,ibatch,iprot) write (iout,'(9helowest_t,10f15.3)') & elowest_t(1,ib),elowest_t(2,ib) enddo write (iout,*) "Processor",me,me1," finished MPI_Reduce: 2" #endif do ib=1,nbeta(ibatch,iprot) elowest(ib,ibatch,iprot)=elowest_t(1,ib) ind_lowest(ib,ibatch,iprot)=elowest_t(2,ib) enddo enddo ! ibatc enddo ! iprot do iprot=1,nprot do ibatch=1,natlike(iprot)+2 do ib=1,nbeta(ibatch,iprot) write (iout,*) "iprot",iprot," ibatch",ibatch," ib",ib, & " elowest",elowest(ib,ibatch,iprot) enddo enddo enddo c c Allgather to provide all energies to all processors c do iprot=1,nprot do i=1,scount(me1,iprot) e_total_(i)=e_total(indstart(me1,iprot)+i-1,iprot) enddo call MPI_Allgatherv(e_total_(1), & scount(me1,iprot),MPI_DOUBLE_PRECISION,e_total(1,iprot), & scount(0,iprot),idispl(0,iprot),MPI_DOUBLE_PRECISION, & Comm1, IERROR) c call MPI_Allgatherv(e_total(indstart(me1,iprot),iprot), c & scount(me1,iprot),MPI_DOUBLE_PRECISION,e_total(1,iprot), c & scount(0,iprot),idispl(0,iprot),MPI_DOUBLE_PRECISION, c & Comm1, IERROR) enddo #endif c c Now determine which conformations will enter the database. c do iprot=1,nprot call restore_molinfo(iprot) c Clear the list of conformations do i=1,min0(ntot(iprot),maxstr) list_conf(i,iprot)=0 enddo do i=1,ntot(iprot) lflag(i)=.false. enddo c Make the list of conformations based on energy cut-off. nn=0 #ifdef DEBUG write (iout,*) "iprot",iprot," ibatch",ibatch," betmin", & betmin(ibatch,iprot) #endif #ifdef DEBUG write (iout,*) "e_lowb",e_lowb(iprot) write (iout,*) "t_lowb",t_lowb(iprot) #endif do ibatch=1,natlike(iprot)+2 do i=1,ntot(iprot) jj = i2ii(i,iprot) #ifdef MPI c write (iout,*) "i",i," ii",ii," indstart", c & indstart(me1,iprot)," indend",indend(me1,iprot) if (i.ge.indstart(me1,iprot).and.i.le.indend(me1,iprot)) & then #endif c write (iout,*) "i",i," ii",ii," kbatch",kbatch(i,iprot), c & " flag",lflag(i) if (.not.lflag(i)) then c if (eini(i,iprot).lt.e_lowb(iprot) .or. c & entfac(i,iprot).lt.t_lowb(iprot)) then c#ifdef DEBUG c write (iout,*) "Conformation",i," eini",eini(i,iprot), c & " entfac",entfac(i,iprot)," e_lowb", c & e_lowb(iprot), c & " t_lowb",t_lowb(iprot) c#endif c lflag(i)=.true. c goto 122 c endif do ib=1,nbeta(ibatch,iprot) dene=betaT(ib,ibatch,iprot)*(e_total(i,iprot) & -elowest(ib,ibatch,iprot))+entfac(i,iprot) #ifdef DEBUG write (iout,*) "beta",betaT(ib,ibatch,iprot), & " i",i," e_total", & e_total(i,iprot), & " elowest",elowest(ib,ibatch,iprot)," dene",dene, & " enecut",enecut(iprot) #endif if (dene.lt.enecut(iprot)) then nn=nn+1 list_conf(nn,iprot)=i lflag(i)=.true. goto 122 endif enddo 122 continue endif #ifdef MPI endif #endif enddo ! j enddo ! ibatch ntot_work(iprot)=nn if (nn.gt.maxstr) then write (iout,*) "Error - after applying cutoff the number", & " of conformations for protein ",i," exceeds MAXSTR:", & nn,maxstr write (iout,*) "The calculation is terminating." call flush(iout) #ifdef MPI call MPI_Finalize(ierror) #endif stop endif call imysort(ntot_work(iprot),list_conf(1,iprot),ipermut) write (iout,*) "Protein",iprot,ntot_work(iprot), & " conformations within scaled energy cut-off=",enecut(iprot), & " found at processor",me #ifdef MPI c c All workers get the complete list of conformations. c call MPI_Allgather(ntot_work(iprot),1,MPI_INTEGER, & scount(0,iprot),1,MPI_INTEGER,Comm1,IERROR) idispl(0,iprot)=0 do i=1,nprocs1-1 idispl(i,iprot)=idispl(i-1,iprot)+scount(i-1,iprot) enddo #ifdef DEBUG write (iout,*) "Protein",iprot," Scount and Idispl" do i=0,nprocs1-1 write (iout,*) i,scount(i,iprot),idispl(i,iprot) enddo write (iout,*) "Protein",i, & " local list of conformations of processor",me do i=1,ntot_work(iprot) write(iout,*) i,list_conf(i,iprot) enddo write (iout,*) "Before REDUCE: ntot_work",ntot_work(iprot) call flush(iout) #endif call MPI_Allreduce(ntot_work(iprot),nn,1,MPI_INTEGER,MPI_SUM, & Comm1,IERROR) ntot_work(iprot)=nn #ifdef DEBUG write (iout,*) "After REDUCE: ntot_work",ntot_work(iprot) call flush(iout) #endif call MPI_Allgatherv(list_conf(1,iprot), & scount(me1,iprot),MPI_INTEGER,list_conf_(1,iprot), & scount(0,iprot),idispl(0,iprot),MPI_INTEGER,Comm1,IERROR) do i=1,ntot_work(iprot) list_conf(i,iprot)=list_conf_(i,iprot) enddo #ifdef DEBUG write (iout,*) "Protein",i, & " global list of conformations of processor",me do i=1,ntot_work(iprot) write(iout,'(2i5,e15.5,33f7.3)')i,list_conf(i,iprot), & e_total(list_conf(i,iprot),iprot) c & ,(nu(k,i,iprot),k=1,natlike(iprot)) enddo call flush(iout) #endif #endif c c Construct the mapping of the new list to the original numbers of c conformations. c do i=1,ntot(iprot) tsil(i,iprot)=0 enddo do i=1,ntot_work(iprot) tsil(list_conf(i,iprot),iprot)=i enddo #ifdef DEBUG write (iout,*) "Protein",i," List-to-conformation mapping" do i=1,ntot(iprot) write(iout,*) i,tsil(i,iprot) enddo #endif enddo ! iprot c c Divide the work again based on the current lists of conformations c call work_partition(.true.) c c If the conformations fit into memory, read them off a DA scratchfile. c do iprot=1,nprot call restore_molinfo(iprot) #ifdef MPI nchunk_conf(iprot)=iroof(scount(me1,iprot),maxstr_proc) #else nchunk_conf(iprot)=iroof(ntot_work(iprot),maxstr_proc) #endif if (nchunk_conf(iprot).eq.1) then write (iout,*) "Protein",iprot, & " in-memory storage of conformations." if (init_ene .or. mod_fourier(nloctyp).gt.0 & .or.mod_fouriertor(nloctyp).gt.0.or.mod_elec.or.mod_scp) then open (icbase,file=bprotfiles(iprot),status="old", & form="unformatted",access="direct",recl=lenrec(iprot)) #ifdef MPI call daread_ccoords(iprot,indstart(me1,iprot), & indend(me1,iprot)) #else call daread_ccoords(iprot,1,ntot_work(iprot)) #endif close(icbase) else open (ientin,file=benefiles(iprot),status="old", & form="unformatted",access="direct",recl=lenrec_ene(iprot)) #ifdef MPI call daread_ene(iprot,indstart(me1,iprot), & indend(me1,iprot)) #else call daread_ene(iprot,1,ntot_work(iprot)) #endif close(ientin) open (icbase,file=bprotfiles(iprot),status="old", & form="unformatted",access="direct",recl=lenrec(iprot)) #ifdef MPI call daread_ccoords(iprot,indstart(me1,iprot), & indend(me1,iprot)) #else call daread_ccoords(iprot,1,ntot_work(iprot)) #endif close(icbase) endif #ifdef DEBUG write (iout,*) "Protein",i, & " global list of conformations of processor",me do i=1,ntot_work(iprot) write(iout,'(2i5,e15.5,33f7.3)')i,list_conf(i,iprot), & e_total(list_conf(i,iprot),iprot) enddo #endif else write (iout,*) "Protein",iprot, & " off-memory storage of conformations; ", & "energy will be evaluated in",nchunk_conf(iprot)," passes." endif enddo do i=1,n_ene ww_orig(i)=ww(i) enddo do iprot=1,nprot #ifdef DEBUG write (iout,*) "E_TOTAL and ETOT_ORIG of protein",iprot #endif do i=1,ntot_work(iprot) etot_orig(i,iprot)=e_total(list_conf(i,iprot),iprot) #ifdef DEBUG write (iout,*) i,list_conf(i,iprot), & e_total(list_conf(i,iprot),iprot),etot_orig(i,iprot) #endif enddo do ibatch=1,natlike(iprot)+2 do ib=1,nbeta(ibatch,iprot) emini = elowest(ib,ibatch,iprot) if (elowest(ib,ibatch,iprot) .lt. emini) & emini = elowest(ib,ibatch,iprot) emin_orig(ib,ibatch,iprot)=emini enddo enddo enddo #ifdef MPI if (me1 .eq. Master) call write_conf_count #else call write_conf_count #endif #ifdef DEBUG write (iout,*) "ELOWEST at the end of MAKE_LIST" do iprot=1,nprot do ibatch=1,natlike(iprot)+2 do ib=1,nbeta(ibatch,iprot) write (iout,*) "iprot",iprot," ibatch",ibatch," ib",ib, & " elowest",elowest(ib,ibatch,iprot) enddo enddo enddo #endif tcpu_fin = tcpu() - tcpu_ini write (iout,*) "Time for creating list of conformations",tcpu_fin call flush(iout) t_func = t_func + tcpu_fin return end c---------------------------------------------------------------------------- subroutine emin_search(iprot) implicit none include "DIMENSIONS" include "DIMENSIONS.ZSCOPT" #ifdef MPI include "mpif.h" integer IERROR,ErrCode,Status(MPI_STATUS_SIZE,10) integer req(10),msg_in(5),msg_out(5),address,size character*1 buffer(8*(2*maxstr_proc*nntyp+8000)) include "COMMON.MPI" #endif include "COMMON.WEIGHTS" include "COMMON.WEIGHTDER" include "COMMON.COMPAR" include "COMMON.ENERGIES" include "COMMON.IOUNITS" include "COMMON.VMCPAR" include "COMMON.NAMES" include "COMMON.INTERACT" include "COMMON.TIME1" include "COMMON.CHAIN" include "COMMON.PROTFILES" include "COMMON.VAR" include "COMMON.GEO" include "COMMON.CLASSES" C Define local variables integer i,ii,iii,kkk,jj,j,k,kk,l,iprot,ib,ibatch,nn integer ipass_conf,istart_conf,iend_conf,Previous,Next double precision energia(0:max_ene) double precision etoti,elowesti,dene double precision tcpu_ini,tcpu_fin,tcpu double precision etot_aux,enepsjk, & emini,elowest_aux(2,maxT) integer iroof,icant external iroof,icant logical lprn do ibatch=1,natlike(iprot)+2 do ib=1,nbeta(ibatch,iprot) elowest(ib,ibatch,iprot)=1.0d20 enddo enddo do ibatch=1,natlike(iprot)+2 do k=1,ntot(iprot) jj=i2ii(k,iprot) #ifdef MPI if (jj.gt.0) then #endif do ib=1,nbeta(ibatch,iprot) etoti=0.0d0 do kk=1,n_ene etoti=etoti+ww(kk)*enetb(jj,kk,iprot) & *escal(kk,ib,ibatch,iprot) enddo if (ib.eq.1 .and. etoti.lt.elowest_ent(1,ibatch,iprot)) & then elowest_ent(1,ibatch,iprot)=etoti elowest_ent(2,ibatch,iprot)=entfac(k,iprot) endif c if (ib.gt.1) c & etoti=etoti+entfac(k,iprot)/betaT(ib,ibatch,iprot) etoti=etoti+entfac(k,iprot)/betaT(ib,ibatch,iprot) if (etoti.lt.elowest(ib,ibatch,iprot)) then elowest(ib,ibatch,iprot)=etoti ind_lowest(ib,ibatch,iprot)=k endif c write (iout,*) ib,betaT(ib,ibatch,iprot),etoti, c & entfac(k,iprot) #ifdef DEBUG write (iout,'(2i5,2e15.5,f8.3)') k,jj, & etoti,elowest(ib,ibatch,iprot),betaT(ib,ibatch,iprot) #endif enddo ! ib #ifdef DEBUG write (iout,'(2i5,20f8.2)') j,jj,(enetb(jj,kk,iprot), & kk=1,n_ene) #endif #ifdef MPI endif #endif enddo ! k enddo ! ibatch return end