2 implicit real*8 (a-h,o-z)
5 include 'COMMON.CONTROL'
9 include 'COMMON.LANGEVIN'
11 include 'COMMON.LANGEVIN.lang0'
13 include 'COMMON.CHAIN'
14 include 'COMMON.DERIV'
16 include 'COMMON.LOCAL'
17 include 'COMMON.INTERACT'
18 include 'COMMON.IOUNITS'
19 include 'COMMON.NAMES'
20 include 'COMMON.TIME1'
22 include 'COMMON.SETUP'
24 include 'COMMON.HAIRPIN'
26 double precision cm(3),L(3),vcm(3)
27 double precision energia(0:n_ene)
28 double precision remd_t_bath(maxprocs)
29 integer iremd_iset(maxprocs)
31 & (maxprocs/4,maxprocs/20,maxprocs/200,maxprocs/200)
32 double precision remd_ene(0:n_ene+4,maxprocs)
33 integer iremd_acc(maxprocs),iremd_tot(maxprocs)
34 integer iremd_acc_usa(maxprocs),iremd_tot_usa(maxprocs)
40 cold integer nup(0:maxprocs),ndown(0:maxprocs)
41 integer rep2i(0:maxprocs),ireqi(maxprocs)
42 integer icache_all(maxprocs)
43 integer status(MPI_STATUS_SIZE),statusi(MPI_STATUS_SIZE,maxprocs)
44 logical synflag,end_of_run,file_exist /.false./,ovrtim
50 if(me.eq.king.or..not.out1file) then
51 write (iout,*) 'MREMD',nodes,'time before',time00-walltime
52 write (iout,*) "NREP=",nrep
56 if (ilen(tmpdir).gt.0 .and. (me.eq.king)) then
57 call copy_to_tmp(pref_orig(:ilen(pref_orig))//"_mremd.rst")
59 mremd_rst_name=prefix(:ilen(prefix))//"_mremd.rst"
61 cd print *,'MREMD',nodes
62 cd print *,'mmm',me,remd_mlist,(remd_m(i),i=1,nrep)
63 cde write (iout,*) "Start MREMD: me",me," t_bath",t_bath
67 do il1=1,max0(mset(il),1)
83 if(me.eq.king.or..not.out1file) then
84 write(iout,*) (i2rep(i),i=0,nodes-1)
85 write(iout,*) (i2set(i),i=0,nodes-1)
90 write(iout,*) i,j,il,il1,i_index(i,j,il,il1)
97 c print *,'i2rep',me,i2rep(me)
98 c print *,'rep2i',(rep2i(i),i=0,nrep)
100 cold if (i2rep(me).eq.nrep) then
103 cold nup(0)=remd_m(i2rep(me)+1)
104 cold k=rep2i(int(i2rep(me)))+1
111 cd print '(i4,a4,100i4)',me,' nup',(nup(i),i=0,nup(0))
113 cold if (i2rep(me).eq.1) then
116 cold ndown(0)=remd_m(i2rep(me)-1)
117 cold k=rep2i(i2rep(me)-2)+1
124 cd print '(i4,a6,100i4)',me,' ndown',(ndown(i),i=0,ndown(0))
127 write (*,*) "Processor",me," rest",rest,"
128 & restart1fie",restart1file
129 if(rest.and.restart1file) then
131 & inquire(file=mremd_rst_name,exist=file_exist)
132 cd write (*,*) me," Before broadcast: file_exist",file_exist
133 call MPI_Bcast(file_exist,1,MPI_LOGICAL,king,CG_COMM,
135 cd write (*,*) me," After broadcast: file_exist",file_exist
137 if(me.eq.king.or..not.out1file)
138 & write (iout,*) 'Master is reading restart1file'
139 call read1restart(i_index)
141 if(me.eq.king.or..not.out1file)
142 & write (iout,*) 'WARNING : no restart1file'
145 if(me.eq.king.or..not.out1file) then
146 write(iout,*) "i2set",(i2set(i),i=0,nodes-1)
147 write(iout,*) "i_index"
152 write(iout,*) i,j,il,il1,i_index(i,j,il,il1)
161 if (rest.and..not.restart1file)
162 & inquire(file=mremd_rst_name,exist=file_exist)
163 if(.not.file_exist.and.rest.and..not.restart1file)
164 & write(iout,*) 'WARNING : no restart file',mremd_rst_name
165 IF (rest.and.file_exist.and..not.restart1file) THEN
166 write (iout,*) 'Master is reading restart file',
168 open(irest2,file=mremd_rst_name,status='unknown')
170 read (irest2,*) (i2rep(i),i=0,nodes-1)
172 read (irest2,*) (ifirst(i),i=1,remd_m(1))
175 read (irest2,*) nupa(0,il),(nupa(i,il),i=1,nupa(0,il))
177 read (irest2,*) ndowna(0,il),
178 & (ndowna(i,il),i=1,ndowna(0,il))
184 read (irest2,*) (mset(i),i=1,nset)
186 read (irest2,*) (i2set(i),i=0,nodes-1)
191 read(irest2,*) (i_index(i,j,il,il1),j=1,remd_m(i))
196 write(iout,*) "i2set",(i2set(i),i=0,nodes-1)
197 write(iout,*) "i_index"
202 write(iout,*) i,j,il,il1,i_index(i,j,il,il1)
211 write (iout,'(a6,1000i5)') "i2rep",(i2rep(i),i=0,nodes-1)
212 write (iout,'(a6,1000i5)') "ifirst",
213 & (ifirst(i),i=1,remd_m(1))
215 write (iout,'(a6,i4,a1,100i4)') "nupa",il,":",
216 & (nupa(i,il),i=1,nupa(0,il))
217 write (iout,'(a6,i4,a1,100i4)') "ndowna",il,":",
218 & (ndowna(i,il),i=1,ndowna(0,il))
220 ELSE IF (.not.(rest.and.file_exist)) THEN
226 if (i2rep(il-1).eq.nrep) then
229 nupa(0,il)=remd_m(i2rep(il-1)+1)
230 k=rep2i(int(i2rep(il-1)))+1
236 if (i2rep(il-1).eq.1) then
239 ndowna(0,il)=remd_m(i2rep(il-1)-1)
240 k=rep2i(i2rep(il-1)-2)+1
248 write (iout,'(a6,100i4)') "ifirst",
249 & (ifirst(i),i=1,remd_m(1))
251 write (iout,'(a6,i4,a1,100i4)') "nupa",il,":",
252 & (nupa(i,il),i=1,nupa(0,il))
253 write (iout,'(a6,i4,a1,100i4)') "ndowna",il,":",
254 & (ndowna(i,il),i=1,ndowna(0,il))
260 c t_bath=retmin+(retmax-retmin)*me/(nodes-1)
261 if(.not.(rest.and.file_exist.and.restart1file)) then
262 if (me .eq. king) then
265 t_bath=retmin+(retmax-retmin)*exp(float(i2rep(me)-nrep))
267 cd print *,'ttt',me,remd_tlist,(remd_t(i),i=1,nrep)
268 if (remd_tlist) t_bath=remd_t(int(i2rep(me)))
273 if(me.eq.king.or..not.out1file)
274 & write(iout,*) me,"iset=",iset,"t_bath=",t_bath
277 stdfp=dsqrt(2*Rb*t_bath/d_time)
279 stdfsc(i)=dsqrt(2*Rb*t_bath/d_time)
282 c print *,'irep',me,t_bath
284 if (me.eq.king .or. .not. out1file)
285 & write (iout,'(a60,f10.5)') "REMD Temperature:",t_bath
286 call rescale_weights(t_bath)
290 c------copy MD--------------
291 c The driver for molecular dynamics subroutines
292 c------------------------------------------------
298 if(me.eq.king.or..not.out1file)
299 & write (iout,'(20(1h=),a20,20(1h=))') "MD calculation started"
305 c Determine the inverse of the inertia matrix.
306 call setup_MD_matrices
307 write(iout,*), "TU", dc(1,0)
310 write(iout,*), "TU2", dc(1,0)
313 if (me.eq.king .or. .not. out1file)
314 & write (iout,'(a60,f10.5)') "REMD restart Temperature:",t_bath
315 stdfp=dsqrt(2*Rb*t_bath/d_time)
317 stdfsc(i)=dsqrt(2*Rb*t_bath/d_time)
319 call rescale_weights(t_bath)
323 t_MDsetup = MPI_Wtime()-tt0
325 t_MDsetup = tcpu()-tt0
328 c Entering the MD loop
334 if (lang.eq.2 .or. lang.eq.3) then
338 call sd_verlet_p_setup
340 call sd_verlet_ciccotti_setup
344 pfric0_mat(i,j,0)=pfric_mat(i,j)
345 afric0_mat(i,j,0)=afric_mat(i,j)
346 vfric0_mat(i,j,0)=vfric_mat(i,j)
347 prand0_mat(i,j,0)=prand_mat(i,j)
348 vrand0_mat1(i,j,0)=vrand_mat1(i,j)
349 vrand0_mat2(i,j,0)=vrand_mat2(i,j)
354 flag_stoch(i)=.false.
358 & "LANG=2 or 3 NOT SUPPORTED. Recompile without -DLANG0"
360 call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
364 else if (lang.eq.1 .or. lang.eq.4) then
368 if (me.eq.king .or. .not. out1file)
369 & write(iout,*) 'Setup time',time00-walltime
372 t_langsetup=MPI_Wtime()-tt0
375 t_langsetup=tcpu()-tt0
380 do while(.not.end_of_run)
382 if(itime.eq.n_timestep.and.me.eq.king) end_of_run=.true.
383 if(mremdsync.and.itime.eq.n_timestep) end_of_run=.true.
385 if (lang.gt.0 .and. surfarea .and.
386 & mod(itime,reset_fricmat).eq.0) then
387 if (lang.eq.2 .or. lang.eq.3) then
391 call sd_verlet_p_setup
393 call sd_verlet_ciccotti_setup
397 pfric0_mat(i,j,0)=pfric_mat(i,j)
398 afric0_mat(i,j,0)=afric_mat(i,j)
399 vfric0_mat(i,j,0)=vfric_mat(i,j)
400 prand0_mat(i,j,0)=prand_mat(i,j)
401 vrand0_mat1(i,j,0)=vrand_mat1(i,j)
402 vrand0_mat2(i,j,0)=vrand_mat2(i,j)
407 flag_stoch(i)=.false.
410 else if (lang.eq.1 .or. lang.eq.4) then
413 write (iout,'(a,i10)')
414 & "Friction matrix reset based on surface area, itime",itime
416 if (reset_vel .and. tbf .and. lang.eq.0
417 & .and. mod(itime,count_reset_vel).eq.0) then
419 if (me.eq.king .or. .not. out1file)
420 & write(iout,'(a,f20.2)')
421 & "Velocities reset to random values, time",totT
424 d_t_old(j,i)=d_t(j,i)
428 if (reset_moment .and. mod(itime,count_reset_moment).eq.0) then
432 d_t(j,0)=d_t(j,0)-vcm(j)
435 kinetic_T=2.0d0/(dimen3*Rb)*EK
436 scalfac=dsqrt(T_bath/kinetic_T)
437 cd write(iout,'(a,f20.2)') "Momenta zeroed out, time",totT
440 d_t_old(j,i)=scalfac*d_t(j,i)
446 c Time-reversible RESPA algorithm
447 c (Tuckerman et al., J. Chem. Phys., 97, 1990, 1992)
448 call RESPA_step(itime)
450 c Variable time step algorithm.
451 call velverlet_step(itime)
455 call brown_step(itime)
457 print *,"Brown dynamics not here!"
459 call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
465 if (mod(itime,ntwe).eq.0) call statout(itime)
467 if (mod(itime,ntwx).eq.0.and..not.traj1file) then
468 write (tytul,'("time",f8.2," temp",f8.1)') totT,t_bath
470 call hairpin(.true.,nharp,iharp)
471 call secondary2(.true.)
472 call pdbout(potE,tytul,ipdb)
477 if (mod(itime,ntwx).eq.0.and.traj1file) then
478 if(ntwx_cache.lt.max_cache_traj_use) then
479 ntwx_cache=ntwx_cache+1
481 if (max_cache_traj_use.ne.1)
482 & print *,itime,"processor ",me," over cache ",ntwx_cache
485 totT_cache(i)=totT_cache(i+1)
486 EK_cache(i)=EK_cache(i+1)
487 potE_cache(i)=potE_cache(i+1)
488 t_bath_cache(i)=t_bath_cache(i+1)
489 Uconst_cache(i)=Uconst_cache(i+1)
490 iset_cache(i)=iset_cache(i+1)
493 qfrag_cache(ii,i)=qfrag_cache(ii,i+1)
496 qpair_cache(ii,i)=qpair_cache(ii,i+1)
499 utheta_cache(ii,i)=utheta_cache(ii,i+1)
500 ugamma_cache(ii,i)=ugamma_cache(ii,i+1)
501 uscdiff_cache(ii,i)=uscdiff_cache(ii,i+1)
507 c_cache(j,ii,i)=c_cache(j,ii,i+1)
513 totT_cache(ntwx_cache)=totT
514 EK_cache(ntwx_cache)=EK
515 potE_cache(ntwx_cache)=potE
516 t_bath_cache(ntwx_cache)=t_bath
517 Uconst_cache(ntwx_cache)=Uconst
518 iset_cache(ntwx_cache)=iset
521 qfrag_cache(i,ntwx_cache)=qfrag(i)
524 qpair_cache(i,ntwx_cache)=qpair(i)
527 utheta_cache(i,ntwx_cache)=utheta(i)
528 ugamma_cache(i,ntwx_cache)=ugamma(i)
529 uscdiff_cache(i,ntwx_cache)=uscdiff(i)
531 C print *,'przed returnbox'
533 C call enerprint(remd_ene(0,i))
536 c_cache(j,i,ntwx_cache)=c(j,i)
541 if ((rstcount.eq.1000.or.itime.eq.n_timestep)
542 & .and..not.restart1file) then
545 open(irest1,file=mremd_rst_name,status='unknown')
546 write (irest1,*) "i2rep"
547 write (irest1,*) (i2rep(i),i=0,nodes-1)
548 write (irest1,*) "ifirst"
549 write (irest1,*) (ifirst(i),i=1,remd_m(1))
551 write (irest1,*) "nupa",il
552 write (irest1,*) nupa(0,il),(nupa(i,il),i=1,nupa(0,il))
553 write (irest1,*) "ndowna",il
554 write (irest1,*) ndowna(0,il),
555 & (ndowna(i,il),i=1,ndowna(0,il))
558 write (irest1,*) "nset"
559 write (irest1,*) nset
560 write (irest1,*) "mset"
561 write (irest1,*) (mset(i),i=1,nset)
562 write (irest1,*) "i2set"
563 write (irest1,*) (i2set(i),i=0,nodes-1)
564 write (irest1,*) "i_index"
568 write(irest1,*) (i_index(i,j,il,il1),j=1,remd_m(i))
576 open(irest2,file=rest2name,status='unknown')
577 write(irest2,*) totT,EK,potE,totE,t_bath
579 write (irest2,'(3e15.5)') (d_t(j,i),j=1,3)
582 write (irest2,'(3e15.5)') (dc(j,i),j=1,3)
585 write (irest2,*) iset
592 c forced synchronization
593 if (mod(itime,i_sync_step).eq.0 .and. me.ne.king
594 & .and. .not. mremdsync) then
596 call mpi_iprobe(0,101,CG_COMM,synflag,status,ierr)
598 call mpi_recv(itime_master, 1, MPI_INTEGER,
599 & 0,101,CG_COMM, status, ierr)
600 call mpi_barrier(CG_COMM, ierr)
601 cdeb if (out1file.or.traj1file) then
602 cdeb call mpi_gather(itime,1,mpi_integer,
603 cdeb & icache_all,1,mpi_integer,king,
606 & call mpi_gather(ntwx_cache,1,mpi_integer,
607 & icache_all,1,mpi_integer,king,
610 & write(iout,*) 'REMD synchro at',itime_master,itime
611 if (itime_master.ge.n_timestep .or. ovrtim())
613 ctime call flush(iout)
618 if ((mod(itime,nstex).eq.0.and.me.eq.king
619 & .or.end_of_run.and.me.eq.king )
620 & .and. .not. mremdsync ) then
623 call mpi_isend(itime,1,MPI_INTEGER,i,101,
624 & CG_COMM, ireqi(i), ierr)
625 cd write(iout,*) 'REMD synchro with',i
628 call mpi_waitall(nodes-1,ireqi,statusi,ierr)
629 call mpi_barrier(CG_COMM, ierr)
631 write(iout,*) 'REMD synchro at',itime,'time=',time01-time00
632 if (out1file.or.traj1file) then
633 cdeb call mpi_gather(itime,1,mpi_integer,
634 cdeb & itime_all,1,mpi_integer,king,
636 cdeb write(iout,'(a19,8000i8)') ' REMD synchro itime',
637 cdeb & (itime_all(i),i=1,nodes)
639 cdeb imin_itime=itime_all(1)
641 cdeb if(itime_all(i).lt.imin_itime) imin_itime=itime_all(i)
643 cdeb ii_write=(imin_itime-imin_itime_old)/ntwx
644 cdeb imin_itime_old=int(imin_itime/ntwx)*ntwx
645 cdeb write(iout,*) imin_itime,imin_itime_old,ii_write
646 call mpi_gather(ntwx_cache,1,mpi_integer,
647 & icache_all,1,mpi_integer,king,
649 c write(iout,'(a19,8000i8)') ' ntwx_cache',
650 c & (icache_all(i),i=1,nodes)
651 ii_write=icache_all(1)
653 if(icache_all(i).lt.ii_write) ii_write=icache_all(i)
655 c write(iout,*) "MIN ii_write=",ii_write
658 ctime call flush(iout)
660 if(mremdsync .and. mod(itime,nstex).eq.0) then
662 if (me.eq.king .or. .not. out1file)
663 & write(iout,*) 'REMD synchro at',itime
666 call mpi_gather(ntwx_cache,1,mpi_integer,
667 & icache_all,1,mpi_integer,king,
670 write(iout,'(a19,8000i8)') ' ntwx_cache',
671 & (icache_all(i),i=1,nodes)
672 ii_write=icache_all(1)
674 if(icache_all(i).lt.ii_write) ii_write=icache_all(i)
676 write(iout,*) "MIN ii_write=",ii_write
682 c Update the time safety limiy
683 if (time001-time00.gt.safety) then
684 safety=time001-time00+600
685 write (iout,*) "****** SAFETY increased to",safety," s"
687 if (ovrtim()) end_of_run=.true.
689 if(synflag.and..not.end_of_run) then
693 write(iout,*) 'REMD before',me,t_bath
695 c call mpi_gather(t_bath,1,mpi_double_precision,
696 c & remd_t_bath,1,mpi_double_precision,king,
698 potEcomp(n_ene+1)=t_bath
700 potEcomp(n_ene+2)=iset
701 if (iset.lt.nset) then
705 potEcomp(n_ene+3)=Uconst
712 potEcomp(n_ene+4)=Uconst
716 call mpi_gather(potEcomp(0),n_ene+5,mpi_double_precision,
717 & remd_ene(0,1),n_ene+5,mpi_double_precision,king,
720 call mpi_gather(elow,1,mpi_double_precision,
721 & elowi,1,mpi_double_precision,king,
723 call mpi_gather(ehigh,1,mpi_double_precision,
724 & ehighi,1,mpi_double_precision,king,
729 if (me.eq.king .or. .not. out1file) then
730 write(iout,*) 'REMD gather times=',time03-time01
734 if (restart1file) call write1rst(i_index)
737 if (me.eq.king .or. .not. out1file) then
738 write(iout,*) 'REMD writing rst time=',time04-time03
741 if (traj1file) call write1traj
743 cdeb call mpi_gather(ntwx_cache,1,mpi_integer,
744 cdeb & icache_all,1,mpi_integer,king,
746 cdeb write(iout,'(a19,8000i8)') ' ntwx_cache after traj1file',
747 cdeb & (icache_all(i),i=1,nodes)
752 if (me.eq.king .or. .not. out1file) then
753 write(iout,*) 'REMD writing traj time=',time05-time04
760 remd_t_bath(i)=remd_ene(n_ene+1,i)
761 iremd_iset(i)=remd_ene(n_ene+2,i)
764 co write(iout,*) 'REMD exchange temp,ene,elow,ehigh'
766 write(iout,'(i4,4f12.5)') i,remd_t_bath(i),remd_ene(0,i),
770 write(iout,*) 'REMD exchange temp,ene'
772 write(iout,'(i4,2f12.5)') i,remd_t_bath(i),remd_ene(0,i)
773 write(iout,'(6f12.5)') (remd_ene(j,i),j=1,n_ene)
776 c-------------------------------------
778 write (iout,*) "Enter exchnge, remd_m",remd_m(1),
781 write (iout,*) "remd_m(1)",remd_m(1)
783 i=ifirst(iran_num(1,remd_m(1)))
789 write (iout,*) "i",i," nupa(0,i)",int(nupa(0,i))
790 if(i.gt.0.and.nupa(0,i).gt.0) then
792 c if (i.eq.1 .and. int(nupa(0,i)).eq.1) then
794 c & "CHUJ ABSOLUTNY!!! No way to sample a distinct replica in MREMD"
796 c call MPI_Abort(MPI_COMM_WORLD,ERRCODE,ierr)
798 c do while (iex.eq.i)
799 c write (iout,*) "upper",nupa(int(nupa(0,i)),i)
800 iex=nupa(iran_num(1,int(nupa(0,i))),i)
802 c write (iout,*) "nupa(0,i)",nupa(0,i)," iex",iex
804 call muca_delta(remd_t_bath,remd_ene,i,iex,delta)
806 c Swap temperatures between conformations i and iex with recalculating the free energies
807 c following temperature changes.
808 ene_iex_iex=remd_ene(0,iex)
809 ene_i_i=remd_ene(0,i)
810 c write (iout,*) "i",i," ene_i_i",ene_i_i,
811 c & " iex",iex," ene_iex_iex",ene_iex_iex
812 c write (iout,*) "rescaling weights with temperature",
815 call rescale_weights(remd_t_bath(i))
817 c write (iout,*) "0,iex",remd_t_bath(i)
818 c call enerprint(remd_ene(0,iex))
820 call sum_energy(remd_ene(0,iex),.false.)
821 ene_iex_i=remd_ene(0,iex)
822 c write (iout,*) "ene_iex_i",remd_ene(0,iex)
824 c write (iout,*) "0,i",remd_t_bath(i)
825 c call enerprint(remd_ene(0,i))
827 call sum_energy(remd_ene(0,i),.false.)
828 c write (iout,*) "ene_i_i",remd_ene(0,i)
830 c write (iout,*) "rescaling weights with temperature",
832 if (real(ene_i_i).ne.real(remd_ene(0,i))) then
833 write (iout,*) "ERROR: inconsistent energies:",i,
834 & ene_i_i,remd_ene(0,i)
836 call rescale_weights(remd_t_bath(iex))
838 c write (iout,*) "0,i",remd_t_bath(iex)
839 call enerprint(remd_ene(0,i))
841 call sum_energy(remd_ene(0,i),.false.)
842 c write (iout,*) "ene_i_iex",remd_ene(0,i)
844 ene_i_iex=remd_ene(0,i)
846 c write (iout,*) "0,iex",remd_t_bath(iex)
847 c call enerprint(remd_ene(0,iex))
849 call sum_energy(remd_ene(0,iex),.false.)
850 if (real(ene_iex_iex).ne.real(remd_ene(0,iex))) then
851 write (iout,*) "ERROR: inconsistent energies:",iex,
852 & ene_iex_iex,remd_ene(0,iex)
854 c write (iout,*) "ene_iex_iex",remd_ene(0,iex)
855 c write (iout,*) "i",i," iex",iex
856 c write (iout,'(4(a,e15.5))') "ene_i_i",ene_i_i,
857 c & " ene_i_iex",ene_i_iex,
858 c & " ene_iex_i",ene_iex_i," ene_iex_iex",ene_iex_iex
860 delta=(ene_iex_iex-ene_i_iex)/(Rb*remd_t_bath(iex))-
861 & (ene_iex_i-ene_i_i)/(Rb*remd_t_bath(i))
863 c write(iout,*) 'delta',delta
864 c delta=(remd_t_bath(i)-remd_t_bath(iex))*
865 c & (remd_ene(i)-remd_ene(iex))/Rb/
866 c & (remd_t_bath(i)*remd_t_bath(iex))
868 if (delta .gt. 50.0d0) then
874 else if (delta.lt.-50.0d0) then
883 iremd_tot(int(i2rep(i-1)))=iremd_tot(int(i2rep(i-1)))+1
884 xxx=ran_number(0.0d0,1.0d0)
885 c write(iout,'(2i4,a6,2f12.5)') i,iex,' delta',delta,xxx
887 if (delta .gt. xxx) then
889 remd_t_bath(i)=remd_t_bath(iex)
891 remd_ene(0,i)=ene_i_iex
892 remd_ene(0,iex)=ene_iex_i
898 ehighi(i)=ehighi(iex)
905 nupa(k,i)=nupa(k,iex)
908 ndowna(k,i)=ndowna(k,iex)
912 if (ifirst(il).eq.i) ifirst(il)=iex
914 if (nupa(k,il).eq.i) then
916 elseif (nupa(k,il).eq.iex) then
921 if (ndowna(k,il).eq.i) then
923 elseif (ndowna(k,il).eq.iex) then
929 iremd_acc(int(i2rep(i-1)))=iremd_acc(int(i2rep(i-1)))+1
931 i2rep(i-1)=i2rep(iex-1)
934 c write(iout,*) 'exchange',i,iex
935 c write (iout,'(a8,100i4)') "@ ifirst",
936 c & (ifirst(k),k=1,remd_m(1))
938 c write (iout,'(a8,i4,a1,100i4)') "@ nupa",il,":",
939 c & (nupa(k,il),k=1,nupa(0,il))
940 c write (iout,'(a8,i4,a1,100i4)') "@ ndowna",il,":",
941 c & (ndowna(k,il),k=1,ndowna(0,il))
946 remd_ene(0,iex)=ene_iex_iex
947 remd_ene(0,i)=ene_i_i
953 cd write (iout,*) "exchange completed"
957 cd write(iout,*) "########",ii
959 i_temp=iran_num(1,nrep)
960 i_mult=iran_num(1,remd_m(i_temp))
961 i_iset=iran_num(1,nset)
962 i_mset=iran_num(1,mset(i_iset))
963 i=i_index(i_temp,i_mult,i_iset,i_mset)
965 cd write(iout,*) "i=",i,i_temp,i_mult,i_iset,i_mset
968 cd write(iout,*) "i_dir=",i_dir
970 if(i_dir.eq.1 .and. remd_m(i_temp+1).gt.0 )then
973 i_mult1=iran_num(1,remd_m(i_temp1))
975 i_mset1=iran_num(1,mset(i_iset1))
976 iex=i_index(i_temp1,i_mult1,i_iset1,i_mset1)
978 elseif(i_dir.eq.2 .and. mset(i_iset+1).gt.0)then
981 i_mult1=iran_num(1,remd_m(i_temp1))
983 i_mset1=iran_num(1,mset(i_iset1))
984 iex=i_index(i_temp1,i_mult1,i_iset1,i_mset1)
985 econstr_temp_i=remd_ene(20,i)
986 econstr_temp_iex=remd_ene(20,iex)
987 remd_ene(20,i)=remd_ene(n_ene+3,i)
988 remd_ene(20,iex)=remd_ene(n_ene+4,iex)
990 elseif(remd_m(i_temp+1).gt.0.and.mset(i_iset+1).gt.0)then
993 i_mult1=iran_num(1,remd_m(i_temp1))
995 i_mset1=iran_num(1,mset(i_iset1))
996 iex=i_index(i_temp1,i_mult1,i_iset1,i_mset1)
997 econstr_temp_i=remd_ene(20,i)
998 econstr_temp_iex=remd_ene(20,iex)
999 remd_ene(20,i)=remd_ene(n_ene+3,i)
1000 remd_ene(20,iex)=remd_ene(n_ene+4,iex)
1006 cd write(iout,*) "iex=",iex,i_temp1,i_mult1,i_iset1,i_mset1
1009 c Swap temperatures between conformations i and iex with recalculating the free energies
1010 c following temperature changes.
1011 ene_iex_iex=remd_ene(0,iex)
1012 ene_i_i=remd_ene(0,i)
1013 co write (iout,*) "rescaling weights with temperature",
1015 call rescale_weights(remd_t_bath(i))
1017 call sum_energy(remd_ene(0,iex),.false.)
1018 ene_iex_i=remd_ene(0,iex)
1019 cd write (iout,*) "ene_iex_i",remd_ene(0,iex)
1020 c call sum_energy(remd_ene(0,i),.false.)
1021 cd write (iout,*) "ene_i_i",remd_ene(0,i)
1022 c write (iout,*) "rescaling weights with temperature",
1023 c & remd_t_bath(iex)
1024 c if (real(ene_i_i).ne.real(remd_ene(0,i))) then
1025 c write (iout,*) "ERROR: inconsistent energies:",i,
1026 c & ene_i_i,remd_ene(0,i)
1028 call rescale_weights(remd_t_bath(iex))
1029 call sum_energy(remd_ene(0,i),.false.)
1030 cd write (iout,*) "ene_i_iex",remd_ene(0,i)
1031 ene_i_iex=remd_ene(0,i)
1032 c call sum_energy(remd_ene(0,iex),.false.)
1033 c if (real(ene_iex_iex).ne.real(remd_ene(0,iex))) then
1034 c write (iout,*) "ERROR: inconsistent energies:",iex,
1035 c & ene_iex_iex,remd_ene(0,iex)
1037 cd write (iout,*) "ene_iex_iex",remd_ene(0,iex)
1038 c write (iout,*) "i",i," iex",iex
1039 cd write (iout,'(4(a,e15.5))') "ene_i_i",ene_i_i,
1040 cd & " ene_i_iex",ene_i_iex,
1041 cd & " ene_iex_i",ene_iex_i," ene_iex_iex",ene_iex_iex
1042 delta=(ene_iex_iex-ene_i_iex)/(Rb*remd_t_bath(iex))-
1043 & (ene_iex_i-ene_i_i)/(Rb*remd_t_bath(i))
1045 cd write(iout,*) 'delta',delta
1046 c delta=(remd_t_bath(i)-remd_t_bath(iex))*
1047 c & (remd_ene(i)-remd_ene(iex))/Rb/
1048 c & (remd_t_bath(i)*remd_t_bath(iex))
1049 if (delta .gt. 50.0d0) then
1054 if (i_dir.eq.1.or.i_dir.eq.3)
1055 & iremd_tot(int(i2rep(i-1)))=iremd_tot(int(i2rep(i-1)))+1
1056 if (i_dir.eq.2.or.i_dir.eq.3)
1057 & iremd_tot_usa(int(i2set(i-1)))=
1058 & iremd_tot_usa(int(i2set(i-1)))+1
1059 xxx=ran_number(0.0d0,1.0d0)
1060 cd write(iout,'(2i4,a6,2f12.5)') i,iex,' delta',delta,xxx
1061 if (delta .gt. xxx) then
1063 remd_t_bath(i)=remd_t_bath(iex)
1064 remd_t_bath(iex)=tmp
1067 iremd_iset(i)=iremd_iset(iex)
1068 iremd_iset(iex)=itmp
1070 remd_ene(0,i)=ene_i_iex
1071 remd_ene(0,iex)=ene_iex_i
1073 if (i_dir.eq.1.or.i_dir.eq.3)
1074 & iremd_acc(int(i2rep(i-1)))=iremd_acc(int(i2rep(i-1)))+1
1077 i2rep(i-1)=i2rep(iex-1)
1080 if (i_dir.eq.2.or.i_dir.eq.3)
1081 & iremd_acc_usa(int(i2set(i-1)))=
1082 & iremd_acc_usa(int(i2set(i-1)))+1
1085 i2set(i-1)=i2set(iex-1)
1088 itmp=i_index(i_temp,i_mult,i_iset,i_mset)
1089 i_index(i_temp,i_mult,i_iset,i_mset)=
1090 & i_index(i_temp1,i_mult1,i_iset1,i_mset1)
1091 i_index(i_temp1,i_mult1,i_iset1,i_mset1)=itmp
1094 remd_ene(0,iex)=ene_iex_iex
1095 remd_ene(0,i)=ene_i_i
1096 remd_ene(20,iex)=econstr_temp_iex
1097 remd_ene(20,i)=econstr_temp_i
1101 cd do il1=1,mset(il)
1104 cd write(iout,*) i,j,il,il1,i_index(i,j,il,il1)
1117 c-------------------------------------
1118 write (iout,*) "NREP",nrep
1120 if(iremd_tot(i).ne.0)
1121 & write(iout,'(a3,i4,2f12.5,i5)') 'ACC',i,remd_t(i)
1122 & ,iremd_acc(i)/(1.0*iremd_tot(i)),iremd_tot(i)
1127 if(iremd_tot_usa(i).ne.0)
1128 & write(iout,'(a10,i4,f12.5,i8)') 'ACC_usampl',i,
1129 & iremd_acc_usa(i)/(1.0*iremd_tot_usa(i)),iremd_tot_usa(i)
1135 cd write (iout,'(a6,100i4)') "ifirst",
1136 cd & (ifirst(i),i=1,remd_m(1))
1138 cd write (iout,'(a5,i4,a1,100i4)') "nup",il,":",
1139 cd & (nupa(i,il),i=1,nupa(0,il))
1140 cd write (iout,'(a5,i4,a1,100i4)') "ndown",il,":",
1141 cd & (ndowna(i,il),i=1,ndowna(0,il))
1146 cd write (iout,*) "Before scatter"
1148 call mpi_scatter(remd_t_bath,1,mpi_double_precision,
1149 & t_bath,1,mpi_double_precision,king,
1151 cd write (iout,*) "After scatter"
1154 & call mpi_scatter(iremd_iset,1,mpi_integer,
1155 & iset,1,mpi_integer,king,
1159 if (me.eq.king .or. .not. out1file) then
1160 write(iout,*) 'REMD scatter time=',time07-time06
1164 call mpi_scatter(elowi,1,mpi_double_precision,
1165 & elow,1,mpi_double_precision,king,
1167 call mpi_scatter(ehighi,1,mpi_double_precision,
1168 & ehigh,1,mpi_double_precision,king,
1171 call rescale_weights(t_bath)
1172 co write (iout,*) "Processor",me,
1173 co & " rescaling weights with temperature",t_bath
1175 stdfp=dsqrt(2*Rb*t_bath/d_time)
1177 stdfsc(i)=dsqrt(2*Rb*t_bath/d_time)
1180 cde write(iout,*) 'REMD after',me,t_bath
1182 if (me.eq.king .or. .not. out1file) then
1183 write(iout,*) 'REMD exchange time=',time08-time00
1189 if (restart1file) then
1190 if (me.eq.king .or. .not. out1file)
1191 & write(iout,*) 'writing restart at the end of run'
1192 call write1rst(i_index)
1195 if (traj1file) call write1traj
1197 cdeb call mpi_gather(ntwx_cache,1,mpi_integer,
1198 cdeb & icache_all,1,mpi_integer,king,
1199 cdeb & CG_COMM,ierr)
1200 cdeb write(iout,'(a40,8000i8)')
1201 cdeb & ' ntwx_cache after traj1file at the end',
1202 cdeb & (icache_all(i),i=1,nodes)
1207 t_MD=MPI_Wtime()-tt0
1211 if (me.eq.king .or. .not. out1file) then
1212 write (iout,'(//35(1h=),a10,35(1h=)/10(/a40,1pe15.5))')
1214 & 'MD calculations setup:',t_MDsetup,
1215 & 'Energy & gradient evaluation:',t_enegrad,
1216 & 'Stochastic MD setup:',t_langsetup,
1217 & 'Stochastic MD step setup:',t_sdsetup,
1219 write (iout,'(/28(1h=),a25,27(1h=))')
1220 & ' End of MD calculation '
1225 c-----------------------------------------------------------------------
1226 subroutine write1rst(i_index)
1227 implicit real*8 (a-h,o-z)
1228 include 'DIMENSIONS'
1231 include 'COMMON.IOUNITS'
1232 include 'COMMON.REMD'
1233 include 'COMMON.SETUP'
1234 include 'COMMON.CHAIN'
1235 include 'COMMON.SBRIDGE'
1236 include 'COMMON.INTERACT'
1238 real d_restart1(3,0:2*(maxres)*maxprocs),r_d(3,0:2*maxres),
1239 & d_restart2(3,0:2*(maxres)*maxprocs)
1243 & (maxprocs/4,maxprocs/20,maxprocs/200,maxprocs/200)
1244 common /przechowalnia/ d_restart1,d_restart2
1249 t5_restart1(4)=t_bath
1250 t5_restart1(5)=Uconst
1252 call mpi_gather(t5_restart1,5,mpi_real,
1253 & t_restart1,5,mpi_real,king,CG_COMM,ierr)
1261 call mpi_gather(r_d,3*2*(nres+1),mpi_real,
1262 & d_restart1,3*2*(nres+1),mpi_real,king,
1265 C write (*,*) dc(j,0),"TU3"
1274 call mpi_gather(r_d,3*2*(nres+1),mpi_real,
1275 & d_restart2,3*2*(nres+1),mpi_real,king,
1280 call xdrfopen_(ixdrf,mremd_rst_name, "w", iret)
1282 call xdrfint_(ixdrf, i2rep(i), iret)
1285 call xdrfint_(ixdrf, ifirst(i), iret)
1289 call xdrfint_(ixdrf, nupa(i,il), iret)
1293 call xdrfint_(ixdrf, ndowna(i,il), iret)
1299 call xdrffloat_(ixdrf, t_restart1(j,il), iret)
1307 & xdrffloat_(ixdrf, d_restart1(j,i+2*(nres+1)*il), iret)
1315 & xdrffloat_(ixdrf, d_restart2(j,i+2*(nres+1)*il), iret)
1321 call xdrfint_(ixdrf, nset, iret)
1323 call xdrfint_(ixdrf,mset(i), iret)
1326 call xdrfint_(ixdrf,i2set(i), iret)
1332 itmp=i_index(i,j,il,il1)
1333 call xdrfint_(ixdrf,itmp, iret)
1340 call xdrfclose_(ixdrf, iret)
1342 call xdrfopen(ixdrf,mremd_rst_name, "w", iret)
1344 call xdrfint(ixdrf, i2rep(i), iret)
1347 call xdrfint(ixdrf, ifirst(i), iret)
1351 call xdrfint(ixdrf, nupa(i,il), iret)
1355 call xdrfint(ixdrf, ndowna(i,il), iret)
1361 call xdrffloat(ixdrf, t_restart1(j,il), iret)
1369 & xdrffloat(ixdrf, d_restart1(j,i+2*(nres+1)*il), iret)
1377 & xdrffloat(ixdrf, d_restart2(j,i+2*(nres+1)*il), iret)
1384 call xdrfint(ixdrf, nset, iret)
1386 call xdrfint(ixdrf,mset(i), iret)
1389 call xdrfint(ixdrf,i2set(i), iret)
1395 itmp=i_index(i,j,il,il1)
1396 call xdrfint(ixdrf,itmp, iret)
1403 call xdrfclose(ixdrf, iret)
1410 subroutine write1traj
1411 implicit real*8 (a-h,o-z)
1412 include 'DIMENSIONS'
1415 include 'COMMON.IOUNITS'
1416 include 'COMMON.REMD'
1417 include 'COMMON.SETUP'
1418 include 'COMMON.CHAIN'
1419 include 'COMMON.SBRIDGE'
1420 include 'COMMON.INTERACT'
1424 real xcoord(3,maxres2+2),prec
1425 real r_qfrag(50),r_qpair(100)
1426 real r_utheta(50),r_ugamma(100),r_uscdiff(100)
1427 real p_qfrag(50*maxprocs),p_qpair(100*maxprocs)
1428 real p_utheta(50*maxprocs),p_ugamma(100*maxprocs),
1429 & p_uscdiff(100*maxprocs)
1430 real p_c(3,(maxres2+2)*maxprocs),r_c(3,maxres2+2)
1431 common /przechowalnia/ p_c
1433 call mpi_bcast(ii_write,1,mpi_integer,
1434 & king,CG_COMM,ierr)
1437 print *,'traj1file',me,ii_write,ntwx_cache
1441 if(me.eq.king) call xdrfopen_(ixdrf,cartname, "a", iret)
1443 if(me.eq.king) call xdrfopen(ixdrf,cartname, "a", iret)
1446 t5_restart1(1)=totT_cache(ii)
1447 t5_restart1(2)=EK_cache(ii)
1448 t5_restart1(3)=potE_cache(ii)
1449 t5_restart1(4)=t_bath_cache(ii)
1450 t5_restart1(5)=Uconst_cache(ii)
1451 call mpi_gather(t5_restart1,5,mpi_real,
1452 & t_restart1,5,mpi_real,king,CG_COMM,ierr)
1454 call mpi_gather(iset_cache(ii),1,mpi_integer,
1455 & iset_restart1,1,mpi_integer,king,CG_COMM,ierr)
1458 r_qfrag(i)=qfrag_cache(i,ii)
1461 r_qpair(i)=qpair_cache(i,ii)
1464 r_utheta(i)=utheta_cache(i,ii)
1465 r_ugamma(i)=ugamma_cache(i,ii)
1466 r_uscdiff(i)=uscdiff_cache(i,ii)
1469 call mpi_gather(r_qfrag,nfrag,mpi_real,
1470 & p_qfrag,nfrag,mpi_real,king,
1472 call mpi_gather(r_qpair,npair,mpi_real,
1473 & p_qpair,npair,mpi_real,king,
1475 call mpi_gather(r_utheta,nfrag_back,mpi_real,
1476 & p_utheta,nfrag_back,mpi_real,king,
1478 call mpi_gather(r_ugamma,nfrag_back,mpi_real,
1479 & p_ugamma,nfrag_back,mpi_real,king,
1481 call mpi_gather(r_uscdiff,nfrag_back,mpi_real,
1482 & p_uscdiff,nfrag_back,mpi_real,king,
1486 write (iout,*) "p_qfrag"
1488 write (iout,*) i,(p_qfrag((i-1)*nfrag+j),j=1,nfrag)
1490 write (iout,*) "p_qpair"
1492 write (iout,*) i,(p_qpair((i-1)*npair+j),j=1,npair)
1498 r_c(j,i)=c_cache(j,i,ii)
1502 call mpi_gather(r_c,3*2*nres,mpi_real,
1503 & p_c,3*2*nres,mpi_real,king,
1509 call xdrffloat_(ixdrf, real(t_restart1(1,il)), iret)
1510 call xdrffloat_(ixdrf, real(t_restart1(3,il)), iret)
1511 call xdrffloat_(ixdrf, real(t_restart1(5,il)), iret)
1512 call xdrffloat_(ixdrf, real(t_restart1(4,il)), iret)
1513 call xdrfint_(ixdrf, nss, iret)
1516 call xdrfint(ixdrf, idssb(j)+nres, iret)
1517 call xdrfint(ixdrf, jdssb(j)+nres, iret)
1519 call xdrfint_(ixdrf, ihpb(j), iret)
1520 call xdrfint_(ixdrf, jhpb(j), iret)
1523 call xdrfint_(ixdrf, nfrag+npair+3*nfrag_back, iret)
1524 call xdrfint_(ixdrf, iset_restart1(il), iret)
1526 call xdrffloat_(ixdrf, p_qfrag(i+(il-1)*nfrag), iret)
1529 call xdrffloat_(ixdrf, p_qpair(i+(il-1)*npair), iret)
1532 call xdrffloat_(ixdrf, p_utheta(i+(il-1)*nfrag_back), iret)
1533 call xdrffloat_(ixdrf, p_ugamma(i+(il-1)*nfrag_back), iret)
1534 call xdrffloat_(ixdrf, p_uscdiff(i+(il-1)*nfrag_back), iret)
1539 xcoord(j,i)=p_c(j,i+(il-1)*nres*2)
1544 xcoord(j,nres+i-nnt+1)=p_c(j,i+nres+(il-1)*nres*2)
1548 call xdrf3dfcoord_(ixdrf, xcoord, itmp, prec, iret)
1552 call xdrffloat(ixdrf, real(t_restart1(1,il)), iret)
1553 call xdrffloat(ixdrf, real(t_restart1(3,il)), iret)
1554 call xdrffloat(ixdrf, real(t_restart1(5,il)), iret)
1555 call xdrffloat(ixdrf, real(t_restart1(4,il)), iret)
1556 call xdrfint(ixdrf, nss, iret)
1559 call xdrfint(ixdrf, idssb(j)+nres, iret)
1560 call xdrfint(ixdrf, jdssb(j)+nres, iret)
1562 call xdrfint(ixdrf, ihpb(j), iret)
1563 call xdrfint(ixdrf, jhpb(j), iret)
1566 call xdrfint(ixdrf, nfrag+npair+3*nfrag_back, iret)
1567 call xdrfint(ixdrf, iset_restart1(il), iret)
1569 call xdrffloat(ixdrf, p_qfrag(i+(il-1)*nfrag), iret)
1572 call xdrffloat(ixdrf, p_qpair(i+(il-1)*npair), iret)
1575 call xdrffloat(ixdrf, p_utheta(i+(il-1)*nfrag_back), iret)
1576 call xdrffloat(ixdrf, p_ugamma(i+(il-1)*nfrag_back), iret)
1577 call xdrffloat(ixdrf, p_uscdiff(i+(il-1)*nfrag_back), iret)
1582 xcoord(j,i)=p_c(j,i+(il-1)*nres*2)
1587 xcoord(j,nres+i-nnt+1)=p_c(j,i+nres+(il-1)*nres*2)
1591 call xdrf3dfcoord(ixdrf, xcoord, itmp, prec, iret)
1597 if(me.eq.king) call xdrfclose_(ixdrf, iret)
1599 if(me.eq.king) call xdrfclose(ixdrf, iret)
1601 do i=1,ntwx_cache-ii_write
1603 totT_cache(i)=totT_cache(ii_write+i)
1604 EK_cache(i)=EK_cache(ii_write+i)
1605 potE_cache(i)=potE_cache(ii_write+i)
1606 t_bath_cache(i)=t_bath_cache(ii_write+i)
1607 Uconst_cache(i)=Uconst_cache(ii_write+i)
1608 iset_cache(i)=iset_cache(ii_write+i)
1611 qfrag_cache(ii,i)=qfrag_cache(ii,ii_write+i)
1614 qpair_cache(ii,i)=qpair_cache(ii,ii_write+i)
1617 utheta_cache(ii,i)=utheta_cache(ii,ii_write+i)
1618 ugamma_cache(ii,i)=ugamma_cache(ii,ii_write+i)
1619 uscdiff_cache(ii,i)=uscdiff_cache(ii,ii_write+i)
1624 c_cache(j,ii,i)=c_cache(j,ii,ii_write+i)
1628 ntwx_cache=ntwx_cache-ii_write
1633 subroutine read1restart(i_index)
1634 implicit real*8 (a-h,o-z)
1635 include 'DIMENSIONS'
1638 include 'COMMON.IOUNITS'
1639 include 'COMMON.REMD'
1640 include 'COMMON.SETUP'
1641 include 'COMMON.CHAIN'
1642 include 'COMMON.SBRIDGE'
1643 include 'COMMON.INTERACT'
1644 real d_restart1(3,0:2*(maxres)*maxprocs),r_d(3,0:2*(maxres)),
1647 & (maxprocs/4,maxprocs/20,maxprocs/200,maxprocs/200)
1648 common /przechowalnia/ d_restart1
1649 write (*,*) "Processor",me," called read1restart"
1652 open(irest2,file=mremd_rst_name,status='unknown')
1653 read(irest2,*,err=334) i
1654 write(iout,*) "Reading old rst in ASCI format"
1656 call read1restart_old
1660 call xdrfopen_(ixdrf,mremd_rst_name, "r", iret)
1663 call xdrfint_(ixdrf, i2rep(i), iret)
1666 call xdrfint_(ixdrf, ifirst(i), iret)
1669 call xdrfint_(ixdrf, nupa(0,il), iret)
1671 call xdrfint_(ixdrf, nupa(i,il), iret)
1674 call xdrfint_(ixdrf, ndowna(0,il), iret)
1676 call xdrfint_(ixdrf, ndowna(i,il), iret)
1681 call xdrffloat_(ixdrf, t_restart1(j,il), iret)
1685 call xdrfopen(ixdrf,mremd_rst_name, "r", iret)
1688 call xdrfint(ixdrf, i2rep(i), iret)
1691 call xdrfint(ixdrf, ifirst(i), iret)
1694 call xdrfint(ixdrf, nupa(0,il), iret)
1696 call xdrfint(ixdrf, nupa(i,il), iret)
1699 call xdrfint(ixdrf, ndowna(0,il), iret)
1701 call xdrfint(ixdrf, ndowna(i,il), iret)
1706 call xdrffloat(ixdrf, t_restart1(j,il), iret)
1711 call mpi_scatter(t_restart1,5,mpi_real,
1712 & t5_restart1,5,mpi_real,king,CG_COMM,ierr)
1716 t_bath=t5_restart1(4)
1721 c read(irest2,'(3e15.5)')
1722 c & (d_restart1(j,i+2*nres*il),j=1,3)
1725 call xdrffloat_(ixdrf, d_restart1(j,i+2*(nres+1)*il), iret)
1727 call xdrffloat(ixdrf, d_restart1(j,i+2*(nres+1)*il), iret)
1733 call mpi_scatter(d_restart1,3*2*(nres+1),mpi_real,
1734 & r_d,3*2*(nres+1),mpi_real,king,CG_COMM,ierr)
1744 c read(irest2,'(3e15.5)')
1745 c & (d_restart1(j,i+2*nres*il),j=1,3)
1749 & xdrffloat_(ixdrf, d_restart1(j,i+2*(nres+1)*il), iret)
1752 & xdrffloat(ixdrf, d_restart1(j,i+2*(nres+1)*il), iret)
1758 call mpi_scatter(d_restart1,3*2*(nres+1),mpi_real,
1759 & r_d,3*2*(nres+1),mpi_real,king,CG_COMM,ierr)
1770 call xdrfint_(ixdrf, nset, iret)
1772 call xdrfint_(ixdrf,mset(i), iret)
1775 call xdrfint_(ixdrf,i2set(i), iret)
1781 call xdrfint_(ixdrf,itmp, iret)
1782 i_index(i,j,il,il1)=itmp
1790 call xdrfint(ixdrf, nset, iret)
1792 call xdrfint(ixdrf,mset(i), iret)
1795 call xdrfint(ixdrf,i2set(i), iret)
1801 call xdrfint(ixdrf,itmp, iret)
1802 i_index(i,j,il,il1)=itmp
1809 call mpi_scatter(i2set,1,mpi_integer,
1810 & iset,1,mpi_integer,king,
1816 if(me.eq.king) close(irest2)
1820 subroutine read1restart_old
1821 implicit real*8 (a-h,o-z)
1822 include 'DIMENSIONS'
1825 include 'COMMON.IOUNITS'
1826 include 'COMMON.REMD'
1827 include 'COMMON.SETUP'
1828 include 'COMMON.CHAIN'
1829 include 'COMMON.SBRIDGE'
1830 include 'COMMON.INTERACT'
1831 real d_restart1(3,0:2*maxres*maxprocs),r_d(3,0:2*maxres),
1833 common /przechowalnia/ d_restart1
1835 open(irest2,file=mremd_rst_name,status='unknown')
1836 read (irest2,*) (i2rep(i),i=0,nodes-1)
1837 read (irest2,*) (ifirst(i),i=1,remd_m(1))
1839 read (irest2,*) nupa(0,il),(nupa(i,il),i=1,nupa(0,il))
1840 read (irest2,*) ndowna(0,il),
1841 & (ndowna(i,il),i=1,ndowna(0,il))
1844 read(irest2,*) (t_restart1(j,il),j=1,4)
1847 call mpi_scatter(t_restart1,5,mpi_real,
1848 & t5_restart1,5,mpi_real,king,CG_COMM,ierr)
1852 t_bath=t5_restart1(4)
1857 read(irest2,'(3e15.5)')
1858 & (d_restart1(j,i+2*(nres+1)*il),j=1,3)
1862 call mpi_scatter(d_restart1,3*2*(nres+1),mpi_real,
1863 & r_d,3*2*(nres+1),mpi_real,king,CG_COMM,ierr)
1873 read(irest2,'(3e15.5)')
1874 & (d_restart1(j,i+2*(nres+1)*il),j=1,3)
1878 call mpi_scatter(d_restart1,3*2*(nres+1),mpi_real,
1879 & r_d,3*2*(nres+1),mpi_real,king,CG_COMM,ierr)
1885 if(me.eq.king) close(irest2)
1888 c------------------------------------------