subroutine MREMD implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'mpif.h' include 'COMMON.CONTROL' include 'COMMON.VAR' include 'COMMON.MD' #ifndef LANG0 include 'COMMON.LANGEVIN' #else include 'COMMON.LANGEVIN.lang0' #endif include 'COMMON.CHAIN' include 'COMMON.DERIV' include 'COMMON.GEO' include 'COMMON.LOCAL' include 'COMMON.INTERACT' include 'COMMON.IOUNITS' include 'COMMON.NAMES' include 'COMMON.TIME1' include 'COMMON.REMD' include 'COMMON.SETUP' include 'COMMON.MUCA' include 'COMMON.HAIRPIN' double precision cm(3),L(3),vcm(3) double precision energia(0:n_ene) double precision remd_t_bath(maxprocs) double precision remd_ene(0:n_ene+1,maxprocs) integer iremd_acc(maxprocs),iremd_tot(maxprocs) integer ilen,rstcount external ilen character*50 tytul common /gucio/ cm integer itime integer nup(0:maxprocs),ndown(0:maxprocs) integer rep2i(0:maxprocs) integer itime_all(maxprocs) integer status(MPI_STATUS_SIZE) logical synflag,end_of_run,file_exist time00=MPI_WTIME() if(me.eq.king.or..not.out1file) & write (iout,*) 'MREMD',nodes,'time before',time00-walltime synflag=.false. mremd_rst_name=prefix(:ilen(prefix))//"_mremd.rst" cd print *,'MREMD',nodes cd print *,'mmm',me,remd_mlist,(remd_m(i),i=1,nrep) cde write (iout,*) "Start MREMD: me",me," t_bath",t_bath k=0 rep2i(k)=-1 do i=1,nrep iremd_acc(i)=0 iremd_tot(i)=0 do j=1,remd_m(i) i2rep(k)=i rep2i(i)=k k=k+1 enddo enddo c print *,'i2rep',me,i2rep(me) c print *,'rep2i',(rep2i(i),i=0,nrep) if (i2rep(me).eq.nrep) then nup(0)=0 else nup(0)=remd_m(i2rep(me)+1) k=rep2i(i2rep(me))+1 do i=1,nup(0) nup(i)=k k=k+1 enddo endif cd print '(i4,a4,100i4)',me,' nup',(nup(i),i=0,nup(0)) if (i2rep(me).eq.1) then ndown(0)=0 else ndown(0)=remd_m(i2rep(me)-1) k=rep2i(i2rep(me)-2)+1 do i=1,ndown(0) ndown(i)=k k=k+1 enddo endif cd print '(i4,a6,100i4)',me,' ndown',(ndown(i),i=0,ndown(0)) if(rest.and.restart1file) then inquire(file=mremd_rst_name,exist=file_exist) if(file_exist) call read1restart endif if(me.eq.king) then if (rest.and..not.restart1file) & inquire(file=mremd_rst_name,exist=file_exist) IF (rest.and.file_exist.and..not.restart1file) THEN open(irest2,file=mremd_rst_name,status='unknown') read (irest2,*) read (irest2,*) (i2rep(i),i=0,nodes-1) read (irest2,*) read (irest2,*) (ifirst(i),i=1,remd_m(1)) do il=1,nodes read (irest2,*) read (irest2,*) nupa(0,il),(nupa(i,il),i=1,nupa(0,il)) read (irest2,*) read (irest2,*) ndowna(0,il), & (ndowna(i,il),i=1,ndowna(0,il)) enddo close(irest2) write (iout,'(a6,1000i5)') "i2rep",(i2rep(i),i=0,nodes-1) write (iout,'(a6,1000i5)') "ifirst", & (ifirst(i),i=1,remd_m(1)) do il=1,nodes write (iout,'(a6,i4,a1,100i4)') "nupa",il,":", & (nupa(i,il),i=1,nupa(0,il)) write (iout,'(a6,i4,a1,100i4)') "ndowna",il,":", & (ndowna(i,il),i=1,ndowna(0,il)) enddo ELSE do il=1,remd_m(1) ifirst(il)=il enddo do il=1,nodes if (i2rep(il-1).eq.nrep) then nupa(0,il)=0 else nupa(0,il)=remd_m(i2rep(il-1)+1) k=rep2i(i2rep(il-1))+1 do i=1,nupa(0,il) nupa(i,il)=k+1 k=k+1 enddo endif if (i2rep(il-1).eq.1) then ndowna(0,il)=0 else ndowna(0,il)=remd_m(i2rep(il-1)-1) k=rep2i(i2rep(il-1)-2)+1 do i=1,ndowna(0,il) ndowna(i,il)=k+1 k=k+1 enddo endif enddo c write (iout,'(a6,100i4)') "ifirst", c & (ifirst(i),i=1,remd_m(1)) c do il=1,nodes c write (iout,'(a6,i4,a1,100i4)') "nupa",il,":", c & (nupa(i,il),i=1,nupa(0,il)) c write (iout,'(a6,i4,a1,100i4)') "ndowna",il,":", c & (ndowna(i,il),i=1,ndowna(0,il)) c enddo ENDIF endif c c t_bath=retmin+(retmax-retmin)*me/(nodes-1) if(.not.(rest.and.file_exist.and.restart1file)) then if (me .eq. king) then t_bath=retmin else t_bath=retmin+(retmax-retmin)*exp(float(i2rep(me)-nrep)) endif cd print *,'ttt',me,remd_tlist,(remd_t(i),i=1,nrep) if (remd_tlist) t_bath=remd_t(i2rep(me)) endif c stdfp=dsqrt(2*Rb*t_bath/d_time) do i=1,ntyp stdfsc(i)=dsqrt(2*Rb*t_bath/d_time) enddo c print *,'irep',me,t_bath if (.not.rest) then if (me.eq.king .or. .not. out1file) & write (iout,'(a60,f10.5)') "REMD Temperature:",t_bath call rescale_weights(t_bath) endif c------copy MD-------------- c The driver for molecular dynamics subroutines c------------------------------------------------ t_MDsetup=0.0d0 t_langsetup=0.0d0 t_MD=0.0d0 t_enegrad=0.0d0 t_sdsetup=0.0d0 if(me.eq.king.or..not.out1file) & write (iout,'(20(1h=),a20,20(1h=))') "MD calculation started" tt0 = tcpu() c Determine the inverse of the inertia matrix. call setup_MD_matrices c Initialize MD call init_MD if (rest) then if (me.eq.king .or. .not. out1file) & write (iout,'(a60,f10.5)') "REMD restart Temperature:",t_bath stdfp=dsqrt(2*Rb*t_bath/d_time) do i=1,ntyp stdfsc(i)=dsqrt(2*Rb*t_bath/d_time) enddo call rescale_weights(t_bath) endif t_MDsetup = tcpu()-tt0 rstcount=0 c Entering the MD loop tt0 = tcpu() if (lang.eq.2 .or. lang.eq.3) then call setup_fricmat if (lang.eq.2) then call sd_verlet_p_setup else call sd_verlet_ciccotti_setup endif #ifndef LANG0 do i=1,dimen do j=1,dimen pfric0_mat(i,j,0)=pfric_mat(i,j) afric0_mat(i,j,0)=afric_mat(i,j) vfric0_mat(i,j,0)=vfric_mat(i,j) prand0_mat(i,j,0)=prand_mat(i,j) vrand0_mat1(i,j,0)=vrand_mat1(i,j) vrand0_mat2(i,j,0)=vrand_mat2(i,j) enddo enddo #endif flag_stoch(0)=.true. do i=1,maxflag_stoch flag_stoch(i)=.false. enddo else if (lang.eq.1 .or. lang.eq.4) then call setup_fricmat endif time00=MPI_WTIME() if (me.eq.king .or. .not. out1file) & write(iout,*) 'Setup time',time00-walltime call flush(iout) t_langsetup=tcpu()-tt0 tt0=tcpu() itime=0 end_of_run=.false. do while(.not.end_of_run) itime=itime+1 if(itime.eq.n_timestep.and.me.eq.king) end_of_run=.true. if(mremdsync.and.itime.eq.n_timestep) end_of_run=.true. rstcount=rstcount+1 if (lang.gt.0 .and. surfarea .and. & mod(itime,reset_fricmat).eq.0) then if (lang.eq.2 .or. lang.eq.3) then call setup_fricmat if (lang.eq.2) then call sd_verlet_p_setup else call sd_verlet_ciccotti_setup endif #ifndef LANG0 do i=1,dimen do j=1,dimen pfric0_mat(i,j,0)=pfric_mat(i,j) afric0_mat(i,j,0)=afric_mat(i,j) vfric0_mat(i,j,0)=vfric_mat(i,j) prand0_mat(i,j,0)=prand_mat(i,j) vrand0_mat1(i,j,0)=vrand_mat1(i,j) vrand0_mat2(i,j,0)=vrand_mat2(i,j) enddo enddo #endif flag_stoch(0)=.true. do i=1,maxflag_stoch flag_stoch(i)=.false. enddo else if (lang.eq.1 .or. lang.eq.4) then call setup_fricmat endif write (iout,'(a,i10)') & "Friction matrix reset based on surface area, itime",itime endif if (reset_vel .and. tbf .and. lang.eq.0 & .and. mod(itime,count_reset_vel).eq.0) then call random_vel if (me.eq.king .or. .not. out1file) & write(iout,'(a,f20.2)') & "Velocities reset to random values, time",totT do i=0,2*nres do j=1,3 d_t_old(j,i)=d_t(j,i) enddo enddo endif if (reset_moment .and. mod(itime,count_reset_moment).eq.0) then call inertia_tensor call vcm_vel(vcm) do j=1,3 d_t(j,0)=d_t(j,0)-vcm(j) enddo call kinetic(EK) kinetic_T=2.0d0/(dimen*Rb)*EK scalfac=dsqrt(T_bath/kinetic_T) cd write(iout,'(a,f20.2)') "Momenta zeroed out, time",totT do i=0,2*nres do j=1,3 d_t_old(j,i)=scalfac*d_t(j,i) enddo enddo endif if (lang.ne.4) then if (RESPA) then c Time-reversible RESPA algorithm c (Tuckerman et al., J. Chem. Phys., 97, 1990, 1992) call RESPA_step(itime) else c Variable time step algorithm. call velverlet_step(itime) endif else call brown_step(itime) endif if(ntwe.ne.0) then if (mod(itime,ntwe).eq.0) call statout(itime) endif if (mod(itime,ntwx).eq.0.and..not.traj1file) then write (tytul,'("time",f8.2," temp",f8.1)') totT,t_bath if(mdpdb) then call hairpin(.true.,nharp,iharp) call secondary2(.true.) call pdbout(potE,tytul,ipdb) else call cartout(totT) endif endif if ((rstcount.eq.1000.or.itime.eq.n_timestep) & .and..not.restart1file) then if(me.eq.king) then open(irest1,file=mremd_rst_name,status='unknown') write (irest1,*) "i2rep" write (irest1,*) (i2rep(i),i=0,nodes-1) write (irest1,*) "ifirst" write (irest1,*) (ifirst(i),i=1,remd_m(1)) do il=1,nodes write (irest1,*) "nupa",il write (irest1,*) nupa(0,il),(nupa(i,il),i=1,nupa(0,il)) write (irest1,*) "ndowna",il write (irest1,*) ndowna(0,il), & (ndowna(i,il),i=1,ndowna(0,il)) enddo close(irest1) endif open(irest2,file=rest2name,status='unknown') write(irest2,*) totT,EK,potE,totE,t_bath do i=1,2*nres write (irest2,'(3e15.5)') (d_t(j,i),j=1,3) enddo do i=1,2*nres write (irest2,'(3e15.5)') (dc(j,i),j=1,3) enddo close(irest2) rstcount=0 endif c REMD - exchange c forced synchronization if (mod(itime,100).eq.0 .and. me.ne.king & .and. .not. mremdsync) then synflag=.false. call mpi_iprobe(0,101,CG_COMM,synflag,status,ierr) if (synflag) then call mpi_recv(itime_master, 1, MPI_INTEGER, & 0,101,CG_COMM, status, ierr) if (.not. out1file) then write(iout,*) 'REMD synchro at',itime_master,itime else call mpi_gather(itime,1,mpi_integer, & itime_all,1,mpi_integer,king, & CG_COMM,ierr) endif if(itime_master.ge.n_timestep) end_of_run=.true. call flush(iout) endif endif c REMD - exchange if ((mod(itime,nstex).eq.0.and.me.eq.king & .or.end_of_run.and.me.eq.king ) & .and. .not. mremdsync ) then synflag=.true. time00=MPI_WTIME() do i=1,nodes-1 call mpi_isend(itime,1,MPI_INTEGER,i,101, & CG_COMM, ireq, ierr) cd write(iout,*) 'REMD synchro with',i cd call flush(iout) enddo time02=MPI_WTIME() write(iout,*) 'REMD synchro at',itime,'time=',time02-time00 if (out1file) then call mpi_gather(itime,1,mpi_integer, & itime_all,1,mpi_integer,king, & CG_COMM,ierr) write(iout,'(a18,8000i8)') 'REMD synchro itime', & (itime_all(i),i=1,nodes) endif call flush(iout) endif if(mremdsync .and. mod(itime,nstex).eq.0) then synflag=.true. if (me.eq.king .or. .not. out1file) & write(iout,*) 'REMD synchro at',itime call flush(iout) endif if(synflag.and..not.end_of_run) then time00=MPI_WTIME() synflag=.false. cd write(iout,*) 'REMD before',me,t_bath c call mpi_gather(t_bath,1,mpi_double_precision, c & remd_t_bath,1,mpi_double_precision,king, c & CG_COMM,ierr) potEcomp(n_ene+1)=t_bath call mpi_gather(potEcomp(0),n_ene+2,mpi_double_precision, & remd_ene(0,1),n_ene+2,mpi_double_precision,king, & CG_COMM,ierr) if(lmuca) then call mpi_gather(elow,1,mpi_double_precision, & elowi,1,mpi_double_precision,king, & CG_COMM,ierr) call mpi_gather(ehigh,1,mpi_double_precision, & ehighi,1,mpi_double_precision,king, & CG_COMM,ierr) endif if (restart1file) call write1rst if (traj1file) call write1traj if (me.eq.king) then do i=1,nodes remd_t_bath(i)=remd_ene(n_ene+1,i) enddo if(lmuca) then co write(iout,*) 'REMD exchange temp,ene,elow,ehigh' do i=1,nodes write(iout,'(i4,4f12.5)') i,remd_t_bath(i),remd_ene(0,i), & elowi(i),ehighi(i) enddo else cd write(iout,*) 'REMD exchange temp,ene' c do i=1,nodes co write(iout,'(i4,2f12.5)') i,remd_t_bath(i),remd_ene(0,i) c write(iout,'(6f12.5)') (remd_ene(j,i),j=1,n_ene) c enddo endif c------------------------------------- do irr=1,remd_m(1) i=ifirst(iran_num(1,remd_m(1))) do ii=1,nodes-1 if(i.gt.0.and.nupa(0,i).gt.0) then iex=nupa(iran_num(1,int(nupa(0,i))),i) if (lmuca) then call muca_delta(remd_t_bath,remd_ene,i,iex,delta) else c Swap temperatures between conformations i and iex with recalculating the free energies c following temperature changes. ene_iex_iex=remd_ene(0,iex) ene_i_i=remd_ene(0,i) co write (iout,*) "rescaling weights with temperature", co & remd_t_bath(i) call rescale_weights(remd_t_bath(i)) call sum_energy(remd_ene(0,iex),.false.) ene_iex_i=remd_ene(0,iex) cd write (iout,*) "ene_iex_i",remd_ene(0,iex) call sum_energy(remd_ene(0,i),.false.) cd write (iout,*) "ene_i_i",remd_ene(0,i) c write (iout,*) "rescaling weights with temperature", c & remd_t_bath(iex) if (real(ene_i_i).ne.real(remd_ene(0,i))) then write (iout,*) "ERROR: inconsistent energies:",i, & ene_i_i,remd_ene(0,i) endif call rescale_weights(remd_t_bath(iex)) call sum_energy(remd_ene(0,i),.false.) c write (iout,*) "ene_i_iex",remd_ene(0,i) ene_i_iex=remd_ene(0,i) call sum_energy(remd_ene(0,iex),.false.) if (real(ene_iex_iex).ne.real(remd_ene(0,iex))) then write (iout,*) "ERROR: inconsistent energies:",iex, & ene_iex_iex,remd_ene(0,iex) endif c write (iout,*) "ene_iex_iex",remd_ene(0,iex) c write (iout,*) "i",i," iex",iex co write (iout,'(4(a,e15.5))') "ene_i_i",ene_i_i, co & " ene_i_iex",ene_i_iex, co & " ene_iex_i",ene_iex_i," ene_iex_iex",ene_iex_iex delta=(ene_iex_iex-ene_i_iex)/(Rb*remd_t_bath(iex))- & (ene_iex_i-ene_i_i)/(Rb*remd_t_bath(i)) delta=-delta c write(iout,*) 'delta',delta c delta=(remd_t_bath(i)-remd_t_bath(iex))* c & (remd_ene(i)-remd_ene(iex))/Rb/ c & (remd_t_bath(i)*remd_t_bath(iex)) endif if (delta .gt. 50.0d0) then delta=0.0d0 else #ifdef OSF if(isnan(delta))then delta=0.0d0 else if (delta.lt.-50.0d0) then delta=dexp(50.0d0) else delta=dexp(-delta) endif #else delta=dexp(-delta) #endif endif iremd_tot(i2rep(i-1))=iremd_tot(i2rep(i-1))+1 xxx=ran_number(0.0d0,1.0d0) co write(iout,'(2i4,a6,2f12.5)') i,iex,' delta',delta,xxx if (delta .gt. xxx) then tmp=remd_t_bath(i) remd_t_bath(i)=remd_t_bath(iex) remd_t_bath(iex)=tmp remd_ene(0,i)=ene_i_iex remd_ene(0,iex)=ene_iex_i if(lmuca) then tmp=elowi(i) elowi(i)=elowi(iex) elowi(iex)=tmp tmp=ehighi(i) ehighi(i)=ehighi(iex) ehighi(iex)=tmp endif do k=0,nodes itmp=nupa(k,i) nupa(k,i)=nupa(k,iex) nupa(k,iex)=itmp itmp=ndowna(k,i) ndowna(k,i)=ndowna(k,iex) ndowna(k,iex)=itmp enddo do il=1,nodes if (ifirst(il).eq.i) ifirst(il)=iex do k=1,nupa(0,il) if (nupa(k,il).eq.i) then nupa(k,il)=iex elseif (nupa(k,il).eq.iex) then nupa(k,il)=i endif enddo do k=1,ndowna(0,il) if (ndowna(k,il).eq.i) then ndowna(k,il)=iex elseif (ndowna(k,il).eq.iex) then ndowna(k,il)=i endif enddo enddo iremd_acc(i2rep(i-1))=iremd_acc(i2rep(i-1))+1 itmp=i2rep(i-1) i2rep(i-1)=i2rep(iex-1) i2rep(iex-1)=itmp cd write(iout,*) 'exchange',i,iex cd write (iout,'(a8,100i4)') "@ ifirst", cd & (ifirst(k),k=1,remd_m(1)) cd do il=1,nodes cd write (iout,'(a8,i4,a1,100i4)') "@ nupa",il,":", cd & (nupa(k,il),k=1,nupa(0,il)) cd write (iout,'(a8,i4,a1,100i4)') "@ ndowna",il,":", cd & (ndowna(k,il),k=1,ndowna(0,il)) cd enddo else remd_ene(0,iex)=ene_iex_iex remd_ene(0,i)=ene_i_i i=iex endif endif enddo enddo c------------------------------------- do i=1,nrep if(iremd_tot(i).ne.0) & write(iout,'(a3,i4,2f12.5)') 'ACC',i,remd_t(i) & ,iremd_acc(i)/(1.0*iremd_tot(i)) enddo cd write (iout,'(a6,100i4)') "ifirst", cd & (ifirst(i),i=1,remd_m(1)) cd do il=1,nodes cd write (iout,'(a5,i4,a1,100i4)') "nup",il,":", cd & (nupa(i,il),i=1,nupa(0,il)) cd write (iout,'(a5,i4,a1,100i4)') "ndown",il,":", cd & (ndowna(i,il),i=1,ndowna(0,il)) cd enddo endif call mpi_scatter(remd_t_bath,1,mpi_double_precision, & t_bath,1,mpi_double_precision,king, & CG_COMM,ierr) if(lmuca) then call mpi_scatter(elowi,1,mpi_double_precision, & elow,1,mpi_double_precision,king, & CG_COMM,ierr) call mpi_scatter(ehighi,1,mpi_double_precision, & ehigh,1,mpi_double_precision,king, & CG_COMM,ierr) endif call rescale_weights(t_bath) co write (iout,*) "Processor",me, co & " rescaling weights with temperature",t_bath stdfp=dsqrt(2*Rb*t_bath/d_time) do i=1,ntyp stdfsc(i)=dsqrt(2*Rb*t_bath/d_time) enddo cde write(iout,*) 'REMD after',me,t_bath time02=MPI_WTIME() if (me.eq.king .or. .not. out1file) then write(iout,*) 'REMD exchange time=',time02-time00 call flush(iout) endif endif enddo if (restart1file) then if (me.eq.king .or. .not. out1file) & write(iout,*) 'writing restart at the end of run' call write1rst endif if (traj1file) call write1traj t_MD=tcpu()-tt0 if (me.eq.king .or. .not. out1file) then write (iout,'(//35(1h=),a10,35(1h=)/10(/a40,1pe15.5))') & ' Timing ', & 'MD calculations setup:',t_MDsetup, & 'Energy & gradient evaluation:',t_enegrad, & 'Stochastic MD setup:',t_langsetup, & 'Stochastic MD step setup:',t_sdsetup, & 'MD steps:',t_MD write (iout,'(/28(1h=),a25,27(1h=))') & ' End of MD calculation ' endif return end subroutine write1rst implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'mpif.h' include 'COMMON.MD' include 'COMMON.IOUNITS' include 'COMMON.REMD' include 'COMMON.SETUP' include 'COMMON.CHAIN' include 'COMMON.SBRIDGE' include 'COMMON.INTERACT' real d_restart1(3,2*maxres*maxprocs),r_d(3,2*maxres), & d_restart2(3,2*maxres*maxprocs) real t5_restart1(5) integer iret,itmp t5_restart1(1)=totT t5_restart1(2)=EK t5_restart1(3)=potE t5_restart1(4)=t_bath t5_restart1(5)=Uconst call mpi_gather(t5_restart1,5,mpi_real, & t_restart1,5,mpi_real,king,CG_COMM,ierr) do i=1,2*nres do j=1,3 r_d(j,i)=d_t(j,i) enddo enddo call mpi_gather(r_d,3*2*nres,mpi_real, & d_restart1,3*2*nres,mpi_real,king, & CG_COMM,ierr) do i=1,2*nres do j=1,3 r_d(j,i)=dc(j,i) enddo enddo call mpi_gather(r_d,3*2*nres,mpi_real, & d_restart2,3*2*nres,mpi_real,king, & CG_COMM,ierr) if(me.eq.king) then open(irest1,file=mremd_rst_name,status='unknown') write (irest1,*) (i2rep(i),i=0,nodes-1) write (irest1,*) (ifirst(i),i=1,remd_m(1)) do il=1,nodes write (irest1,*) nupa(0,il),(nupa(i,il),i=1,nupa(0,il)) write (irest1,*) ndowna(0,il), & (ndowna(i,il),i=1,ndowna(0,il)) enddo do il=1,nodes write(irest1,*) (t_restart1(j,il),j=1,4) enddo do il=0,nodes-1 do i=1,2*nres write(irest1,'(3e15.5)') (d_restart1(j,i+2*nres*il),j=1,3) enddo enddo do il=0,nodes-1 do i=1,2*nres write(irest1,'(3e15.5)') (d_restart2(j,i+2*nres*il),j=1,3) enddo enddo close(irest1) endif return end subroutine write1traj implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'mpif.h' include 'COMMON.MD' include 'COMMON.IOUNITS' include 'COMMON.REMD' include 'COMMON.SETUP' include 'COMMON.CHAIN' include 'COMMON.SBRIDGE' include 'COMMON.INTERACT' real t5_restart1(5) integer iret,itmp real xcoord(3,maxres2+2),prec real r_qfrag(50),r_qpair(100) real p_qfrag(50*maxprocs),p_qpair(100*maxprocs) real p_c(3,(maxres2+2)*maxprocs),r_c(3,maxres2+2) if(.not.restart1file) then t5_restart1(1)=totT t5_restart1(2)=EK t5_restart1(3)=potE t5_restart1(4)=t_bath t5_restart1(5)=Uconst call mpi_gather(t5_restart1,5,mpi_real, & t_restart1,5,mpi_real,king,CG_COMM,ierr) endif do i=1,nfrag r_qfrag(i)=qfrag(i) enddo do i=1,npair r_qpair(i)=qpair(i) enddo call mpi_gather(r_qfrag,nfrag,mpi_real, & p_qfrag,nfrag,mpi_real,king, & CG_COMM,ierr) call mpi_gather(qpair,nfrag,mpi_real, & p_qpair,nfrag,mpi_real,king, & CG_COMM,ierr) do i=1,nres*2 do j=1,3 r_c(j,i)=c(j,i) enddo enddo call mpi_gather(r_c,3*2*nres,mpi_real, & p_c,3*2*nres,mpi_real,king, & CG_COMM,ierr) if(me.eq.king) then call xdrfopen(ixdrf,cartname, "a", iret) do il=1,nodes call xdrffloat(ixdrf, real(t_restart1(1,il)), iret) call xdrffloat(ixdrf, real(t_restart1(3,il)), iret) call xdrffloat(ixdrf, real(t_restart1(5,il)), iret) call xdrffloat(ixdrf, real(t_restart1(4,il)), iret) call xdrfint(ixdrf, nss, iret) do j=1,nss call xdrfint(ixdrf, ihpb(j), iret) call xdrfint(ixdrf, jhpb(j), iret) enddo call xdrfint(ixdrf, nfrag+npair, iret) do i=1,nfrag call xdrffloat(ixdrf, p_qfrag(i+(il-1)*nfrag), iret) enddo do i=1,npair call xdrffloat(ixdrf, p_qpair(i+(il-1)*npair), iret) enddo prec=10000.0 do i=1,nres do j=1,3 xcoord(j,i)=p_c(j,i+(il-1)*nres*2) enddo enddo do i=nnt,nct do j=1,3 xcoord(j,nres+i-nnt+1)=p_c(j,i+nres+(il-1)*nres*2) enddo enddo itmp=nres+nct-nnt+1 call xdrf3dfcoord(ixdrf, xcoord, itmp, prec, iret) enddo call xdrfclose(ixdrf, iret) endif return end subroutine read1restart implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'mpif.h' include 'COMMON.MD' include 'COMMON.IOUNITS' include 'COMMON.REMD' include 'COMMON.SETUP' include 'COMMON.CHAIN' include 'COMMON.SBRIDGE' include 'COMMON.INTERACT' real d_restart1(3,2*maxres*maxprocs),r_d(3,2*maxres), & t5_restart1(5) if(me.eq.king)then open(irest2,file=mremd_rst_name,status='unknown') read (irest2,*) (i2rep(i),i=0,nodes-1) read (irest2,*) (ifirst(i),i=1,remd_m(1)) do il=1,nodes read (irest2,*) nupa(0,il),(nupa(i,il),i=1,nupa(0,il)) read (irest2,*) ndowna(0,il), & (ndowna(i,il),i=1,ndowna(0,il)) enddo do il=1,nodes read(irest2,*) (t_restart1(j,il),j=1,4) enddo endif call mpi_scatter(t_restart1,5,mpi_real, & t5_restart1,5,mpi_real,king,CG_COMM,ierr) totT=t5_restart1(1) EK=t5_restart1(2) potE=t5_restart1(3) t_bath=t5_restart1(4) if(me.eq.king)then do il=0,nodes-1 do i=1,2*nres read(irest2,'(3e15.5)') & (d_restart1(j,i+2*nres*il),j=1,3) enddo enddo endif call mpi_scatter(d_restart1,3*2*nres,mpi_real, & r_d,3*2*nres,mpi_real,king,CG_COMM,ierr) do i=1,2*nres do j=1,3 d_t(j,i)=r_d(j,i) enddo enddo if(me.eq.king)then do il=0,nodes-1 do i=1,2*nres read(irest2,'(3e15.5)') & (d_restart1(j,i+2*nres*il),j=1,3) enddo enddo endif call mpi_scatter(d_restart1,3*2*nres,mpi_real, & r_d,3*2*nres,mpi_real,king,CG_COMM,ierr) do i=1,2*nres do j=1,3 dc(j,i)=r_d(j,i) enddo enddo if(me.eq.king) close(irest2) return end