correction to the last commit - cleaning of homo restraints
[unres.git] / source / cluster / wham / src-M / read_coords.F
1       subroutine read_coords(ncon,*)
2       implicit none
3       include "DIMENSIONS"
4       include "sizesclu.dat"
5 #ifdef MPI
6       include "mpif.h"
7       integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
8       include "COMMON.MPI"
9 #endif
10       include "COMMON.CONTROL"
11       include "COMMON.CHAIN"
12       include "COMMON.INTERACT"
13       include "COMMON.IOUNITS"
14       include "COMMON.VAR"
15       include "COMMON.SBRIDGE"
16       include "COMMON.GEO"
17       include "COMMON.CLUSTER"
18       character*3 liczba
19       integer ncon
20       integer i,j,jj,jjj,jj_old,icount,k,kk,l,ii,if,ib,
21      &  nn,nn1,inan
22       integer ixdrf,iret,itmp
23       real*4 prec,reini,refree,rmsdev
24       integer nrec,nlines,iscor,lenrec,lenrec_in
25       double precision energ,t_acq,tcpu
26       integer ilen,iroof
27       external ilen,iroof
28       double precision rjunk
29       integer ntot_all(0:maxprocs-1)
30       logical lerr
31       double precision energia(0:max_ene),etot
32       real*4 csingle(3,maxres2+2)
33       integer Previous,Next
34       character*256 bprotfiles
35 c      print *,"Processor",me," calls read_protein_data"
36 #ifdef MPI
37       if (me.eq.master) then
38         Previous=MPI_PROC_NULL
39       else
40         Previous=me-1
41       endif
42       if (me.eq.nprocs-1) then
43         Next=MPI_PROC_NULL
44       else
45         Next=me+1
46       endif
47 c Set the scratchfile names
48       write (liczba,'(bz,i3.3)') me
49 #endif
50 c 1/27/05 AL Change stored coordinates to single precision and don't store 
51 c         energy components in the binary databases.
52       lenrec=12*(nres+nct-nnt+1)+4*(2*nss+2)+16
53       lenrec_in=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
54 #ifdef DEBUG
55       write (iout,*) "nres",nres," nnt",nnt," nct",nct," nss", nss
56       write (iout,*) "lenrec_in",lenrec_in
57 #endif
58       bprotfiles=scratchdir(:ilen(scratchdir))//
59      &       "/"//prefix(:ilen(prefix))//liczba//".xbin"
60
61 #ifdef CHUJ
62       ICON=1
63   123 continue
64       if (from_cart .and. .not. from_bx .and. .not. from_cx) then
65         if (efree) then
66         read (intin,*,end=13,err=11) energy(icon),totfree(icon),
67      &    rmstb(icon),
68      &    nss_all(icon),(ihpb_all(ii,icon),jhpb_all(i,icon),
69      &    i=1,nss_all(icon)),iscore(icon)
70         else
71         read (intin,*,end=13,err=11) energy(icon),rmstb(icon),
72      &    nss_all(icon),(ihpb_all(ii,icon),jhpb_all(i,icon),
73      &    i=1,nss_all(icon)),iscore(icon)
74         endif
75         read (intin,'(8f10.5)',end=13,err=10) 
76      &    ((allcart(j,i,icon),j=1,3),i=1,nres),
77      &    ((allcart(j,i+nres,icon),j=1,3),i=nnt,nct)
78         print *,icon,energy(icon),nss_all(icon),rmstb(icon)
79       else 
80         read(intin,'(a80)',end=13,err=12) lineh
81         read(lineh(:5),*,err=8) ic
82         if (efree) then
83         read(lineh(6:),*,err=8) energy(icon)
84         else
85         read(lineh(6:),*,err=8) energy(icon)
86         endif
87         goto 9
88     8   ic=1
89         print *,'error, assuming e=1d10',lineh
90         energy(icon)=1d10
91         nss=0
92     9   continue
93 cold        read(lineh(18:),*,end=13,err=11) nss_all(icon)
94         ii = index(lineh(15:)," ")+15
95         read(lineh(ii:),*,end=13,err=11) nss_all(icon)
96         IF (NSS_all(icon).LT.9) THEN
97           read (lineh(20:),*,end=102)
98      &    (IHPB_all(I,icon),JHPB_all(I,icon),I=1,NSS_all(icon)),
99      &    iscore(icon)
100         ELSE
101           read (lineh(20:),*,end=102) 
102      &           (IHPB_all(I,icon),JHPB_all(I,icon),I=1,8)
103           read (intin,*) (IHPB_all(I,icon),JHPB_all(I,icon),
104      &      I=9,NSS_all(icon)),iscore(icon)
105         ENDIF
106
107   102   continue  
108
109         PRINT *,'IC:',IC,' ENERGY:',ENERGY(ICON)
110         call read_angles(intin,*13)
111         do i=1,nres
112           phiall(i,icon)=phi(i)
113           thetall(i,icon)=theta(i)
114           alphall(i,icon)=alph(i)
115           omall(i,icon)=omeg(i)
116         enddo
117       endif
118       ICON=ICON+1
119       GOTO 123
120 C
121 C CALCULATE DISTANCES
122 C
123    10 print *,'something wrong with angles'
124       goto 13
125    11 print *,'something wrong with NSS',nss
126       goto 13
127    12 print *,'something wrong with header'
128
129    13 NCON=ICON-1
130
131 #endif
132       call flush(iout)
133       jj_old=1
134       open (icbase,file=bprotfiles,status="unknown",
135      &   form="unformatted",access="direct",recl=lenrec)
136 c Read conformations from binary DA files (one per batch) and write them to 
137 c a binary DA scratchfile.
138       jj=0
139       jjj=0
140 #ifdef MPI
141       write (liczba,'(bz,i3.3)') me
142       IF (ME.EQ.MASTER) THEN
143 c Only the master reads the database; it'll send it to the other procs
144 c through a ring.
145 #endif
146         t_acq = tcpu()
147         icount=0
148
149         if (from_bx) then
150
151           open (intin,file=intinname,status="old",form="unformatted",
152      &            access="direct",recl=lenrec_in)
153
154         else if (from_cx) then
155 #if (defined(AIX) && !defined(JUBL))
156           call xdrfopen_(ixdrf,intinname, "r", iret)
157 #else
158           call xdrfopen(ixdrf,intinname, "r", iret)
159 #endif
160           prec=10000.0
161           write (iout,*) "xdrfopen: iret",iret
162           if (iret.eq.0) then
163             write (iout,*) "Error: coordinate file ",
164      &       intinname(:ilen(intinname))," does not exist."
165             call flush(iout)
166 #ifdef MPI
167             call MPI_ABORT(MPI_COMM_WORLD,IERROR,ERRCODE)
168 #endif
169             stop
170           endif
171         else
172           write (iout,*) "Error: coordinate format not specified"
173           call flush(iout)
174 #ifdef MPI
175           call MPI_ABORT(MPI_COMM_WORLD,IERROR,ERRCODE)
176 #else
177           stop
178 #endif
179         endif
180
181 C#define DEBUG
182 #ifdef DEBUG
183         write (iout,*) "Opening file ",intinname(:ilen(intinname))
184         write (iout,*) "lenrec",lenrec_in
185         call flush(iout)
186 #endif
187 C#undef DEBUG
188 c        write (iout,*) "maxconf",maxconf
189         i=0
190         do while (.true.)
191            i=i+1
192            if (i.gt.maxconf) then
193              write (iout,*) "Error: too many conformations ",
194      &        "(",maxconf,") maximum."
195 #ifdef MPI
196              call MPI_Abort(MPI_COMM_WORLD,errcode,ierror)
197 #endif
198              stop
199            endif
200 c          write (iout,*) "i",i
201 c          call flush(iout)
202           if (from_bx) then
203             read(intin,err=101,end=101) 
204      &       ((csingle(l,k),l=1,3),k=1,nres),
205      &       ((csingle(l,k+nres),l=1,3),k=nnt,nct),
206      &       nss,(ihpb(k),jhpb(k),k=1,nss),
207      &       energy(jj+1),
208      &       entfac(jj+1),rmstb(jj+1),iscor
209              do j=1,2*nres
210                do k=1,3
211                  c(k,j)=csingle(k,j)
212                enddo
213              enddo
214           else
215 #if (defined(AIX) && !defined(JUBL))
216             call xdrf3dfcoord_(ixdrf, csingle, itmp, prec, iret)
217             if (iret.eq.0) goto 101
218             call xdrfint_(ixdrf, nss, iret)
219             if (iret.eq.0) goto 101
220             do j=1,nss
221            if (dyn_ss) then
222             call xdrfint(ixdrf, idssb(j), iret)
223             call xdrfint(ixdrf, jdssb(j), iret)
224         idssb(j)=idssb(j)-nres
225         jdssb(j)=jdssb(j)-nres
226            else
227               call xdrfint_(ixdrf, ihpb(j), iret)
228               if (iret.eq.0) goto 101
229               call xdrfint_(ixdrf, jhpb(j), iret)
230               if (iret.eq.0) goto 101
231            endif
232             enddo
233             call xdrffloat_(ixdrf,reini,iret)
234             if (iret.eq.0) goto 101
235             call xdrffloat_(ixdrf,refree,iret)
236             if (iret.eq.0) goto 101
237             call xdrffloat_(ixdrf,rmsdev,iret)
238             if (iret.eq.0) goto 101
239             call xdrfint_(ixdrf,iscor,iret)
240             if (iret.eq.0) goto 101
241 #else
242 c            write (iout,*) "calling xdrf3dfcoord"
243             call xdrf3dfcoord(ixdrf, csingle, itmp, prec, iret)
244 c            write (iout,*) "iret",iret
245 c            call flush(iout)
246             if (iret.eq.0) goto 101
247             call xdrfint(ixdrf, nss, iret)
248 c            write (iout,*) "iret",iret
249 c            write (iout,*) "nss",nss
250             call flush(iout)
251             if (iret.eq.0) goto 101
252             do k=1,nss
253            if (dyn_ss) then
254             call xdrfint(ixdrf, idssb(k), iret)
255             call xdrfint(ixdrf, jdssb(k), iret)
256             else
257               call xdrfint(ixdrf, ihpb(k), iret)
258               if (iret.eq.0) goto 101
259               call xdrfint(ixdrf, jhpb(k), iret)
260               if (iret.eq.0) goto 101
261             endif
262             enddo
263             call xdrffloat(ixdrf,reini,iret)
264             if (iret.eq.0) goto 101
265             call xdrffloat(ixdrf,refree,iret)
266             if (iret.eq.0) goto 101
267             call xdrffloat(ixdrf,rmsdev,iret)
268             if (iret.eq.0) goto 101
269             call xdrfint(ixdrf,iscor,iret)
270             if (iret.eq.0) goto 101
271 #endif
272             energy(jj+1)=reini
273             entfac(jj+1)=refree
274             rmstb(jj+1)=rmsdev
275             do k=1,nres
276               do l=1,3
277                 c(l,k)=csingle(l,k)
278               enddo
279             enddo
280             do k=nnt,nct
281               do l=1,3
282                 c(l,nres+k)=csingle(l,nres+k-nnt+1)
283               enddo
284             enddo
285           endif
286 C#define DEBUG
287 #ifdef DEBUG
288           write (iout,'(5hREAD ,i5,3f15.4,i10)') 
289      &     jj+1,energy(jj+1),entfac(jj+1),
290      &     rmstb(jj+1),iscor
291           write (iout,*) "Conformation",jjj+1,jj+1
292           write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
293           write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
294           call flush(iout)
295 #endif
296 C#undef DEBUG
297           call add_new_cconf(jjj,jj,jj_old,icount,Next)
298         enddo
299   101   continue
300         write (iout,*) i-1," conformations read from DA file ",
301      &    intinname(:ilen(intinname))
302         write (iout,*) jj," conformations read so far"
303         if (from_bx) then
304           close(intin)
305         else
306 #if (defined(AIX) && !defined(JUBL))
307           call xdrfclose_(ixdrf, iret)
308 #else
309           call xdrfclose(ixdrf, iret)
310 #endif
311         endif
312 #ifdef MPI
313 c#ifdef DEBUG   
314         write (iout,*) "jj_old",jj_old," jj",jj
315 c#endif
316         call write_and_send_cconf(icount,jj_old,jj,Next)
317         call MPI_Send(0,1,MPI_INTEGER,Next,570,
318      &             MPI_COMM_WORLD,IERROR)
319         jj_old=jj+1
320 #else
321         call write_and_send_cconf(icount,jj_old,jj,Next)
322 #endif
323         t_acq = tcpu() - t_acq
324 #ifdef MPI
325         write (iout,*) "Processor",me,
326      &    " time for conformation read/send",t_acq
327       ELSE
328 c A worker gets the confs from the master and sends them to its neighbor
329         t_acq = tcpu()
330         call receive_and_pass_cconf(icount,jj_old,jj,
331      &    Previous,Next)
332         t_acq = tcpu() - t_acq
333       ENDIF
334 #endif
335       ncon=jj
336 c      close(icbase)
337       close(intin)
338
339       write(iout,*)"A total of",ncon," conformations read."
340
341 #ifdef MPI
342 c Check if everyone has the same number of conformations
343       call MPI_Allgather(ncon,1,MPI_INTEGER,
344      &  ntot_all(0),1,MPI_INTEGER,MPI_Comm_World,IERROR)
345       lerr=.false.
346       do i=0,nprocs-1
347         if (i.ne.me) then
348             if (ncon.ne.ntot_all(i)) then
349               write (iout,*) "Number of conformations at processor",i,
350      &         " differs from that at processor",me,
351      &         ncon,ntot_all(i)
352               lerr = .true.
353             endif
354         endif
355       enddo
356       if (lerr) then
357         write (iout,*)
358         write (iout,*) "Number of conformations read by processors"
359         write (iout,*)
360         do i=0,nprocs-1
361           write (iout,'(8i10)') i,ntot_all(i)
362         enddo
363         write (iout,*) "Calculation terminated."
364         call flush(iout)
365         return1
366       endif
367       return
368 #endif
369  1111 write(iout,*) "Error opening coordinate file ",
370      & intinname(:ilen(intinname))
371       call flush(iout)
372       return1
373       end
374 c------------------------------------------------------------------------------
375       subroutine add_new_cconf(jjj,jj,jj_old,icount,Next)
376       implicit none
377       include "DIMENSIONS"
378       include "sizesclu.dat"
379       include "COMMON.CLUSTER"
380       include "COMMON.CONTROL"
381       include "COMMON.CHAIN"
382       include "COMMON.INTERACT"
383       include "COMMON.LOCAL"
384       include "COMMON.IOUNITS"
385       include "COMMON.NAMES"
386       include "COMMON.VAR"
387       include "COMMON.SBRIDGE"
388       include "COMMON.GEO"
389       integer i,j,jj,jjj,jj_old,icount,k,kk,l,ii,ib
390      &  nn,nn1,inan,Next,itj,chalen
391       double precision etot,energia(0:max_ene)
392       jjj=jjj+1
393       chalen=int((nct-nnt+2)/symetr)
394       call int_from_cart1(.false.)
395       do j=nnt+1,nct
396         if ((vbld(j).lt.2.0d0 .or. vbld(j).gt.6.0d0)
397      &      .and.(itype(j).ne.ntyp1)) then
398          if (j.gt.2) then
399           if (itel(j).ne.0 .and. itel(j-1).ne.0) then
400           write (iout,*) "Conformation",jjj,jj+1
401           write (iout,*) "Bad CA-CA bond length",j," ",vbld(j),itel(j),
402      & chalen
403           write (iout,*) "The Cartesian geometry is:"
404           write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
405           write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
406           write (iout,*) "The internal geometry is:"
407           write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
408           write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
409           write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
410           write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
411           write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
412           write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
413           write (iout,*) 
414      &      "This conformation WILL NOT be added to the database."
415           return
416           endif
417          endif
418         endif
419       enddo
420       do j=nnt,nct
421         itj=itype(j)
422         if (itype(j).ne.10 .and. (vbld(nres+j)-dsc(iabs(itj))).gt.5.0d0
423      &  .and. itype(j).ne.ntyp1) then
424           write (iout,*) "Conformation",jjj,jj+1
425           write (iout,*) "Bad CA-SC bond length",j," ",vbld(nres+j)
426           write (iout,*) "The Cartesian geometry is:"
427           write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
428           write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
429           write (iout,*) "The internal geometry is:"
430           write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
431           write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
432           write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
433           write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
434           write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
435           write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
436           write (iout,*) 
437      &      "This conformation WILL NOT be added to the database."
438           return
439         endif
440       enddo
441       do j=3,nres
442         if (theta(j).le.0.0d0) then
443           write (iout,*) 
444      &      "Zero theta angle(s) in conformation",jjj,jj+1
445           write (iout,*) "The Cartesian geometry is:"
446           write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
447           write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
448           write (iout,*) "The internal geometry is:"
449           write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
450           write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
451           write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
452           write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
453           write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
454           write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
455           write (iout,*)
456      &      "This conformation WILL NOT be added to the database."
457           return
458         endif
459         if (theta(j).gt.179.97*deg2rad) theta(j)=179.97*deg2rad
460       enddo
461       jj=jj+1
462 #ifdef DEBUG
463       write (iout,*) "Conformation",jjj,jj
464       write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
465       write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
466       write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
467       write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
468       write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
469       write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
470       write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct)
471       write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
472       write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
473       write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
474       write (iout,'(e15.5,16i5)') entfac(icount+1)
475 c     &        iscore(icount+1,0)
476 #endif
477       icount=icount+1
478       call store_cconf_from_file(jj,icount)
479       if (icount.eq.maxstr_proc) then
480 #ifdef DEBUG
481         write (iout,* ) "jj_old",jj_old," jj",jj
482 #endif
483         call write_and_send_cconf(icount,jj_old,jj,Next)
484         jj_old=jj+1
485         icount=0
486       endif
487       return
488       end
489 c------------------------------------------------------------------------------
490       subroutine store_cconf_from_file(jj,icount)
491       implicit none
492       include "DIMENSIONS"
493       include "sizesclu.dat"
494       include "COMMON.CLUSTER"
495       include "COMMON.CHAIN"
496       include "COMMON.SBRIDGE"
497       include "COMMON.INTERACT"
498       include "COMMON.IOUNITS"
499       include "COMMON.VAR"
500       integer i,j,jj,icount
501 c Store the conformation that has been read in
502       do i=1,2*nres
503         do j=1,3
504           allcart(j,i,icount)=c(j,i)
505         enddo
506       enddo
507       nss_all(icount)=nss
508       do i=1,nss
509         ihpb_all(i,icount)=ihpb(i)
510         jhpb_all(i,icount)=jhpb(i)
511       enddo
512       return
513       end
514 c------------------------------------------------------------------------------
515       subroutine write_and_send_cconf(icount,jj_old,jj,Next)
516       implicit none
517       include "DIMENSIONS"
518       include "sizesclu.dat"
519 #ifdef MPI
520       include "mpif.h"
521       integer IERROR
522       include "COMMON.MPI"
523 #endif
524       include "COMMON.CHAIN"
525       include "COMMON.SBRIDGE"
526       include "COMMON.INTERACT"
527       include "COMMON.IOUNITS"
528       include "COMMON.CLUSTER"
529       include "COMMON.VAR"
530       integer icount,jj_old,jj,Next
531 c Write the structures to a scratch file
532 #ifdef MPI
533 c Master sends the portion of conformations that have been read in to the neighbor
534 #ifdef DEBUG
535       write (iout,*) "Processor",me," entered WRITE_AND_SEND_CONF"
536       call flush(iout)
537 #endif
538       call MPI_Send(icount,1,MPI_INTEGER,Next,570,MPI_COMM_WORLD,IERROR)
539       call MPI_Send(nss_all(1),icount,MPI_INTEGER,
540      &    Next,571,MPI_COMM_WORLD,IERROR)
541       call MPI_Send(ihpb_all(1,1),icount,MPI_INTEGER,
542      &    Next,572,MPI_COMM_WORLD,IERROR)
543       call MPI_Send(jhpb_all(1,1),icount,MPI_INTEGER,
544      &    Next,573,MPI_COMM_WORLD,IERROR)
545       call MPI_Send(rmstb(jj_old),icount,MPI_DOUBLE_PRECISION,
546      &    Next,577,MPI_COMM_WORLD,IERROR)
547       call MPI_Send(entfac(jj_old),icount,MPI_DOUBLE_PRECISION,
548      &    Next,579,MPI_COMM_WORLD,IERROR)
549       call MPI_Send(allcart(1,1,1),3*icount*maxres2,
550      &    MPI_REAL,Next,580,MPI_COMM_WORLD,IERROR)
551 #endif
552       call dawrite_ccoords(jj_old,jj,icbase)
553       return
554       end
555 c------------------------------------------------------------------------------
556 #ifdef MPI
557       subroutine receive_and_pass_cconf(icount,jj_old,jj,Previous,
558      &  Next)
559       implicit none
560       include "DIMENSIONS"
561       include "sizesclu.dat"
562       include "mpif.h"
563       integer IERROR,STATUS(MPI_STATUS_SIZE)
564       include "COMMON.MPI"
565       include "COMMON.CHAIN"
566       include "COMMON.SBRIDGE"
567       include "COMMON.INTERACT"
568       include "COMMON.IOUNITS"
569       include "COMMON.VAR"
570       include "COMMON.GEO"
571       include "COMMON.CLUSTER"
572       integer i,j,k,icount,jj_old,jj,Previous,Next
573       icount=1
574 #ifdef DEBUG
575       write (iout,*) "Processor",me," entered RECEIVE_AND_PASS_CONF"
576       call flush(iout)
577 #endif
578       do while (icount.gt.0) 
579       call MPI_Recv(icount,1,MPI_INTEGER,Previous,570,MPI_COMM_WORLD,
580      &     STATUS,IERROR)
581       call MPI_Send(icount,1,MPI_INTEGER,Next,570,MPI_COMM_WORLD,
582      &     IERROR)
583 #ifdef DEBUG
584       write (iout,*) "Processor",me," icount",icount
585 #endif
586       if (icount.eq.0) return
587       call MPI_Recv(nss_all(1),icount,MPI_INTEGER,
588      &    Previous,571,MPI_COMM_WORLD,STATUS,IERROR)
589       call MPI_Send(nss_all(1),icount,MPI_INTEGER,
590      &  Next,571,MPI_COMM_WORLD,IERROR)
591       call MPI_Recv(ihpb_all(1,1),icount,MPI_INTEGER,
592      &    Previous,572,MPI_COMM_WORLD,STATUS,IERROR)
593       call MPI_Send(ihpb_all(1,1),icount,MPI_INTEGER,
594      &  Next,572,MPI_COMM_WORLD,IERROR)
595       call MPI_Recv(jhpb_all(1,1),icount,MPI_INTEGER,
596      &    Previous,573,MPI_COMM_WORLD,STATUS,IERROR)
597       call MPI_Send(jhpb_all(1,1),icount,MPI_INTEGER,
598      &  Next,573,MPI_COMM_WORLD,IERROR)
599       call MPI_Recv(rmstb(jj_old),icount,MPI_DOUBLE_PRECISION,
600      &  Previous,577,MPI_COMM_WORLD,STATUS,IERROR)
601       call MPI_Send(rmstb(jj_old),icount,MPI_DOUBLE_PRECISION,
602      &  Next,577,MPI_COMM_WORLD,IERROR)
603       call MPI_Recv(entfac(jj_old),icount,MPI_DOUBLE_PRECISION,
604      &  Previous,579,MPI_COMM_WORLD,STATUS,IERROR)
605       call MPI_Send(entfac(jj_old),icount,MPI_DOUBLE_PRECISION,
606      &  Next,579,MPI_COMM_WORLD,IERROR)
607       call MPI_Recv(allcart(1,1,1),3*icount*maxres2,
608      &  MPI_REAL,Previous,580,MPI_COMM_WORLD,STATUS,IERROR)
609       call MPI_Send(allcart(1,1,1),3*icount*maxres2,
610      &  MPI_REAL,Next,580,MPI_COMM_WORLD,IERROR)
611       jj=jj_old+icount-1
612       call dawrite_ccoords(jj_old,jj,icbase)
613       jj_old=jj+1
614 #ifdef DEBUG
615       write (iout,*) "Processor",me," received",icount," conformations"
616       do i=1,icount
617         write (iout,'(8f10.4)') (allcart(l,k,i),l=1,3,k=1,nres)
618         write (iout,'(8f10.4)')((allcart(l,k,i+nres),l=1,3,k=nnt,nct)
619         write (iout,'(e15.5,16i5)') entfac(i)
620       enddo
621 #endif
622       enddo
623       return
624       end
625 #endif
626 c------------------------------------------------------------------------------
627       subroutine daread_ccoords(istart_conf,iend_conf)
628       implicit none
629       include "DIMENSIONS"
630       include "sizesclu.dat"
631 #ifdef MPI
632       include "mpif.h"
633       include "COMMON.MPI"
634 #endif
635       include "COMMON.CHAIN"
636       include "COMMON.CLUSTER"
637       include "COMMON.IOUNITS"
638       include "COMMON.INTERACT"
639       include "COMMON.VAR"
640       include "COMMON.SBRIDGE"
641       include "COMMON.GEO"
642       integer istart_conf,iend_conf
643       integer i,j,ij,ii,iii
644       integer len
645       character*16 form,acc
646       character*80 nam
647 c
648 c Read conformations off a DA scratchfile.
649 c
650 C#define DEBUG
651 #ifdef DEBUG
652       write (iout,*) "DAREAD_COORDS"
653       write (iout,*) "istart_conf",istart_conf," iend_conf",iend_conf
654       inquire(unit=icbase,name=nam,recl=len,form=form,access=acc)
655       write (iout,*) "len=",len," form=",form," acc=",acc
656       write (iout,*) "nam=",nam
657       call flush(iout)
658 #endif
659       do ii=istart_conf,iend_conf
660         ij = ii - istart_conf + 1
661         iii=list_conf(ii)
662 #ifdef DEBUG
663         write (iout,*) "Reading binary file, record",iii," ii",ii
664         call flush(iout)
665 #endif
666         if (dyn_ss) then
667         read(icbase,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres),
668      &    ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),
669 c     &    nss_all(ij),(ihpb_all(i,ij),jhpb_all(i,ij),i=1,nss),
670      &    entfac(ii),rmstb(ii)
671         else
672         read(icbase,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres),
673      &    ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),
674      &    nss_all(ij),(ihpb_all(i,ij),jhpb_all(i,ij),i=1,nss),
675      &    entfac(ii),rmstb(ii)
676          endif
677 #ifdef DEBUG
678         write (iout,*) ii,iii,ij,entfac(ii)
679         write (iout,'(8f10.5)') ((allcart(j,i,ij),j=1,3),i=1,nres)
680         write (iout,'(8f10.4)') ((allcart(j,i,ij),j=1,3),
681      &    i=nnt+nres,nct+nres)
682         write (iout,'(2e15.5)') entfac(ij)
683         write (iout,'(16i5)') nss_all(ij),(ihpb_all(i,ij),
684      &    jhpb_all(i,ij),i=1,nss)
685         call flush(iout)
686 #endif
687 C#undef DEBUG
688       enddo
689       write (iout,*) "just before leave"
690       call flush(iout)
691       return
692       end
693 c------------------------------------------------------------------------------
694       subroutine dawrite_ccoords(istart_conf,iend_conf,unit_out)
695       implicit none
696       include "DIMENSIONS"
697       include "sizesclu.dat"
698 #ifdef MPI
699       include "mpif.h"
700       include "COMMON.MPI"
701 #endif
702       include "COMMON.CHAIN"
703       include "COMMON.INTERACT"
704       include "COMMON.IOUNITS"
705       include "COMMON.VAR"
706       include "COMMON.SBRIDGE"
707       include "COMMON.GEO"
708       include "COMMON.CLUSTER"
709       integer istart_conf,iend_conf
710       integer i,j,ii,ij,iii,unit_out
711       integer len
712       character*16 form,acc
713       character*32 nam
714 c
715 c Write conformations to a DA scratchfile.
716 c
717 #ifdef DEBUG
718       write (iout,*) "DAWRITE_COORDS"
719       write (iout,*) "istart_conf",istart_conf," iend_conf",iend_conf
720       write (iout,*) "lenrec",lenrec
721       inquire(unit=unit_out,name=nam,recl=len,form=form,access=acc)
722       write (iout,*) "len=",len," form=",form," acc=",acc
723       write (iout,*) "nam=",nam
724       call flush(iout)
725 #endif
726       do ii=istart_conf,iend_conf
727         iii=list_conf(ii)
728         ij = ii - istart_conf + 1
729 #ifdef DEBUG
730         write (iout,*) "Writing binary file, record",iii," ii",ii
731         call flush(iout)
732 #endif
733        if (dyn_ss) then
734         write(unit_out,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres),
735      &    ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),
736 c     &    nss_all(ij),(ihpb_all(i,ij),jhpb_all(i,ij),i=1,nss_all(ij))
737      &    entfac(ii),rmstb(ii)
738         else
739         write(unit_out,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres),
740      &    ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),
741      &    nss_all(ij),(ihpb_all(i,ij),jhpb_all(i,ij),i=1,nss_all(ij)),
742      &    entfac(ii),rmstb(ii)
743        endif
744 #ifdef DEBUG
745         write (iout,'(8f10.5)') ((allcart(j,i,ij),j=1,3),i=1,nres)
746         write (iout,'(8f10.4)') ((allcart(j,i,ij),j=1,3),i=nnt+nres,
747      &   nct+nres)
748         write (iout,'(2e15.5)') entfac(ij)
749         write (iout,'(16i5)') nss_all(ij),(ihpb(i,ij),jhpb(i,ij),i=1,
750      &   nss_all(ij))
751         call flush(iout)
752 #endif
753       enddo
754       return
755       end