ctest wham & cluster src-M
[unres.git] / source / unres / src_CSA / together.F
1 #ifdef MPI
2       Subroutine together
3 c  feeds tasks for parallel processing
4       implicit real*8 (a-h,o-z)
5       include 'DIMENSIONS'
6       include 'mpif.h'
7       real ran1,ran2
8       include 'COMMON.CSA'
9       include 'COMMON.BANK'
10       include 'COMMON.IOUNITS'
11       include 'COMMON.CHAIN'
12       include 'COMMON.TIME1'
13       include 'COMMON.SETUP'
14       include 'COMMON.VAR'
15       include 'COMMON.GEO'
16       include 'COMMON.CONTROL'
17       include 'COMMON.SBRIDGE'
18       real tcpu
19       double precision time_start,time_start_c,time0f,time0i
20       logical ovrtim,sync_iter,timeout,flag,timeout1
21       dimension muster(mpi_status_size)
22       dimension t100(0:100),indx(mxio)
23       dimension xout(maxvar),eout(mxch*(mxch+1)/2+1),ind(9)
24       dimension cout(2)
25       parameter (rad=1.745329252d-2)
26
27 cccccccccccccccccccccccccccccccccccccccccccccccc
28       IF (ME.EQ.KING) THEN
29
30        time0f=MPI_WTIME()
31        ilastnstep=1
32        sync_iter=.false.
33        numch=1
34        nrmsdb=0
35        nrmsdb1=0
36        rmsdbc1c=rmsdbc1
37        nstep=0
38        call csa_read
39        call make_array
40
41        if(iref.ne.0) call from_int(1,0,idum)
42
43 c To minimize input conformation (bank conformation)
44 c Output to $mol.reminimized
45        if (irestart.lt.0) then
46         call read_bank(0,nft,cutdifr)
47         if (irestart.lt.-10) then
48          p_cut=nres*4.d0
49          call prune_bank(p_cut)
50          return
51         endif
52         call reminimize(jlee)
53         return
54        endif
55
56        if (irestart.eq.0) then
57         call initial_write
58         nbank=nconf
59         ntbank=nconf
60         if (ntbankm.eq.0) ntbank=0
61         nstep=0
62         nft=0
63         do i=1,mxio
64          ibank(i)=0
65          jbank(i)=0
66         enddo
67        else
68         call restart_write
69 c!bankt        call read_bankt(jlee,nft,cutdifr)
70         call read_bank(jlee,nft,cutdifr)
71         call read_rbank(jlee,adif)
72         if(iref.ne.0) call from_int(1,0,idum)
73        endif
74
75        nstmax=nstmax+nstep
76        ntrial=n1+n2+n3+n4+n5+n6+n7+n8
77        ntry=ntrial+1
78        ntry=ntry*nseed
79
80 c ntrial : number of trial conformations per seed.
81 c ntry : total number of trial conformations including seed conformations.
82
83        idum2=-123
84 #ifdef G77
85        imax=2**30-1
86 #else
87        imax=2**31-1
88 #endif
89        ENDIF
90
91        call mpi_bcast(jend,1,mpi_integer,0,CG_COMM,ierr)
92 cccccccccccccccccccccccccccccccccccccccc
93        do 300 jlee=1,jend
94 cccccccccccccccccccccccccccccccccccccccc
95   331   continue 
96         IF (ME.EQ.KING) THEN
97         if(sync_iter) goto 333
98         idum=-  ran2(idum2)*imax
99         if(jlee.lt.jstart) goto 300
100
101 C Restart the random number generator for conformation generation
102
103         if(irestart.gt.0) then
104          idum2=idum2+nstep
105          if(idum2.le.0) idum2=-idum2+1
106          idum=-  ran2(idum2)*imax
107         endif
108
109         idumm=idum
110         call vrndst(idumm)
111
112         open(icsa_seed,file=csa_seed,status="old")
113         write(icsa_seed,*) "jlee : ",jlee
114         close(icsa_seed)
115
116       call history_append
117       write(icsa_history,*) "number of procs is ",nodes
118       write(icsa_history,*) jlee,idum,idum2
119       close(icsa_history)
120
121 cccccccccccccccccccccccccccccccccccccccccccccccc
122   333 icycle=0
123
124        call history_append
125         write(icsa_history,*) "nbank is ",nbank
126        close(icsa_history)
127
128       if(irestart.eq.1) goto 111
129       if(irestart.eq.2) then
130        icycle=0
131        do i=1,nbank
132         ibank(i)=1
133        enddo
134        do i=nbank+1,nbank+nconf
135         ibank(i)=0
136        enddo
137       endif
138
139 c  start energy minimization
140       nconfr=max0(nconf+nadd,nodes-1)
141       if (sync_iter) nconf_in=0
142 c  king-emperor - feed input and sort output
143        write (iout,*) "NCONF_IN",nconf_in
144        m=0
145        if (nconf_in.gt.0) then
146 c al 7/2/00 - added possibility to read in some of the initial conformations
147         do m=1,nconf_in
148           read (intin,'(i5)',end=11,err=12) iconf
149    12     continue
150           write (iout,*) "write READ_ANGLES",iconf,m
151           call read_angles(intin,*11)
152           if (iref.eq.0) then
153             mm=m
154           else
155             mm=m+1
156           endif
157           do j=2,nres-1
158             dihang_in(1,j,1,mm)=theta(j+1)
159             dihang_in(2,j,1,mm)=phi(j+2)
160             dihang_in(3,j,1,mm)=alph(j)
161             dihang_in(4,j,1,mm)=omeg(j)
162           enddo
163         enddo ! m
164         goto 13
165    11   write (iout,*) nconf_in," conformations requested, but only",
166      &   m-1," found in the angle file."
167         nconf_in=m-1
168    13   continue
169         m=nconf_in
170         write (iout,*) nconf_in,
171      &    " initial conformations have been read in."
172        endif
173        if (iref.eq.0) then
174         if (nconfr.gt.nconf_in) then
175           call make_ranvar(nconfr,m,idum)
176           write (iout,*) nconfr-nconf_in,
177      &     " conformations have been generated randomly."
178         endif
179        else
180         nconfr=nconfr*2
181         call from_int(nconfr,m,idum)
182 c       call from_pdb(nconfr,idum)
183        endif
184        write (iout,*) 'Exitted from make_ranvar nconfr=',nconfr
185        write (*,*) 'Exitted from make_ranvar nconfr=',nconfr
186        do m=1,nconfr
187           write (iout,*) 'Initial conformation',m
188           write(iout,'(8f10.4)') (rad2deg*dihang_in(1,j,1,m),j=2,nres-1)
189           write(iout,'(8f10.4)') (rad2deg*dihang_in(2,j,1,m),j=2,nres-1)
190           write(iout,'(8f10.4)') (rad2deg*dihang_in(3,j,1,m),j=2,nres-1)
191           write(iout,'(8f10.4)') (rad2deg*dihang_in(4,j,1,m),j=2,nres-1)
192        enddo 
193        write(iout,*)'Calling FEEDIN NCONF',nconfr
194        time1i=MPI_WTIME()
195        call feedin(nconfr,nft)
196        write (iout,*) ' Time for first bank min.',MPI_WTIME()-time1i
197        call  history_append
198         write(icsa_history,*) jlee,nft,nbank
199         write(icsa_history,851) (etot(i),i=1,nconfr)
200         write(icsa_history,850) (rmsn(i),i=1,nconfr)
201         write(icsa_history,850) (pncn(i),i=1,nconfr)
202         write(icsa_history,*)
203        close(icsa_history)
204       ELSE
205 c To minimize input conformation (bank conformation)
206 c Output to $mol.reminimized   
207        if (irestart.lt.0) then 
208         call reminimize(jlee)
209         return
210        endif
211        if (irestart.eq.1) goto 111
212 c  soldier - perform energy minimization
213  334   call minim_jlee
214
215
216       ENDIF
217
218 ccccccccccccccccccccccccccccccccccc
219 c need to syncronize all procs
220       call mpi_barrier(CG_COMM,ierr)
221       if (ierr.ne.0) then
222        print *, ' cannot synchronize MPI'
223        stop
224       endif
225 ccccccccccccccccccccccccccccccccccc
226
227       IF (ME.EQ.KING) THEN
228
229 c      print *,"ok after minim"
230       nstep=nstep+nconf
231       if(irestart.eq.2) then
232        nbank=nbank+nconf
233 c      ntbank=ntbank+nconf
234        if(ntbank.gt.ntbankm) ntbank=ntbankm
235       endif
236 c      print *,"ok before indexx"
237       if(iref.eq.0) then
238        call indexx(nconfr,etot,indx)
239       else
240 c cc/al 7/6/00
241        do k=1,nconfr
242          indx(k)=k
243        enddo
244        call indexx(nconfr-nconf_in,rmsn(nconf_in+1),indx(nconf_in+1))
245        do k=nconf_in+1,nconfr
246          indx(k)=indx(k)+nconf_in
247        enddo
248 c cc/al
249 c       call indexx(nconfr,rmsn,indx)
250       endif
251 c      print *,"ok after indexx"
252       do im=1,nconf
253        m=indx(im)
254        if (m.gt.mxio .or. m.lt.1) then
255          write (iout,*) 'Dimension ERROR in TOGEHER: IM',im,' M',m
256          call mpi_abort(mpi_comm_world,ierror,ierrcode)
257        endif
258        jbank(im+nbank-nconf)=0
259        bene(im+nbank-nconf)=etot(m)
260        rene(im+nbank-nconf)=etot(m)
261 c!bankt       btene(im)=etot(m)
262 c
263        brmsn(im+nbank-nconf)=rmsn(m)
264        bpncn(im+nbank-nconf)=pncn(m)
265        rrmsn(im+nbank-nconf)=rmsn(m)
266        rpncn(im+nbank-nconf)=pncn(m)
267        if (im+nbank-nconf.gt.mxio .or. im+nbank-nconf.lt.1) then
268          write (iout,*) 'Dimension ERROR in TOGEHER: IM',im,
269      &   ' NBANK',nbank,' NCONF',nconf,' IM+NBANK-NCONF',im+nbank-nconf
270          call mpi_abort(mpi_comm_world,ierror,ierrcode)
271        endif
272        do k=1,numch
273         do j=2,nres-1
274          do i=1,4
275           bvar(i,j,k,im+nbank-nconf)=dihang(i,j,k,m)
276           rvar(i,j,k,im+nbank-nconf)=dihang(i,j,k,m)
277 c!bankt          btvar(i,j,k,im)=dihang(i,j,k,m)
278 c
279          enddo
280         enddo
281        enddo
282        if(iref.eq.1) then
283         if(brmsn(im+nbank-nconf).gt.rmscut.or.
284      &     bpncn(im+nbank-nconf).lt.pnccut) ibank(im+nbank-nconf)=9
285        endif
286        if(vdisulf) then
287            bvar_ns(im+nbank-nconf)=ns-2*nss
288            k=0
289            do i=1,ns
290              j=1
291              do while( iss(i).ne.ihpb(j)-nres .and. 
292      &                 iss(i).ne.jhpb(j)-nres .and. j.le.nss)
293               j=j+1 
294              enddo
295              if (j.gt.nss) then            
296                k=k+1
297                bvar_s(k,im+nbank-nconf)=iss(i)
298              endif
299            enddo
300        endif
301        bvar_nss(im+nbank-nconf)=nss
302        do i=1,nss
303            bvar_ss(1,i,im+nbank-nconf)=ihpb(i)
304            bvar_ss(2,i,im+nbank-nconf)=jhpb(i)
305        enddo
306       enddo
307       ENDIF
308
309   111 continue
310
311       IF (ME.EQ.KING) THEN
312
313       call find_max
314       call find_min
315  
316       if (tm_score) then
317        call get_diff_p
318       else
319        call get_diff
320       endif
321       if(nbank.eq.nconf.and.irestart.eq.0) then
322        adif=avedif
323       endif
324
325       write (iout,*) "AVEDIF",avedif
326       cutdif=adif/cut1
327       ctdif1=adif/cut2
328
329 cd      print *,"adif,xctdif,cutdifr"
330 cd      print *,adif,xctdif,cutdifr
331        nst=ntotal/ntrial/nseed
332        xctdif=(cutdif/ctdif1)**(-1.0/nst)
333        if(irestart.ge.1) call estimate_cutdif(adif,xctdif,cutdifr)
334 c       print *,"ok after estimate"
335
336       irestart=0
337
338        call write_rbank(jlee,adif,nft)
339        call write_bank(jlee,nft)
340 c!bankt       call write_bankt(jlee,nft)
341 c       call write_bank1(jlee)
342        call  history_append
343         write(icsa_history,*) "xctdif: ", xctdif,nst,adif/cut1,ctdif1
344         write(icsa_history,851) (bene(i),i=1,nbank)
345         write(icsa_history,850) (brmsn(i),i=1,nbank)
346         write(icsa_history,850) (bpncn(i),i=1,nbank)
347         close(icsa_history)
348   850 format(10f8.3)
349   851 format(5e15.6)
350
351       ifar=nseed/4*3+1
352       ifar=nseed+1
353       ENDIF
354     
355
356       finished=.false.
357       iter = 0
358       irecv = 0
359       isent =0
360       ifrom= 0
361       time0i=MPI_WTIME()
362       time1i=time0i
363       time_start_c=time0i
364       if (.not.sync_iter) then 
365         time_start=time0i
366         nft00=nft
367       else
368         sync_iter=.false.
369       endif
370       nft00_c=nft
371       nft0i=nft
372 ccccccccccccccccccccccccccccccccccccccc
373       do while (.not. finished)
374 ccccccccccccccccccccccccccccccccccccccc
375 crc      print *,"iter ", iter,' isent=',isent
376
377       IF (ME.EQ.KING) THEN
378 c  start energy minimization
379
380        if (isent.eq.0) then
381 c  king-emperor - select seeds & make var & feed input
382 cd        print *,'generating new conf',ntrial,MPI_WTIME()
383         call select_is(nseed,ifar,idum)
384
385         open(icsa_seed,file=csa_seed,status="old")
386         write(icsa_seed,39) 
387      &    jlee,icycle,nstep,(is(i),bene(is(i)),i=1,nseed)
388         close(icsa_seed)
389         call  history_append
390         write(icsa_history,40) jlee,icycle,nstep,cutdif,ibmin,ibmax,
391      *   ebmin,ebmax,nft,iuse,nbank,ntbank
392         close(icsa_history)
393
394          
395
396         call make_var(ntry,idum,iter)
397 cd        print *,'new trial generated',ntrial,MPI_WTIME()
398            time2i=MPI_WTIME()
399            write (iout,'(a20,i4,f12.2)') 
400      &       'Time for make trial',iter+1,time2i-time1i
401        endif
402
403 crc        write(iout,*)'1:Calling FEEDIN NTRY',NTRY,' ntrial',ntrial
404 crc        call feedin(ntry,nft)
405
406        isent=isent+1
407        if (isent.ge.nodes.or.iter.gt.0)  then
408 ct            print *,'waiting ',MPI_WTIME()
409             irecv=irecv+1
410             call recv(0,ifrom,xout,eout,ind,timeout)
411 ct            print *,'   ',irecv,' received from',ifrom,MPI_WTIME()
412
413             if(tm_score) then
414               nft=nft+ind(3)
415               movernx(irecv)=iabs(ind(5))
416               call getx(ind,xout,eout,cout,rad,iw_pdb,irecv)
417               if(vdisulf) then
418                nss_out(irecv)=nss
419                do i=1,nss
420                 iss_out(i,irecv)=ihpb(i)
421                 jss_out(i,irecv)=jhpb(i)  
422                enddo
423               endif
424               if(iw_pdb.gt.0) 
425      &          call write_csa_pdb(xout,eout,nft,irecv,iw_pdb)
426             endif
427
428             if(tm_score.and.eout(1).lt.ebmax) then
429              if(iref.eq.0   .or.
430      &         (rmsn(irecv).le.rmscut.and.pncn(irecv).ge.pnccut))
431      &         call refresh_bank_master_tmscore(ifrom,eout(1),irecv)
432             endif
433        else
434             ifrom=ifrom+1
435        endif
436
437 ct            print *,'sending to',ifrom,MPI_WTIME()
438        call send(isent,ifrom,iter)
439 ct            print *,isent,' sent ',MPI_WTIME()
440
441 c store results -----------------------------------------------
442        if ((isent.ge.nodes.or.iter.gt.0).and..not.tm_score)  then
443          nft=nft+ind(3)
444          movernx(irecv)=iabs(ind(5))
445          call getx(ind,xout,eout,cout,rad,iw_pdb,irecv)
446          if(vdisulf) then
447              nss_out(irecv)=nss
448              do i=1,nss
449                iss_out(i,irecv)=ihpb(i)
450                jss_out(i,irecv)=jhpb(i)  
451              enddo
452          endif
453          if(iw_pdb.gt.0) 
454      &          call write_csa_pdb(xout,eout,nft,irecv,iw_pdb)
455        endif
456 c--------------------------------------------------------------
457        if (isent.eq.ntry) then
458            time1i=MPI_WTIME()
459            write (iout,'(a18,f12.2,a14,f10.2)') 
460      &       'Nonsetup time     ',time1i-time_start_c,
461      &       ' sec, Eval/s =',(nft-nft00_c)/(time1i-time_start_c)
462            write (iout,'(a14,i4,f12.2,a14,f10.2)') 
463      &       'Time for iter ',iter+1,time1i-time0i,
464      &       ' sec, Eval/s =',(nft-nft0i)/(time1i-time0i)
465            time0i=time1i
466            nft0i=nft
467            cutdif=cutdif*xctdif
468            if(cutdif.lt.ctdif1) cutdif=ctdif1
469            if (iter.eq.0) then
470               print *,'UPDATING ',ntry-nodes+1,irecv
471               write(iout,*) 'UPDATING ',ntry-nodes+1
472               iter=iter+1
473 c----------------- call update(ntry-nodes+1) -------------------
474               nstep=nstep+ntry-nseed-(nodes-1)
475               if (tm_score) then
476 ctm               call refresh_bank(ntry)
477                call print_mv_stat
478                   do i=0,mxmv
479                    do j=1,3
480                     nstatnx_tot(i,j)=nstatnx_tot(i,j)+nstatnx(i,j)
481                     nstatnx(i,j)=0
482                    enddo
483                   enddo
484               else
485                call refresh_bank(ntry-nodes+1)
486               endif
487 c!bankt              call refresh_bankt(ntry-nodes+1)
488            else
489 c----------------- call update(ntry) ---------------------------
490               iter=iter+1
491               print *,'UPDATING ',ntry,irecv
492               write(iout,*) 'UPDATING ',ntry
493               nstep=nstep+ntry-nseed
494               if (tm_score) then
495 ctm               call refresh_bank(ntry)
496                call print_mv_stat
497                   do i=0,mxmv
498                    do j=1,3
499                     nstatnx_tot(i,j)=nstatnx_tot(i,j)+nstatnx(i,j)
500                     nstatnx(i,j)=0
501                    enddo
502                   enddo
503               else
504                call refresh_bank(ntry)
505               endif
506 c!bankt              call refresh_bankt(ntry)
507            endif         
508 c----------------------------------------------------------------- 
509
510            call write_bank(jlee,nft)
511 c!bankt           call write_bankt(jlee,nft)
512            call find_min
513
514            time1i=MPI_WTIME()
515            write (iout,'(a20,i4,f12.2)') 
516      &       'Time for refresh ',iter,time1i-time0i
517
518            if(ebmin.lt.estop) finished=.true.
519            if(icycle.gt.icmax) then
520                call write_bank1(jlee)
521                do i=1,nbank
522 c                 ibank(i)=2
523                  ibank(i)=1
524                enddo
525                nbank=nbank+nconf
526                if(nbank.gt.nbankm) then 
527                    nbank=nbank-nconf
528                    finished=.true.
529                else
530 crc                   goto 333
531                    sync_iter=.true.
532                endif
533            endif
534            if(nstep.gt.nstmax) finished=.true.
535
536            if(finished.or.sync_iter) then
537             do ij=1,nodes-1
538               call recv(1,ifrom,xout,eout,ind,timeout)
539               if (timeout) then
540                 nstep=nstep+ij-1
541                 print *,'ERROR worker is not responding'
542                 write(iout,*) 'ERROR worker is not responding'
543                 time1i=MPI_WTIME()-time_start_c
544                 print *,'End of cycle, master time for ',iter,' iters ',
545      &             time1i,'sec, Eval/s ',(nft-nft00_c)/time1i
546                 write (iout,*) 'End of cycle, master time for ',iter,
547      &             ' iters ',time1i,' sec'
548                 write (iout,*) 'Total eval/s ',(nft-nft00_c)/time1i
549                 print *,'UPDATING ',ij-1
550                 write(iout,*) 'UPDATING ',ij-1
551                 call flush(iout)
552                 call refresh_bank(ij-1)
553 c!bankt                call refresh_bankt(ij-1)
554                 goto 1002
555               endif
556 c              print *,'node ',ifrom,' finished ',ij,nft
557               write(iout,*) 'node ',ifrom,' finished ',ij,nft
558               call flush(iout)
559               nft=nft+ind(3)
560               movernx(ij)=iabs(ind(5))
561               call getx(ind,xout,eout,cout,rad,iw_pdb,ij)
562               if(vdisulf) then
563                nss_out(ij)=nss
564                do i=1,nss
565                  iss_out(i,ij)=ihpb(i)
566                  jss_out(i,ij)=jhpb(i)  
567                enddo
568               endif
569               if(iw_pdb.gt.0) 
570      &          call write_csa_pdb(xout,eout,nft,ij,iw_pdb)
571             enddo
572             nstep=nstep+nodes-1
573 crc            print *,'---------bcast finished--------',finished
574             time1i=MPI_WTIME()-time_start_c
575             print *,'End of cycle, master time for ',iter,' iters ',
576      &             time1i,'sec, Eval/s ',(nft-nft00_c)/time1i
577             write (iout,*) 'End of cycle, master time for ',iter,
578      &             ' iters ',time1i,' sec'
579             write (iout,*) 'Total eval/s ',(nft-nft00_c)/time1i
580
581 ctimeout            call mpi_bcast(finished,1,mpi_logical,0,CG_COMM,ierr)
582 ctimeout            call mpi_bcast(sync_iter,1,mpi_logical,0,
583 ctimeout     &                                              CG_COMM,ierr)
584             do ij=1,nodes-1 
585                tstart=MPI_WTIME()
586                call mpi_issend(finished,1,mpi_logical,ij,idchar,
587      &             CG_COMM,ireq,ierr)
588                call mpi_issend(sync_iter,1,mpi_logical,ij,idchar,
589      &             CG_COMM,ireq2,ierr)
590                flag=.false.  
591                timeout1=.false.
592                do while(.not. (flag .or. timeout1))
593                  call MPI_TEST(ireq2,flag,muster,ierr)
594                  tend1=MPI_WTIME()
595                  if(tend1-tstart.gt.60) then 
596                   print *,'ERROR worker ',ij,' is not responding'
597                   write(iout,*) 'ERROR worker ',ij,' is not responding'
598                   timeout1=.true.
599                  endif
600                enddo
601                if(timeout1) then
602                 write(iout,*) 'worker ',ij,' NOT OK ',tend1-tstart
603                 timeout=.true.
604                else
605                 write(iout,*) 'worker ',ij,' OK ',tend1-tstart
606                endif
607             enddo
608             print *,'UPDATING ',nodes-1,ij
609             write(iout,*) 'UPDATING ',nodes-1
610             call refresh_bank(nodes-1)
611 c!bankt            call refresh_bankt(nodes-1)
612  1002       continue
613             call write_bank(jlee,nft)
614 c!bankt            call write_bankt(jlee,nft)
615             call find_min
616
617             do i=0,mxmv  
618               do j=1,3
619                 nstatnx_tot(i,j)=nstatnx_tot(i,j)+nstatnx(i,j)
620                 nstatnx(i,j)=0
621               enddo
622             enddo
623
624             write(iout,*)'### Total stats:'
625             do i=0,mxmv
626              if(nstatnx_tot(i,1).ne.0) then
627               if (i.le.9) then
628               write(iout,'(a5,i1,a7,i6,a7,i4,a5,i4,a5,f5.1)') 
629      &        '### N',i,' total=',nstatnx_tot(i,1),
630      &      ' close=',nstatnx_tot(i,2),' far=',nstatnx_tot(i,3),'%acc',
631      &       (nstatnx_tot(i,2)+nstatnx_tot(i,3))*100.0/nstatnx_tot(i,1)
632               else
633               write(iout,'(a4,i2,a7,i6,a7,i4,a5,i4,a5,f5.1)') 
634      &        '###N',i,' total=',nstatnx_tot(i,1),
635      &      ' close=',nstatnx_tot(i,2),' far=',nstatnx_tot(i,3),'%acc',
636      &       (nstatnx_tot(i,2)+nstatnx_tot(i,3))*100.0/nstatnx_tot(i,1)
637               endif
638              else
639               if (i.le.9) then
640               write(iout,'(a5,i1,a7,i6,a7,i4,a5,i4,a5,f5.1)') 
641      &          '### N',i,' total=',nstatnx_tot(i,1),
642      &          ' close=',nstatnx_tot(i,2),' far=',nstatnx_tot(i,3),
643      &          ' %acc',0.0
644               else
645               write(iout,'(a4,i2,a7,i6,a7,i4,a5,i4,a5,f5.1)') 
646      &          '###N',i,' total=',nstatnx_tot(i,1),
647      &          ' close=',nstatnx_tot(i,2),' far=',nstatnx_tot(i,3),
648      &          ' %acc',0.0
649               endif
650              endif
651             enddo
652
653            endif
654            if(sync_iter) goto 331
655
656    39 format(2i3,i7,5(i4,f8.3,1x),19(/,13x,5(i4,f8.3,1x)))
657    40 format(2i2,i8,f8.1,2i4,2(1pe14.5),i10,3i4)
658    43 format(10i8)
659    44 format('jlee =',i3,':',4f10.1,' E =',f15.5,i7,i10)
660
661            isent=0
662            irecv=0
663        endif
664       ELSE
665        if (tm_score) then
666         call get_diff_p
667        endif
668 c  soldier - perform energy minimization
669         call minim_jlee
670         print *,'End of minim, proc',me,'time ',MPI_WTIME()-time_start
671         write (iout,*) 'End of minim, proc',me,'time ',
672      &                  MPI_WTIME()-time_start
673         call flush(iout)
674 ctimeout        call mpi_bcast(finished,1,mpi_logical,0,CG_COMM,ierr)
675 ctimeout        call mpi_bcast(sync_iter,1,mpi_logical,0,CG_COMM,ierr)
676          call mpi_recv(finished,1,mpi_logical,0,idchar,
677      *                 CG_COMM,muster,ierr)             
678          call mpi_recv(sync_iter,1,mpi_logical,0,idchar,
679      *                 CG_COMM,muster,ierr)             
680         if(sync_iter) goto 331
681       ENDIF
682
683 ccccccccccccccccccccccccccccccccccccccc
684       enddo
685 ccccccccccccccccccccccccccccccccccccccc
686
687       IF (ME.EQ.KING) THEN
688         call  history_append
689         write(icsa_history,40) jlee,icycle,nstep,cutdif,ibmin,ibmax,
690      *  ebmin,ebmax,nft,iuse,nbank,ntbank
691
692         write(icsa_history,44) jlee,0.0,0.0,0.0,
693      &   0.0,ebmin,nstep,nft
694         write(icsa_history,*)
695        close(icsa_history)
696
697        time1i=MPI_WTIME()-time_start
698        print *,'End of RUN, master time ',
699      &             time1i,'sec, Eval/s ',(nft-nft00)/time1i
700        write (iout,*) 'End of RUN, master time  ',
701      &             time1i,' sec'
702        write (iout,*) 'Total eval/s ',(nft-nft00)/time1i
703
704        if(timeout) then 
705         write(iout,*) '!!!! ERROR worker was not responding'
706         write(iout,*) '!!!! cannot finish work normally'
707         write(iout,*) 'Processor0 is calling MPI_ABORT'
708         print *,'!!!! ERROR worker was not responding'
709         print *,'!!!! cannot finish work normally'
710         print *,'Processor0 is calling MPI_ABORT'
711         call flush(iout)
712         call mpi_abort(mpi_comm_world, 111, ierr)
713        endif
714       ENDIF
715
716 cccccccccccccccccccccccccccccc
717   300 continue
718 cccccccccccccccccccccccccccccc
719
720       return
721       end
722 #else
723       Subroutine together
724 c  feeds tasks for parallel processing
725       implicit real*8 (a-h,o-z)
726       include 'DIMENSIONS'
727       include 'COMMON.IOUNITS'
728       write (iout,*) "Unsupported option for the serial version"
729       return
730       end
731 #endif
732 #ifdef MPI
733 c-------------------------------------------------
734       subroutine feedin(nconf,nft)
735 c  sends out starting conformations and receives results of energy minimization
736       implicit real*8 (a-h,o-z)
737       include 'DIMENSIONS'
738       include 'COMMON.VAR'
739       include 'COMMON.IOUNITS'
740       include 'COMMON.CONTROL'
741       include 'mpif.h'
742       dimension xin(maxvar),xout(maxvar),eout(mxch*(mxch+1)/2+1),
743      *          cout(2),ind(9),info(12)
744       dimension muster(mpi_status_size)
745       include 'COMMON.SETUP'
746       parameter (rad=1.745329252d-2)
747
748       print *,'FEEDIN: NCONF=',nconf
749       mm=0
750 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
751       if (nconf .lt. nodes-1) then
752         write (*,*) 'FATAL ERROR in FEEDIN, nconf < nodes -1',
753      &   nconf,nodes-1 
754         write (iout,*) 'FATAL ERROR in FEEDIN, nconf < nodes -1',
755      &   nconf,nodes-1 
756         call mpi_abort(mpi_comm_world,ierror,ierrcode)
757       endif
758       do n=1,nconf
759 c  pull out external and internal variables for next start
760         call putx(xin,n,rad)
761 !        write (iout,*) 'XIN from FEEDIN N=',n
762 !        write(iout,'(8f10.4)') (xin(j),j=1,nvar)
763         mm=mm+1
764         if (mm.lt.nodes) then
765 c  feed task to soldier
766 !       print *, ' sending input for start # ',n
767          info(1)=n
768          info(2)=-1
769          info(3)=0
770          info(4)=0
771          info(5)=0
772          info(6)=0
773          call mpi_send(info,12,mpi_integer,mm,idint,CG_COMM,
774      *                  ierr)
775          call mpi_send(xin,nvar,mpi_double_precision,mm,
776      *                  idreal,CG_COMM,ierr)
777         else
778 c  find an available soldier
779          call mpi_recv(ind,9,mpi_integer,mpi_any_source,idint,
780      *                 CG_COMM,muster,ierr)
781 !        print *, ' receiving output from start # ',ind(1)
782          man=muster(mpi_source)
783 c  receive final energies and variables
784          nft=nft+ind(3)
785          call mpi_recv(eout,1,mpi_double_precision,
786      *                 man,idreal,CG_COMM,muster,ierr)
787 !         print *,eout 
788 #ifdef CO_BIAS
789          call mpi_recv(co,1,mpi_double_precision,
790      *                 man,idreal,CG_COMM,muster,ierr)
791          write (iout,'(a15,f3.2,$)') ' BIAS by contact order*100 ',co
792 #endif
793          call mpi_recv(xout,nvar,mpi_double_precision,
794      *                 man,idreal,CG_COMM,muster,ierr)
795 !         print *,nvar , ierr
796 c  feed next task to soldier
797 !        print *, ' sending input for start # ',n
798          info(1)=n
799          info(2)=-1
800          info(3)=0
801          info(4)=0
802          info(5)=0
803          info(6)=0
804          info(7)=0
805          info(8)=0
806          info(9)=0
807          call mpi_send(info,12,mpi_integer,man,idint,CG_COMM,
808      *                  ierr)
809          call mpi_send(xin,nvar,mpi_double_precision,man,
810      *                  idreal,CG_COMM,ierr)
811 c  retrieve latest results
812          call getx(ind,xout,eout,cout,rad,iw_pdb,ind(1))
813          if(iw_pdb.gt.0) 
814      &        call write_csa_pdb(xout,eout,nft,ind(1),iw_pdb)
815         endif
816       enddo
817 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
818 c  no more input
819 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
820       do j=1,nodes-1
821 c  wait for a soldier
822        call mpi_recv(ind,9,mpi_integer,mpi_any_source,idint,
823      *               CG_COMM,muster,ierr)
824 crc       if (ierr.ne.0) go to 30
825 !      print *, ' receiving output from start # ',ind(1)
826        man=muster(mpi_source)
827 c  receive final energies and variables
828        nft=nft+ind(3)
829        call mpi_recv(eout,1,
830      *               mpi_double_precision,man,idreal,
831      *               CG_COMM,muster,ierr)
832 !       print *,eout
833 #ifdef CO_BIAS
834          call mpi_recv(co,1,mpi_double_precision,
835      *                 man,idreal,CG_COMM,muster,ierr)
836          write (iout,'(a15,f3.2,$)') ' BIAS by contact order*100 ',co
837 #endif
838 crc       if (ierr.ne.0) go to 30
839        call mpi_recv(xout,nvar,mpi_double_precision,
840      *               man,idreal,CG_COMM,muster,ierr)
841 !       print *,nvar , ierr
842 crc       if (ierr.ne.0) go to 30
843 c  halt soldier
844        info(1)=0
845        info(2)=-1
846        info(3)=0 
847        info(4)=0
848        info(5)=0
849        info(6)=0
850        call mpi_send(info,12,mpi_integer,man,idint,CG_COMM,
851      *                ierr)
852 c  retrieve results
853        call getx(ind,xout,eout,cout,rad,iw_pdb,ind(1))
854        if(iw_pdb.gt.0) 
855      &          call write_csa_pdb(xout,eout,nft,ind(1),iw_pdb)
856       enddo
857 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
858       return
859    10 print *, ' dispatching error'
860       call mpi_abort(mpi_comm_world,ierror,ierrcode)
861       return
862    20 print *, ' communication error'
863       call mpi_abort(mpi_comm_world,ierror,ierrcode)
864       return
865    30 print *, ' receiving error'
866       call mpi_abort(mpi_comm_world,ierror,ierrcode)
867       return
868       end
869 cccccccccccccccccccccccccccccccccccccccccccccccccc
870       subroutine getx(ind,xout,eout,cout,rad,iw_pdb,k)
871 c  receives and stores data from soldiers
872       implicit real*8 (a-h,o-z)
873       include 'DIMENSIONS'
874       include 'COMMON.IOUNITS'
875       include 'COMMON.CSA'
876       include 'COMMON.BANK'
877       include 'COMMON.VAR'
878       include 'COMMON.CHAIN'
879       include 'COMMON.CONTACTS'
880       dimension ind(9),xout(maxvar),eout(mxch*(mxch+1)/2+1)
881       dimension cout(2)
882 cjlee
883       double precision przes(3),obr(3,3)
884       logical non_conv
885 cjlee
886       iw_pdb=2
887       if (k.gt.mxio .or. k.lt.1) then 
888         write (iout,*) 
889      &   'ERROR - dimensions of ANGMIN have been exceeded K=',k
890         call mpi_abort(mpi_comm_world,ierror,ierrcode)
891       endif
892 c  store ind()
893       do j=1,9
894        indb(k,j)=ind(j)
895       enddo
896 c  store energies
897       etot(k)=eout(1)
898 c  retrieve dihedral angles etc
899       call var_to_geom(nvar,xout)
900       do j=2,nres-1
901         dihang(1,j,1,k)=theta(j+1)
902         dihang(2,j,1,k)=phi(j+2)
903         dihang(3,j,1,k)=alph(j)
904         dihang(4,j,1,k)=omeg(j)
905       enddo
906       dihang(2,nres-1,1,k)=0.0d0
907 cjlee
908       if(iref.eq.0) then 
909        iw_pdb=1
910 cd       write(iout,'(i3,a3,i4,i5,a6,1pe12.4,a4,i3,i4)') 
911 cd     &      ind(2),' e ',ind(3),ind(1),' etot ',etot(k),' mv ',
912 cd     &      ind(5),ind(4)
913        return
914       endif
915       call chainbuild
916 c     call dihang_to_c(dihang(1,1,1,k))
917 c     call fitsq(rms,c(1,1),crefjlee(1,1),nres,przes,obr,non_conv)
918 c     call fitsq(rms,c(1,2),crefjlee(1,2),nres-1,przes,obr,non_conv)
919 c           call fitsq(rms,c(1,nstart_seq),crefjlee(1,nstart_sup),
920 c    &                 nsup,przes,obr,non_conv)
921 c      rmsn(k)=dsqrt(rms)
922
923        call rmsd_csa(rmsn(k))
924        call contact(.false.,ncont,icont,co)
925        pncn(k)=contact_fract(ncont,ncont_ref,icont,icont_ref)     
926
927 cd       write(iout,'(i3,a3,i4,i5,a6,1pe12.4,a5
928 cd     &      ,0pf5.2,a5,f5.1,a,f6.3,a4,i3,i4)') 
929 cd     &    ind(2),' e ',ind(3),ind(1),' etot ',etot(k),' rms ',
930 cd     &    rmsn(k),' %NC ',pncn(k)*100,' cont.order',co,' mv ',
931 cd     &    ind(5),ind(4)
932
933  
934       if (rmsn(k).gt.rmscut.or.pncn(k).lt.pnccut) iw_pdb=0
935       return
936       end
937 cccccccccccccccccccccccccccccccccccccccccccccccccc
938       subroutine putx(xin,n,rad)
939 c  gets starting variables
940       implicit real*8 (a-h,o-z)
941       include 'DIMENSIONS'
942       include 'COMMON.CSA'
943       include 'COMMON.BANK'
944       include 'COMMON.VAR'
945       include 'COMMON.CHAIN'
946       include 'COMMON.IOUNITS'
947       dimension xin(maxvar)
948
949 c  pull out starting values for variables
950 !       write (iout,*)'PUTX: N=',n
951       do m=1,numch
952 !        write (iout,'(8f10.4)') (dihang_in(1,j,m,n),j=2,nres-1)
953 !        write (iout,'(8f10.4)') (dihang_in(2,j,m,n),j=2,nres-1)
954 !        write (iout,'(8f10.4)') (dihang_in(3,j,m,n),j=2,nres-1)
955 !        write (iout,'(8f10.4)') (dihang_in(4,j,m,n),j=2,nres-1)
956        do j=2,nres-1
957         theta(j+1)=dihang_in(1,j,m,n)
958         phi(j+2)=dihang_in(2,j,m,n)
959         alph(j)=dihang_in(3,j,m,n)
960         omeg(j)=dihang_in(4,j,m,n)
961        enddo
962       enddo
963 c  set up array of variables
964       call geom_to_var(nvar,xin)
965 !       write (iout,*) 'xin in PUTX N=',n 
966 !       call intout
967 !       write (iout,'(8f10.4)') (xin(i),i=1,nvar) 
968       return
969       end
970 c--------------------------------------------------------
971       subroutine putx2(xin,iff,n)
972 c  gets starting variables
973       implicit real*8 (a-h,o-z)
974       include 'DIMENSIONS'
975       include 'COMMON.CSA'
976       include 'COMMON.BANK'
977       include 'COMMON.VAR'
978       include 'COMMON.CHAIN'
979       include 'COMMON.IOUNITS'
980       dimension xin(maxvar),iff(maxres)
981
982 c  pull out starting values for variables
983       do m=1,numch
984        do j=2,nres-1
985         theta(j+1)=dihang_in2(1,j,m,n)
986         phi(j+2)=dihang_in2(2,j,m,n)
987         alph(j)=dihang_in2(3,j,m,n)
988         omeg(j)=dihang_in2(4,j,m,n)
989        enddo
990       enddo
991 c  set up array of variables
992       call geom_to_var(nvar,xin)
993        
994       do i=1,nres
995         iff(i)=iff_in(i,n)
996       enddo
997       return
998       end
999
1000 c-------------------------------------------------------
1001       subroutine prune_bank(p_cut)
1002       implicit real*8 (a-h,o-z)
1003       include 'DIMENSIONS'
1004       include 'mpif.h'
1005       include 'COMMON.CSA'
1006       include 'COMMON.BANK'
1007       include 'COMMON.IOUNITS'
1008       include 'COMMON.CHAIN'
1009       include 'COMMON.TIME1'
1010       include 'COMMON.SETUP'
1011 c---------------------------
1012 c This subroutine prunes bank conformations using p_cut
1013 c---------------------------
1014
1015       nprune=0
1016       nprune=nprune+1
1017       m=1
1018       do k=1,numch
1019        do j=2,nres-1
1020         do i=1,4
1021          dihang(i,j,k,nprune)=bvar(i,j,k,m)
1022         enddo
1023        enddo
1024       enddo
1025       bene(nprune)=bene(m)
1026       brmsn(nprune)=brmsn(m)
1027       bpncn(nprune)=bpncn(m) 
1028
1029       do m=2,nbank
1030        ddmin=9.d190
1031        do ip=1,nprune
1032         call get_diff12(dihang(1,1,1,ip),bvar(1,1,1,m),diff) 
1033         if(diff.lt.p_cut) goto 100
1034         if(diff.lt.ddmin) ddmin=diff
1035        enddo
1036        nprune=nprune+1
1037        do k=1,numch
1038         do j=2,nres-1
1039          do i=1,4
1040           dihang(i,j,k,nprune)=bvar(i,j,k,m)
1041          enddo
1042         enddo
1043        enddo
1044        bene(nprune)=bene(m)
1045        brmsn(nprune)=brmsn(m)
1046        bpncn(nprune)=bpncn(m)
1047   100  continue
1048        write (iout,*) 'Pruning :',m,nprune,p_cut,ddmin
1049       enddo
1050       nbank=nprune
1051       print *, 'Pruning :',m,nprune,p_cut
1052       call write_bank(0,0)
1053
1054       return
1055       end
1056 c-------------------------------------------------------
1057
1058       subroutine reminimize(jlee)
1059       implicit real*8 (a-h,o-z)
1060       include 'DIMENSIONS'
1061       include 'mpif.h'
1062       include 'COMMON.CSA'
1063       include 'COMMON.BANK'
1064       include 'COMMON.IOUNITS'
1065       include 'COMMON.CHAIN'
1066       include 'COMMON.TIME1'
1067       include 'COMMON.SETUP'
1068 c---------------------------
1069 c This subroutine re-minimizes bank conformations:
1070 c---------------------------
1071
1072        ntry=nbank
1073
1074        call find_max
1075        call find_min
1076
1077        if (me.eq.king) then
1078         open(icsa_history,file=csa_history,status="old")
1079          write(icsa_history,*) "Re-minimization",nodes,"nodes"
1080          write(icsa_history,851) (bene(i),i=1,nbank)
1081          write(icsa_history,40) jlee,icycle,nstep,cutdif,ibmin,ibmax,
1082      *   ebmin,ebmax,nft,iuse,nbank,ntbank
1083         close(icsa_history)
1084         do index=1,ntry
1085          do k=1,numch
1086           do j=2,nres-1
1087            do i=1,4
1088             dihang_in(i,j,k,index)=bvar(i,j,k,index)
1089            enddo
1090           enddo
1091          enddo
1092         enddo
1093         nft=0
1094         call feedin(ntry,nft)
1095        else
1096         call minim_jlee
1097        endif
1098
1099        call find_max
1100        call find_min
1101
1102        if (me.eq.king) then
1103         do i=1,ntry
1104          call replace_bvar(i,i)
1105         enddo
1106         open(icsa_history,file=csa_history,status="old")
1107          write(icsa_history,40) jlee,icycle,nstep,cutdif,ibmin,ibmax,
1108      *   ebmin,ebmax,nft,iuse,nbank,ntbank
1109          write(icsa_history,851) (bene(i),i=1,nbank)
1110         close(icsa_history)
1111         call write_bank_reminimized(jlee,nft)
1112        endif
1113
1114    40 format(2i2,i8,f8.1,2i4,2(1pe14.5),i10,3i4)
1115   851 format(5e15.6)
1116   850 format(5e15.10)
1117 c  850 format(10f8.3)
1118
1119       return
1120       end
1121 c-------------------------------------------------------
1122       subroutine send(n,mm,it)
1123 c  sends out starting conformation for minimization
1124       implicit real*8 (a-h,o-z)
1125       include 'DIMENSIONS'
1126       include 'COMMON.VAR'
1127       include 'COMMON.IOUNITS'
1128       include 'COMMON.CONTROL'
1129       include 'COMMON.BANK'
1130       include 'COMMON.CHAIN'
1131       include 'mpif.h'
1132       dimension xin(maxvar),xout(maxvar),eout(mxch*(mxch+1)/2+1),
1133      *          cout(2),ind(8),xin2(maxvar),iff(maxres),info(12)
1134       dimension muster(mpi_status_size)
1135       include 'COMMON.SETUP'
1136       parameter (rad=1.745329252d-2)
1137
1138       if (isend2(n).eq.0) then
1139 c  pull out external and internal variables for next start
1140         call putx(xin,n,rad)
1141         info(1)=n
1142         info(2)=it
1143         info(3)=movenx(n)
1144         info(4)=nss_in(n)
1145         info(5)=parent(1,n)
1146         info(6)=parent(2,n)
1147
1148         if (movenx(n).eq.14.or.movenx(n).eq.17) then
1149           info(7)=idata(1,n)
1150           info(8)=idata(2,n)
1151         else if (movenx(n).eq.16) then
1152           info(7)=idata(1,n)
1153           info(8)=idata(2,n)
1154           info(10)=idata(3,n)
1155           info(11)=idata(4,n)
1156           info(12)=idata(5,n)
1157         else
1158          info(7)=0
1159          info(8)=0
1160          info(10)=0
1161          info(11)=0
1162          info(12)=0
1163         endif
1164
1165         if (movenx(n).eq.15) then
1166          info(9)=parent(3,n)
1167         else
1168          info(9)=0
1169         endif
1170         call mpi_send(info,12,mpi_integer,mm,idint,CG_COMM,
1171      *                  ierr)
1172         call mpi_send(xin,nvar,mpi_double_precision,mm,
1173      *                  idreal,CG_COMM,ierr)
1174       else
1175 c  distfit & minimization for n7 move
1176         info(1)=-n
1177         info(2)=it
1178         info(3)=movenx(n)
1179         info(4)=nss_in(n)
1180         info(5)=parent(1,n)
1181         info(6)=parent(2,n)
1182         info(7)=0
1183         info(8)=0
1184         info(9)=0
1185         call mpi_send(info,12,mpi_integer,mm,idint,CG_COMM,
1186      *                  ierr)
1187         call putx2(xin,iff,isend2(n))
1188         call mpi_send(xin,nvar,mpi_double_precision,mm,
1189      *                  idreal,CG_COMM,ierr)
1190         call mpi_send(iff,nres,mpi_integer,mm,
1191      *                  idint,CG_COMM,ierr)
1192         call putx(xin2,n,rad)
1193         call mpi_send(xin2,nvar,mpi_double_precision,mm,
1194      *                  idreal,CG_COMM,ierr)
1195       endif 
1196       if (vdisulf.and.nss_in(n).ne.0) then
1197         call mpi_send(iss_in(1,n),nss_in(n),mpi_integer,mm,
1198      *                  idint,CG_COMM,ierr)
1199         call mpi_send(jss_in(1,n),nss_in(n),mpi_integer,mm,
1200      *                  idint,CG_COMM,ierr)
1201       endif
1202       return
1203       end
1204 c-------------------------------------------------
1205
1206       subroutine recv(ihalt,man,xout,eout,ind,tout)
1207 c  receives results of energy minimization
1208       implicit real*8 (a-h,o-z)
1209       include 'DIMENSIONS'
1210       include 'COMMON.VAR'
1211       include 'COMMON.IOUNITS'
1212       include 'COMMON.CONTROL'
1213       include 'COMMON.SBRIDGE'
1214       include 'COMMON.BANK'
1215       include 'COMMON.CHAIN'
1216       include 'mpif.h'
1217       dimension xin(maxvar),xout(maxvar),eout(mxch*(mxch+1)/2+1),
1218      *          cout(2),ind(9),info(12)
1219       dimension muster(mpi_status_size)
1220       include 'COMMON.SETUP'
1221       logical tout,flag
1222       double precision twait,tstart,tend1
1223       parameter(twait=600.0d0)
1224
1225 c  find an available soldier
1226        tout=.false.
1227        flag=.false.
1228        tstart=MPI_WTIME()
1229        do while(.not. (flag .or. tout))
1230          call MPI_IPROBE(mpi_any_source,idint,CG_COMM,flag, 
1231      *            muster,ierr)
1232          tend1=MPI_WTIME()
1233          if(tend1-tstart.gt.twait .and. ihalt.eq.1) tout=.true.
1234 c_error         if(tend1-tstart.gt.twait) tout=.true.
1235        enddo
1236        if (tout) then 
1237          write(iout,*) 'ERROR = timeout for recv ',tend1-tstart
1238          call flush(iout)
1239          return
1240        endif
1241        man=muster(mpi_source)        
1242
1243 ctimeout         call mpi_recv(ind,9,mpi_integer,mpi_any_source,idint,
1244 ctimeout     *                 CG_COMM,muster,ierr)
1245 !        print *, ' receiving output from start # ',ind(1)
1246 ct         print *,'receiving ',MPI_WTIME()
1247 ctimeout         man=muster(mpi_source)
1248          call mpi_recv(ind,9,mpi_integer,man,idint,
1249      *                 CG_COMM,muster,ierr)
1250 ctimeout
1251 c  receive final energies and variables
1252          call mpi_recv(eout,1,mpi_double_precision,
1253      *                 man,idreal,CG_COMM,muster,ierr)
1254 !         print *,eout 
1255 #ifdef CO_BIAS
1256          call mpi_recv(co,1,mpi_double_precision,
1257      *                 man,idreal,CG_COMM,muster,ierr)
1258          write (iout,'(a15,f3.2,$)') ' BIAS by contact order*100 ',co
1259 #endif
1260          call mpi_recv(xout,nvar,mpi_double_precision,
1261      *                 man,idreal,CG_COMM,muster,ierr)
1262 !         print *,nvar , ierr
1263          if(vdisulf) nss=ind(6)
1264          if(vdisulf.and.nss.ne.0) then
1265           call mpi_recv(ihpb,nss,mpi_integer,
1266      *                 man,idint,CG_COMM,muster,ierr)         
1267           call mpi_recv(jhpb,nss,mpi_integer,
1268      *                 man,idint,CG_COMM,muster,ierr)
1269          endif
1270 c  halt soldier
1271        if(ihalt.eq.1) then 
1272 c        print *,'sending halt to ',man
1273         write(iout,*) 'sending halt to ',man
1274         info(1)=0
1275         info(2)=0
1276         call mpi_send(info,12,mpi_integer,man,idint,CG_COMM,ierr)
1277        endif
1278       return
1279       end
1280
1281 c---------------------------------------------------------- 
1282       subroutine history_append 
1283       implicit real*8 (a-h,o-z)
1284       include 'DIMENSIONS'
1285       include 'COMMON.IOUNITS'
1286
1287 #if defined(AIX) || defined(PGI)
1288        open(icsa_history,file=csa_history,position="append")
1289 #else
1290        open(icsa_history,file=csa_history,access="append")
1291 #endif
1292        return
1293        end
1294 #endif