2 !-----------------------------------------------------------------------------
6 use io_base !, only: ilen
7 use geometry_data, only: nres,c,nres_molec,boxxsize,boxysize,boxzsize
8 use energy_data, only: nnt,nct,nss
9 use control_data, only: lside
11 !-----------------------------------------------------------------------------
14 !-----------------------------------------------------------------------------
16 !-----------------------------------------------------------------------------
18 !-----------------------------------------------------------------------------
19 SUBROUTINE WRTCLUST(NCON,ICUT,PRINTANG,PRINTPDB,printmol2,ib)
21 use hc_, only: ioffset
22 use control_data, only: lprint_cart,lprint_int,titel
23 use geometry, only: int_from_cart1,nres
24 ! implicit real*8 (a-h,o-z)
25 ! include 'DIMENSIONS'
26 ! include 'sizesclu.dat'
27 integer,parameter :: num_in_line=5
28 LOGICAL :: PRINTANG(max_cut)
29 integer :: PRINTPDB(max_cut),printmol2(max_cut)
30 ! include 'COMMON.CONTROL'
31 ! include 'COMMON.HEADER'
32 ! include 'COMMON.CHAIN'
33 ! include 'COMMON.VAR'
34 ! include 'COMMON.CLUSTER'
35 ! include 'COMMON.IOUNITS'
36 ! include 'COMMON.GEO'
37 ! include 'COMMON.FREE'
38 ! include 'COMMON.TEMPFAC'
39 real(kind=8) :: rmsave(maxgr)
40 CHARACTER(len=64) :: prefixp,NUMM,MUMM,EXTEN,extmol
41 character(len=80) :: cfname
42 character(len=8) :: ctemper
43 DATA EXTEN /'.pdb'/,extmol /'.mol2'/,NUMM /'000'/,&
46 integer :: ncon,icut,ib
47 integer :: i,ii,ii1,ii2,igr,ind1,ind2,ico,icon,&
48 irecord,nrecord,j,k,jj,ind,ncon_lim,ncon_out
49 real(kind=8) :: temper,curr_dist,emin,qpart,boltz,&
50 ave_dim,amax_dim,emin1
53 allocate(tempfac(2,nres))
58 ! print *,"calling WRTCLUST",ncon
59 ! write (iout,*) "ICUT",icut," PRINTPDB ",PRINTPDB(icut)
62 temper=1.0d0/(beta_h(ib)*1.987d-3)
63 if (temper.lt.100.0d0) then
64 write(ctemper,'(f3.0)') temper
66 else if (temper.lt.1000.0) then
67 write (ctemper,'(f4.0)') temper
70 write (ctemper,'(f5.0)') temper
74 do i=1,ncon*(ncon-1)/2
77 close(80,status='delete')
79 ! PRINT OUT THE RESULTS OF CLUSTER ANALYSIS
81 ii1= index(intinname,'/')
86 ii2=index(intinname(ii1:),'/')
88 ii = ii1+index(intinname(ii1:),'.')-1
94 prefixp=intinname(ii1:ii)
95 !d print *,icut,printang(icut),printpdb(icut),printmol2(icut)
96 !d print *,'ecut=',ecut
99 WRITE (iout,200) IGR,totfree_gr(igr)/beta_h(ib),LICZ(IGR)
100 NRECORD=LICZ(IGR)/num_in_line
102 DO 63 IRECORD=1,NRECORD
103 IND2=IND1+num_in_line-1
104 WRITE (iout,300) (list_conf(NCONF(IGR,ICO)),&
105 totfree(NCONF(IGR,ICO))/beta_h(ib),ICO=IND1,IND2)
108 WRITE (iout,300) (list_conf(NCONF(IGR,ICO)),&
109 totfree(NCONF(IGR,ICO))/beta_h(ib),ICO=IND1,LICZ(IGR))
111 ICON=list_conf(NCONF(IGR,1))
112 ! WRITE (iout,'(16F5.0)') (rad2deg*phiall(IND,icon),IND=4,nphi+3)
113 ! 12/8/93 Estimation of "diameters" of the subsequent families.
116 ! write (iout,*) "ecut",ecut
119 if (totfree(ii)-emin .gt. ecut) goto 10
124 ind=ioffset(ncon,ii,jj)
126 ind=ioffset(ncon,jj,ii)
128 ! write (iout,*) " ncon",ncon,"i",i," j",j," ii",ii," jj",jj,
131 curr_dist=dabs(diss(ind)+0.0d0)
132 ! write(iout,'(i10,4i4,f12.4)') ind,ii,jj,list_conf(ii),
133 ! & list_conf(jj),curr_dist
134 if (curr_dist .gt. amax_dim) amax_dim=curr_dist
135 ave_dim=ave_dim+curr_dist**2
138 10 if (licz(igr) .gt. 1) &
139 ave_dim=sqrt(ave_dim/(licz(igr)*(licz(igr)-1)/2))
140 write (iout,'(/A,F8.1,A,F8.1)') &
141 'Max. distance in the family:',amax_dim,&
142 '; average distance in the family:',ave_dim
147 boltz=dexp(-totfree(icon))
148 rmsave(igr)=rmsave(igr)+boltz*rmstb(icon)
151 rmsave(igr)=rmsave(igr)/qpart
152 write (iout,'(a,f5.2,a)') "Average RMSD",rmsave(igr)," A"
155 WRITE (iout,500) (list_conf(I),IASS(I),I=1,NCON)
156 ! print *,icut,printang(icut)
157 IF (PRINTANG(ICUT) .and. (lprint_cart .or. lprint_int)) then
159 ! print *,'emin',emin,' ngr',ngr
160 if (lprint_cart) then
161 cfname=prefixp(:ilen(prefixp))//"_T"//ctemper(:ilen(ctemper))&
164 cfname=prefixp(:ilen(prefixp))//"_T"//ctemper(:ilen(ctemper))&
169 if (totfree_gr(igr)-emin.le.ecut) then
170 if (lprint_cart) then
171 call cartout(igr,icon,totfree(icon)/beta_h(ib),&
172 totfree_gr(igr)/beta_h(ib),&
175 ! print '(a)','calling briefout'
178 c(j,i)=allcart(j,i,icon)
181 call int_from_cart1(.false.)
182 !el call briefout(igr,iscore(icon),totfree(icon)/beta_h(ib),&
183 !el totfree_gr(igr),nss_all(icon),ihpb_all(1,icon),&
184 !el jhpb_all(1,icon),cfname)
185 call briefout(igr,totfree(icon)/beta_h(ib),&
187 ! print '(a)','exit briefout'
193 IF (PRINTPDB(ICUT).gt.0) THEN
194 ! Write out a number of conformations from each family in PDB format and
195 ! create InsightII command file for their displaying in different colors
196 cfname=prefixp(:ilen(prefixp))//"_T"//ctemper(:ilen(ctemper))&
198 write (iout,*) "cfname",cfname
199 OPEN(ipdb,FILE=CFNAME,STATUS='UNKNOWN',FORM='FORMATTED')
200 write (ipdb,'(a,f8.2)') &
201 "REMAR AVERAGE CONFORMATIONS AT TEMPERATURE",temper
207 DO WHILE(I.LE.NGR .AND. totfree_gr(i)-EMIN.LE.ECUT)
208 ! write (iout,*) "i",i," ngr",ngr,totfree_gr(I),EMIN,ecut
209 write (NUMM,'(bz,i4.4)') i
210 ncon_lim=min0(licz(i),printpdb(icut))
211 cfname=prefixp(:ilen(prefixp))//"_T"//ctemper(:ilen(ctemper))&
212 //"K_"//numm(:ilen(numm))//exten
213 OPEN(ipdb,FILE=CFNAME,STATUS='UNKNOWN',FORM='FORMATTED')
214 write (ipdb,'("REMARK CLUSTER",i5," FREE ENERGY",1pe14.5," AVE RMSD",0pf5.2)')&
215 i,totfree_gr(i)/beta_h(ib),rmsave(i) !'
216 ! Write conformations of the family i to PDB files
218 do while (ncon_out.lt.printpdb(icut) .and. &
219 ncon_out.lt.licz(i).and. &
220 totfree(nconf(i,ncon_out+1))-EMIN1.LE.ECUT)
222 ! write (iout,*) i,ncon_out,nconf(i,ncon_out),
223 ! & totfree(nconf(i,ncon_out)),emin1,ecut
225 write (iout,*) "ncon_out",ncon_out
235 c(k,ii)=allcart(k,ii,icon)
238 call pdboutC(totfree(icon)/beta_h(ib),rmstb(icon),titel)
239 write (ipdb,'("TER")')
242 ! Average structures and structures closest to average
243 cfname=prefixp(:ilen(prefixp))//"_T"//ctemper(:ilen(ctemper))&
245 OPEN(ipdb,FILE=CFNAME,STATUS='UNKNOWN',FORM='FORMATTED',&
248 write (ipdb,'(a,i5)') "REMARK CLUSTER",i
249 call pdboutC(totfree_gr(i)/beta_h(ib),rmsave(i),titel)
250 write (ipdb,'("TER")')
251 call closest_coord(i)
252 call pdboutC(totfree_gr(i)/beta_h(ib),rmsave(i),titel)
253 write (ipdb,'("TER")')
260 IF (printmol2(icut).gt.0) THEN
261 ! Write out a number of conformations from each family in PDB format and
262 ! create InsightII command file for their displaying in different colors
267 DO WHILE(I.LE.NGR .AND. totfree_gr(i)-EMIN.LE.ECUT)
268 write (NUMM,'(bz,i4.4)') i
269 cfname=prefixp(:ilen(prefixp))//"_T"//ctemper(:ilen(ctemper))&
270 //"K_"//numm(:ilen(numm))//extmol
271 OPEN(imol2,FILE=CFNAME,STATUS='UNKNOWN',FORM='FORMATTED')
273 do while (ncon_out.lt.printmol2(icut) .and. &
274 ncon_out.lt.licz(i).and. &
275 totfree(nconf(i,ncon_out+1))-EMIN1.LE.ECUT)
282 c(k,ii)=allcart(k,ii,icon)
285 CALL MOL2OUT(totfree(icon)/beta_h(ib),'STRUCTURE'//numm)
293 100 FORMAT (//'THERE ARE ',I4,' FAMILIES OF CONFORMATIONS')
294 200 FORMAT (/'FAMILY ',I4,' WITH TOTAL FREE ENERGY',1pE15.5,&
295 ' CONTAINS ',I4,' CONFORMATION(S): ')
296 ! 300 FORMAT ( 8(I4,F6.1))
297 300 FORMAT (5(I4,1pe12.3))
298 400 FORMAT (//'ASSIGNMENT OF CONSECUTIVE CONFORMATIONS TO FAMILIES:')
299 500 FORMAT (8(2I4,2X))
300 600 FORMAT ('REMARK FAMILY',I4,' CONFORMATION',I4,' ENERGY ',E15.6)
302 END SUBROUTINE WRTCLUST
303 !------------------------------------------------------------------------------
304 subroutine ave_coord(igr)
306 use control_data, only:lside
307 use regularize_, only:fitsq,matvec
309 ! include 'DIMENSIONS'
310 ! include 'sizesclu.dat'
311 ! include 'COMMON.CONTROL'
312 ! include 'COMMON.CLUSTER'
313 ! include 'COMMON.CHAIN'
314 ! include 'COMMON.INTERACT'
315 ! include 'COMMON.VAR'
316 ! include 'COMMON.TEMPFAC'
317 ! include 'COMMON.IOUNITS'
319 real(kind=8) :: przes(3),obrot(3,3)
320 real(kind=8) :: xx(3,2*nres),yy(3,2*nres),csq(3,2*nres) !(3,maxres2)
322 integer :: i,ii,j,k,icon,jcon,igr
323 real(kind=8) :: rms,boltz,qpart,cwork(3,2*nres),cref1(3,2*nres) !(3,maxres2)
324 ! write (iout,*) "AVE_COORD: igr",igr
327 boltz = dexp(-totfree(jcon)+eref)
331 c(j,i)=allcart(j,i,jcon)*boltz
332 cref1(j,i)=allcart(j,i,jcon)
333 csq(j,i)=allcart(j,i,jcon)**2*boltz
343 xx(j,ii)=allcart(j,i,jcon)
348 ! if (itype(i).ne.10) then
351 xx(j,ii)=allcart(j,i+nres,jcon)
352 yy(j,ii)=cref1(j,i+nres)
356 call fitsq(rms,xx(1,1),yy(1,1),ii,przes,obrot,non_conv)
360 cwork(j,i)=allcart(j,i,jcon)
363 call fitsq(rms,cwork(1,nnt),cref1(1,nnt),nct-nnt+1,przes,obrot &
366 ! write (iout,*) "rms",rms
368 ! write (iout,'(i3,f10.5,5x,3f10.5)')i,przes(i),(obrot(i,j),j=1,3)
371 print *,'error, rms^2 = ',rms,icon,jcon
374 if (non_conv) print *,non_conv,icon,jcon
375 boltz=dexp(-totfree(jcon)+eref)
376 qpart = qpart + boltz
379 xx(j,i)=allcart(j,i,jcon)
382 call matvec(cwork,obrot,xx,2*nres)
384 ! write (iout,'(i5,2(3f10.5,5x))') i,(cwork(j,i),j=1,3),
385 ! & (allcart(j,i,jcon),j=1,3)
387 cwork(j,i)=cwork(j,i)+przes(j)
388 c(j,i)=c(j,i)+cwork(j,i)*boltz
389 csq(j,i)=csq(j,i)+cwork(j,i)**2*boltz
396 csq(j,i)=csq(j,i)/qpart-c(j,i)**2
398 ! write (iout,'(i5,3f10.5)') i,(csq(j,i),j=1,3)
404 tempfac(1,i)=tempfac(1,i)+csq(j,i)
405 tempfac(2,i)=tempfac(2,i)+csq(j,i+nres)
407 tempfac(1,i)=dsqrt(tempfac(1,i))
408 tempfac(2,i)=dsqrt(tempfac(2,i))
411 end subroutine ave_coord
412 !------------------------------------------------------------------------------
413 subroutine closest_coord(igr)
415 use regularize_, only:fitsq
417 ! include 'DIMENSIONS'
418 ! include 'sizesclu.dat'
419 ! include 'COMMON.IOUNITS'
420 ! include 'COMMON.CONTROL'
421 ! include 'COMMON.CLUSTER'
422 ! include 'COMMON.CHAIN'
423 ! include 'COMMON.INTERACT'
424 ! include 'COMMON.VAR'
426 real(kind=8) :: przes(3),obrot(3,3)
427 real(kind=8) :: xx(3,2*nres),yy(3,2*nres) !(3,maxres2)
428 integer :: i,ii,j,k,icon,jcon,jconmin,igr
429 real(kind=8) :: rms,rmsmin,cwork(3,2*nres)
439 xx(j,ii)=allcart(j,i,jcon)
444 ! if (itype(i).ne.10) then
447 xx(j,ii)=allcart(j,i+nres,jcon)
452 call fitsq(rms,xx(1,1),yy(1,1),ii,przes,obrot,non_conv)
456 cwork(j,i)=allcart(j,i,jcon)
459 call fitsq(rms,cwork(1,nnt),c(1,nnt),nct-nnt+1,przes,obrot,&
463 print *,'error, rms^2 = ',rms,icon,jcon
466 ! write (iout,*) "jcon",jcon," rms",rms," rmsmin",rmsmin
467 if (non_conv) print *,non_conv,icon,jcon
468 if (rms.lt.rmsmin) then
473 ! write (iout,*) "rmsmin",rmsmin," rms",rms
477 c(j,i)=allcart(j,i,jconmin)
481 end subroutine closest_coord
482 !-----------------------------------------------------------------------------
484 !-----------------------------------------------------------------------------
485 subroutine read_coords(ncon,*)
487 use energy_data, only: ihpb,jhpb,max_ene
488 use control_data, only: from_bx,from_cx
489 use control, only: tcpu
491 ! include "DIMENSIONS"
492 ! include "sizesclu.dat"
496 integer :: IERROR,ERRCODE !,STATUS(MPI_STATUS_SIZE)
497 ! include "COMMON.MPI"
499 use MPI_data, only: nprocs
501 ! include "COMMON.CONTROL"
502 ! include "COMMON.CHAIN"
503 ! include "COMMON.INTERACT"
504 ! include "COMMON.IOUNITS"
505 ! include "COMMON.VAR"
506 ! include "COMMON.SBRIDGE"
507 ! include "COMMON.GEO"
508 ! include "COMMON.CLUSTER"
509 character(len=3) :: liczba
511 integer :: i,j,jj,jjj,jj_old,icount,k,kk,l,ii,if,ib,&
513 integer :: ixdrf,iret,itmp
514 real(kind=4) :: prec,reini,refree,rmsdev
515 integer :: nrec,nlines,iscor,lenrec,lenrec_in
516 real(kind=8) :: energ,t_acq !,tcpu
517 !el integer ilen,iroof
518 !el external ilen,iroof
519 real(kind=8) :: rjunk
520 integer :: ntot_all(0:nprocs-1) !(0:maxprocs-1)
522 real(kind=8) :: energia(0:max_ene),etot
523 real(kind=4) :: csingle(3,2*nres+2)
524 integer :: Previous,Next
525 character(len=256) :: bprotfiles
526 ! print *,"Processor",me," calls read_protein_data"
528 if (me.eq.master) then
529 Previous=MPI_PROC_NULL
533 if (me.eq.nprocs-1) then
538 ! Set the scratchfile names
539 write (liczba,'(bz,i3.3)') me
541 allocate(STATUS(MPI_STATUS_SIZE))
543 ! 1/27/05 AL Change stored coordinates to single precision and don't store
544 ! energy components in the binary databases.
545 lenrec=12*(nres+nct-nnt+1)+4*(2*nss+2)+16
546 lenrec_in=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
548 write (iout,*) "nres",nres," nnt",nnt," nct",nct," nss", nss
549 write (iout,*) "lenrec_in",lenrec_in
551 bprotfiles=scratchdir(:ilen(scratchdir))// &
552 "/"//prefix(:ilen(prefix))//liczba//".xbin"
554 ! allocate cluster arrays
555 allocate(energy(0:maxconf),totfree(0:maxconf)) !(0:maxconf)
556 allocate(entfac(maxconf)) !(maxconf)
557 allocate(rmstb(maxconf)) !(maxconf)
558 allocate(allcart(3,2*nres,maxstr_proc)) !(3,maxres2,maxstr_proc)
559 allocate(nss_all(maxstr_proc)) !(maxstr_proc)
560 allocate(ihpb_all(maxss,maxstr_proc),jhpb_all(maxss,maxstr_proc))!(maxss,maxstr_proc)
561 allocate(iscore(maxconf)) !(maxconf)
569 if (from_cart .and. .not. from_bx .and. .not. from_cx) then
571 read (intin,*,end=13,err=11) energy(icon),totfree(icon),&
573 nss_all(icon),(ihpb_all(ii,icon),jhpb_all(i,icon),&
574 i=1,nss_all(icon)),iscore(icon)
576 read (intin,*,end=13,err=11) energy(icon),rmstb(icon),&
577 nss_all(icon),(ihpb_all(ii,icon),jhpb_all(i,icon),&
578 i=1,nss_all(icon)),iscore(icon)
580 read (intin,'(8f10.5)',end=13,err=10) &
581 ((allcart(j,i,icon),j=1,3),i=1,nres),&
582 ((allcart(j,i+nres,icon),j=1,3),i=nnt,nct)
583 print *,icon,energy(icon),nss_all(icon),rmstb(icon)
585 read(intin,'(a80)',end=13,err=12) lineh
586 read(lineh(:5),*,err=8) ic
588 read(lineh(6:),*,err=8) energy(icon)
590 read(lineh(6:),*,err=8) energy(icon)
594 print *,'error, assuming e=1d10',lineh
598 !old read(lineh(18:),*,end=13,err=11) nss_all(icon)
599 ii = index(lineh(15:)," ")+15
600 read(lineh(ii:),*,end=13,err=11) nss_all(icon)
601 IF (NSS_all(icon).LT.9) THEN
602 read (lineh(20:),*,end=102) &
603 (IHPB_all(I,icon),JHPB_all(I,icon),I=1,NSS_all(icon)),&
606 read (lineh(20:),*,end=102) &
607 (IHPB_all(I,icon),JHPB_all(I,icon),I=1,8)
608 read (intin,*) (IHPB_all(I,icon),JHPB_all(I,icon),&
609 I=9,NSS_all(icon)),iscore(icon)
614 PRINT *,'IC:',IC,' ENERGY:',ENERGY(ICON)
615 call read_angles(intin,*13)
617 phiall(i,icon)=phi(i)
618 thetall(i,icon)=theta(i)
619 alphall(i,icon)=alph(i)
620 omall(i,icon)=omeg(i)
626 ! CALCULATE DISTANCES
628 10 print *,'something wrong with angles'
630 11 print *,'something wrong with NSS',nss
632 12 print *,'something wrong with header'
639 open (icbase,file=bprotfiles,status="unknown",&
640 form="unformatted",access="direct",recl=lenrec)
641 ! Read conformations from binary DA files (one per batch) and write them to
642 ! a binary DA scratchfile.
646 write (liczba,'(bz,i3.3)') me
647 IF (ME.EQ.MASTER) THEN
648 ! Only the master reads the database; it'll send it to the other procs
656 open (intin,file=intinname,status="old",form="unformatted",&
657 access="direct",recl=lenrec_in)
659 else if (from_cx) then
660 #if (defined(AIX) && !defined(JUBL))
661 call xdrfopen_(ixdrf,intinname, "r", iret)
663 call xdrfopen(ixdrf,intinname, "r", iret)
666 write (iout,*) "xdrfopen: iret",iret
668 write (iout,*) "Error: coordinate file ",&
669 intinname(:ilen(intinname))," does not exist."
672 call MPI_ABORT(MPI_COMM_WORLD,IERROR,ERRCODE)
677 write (iout,*) "Error: coordinate format not specified"
680 call MPI_ABORT(MPI_COMM_WORLD,IERROR,ERRCODE)
688 write (iout,*) "Opening file ",intinname(:ilen(intinname))
689 write (iout,*) "lenrec",lenrec_in
693 ! write (iout,*) "maxconf",maxconf
697 !el if (i.gt.maxconf) then
698 !el write (iout,*) "Error: too many conformations ",&
699 !el "(",maxconf,") maximum."
701 !el call MPI_Abort(MPI_COMM_WORLD,errcode,ierror)
705 ! write (iout,*) "i",i
708 read(intin,err=101,end=101) &
709 ((csingle(l,k),l=1,3),k=1,nres),&
710 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
711 nss,(ihpb(k),jhpb(k),k=1,nss),&
713 entfac(jj+1),rmstb(jj+1),iscor
720 #if (defined(AIX) && !defined(JUBL))
721 call xdrf3dfcoord_(ixdrf, csingle, itmp, prec, iret)
722 if (iret.eq.0) goto 101
723 call xdrfint_(ixdrf, nss, iret)
724 if (iret.eq.0) goto 101
726 call xdrfint_(ixdrf, ihpb(j), iret)
727 if (iret.eq.0) goto 101
728 call xdrfint_(ixdrf, jhpb(j), iret)
729 if (iret.eq.0) goto 101
731 call xdrffloat_(ixdrf,reini,iret)
732 if (iret.eq.0) goto 101
733 call xdrffloat_(ixdrf,refree,iret)
734 if (iret.eq.0) goto 101
735 call xdrffloat_(ixdrf,rmsdev,iret)
736 if (iret.eq.0) goto 101
737 call xdrfint_(ixdrf,iscor,iret)
738 if (iret.eq.0) goto 101
740 ! write (iout,*) "calling xdrf3dfcoord"
741 call xdrf3dfcoord(ixdrf, csingle, itmp, prec, iret)
742 ! write (iout,*) "iret",iret
744 if (iret.eq.0) goto 101
745 call xdrfint(ixdrf, nss, iret)
746 write (iout,*) "iret",iret
747 ! write (iout,*) "nss",nss,i,"TUTU"
749 if (iret.eq.0) goto 101
751 call xdrfint(ixdrf, ihpb(k), iret)
752 if (iret.eq.0) goto 101
753 call xdrfint(ixdrf, jhpb(k), iret)
754 if (iret.eq.0) goto 101
755 ! write(iout,*), ihpb(k),jhpb(k),"TUTU"
757 call xdrffloat(ixdrf,reini,iret)
758 if (iret.eq.0) goto 101
759 ! write(iout,*) "TUTU", reini
760 call xdrffloat(ixdrf,refree,iret)
761 ! write(iout,*) "TUTU", refree
762 if (iret.eq.0) goto 101
763 call xdrffloat(ixdrf,rmsdev,iret)
764 if (iret.eq.0) goto 101
765 call xdrfint(ixdrf,iscor,iret)
766 if (iret.eq.0) goto 101
778 c(l,nres+k)=csingle(l,nres+k-nnt+1)
783 write (iout,'(5hREAD ,i5,3f15.4,i10)') &
784 jj+1,energy(jj+1),entfac(jj+1),&
786 write (iout,*) "Conformation",jjj+1,jj+1
787 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
788 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
791 call add_new_cconf(jjj,jj,jj_old,icount,Next)
794 write (iout,*) i-1," conformations read from DA file ",&
795 intinname(:ilen(intinname))
796 write (iout,*) jj," conformations read so far"
800 #if (defined(AIX) && !defined(JUBL))
801 call xdrfclose_(ixdrf, iret)
803 call xdrfclose(ixdrf, iret)
808 write (iout,*) "jj_old",jj_old," jj",jj
810 call write_and_send_cconf(icount,jj_old,jj,Next)
811 call MPI_Send(0,1,MPI_INTEGER,Next,570,&
812 MPI_COMM_WORLD,IERROR)
815 call write_and_send_cconf(icount,jj_old,jj,Next)
817 t_acq = tcpu() - t_acq
819 write (iout,*) "Processor",me,&
820 " time for conformation read/send",t_acq
822 ! A worker gets the confs from the master and sends them to its neighbor
824 call receive_and_pass_cconf(icount,jj_old,jj,&
826 t_acq = tcpu() - t_acq
833 write(iout,*)"A total of",ncon," conformations read."
835 allocate(enetb(1:max_ene,ncon)) !(max_ene,maxstr_proc)
837 ! Check if everyone has the same number of conformations
838 call MPI_Allgather(ncon,1,MPI_INTEGER,&
839 ntot_all(0),1,MPI_INTEGER,MPI_Comm_World,IERROR)
843 if (ncon.ne.ntot_all(i)) then
844 write (iout,*) "Number of conformations at processor",i,&
845 " differs from that at processor",me,&
853 write (iout,*) "Number of conformations read by processors"
856 write (iout,'(8i10)') i,ntot_all(i)
858 write (iout,*) "Calculation terminated."
864 1111 write(iout,*) "Error opening coordinate file ",&
865 intinname(:ilen(intinname))
868 end subroutine read_coords
869 !------------------------------------------------------------------------------
870 subroutine add_new_cconf(jjj,jj,jj_old,icount,Next)
872 use geometry_data, only: vbld,rad2deg,theta,phi,alph,omeg,deg2rad
873 use energy_data, only: itel,itype,dsc,max_ene,molnum
874 use control_data, only: symetr
875 use geometry, only: int_from_cart1
877 ! include "DIMENSIONS"
878 ! include "sizesclu.dat"
879 ! include "COMMON.CLUSTER"
880 ! include "COMMON.CONTROL"
881 ! include "COMMON.CHAIN"
882 ! include "COMMON.INTERACT"
883 ! include "COMMON.LOCAL"
884 ! include "COMMON.IOUNITS"
885 ! include "COMMON.NAMES"
886 ! include "COMMON.VAR"
887 ! include "COMMON.SBRIDGE"
888 ! include "COMMON.GEO"
889 integer :: i,j,jj,jjj,jj_old,icount,k,kk,l,ii,ib,&
890 nn,nn1,inan,Next,itj,chalen,mnum
891 real(kind=8) :: etot,energia(0:max_ene)
893 chalen=int((nct-nnt+2)/symetr)
894 call int_from_cart1(.false.)
897 ! write (iout,*) "Check atom",j
899 if (itype(j,mnum).eq.ntyp1_molec(mnum)) cycle
900 if (itype(j+1,molnum(j+1)).eq.ntyp1_molec(molnum(j+1))) cycle
902 if ((vbld(j).lt.2.0d0 .or. vbld(j).gt.5.0d0).and.(mnum.ne.5)) then
904 if (itel(j).ne.0 .and. itel(j-1).ne.0) then
905 write (iout,*) "Conformation",jjj,jj+1
906 write (iout,*) "Bad CA-CA bond length",j," ",vbld(j),itel(j),&
908 write (iout,*) "The Cartesian geometry is:"
909 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
910 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
911 write (iout,*) "The internal geometry is:"
912 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
913 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
914 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
915 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
916 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
917 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
919 "This conformation WILL NOT be added to the database."
929 if (itype(j,1).ne.10 .and. (itype(j,mnum).ne.ntyp1_molec(mnum)) &
930 .and. (vbld(nres+j)-dsc(iabs(itj))) &
932 write (iout,*) "Conformation",jjj,jj+1
933 write (iout,*) "Bad CA-SC bond length",j," ",vbld(nres+j)
934 write (iout,*) "The Cartesian geometry is:"
935 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
936 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
937 write (iout,*) "The internal geometry is:"
938 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
939 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
940 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
941 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
942 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
943 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
945 "This conformation WILL NOT be added to the database."
950 if (theta(j).le.0.0d0) then
952 "Zero theta angle(s) in conformation",jjj,jj+1
953 write (iout,*) "The Cartesian geometry is:"
954 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
955 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
956 write (iout,*) "The internal geometry is:"
957 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
958 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
959 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
960 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
961 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
962 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
964 "This conformation WILL NOT be added to the database."
967 if (theta(j).gt.179.97*deg2rad) theta(j)=179.97*deg2rad
971 write (iout,*) "Conformation",jjj,jj
972 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
973 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
974 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
975 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
976 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
977 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
978 write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct)
979 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
980 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
981 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
982 write (iout,'(e15.5,16i5)') entfac(icount+1)
983 ! & iscore(icount+1,0)
986 call store_cconf_from_file(jj,icount)
987 if (icount.eq.maxstr_proc) then
989 write (iout,* ) "jj_old",jj_old," jj",jj
991 call write_and_send_cconf(icount,jj_old,jj,Next)
996 end subroutine add_new_cconf
997 !------------------------------------------------------------------------------
998 subroutine store_cconf_from_file(jj,icount)
1000 use energy_data, only: ihpb,jhpb
1002 ! include "DIMENSIONS"
1003 ! include "sizesclu.dat"
1004 ! include "COMMON.CLUSTER"
1005 ! include "COMMON.CHAIN"
1006 ! include "COMMON.SBRIDGE"
1007 ! include "COMMON.INTERACT"
1008 ! include "COMMON.IOUNITS"
1009 ! include "COMMON.VAR"
1010 integer :: i,j,jj,icount
1011 ! Store the conformation that has been read in
1014 allcart(j,i,icount)=c(j,i)
1019 ihpb_all(i,icount)=ihpb(i)
1020 jhpb_all(i,icount)=jhpb(i)
1023 end subroutine store_cconf_from_file
1024 !------------------------------------------------------------------------------
1025 subroutine write_and_send_cconf(icount,jj_old,jj,Next)
1028 ! include "DIMENSIONS"
1029 ! include "sizesclu.dat"
1034 ! include "COMMON.MPI"
1036 ! include "COMMON.CHAIN"
1037 ! include "COMMON.SBRIDGE"
1038 ! include "COMMON.INTERACT"
1039 ! include "COMMON.IOUNITS"
1040 ! include "COMMON.CLUSTER"
1041 ! include "COMMON.VAR"
1042 integer :: icount,jj_old,jj,Next
1043 ! Write the structures to a scratch file
1045 ! Master sends the portion of conformations that have been read in to the neighbor
1047 write (iout,*) "Processor",me," entered WRITE_AND_SEND_CONF"
1050 call MPI_Send(icount,1,MPI_INTEGER,Next,570,MPI_COMM_WORLD,IERROR)
1051 call MPI_Send(nss_all(1),icount,MPI_INTEGER,&
1052 Next,571,MPI_COMM_WORLD,IERROR)
1053 call MPI_Send(ihpb_all(1,1),icount,MPI_INTEGER,&
1054 Next,572,MPI_COMM_WORLD,IERROR)
1055 call MPI_Send(jhpb_all(1,1),icount,MPI_INTEGER,&
1056 Next,573,MPI_COMM_WORLD,IERROR)
1057 call MPI_Send(rmstb(jj_old),icount,MPI_DOUBLE_PRECISION,&
1058 Next,577,MPI_COMM_WORLD,IERROR)
1059 call MPI_Send(entfac(jj_old),icount,MPI_DOUBLE_PRECISION,&
1060 Next,579,MPI_COMM_WORLD,IERROR)
1061 call MPI_Send(allcart(1,1,1),3*icount*2*nres,&
1062 MPI_REAL,Next,580,MPI_COMM_WORLD,IERROR)
1064 call dawrite_ccoords(jj_old,jj,icbase)
1066 end subroutine write_and_send_cconf
1067 !------------------------------------------------------------------------------
1069 subroutine receive_and_pass_cconf(icount,jj_old,jj,Previous,Next)
1073 ! include "DIMENSIONS"
1074 ! include "sizesclu.dat"
1076 integer :: IERROR !,STATUS(MPI_STATUS_SIZE)
1077 ! include "COMMON.MPI"
1078 ! include "COMMON.CHAIN"
1079 ! include "COMMON.SBRIDGE"
1080 ! include "COMMON.INTERACT"
1081 ! include "COMMON.IOUNITS"
1082 ! include "COMMON.VAR"
1083 ! include "COMMON.GEO"
1084 ! include "COMMON.CLUSTER"
1085 integer :: i,j,k,icount,jj_old,jj,Previous,Next
1088 write (iout,*) "Processor",me," entered RECEIVE_AND_PASS_CONF"
1091 do while (icount.gt.0)
1092 call MPI_Recv(icount,1,MPI_INTEGER,Previous,570,MPI_COMM_WORLD,&
1094 call MPI_Send(icount,1,MPI_INTEGER,Next,570,MPI_COMM_WORLD,&
1097 write (iout,*) "Processor",me," icount",icount
1099 if (icount.eq.0) return
1100 call MPI_Recv(nss_all(1),icount,MPI_INTEGER,&
1101 Previous,571,MPI_COMM_WORLD,STATUS,IERROR)
1102 call MPI_Send(nss_all(1),icount,MPI_INTEGER,&
1103 Next,571,MPI_COMM_WORLD,IERROR)
1104 call MPI_Recv(ihpb_all(1,1),icount,MPI_INTEGER,&
1105 Previous,572,MPI_COMM_WORLD,STATUS,IERROR)
1106 call MPI_Send(ihpb_all(1,1),icount,MPI_INTEGER,&
1107 Next,572,MPI_COMM_WORLD,IERROR)
1108 call MPI_Recv(jhpb_all(1,1),icount,MPI_INTEGER,&
1109 Previous,573,MPI_COMM_WORLD,STATUS,IERROR)
1110 call MPI_Send(jhpb_all(1,1),icount,MPI_INTEGER,&
1111 Next,573,MPI_COMM_WORLD,IERROR)
1112 call MPI_Recv(rmstb(jj_old),icount,MPI_DOUBLE_PRECISION,&
1113 Previous,577,MPI_COMM_WORLD,STATUS,IERROR)
1114 call MPI_Send(rmstb(jj_old),icount,MPI_DOUBLE_PRECISION,&
1115 Next,577,MPI_COMM_WORLD,IERROR)
1116 call MPI_Recv(entfac(jj_old),icount,MPI_DOUBLE_PRECISION,&
1117 Previous,579,MPI_COMM_WORLD,STATUS,IERROR)
1118 call MPI_Send(entfac(jj_old),icount,MPI_DOUBLE_PRECISION,&
1119 Next,579,MPI_COMM_WORLD,IERROR)
1120 call MPI_Recv(allcart(1,1,1),3*icount*2*nres,&
1121 MPI_REAL,Previous,580,MPI_COMM_WORLD,STATUS,IERROR)
1122 call MPI_Send(allcart(1,1,1),3*icount*2*nres,&
1123 MPI_REAL,Next,580,MPI_COMM_WORLD,IERROR)
1125 call dawrite_ccoords(jj_old,jj,icbase)
1128 write (iout,*) "Processor",me," received",icount," conformations"
1130 write (iout,'(8f10.4)') (allcart(l,k,i),l=1,3,k=1,nres)
1131 write (iout,'(8f10.4)')((allcart(l,k,i+nres),l=1,3,k=nnt,nct)
1132 write (iout,'(e15.5,16i5)') entfac(i)
1137 end subroutine receive_and_pass_cconf
1139 !------------------------------------------------------------------------------
1140 subroutine daread_ccoords(istart_conf,iend_conf)
1142 use energy_data, only: dyn_ss
1145 ! include "DIMENSIONS"
1146 ! include "sizesclu.dat"
1150 ! include "COMMON.MPI"
1152 ! include "COMMON.CHAIN"
1153 ! include "COMMON.CLUSTER"
1154 ! include "COMMON.IOUNITS"
1155 ! include "COMMON.INTERACT"
1156 ! include "COMMON.VAR"
1157 ! include "COMMON.SBRIDGE"
1158 ! include "COMMON.GEO"
1159 integer :: istart_conf,iend_conf
1160 integer :: i,j,ij,ii,iii
1162 character(len=16) :: form,acc
1163 character(len=32) :: nam
1165 ! Read conformations off a DA scratchfile.
1168 write (iout,*) "DAREAD_COORDS"
1169 write (iout,*) "istart_conf",istart_conf," iend_conf",iend_conf
1170 inquire(unit=icbase,name=nam,recl=len,form=form,access=acc)
1171 write (iout,*) "len=",len," form=",form," acc=",acc
1172 write (iout,*) "nam=",nam
1175 do ii=istart_conf,iend_conf
1176 ij = ii - istart_conf + 1
1179 write (iout,*) "Reading binary file, record",iii," ii",ii
1183 read(icbase,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres), &
1184 ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),&
1185 entfac(ii),rmstb(ii)
1187 read(icbase,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres),&
1188 ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),&
1189 nss_all(ij),(ihpb_all(i,ij),jhpb_all(i,ij),i=1,nss),&
1190 entfac(ii),rmstb(ii)
1194 write (iout,*) ii,iii,ij,entfac(ii)
1195 write (iout,'(8f10.5)') ((allcart(j,i,ij),j=1,3),i=1,nres)
1196 write (iout,'(8f10.4)') ((allcart(j,i,ij),j=1,3),&
1197 i=nnt+nres,nct+nres)
1198 write (iout,'(2e15.5)') entfac(ij)
1199 write (iout,'(16i5)') nss_all(ij),(ihpb_all(i,ij),&
1200 jhpb_all(i,ij),i=1,nss)
1205 end subroutine daread_ccoords
1206 !------------------------------------------------------------------------------
1207 subroutine dawrite_ccoords(istart_conf,iend_conf,unit_out)
1210 ! include "DIMENSIONS"
1211 ! include "sizesclu.dat"
1212 use energy_data, only: dyn_ss
1217 ! include "COMMON.MPI"
1219 ! include "COMMON.CHAIN"
1220 ! include "COMMON.INTERACT"
1221 ! include "COMMON.IOUNITS"
1222 ! include "COMMON.VAR"
1223 ! include "COMMON.SBRIDGE"
1224 ! include "COMMON.GEO"
1225 ! include "COMMON.CLUSTER"
1226 integer :: istart_conf,iend_conf
1227 integer :: i,j,ii,ij,iii,unit_out
1229 character(len=16) :: form,acc
1230 character(len=32) :: nam
1232 ! Write conformations to a DA scratchfile.
1235 write (iout,*) "DAWRITE_COORDS"
1236 write (iout,*) "istart_conf",istart_conf," iend_conf",iend_conf
1237 write (iout,*) "lenrec",lenrec
1238 inquire(unit=unit_out,name=nam,recl=len,form=form,access=acc)
1239 write (iout,*) "len=",len," form=",form," acc=",acc
1240 write (iout,*) "nam=",nam
1243 do ii=istart_conf,iend_conf
1245 ij = ii - istart_conf + 1
1247 write (iout,*) "Writing binary file, record",iii," ii",ii
1251 write(unit_out,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres),&
1252 ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),&
1253 entfac(ii),rmstb(ii)
1255 write(unit_out,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres), &
1256 ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),&
1257 nss_all(ij),(ihpb_all(i,ij),jhpb_all(i,ij),i=1,nss_all(ij)),&
1258 entfac(ii),rmstb(ii)
1262 write (iout,'(8f10.5)') ((allcart(j,i,ij),j=1,3),i=1,nres)
1263 write (iout,'(8f10.4)') ((allcart(j,i,ij),j=1,3),i=nnt+nres,&
1265 write (iout,'(2e15.5)') entfac(ij)
1266 write (iout,'(16i5)') nss_all(ij),(ihpb(i,ij),jhpb(i,ij),i=1,&
1272 end subroutine dawrite_ccoords
1273 !-----------------------------------------------------------------------------
1275 !-----------------------------------------------------------------------------
1276 subroutine read_control
1278 ! Read molecular data
1280 use energy_data, only: rescale_mode,distchainmax,ipot !,temp0
1281 use control_data, only: titel,outpdb,outmol2,refstr,pdbref,&
1282 iscode,symetr,punch_dist,print_dist,nstart,nend,&
1283 caonly,iopt,efree,lprint_cart,lprint_int,rlamb_ele,&
1284 r_cut_ele,nclust,tor_mode,scelemode
1286 ! include 'DIMENSIONS'
1287 ! include 'sizesclu.dat'
1288 ! include 'COMMON.IOUNITS'
1289 ! include 'COMMON.TIME1'
1290 ! include 'COMMON.SBRIDGE'
1291 ! include 'COMMON.CONTROL'
1292 ! include 'COMMON.CLUSTER'
1293 ! include 'COMMON.CHAIN'
1294 ! include 'COMMON.HEADER'
1295 ! include 'COMMON.FFIELD'
1296 ! include 'COMMON.FREE'
1297 ! include 'COMMON.INTERACT'
1298 character(len=320) :: controlcard !,ucase
1300 ! include 'COMMON.INFO'
1304 read (INP,'(a80)') titel
1305 call card_concat(controlcard,.true.)
1307 call readi(controlcard,'NRES',nres_molec(1),0)
1308 call readi(controlcard,'NRES_NUCL',nres_molec(2),0)
1309 call readi(controlcard,'NRES_CAT',nres_molec(5),0)
1312 nres=nres_molec(i)+nres
1314 print *,"TU",nres_molec(:)
1316 ! call alloc_clust_arrays
1317 allocate(rcutoff(max_cut+1)) !(max_cut+1)
1318 allocate(beta_h(maxT)) !(maxT)
1319 allocate(mult(nres)) !(maxres)
1322 call readi(controlcard,'RESCALE',rescale_mode,2)
1323 call reada(controlcard,'DISTCHAINMAX',distchainmax,50.0d0)
1324 write (iout,*) "DISTCHAINMAX",distchainmax
1325 call readi(controlcard,'PDBOUT',outpdb,0)
1326 call readi(controlcard,'MOL2OUT',outmol2,0)
1327 refstr=(index(controlcard,'REFSTR').gt.0)
1328 write (iout,*) "REFSTR",refstr
1329 pdbref=(index(controlcard,'PDBREF').gt.0)
1330 iscode=index(controlcard,'ONE_LETTER')
1331 tree=(index(controlcard,'MAKE_TREE').gt.0)
1332 min_var=(index(controlcard,'MINVAR').gt.0)
1333 plot_tree=(index(controlcard,'PLOT_TREE').gt.0)
1334 punch_dist=(index(controlcard,'PUNCH_DIST').gt.0)
1335 call readi(controlcard,'NCUT',ncut,1)
1336 call readi(controlcard,'SYM',symetr,1)
1337 write (iout,*) 'sym', symetr
1338 call readi(controlcard,'NSTART',nstart,0)
1339 call reada(controlcard,'BOXX',boxxsize,100.0d0)
1340 call reada(controlcard,'BOXY',boxysize,100.0d0)
1341 call reada(controlcard,'BOXZ',boxzsize,100.0d0)
1342 call readi(controlcard,'NCLUST',nclust,5)
1343 ! ions=index(controlcard,"IONS").gt.0
1344 call readi(controlcard,"SCELEMODE",scelemode,0)
1345 print *,"SCELE",scelemode
1346 call readi(controlcard,'TORMODE',tor_mode,0)
1347 !C if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0) then
1348 write(iout,*) "torsional and valence angle mode",tor_mode
1349 call reada(controlcard,"R_CUT_ELE",r_cut_ele,25.0d0)
1350 call reada(controlcard,"LAMBDA_ELE",rlamb_ele,0.3d0)
1351 write(iout,*) "R_CUT_ELE=",r_cut_ele
1352 call readi(controlcard,'NEND',nend,0)
1353 call reada(controlcard,'ECUT',ecut,10.0d0)
1354 call reada(controlcard,'PROB',prob_limit,0.99d0)
1355 write (iout,*) "Probability limit",prob_limit
1356 lgrp=(index(controlcard,'LGRP').gt.0)
1357 caonly=(index(controlcard,'CA_ONLY').gt.0)
1358 print_dist=(index(controlcard,'PRINT_DIST').gt.0)
1359 call multreada(controlcard,'CUTOFF',rcutoff,ncut,-1.0d0)
1360 call readi(controlcard,'IOPT',iopt,2)
1361 lside = index(controlcard,"SIDE").gt.0
1362 efree = index(controlcard,"EFREE").gt.0
1363 call readi(controlcard,'NTEMP',nT,1)
1364 write (iout,*) "nT",nT
1365 !el call reada(controlcard,'TEMP0',temp0,300.0d0) !el
1366 call multreada(controlcard,'TEMPER',beta_h,nT,300.0d0)
1367 write (iout,*) "nT",nT
1368 write (iout,*) 'beta_h',(beta_h(i),i=1,nT)
1370 beta_h(i)=1.0d0/(1.987D-3*beta_h(i))
1372 write (iout,*) 'beta_h',(beta_h(i),i=1,nT)
1373 lprint_cart=index(controlcard,"PRINT_CART") .gt.0
1374 lprint_int=index(controlcard,"PRINT_INT") .gt.0
1377 end subroutine read_control
1378 !-----------------------------------------------------------------------------
1381 ! Read molecular data.
1383 use geometry_data, only: nsup,cref,nres0,nstart_sup,nstart_seq,dc
1384 use energy_data!, only: wsc,wscp,welec,wbond,wstrain,wtor,wtor_d,&
1385 ! wang,wscloc,wcorr,wcorr4,wcorr5,wcorr6,wsccor,wel_loc,&
1386 ! wturn3,wturn4,wturn6,wvdwpp,weights
1387 use control_data, only: titel,nstart,nend,pdbref,refstr,iscode,&
1389 use geometry, only: chainbuild,alloc_geo_arrays
1390 use energy, only: alloc_ener_arrays
1391 use io_base, only: rescode
1392 use control, only: setup_var,init_int_table
1393 use conform_compar, only: contact
1395 ! include 'DIMENSIONS'
1396 ! include 'COMMON.IOUNITS'
1397 ! include 'COMMON.GEO'
1398 ! include 'COMMON.VAR'
1399 ! include 'COMMON.INTERACT'
1400 ! include 'COMMON.LOCAL'
1401 ! include 'COMMON.NAMES'
1402 ! include 'COMMON.CHAIN'
1403 ! include 'COMMON.FFIELD'
1404 ! include 'COMMON.SBRIDGE'
1405 ! include 'COMMON.HEADER'
1406 ! include 'COMMON.CONTROL'
1407 ! include 'COMMON.CONTACTS'
1408 ! include 'COMMON.TIME1'
1410 ! include 'COMMON.INFO'
1412 character(len=4) :: sequence(nres) !(maxres)
1413 character(len=800) :: weightcard
1415 real(kind=8) :: x(6*nres) !(maxvar)
1416 real(kind=8) :: ssscale
1417 integer :: itype_pdb(nres) !(maxres)
1419 integer :: i,j,kkk,mnum
1423 !el allocate(weights(n_ene))
1424 allocate(weights(max_ene))
1425 call alloc_geo_arrays
1426 call alloc_ener_arrays
1427 !-----------------------------
1428 allocate(c(3,2*nres+2)) !(3,maxres2+2) maxres2=2*maxres
1429 allocate(dc(3,0:2*nres+2)) !(3,0:maxres2)
1430 allocate(itype(nres+2,5)) !(maxres)
1431 allocate(molnum(nres+1))
1432 allocate(itel(nres+2))
1446 !--------------------------
1447 ! Read weights of the subsequent energy terms.
1448 call card_concat(weightcard,.true.)
1449 call reada(weightcard,'WSC',wsc,1.0d0)
1450 call reada(weightcard,'WLONG',wsc,wsc)
1451 call reada(weightcard,'WSCP',wscp,1.0d0)
1452 call reada(weightcard,'WELEC',welec,1.0D0)
1453 call reada(weightcard,'WVDWPP',wvdwpp,welec)
1454 call reada(weightcard,'WEL_LOC',wel_loc,1.0D0)
1455 call reada(weightcard,'WCORR4',wcorr4,0.0D0)
1456 call reada(weightcard,'WCORR5',wcorr5,0.0D0)
1457 call reada(weightcard,'WCORR6',wcorr6,0.0D0)
1458 call reada(weightcard,'WTURN3',wturn3,1.0D0)
1459 call reada(weightcard,'WTURN4',wturn4,1.0D0)
1460 call reada(weightcard,'WTURN6',wturn6,1.0D0)
1461 call reada(weightcard,'WSTRAIN',wstrain,1.0D0)
1462 call reada(weightcard,'WSCCOR',wsccor,1.0D0)
1463 call reada(weightcard,'WBOND',wbond,1.0D0)
1464 call reada(weightcard,'WTOR',wtor,1.0D0)
1465 call reada(weightcard,'WTORD',wtor_d,1.0D0)
1466 call reada(weightcard,'WANG',wang,1.0D0)
1467 call reada(weightcard,'WSCLOC',wscloc,1.0D0)
1468 call reada(weightcard,'SCAL14',scal14,0.4D0)
1469 call reada(weightcard,'SCALSCP',scalscp,1.0d0)
1470 call reada(weightcard,'CUTOFF',cutoff_corr,7.0d0)
1471 call reada(weightcard,'DELT_CORR',delt_corr,0.5d0)
1472 call reada(weightcard,'WCATCAT',wcatcat,0.0d0)
1473 call reada(weightcard,'WCATPROT',wcatprot,0.0d0)
1474 call reada(weightcard,'TEMP0',temp0,300.0d0) !!! el
1475 if (index(weightcard,'SOFT').gt.0) ipot=6
1476 call reada(weightcard,'WVDWPP_NUCL',wvdwpp_nucl,0.0D0)
1477 call reada(weightcard,'WELPP',welpp,0.0d0)
1478 call reada(weightcard,'WVDWPSB',wvdwpsb,0.0d0)
1479 call reada(weightcard,'WELPSB',welpsb,0.0D0)
1480 call reada(weightcard,'WVDWSB',wvdwsb,0.0d0)
1481 call reada(weightcard,'WELSB',welsb,0.0D0)
1482 call reada(weightcard,'WBOND_NUCL',wbond_nucl,0.0D0)
1483 call reada(weightcard,'WANG_NUCL',wang_nucl,0.0D0)
1484 call reada(weightcard,'WSBLOC',wsbloc,0.0D0)
1485 call reada(weightcard,'WTOR_NUCL',wtor_nucl,0.0D0)
1486 ! print *,"KUR..",wtor_nucl
1487 call reada(weightcard,'WTORD_NUCL',wtor_d_nucl,0.0D0)
1488 call reada(weightcard,'WCORR_NUCL',wcorr_nucl,0.0D0)
1489 call reada(weightcard,'WCORR3_NUC',wcorr3_nucl,0.0D0)
1490 ! 12/1/95 Added weight for the multi-body term WCORR
1491 call reada(weightcard,'WCORRH',wcorr,1.0D0)
1492 call reada(weightcard,'WSCBASE',wscbase,0.0D0)
1493 call reada(weightcard,'WPEPBASE',wpepbase,1.0d0)
1494 call reada(weightcard,'WSCPHO',wscpho,0.0d0)
1495 call reada(weightcard,'WPEPPHO',wpeppho,0.0d0)
1497 call reada(weightcard,"D0CM",d0cm,3.78d0)
1498 call reada(weightcard,"AKCM",akcm,15.1d0)
1499 call reada(weightcard,"AKTH",akth,11.0d0)
1500 call reada(weightcard,"AKCT",akct,12.0d0)
1501 call reada(weightcard,"V1SS",v1ss,-1.08d0)
1502 call reada(weightcard,"V2SS",v2ss,7.61d0)
1503 call reada(weightcard,"V3SS",v3ss,13.7d0)
1504 call reada(weightcard,"EBR",ebr,-5.50D0)
1505 call reada(weightcard,"ATRISS",atriss,0.301D0)
1506 call reada(weightcard,"BTRISS",btriss,0.021D0)
1507 call reada(weightcard,"CTRISS",ctriss,1.001D0)
1508 call reada(weightcard,"DTRISS",dtriss,1.001D0)
1509 call reada(weightcard,"SSSCALE",ssscale,1.0D0)
1510 dyn_ss=(index(weightcard,'DYN_SS').gt.0)
1512 call reada(weightcard,"HT",Ht,0.0D0)
1514 ss_depth=(ebr/wsc-0.25*eps(1,1))*ssscale
1515 Ht=(Ht/wsc-0.25*eps(1,1))*ssscale
1516 akcm=akcm/wsc*ssscale
1517 akth=akth/wsc*ssscale
1518 akct=akct/wsc*ssscale
1519 v1ss=v1ss/wsc*ssscale
1520 v2ss=v2ss/wsc*ssscale
1521 v3ss=v3ss/wsc*ssscale
1523 ss_depth=ebr/wstrain-0.25*eps(1,1)*wsc/wstrain
1526 if (wcorr4.gt.0.0d0) wcorr=wcorr4
1545 !el weights(19)=wsccor !!!!!!!!!!!!!!!!
1547 weights(26)=wvdwpp_nucl
1553 weights(32)=wbond_nucl
1554 weights(33)=wang_nucl
1556 weights(35)=wtor_nucl
1557 weights(36)=wtor_d_nucl
1558 weights(37)=wcorr_nucl
1559 weights(38)=wcorr3_nucl
1561 weights(42)=wcatprot
1567 write (iout,10) wsc,wscp,welec,wvdwpp,wbond,wang,wscloc,wtor,&
1568 wtor_d,wstrain,wel_loc,wcorr,wcorr5,wcorr6,wturn3,&
1569 wturn4,wturn6,wsccor
1570 10 format (/'Energy-term weights (unscaled):'// &
1571 'WSCC= ',f10.6,' (SC-SC)'/ &
1572 'WSCP= ',f10.6,' (SC-p)'/ &
1573 'WELEC= ',f10.6,' (p-p electr)'/ &
1574 'WVDWPP= ',f10.6,' (p-p VDW)'/ &
1575 'WBOND= ',f10.6,' (stretching)'/ &
1576 'WANG= ',f10.6,' (bending)'/ &
1577 'WSCLOC= ',f10.6,' (SC local)'/ &
1578 'WTOR= ',f10.6,' (torsional)'/ &
1579 'WTORD= ',f10.6,' (double torsional)'/ &
1580 'WSTRAIN=',f10.6,' (SS bridges & dist. cnstr.)'/ &
1581 'WEL_LOC=',f10.6,' (multi-body 3-rd order)'/ &
1582 'WCORR4= ',f10.6,' (multi-body 4th order)'/ &
1583 'WCORR5= ',f10.6,' (multi-body 5th order)'/ &
1584 'WCORR6= ',f10.6,' (multi-body 6th order)'/ &
1585 'WTURN3= ',f10.6,' (turns, 3rd order)'/ &
1586 'WTURN4= ',f10.6,' (turns, 4th order)'/ &
1587 'WTURN6= ',f10.6,' (turns, 6th order)'/ &
1588 'WSCCOR= ',f10.6,' (SC-backbone torsinal correalations)')
1590 if (wcorr4.gt.0.0d0) then
1591 write (iout,'(/2a/)') 'Local-electrostatic type correlation ',&
1592 'between contact pairs of peptide groups'
1593 write (iout,'(2(a,f5.3/))') &
1594 'Cutoff on 4-6th order correlation terms: ',cutoff_corr,&
1595 'Range of quenching the correlation terms:',2*delt_corr
1596 else if (wcorr.gt.0.0d0) then
1597 write (iout,'(/2a/)') 'Hydrogen-bonding correlation ',&
1598 'between contact pairs of peptide groups'
1600 write (iout,'(a,f8.3)') &
1601 'Scaling factor of 1,4 SC-p interactions:',scal14
1602 write (iout,'(a,f8.3)') &
1603 'General scaling factor of SC-p interactions:',scalscp
1604 r0_corr=cutoff_corr-delt_corr
1606 aad(i,1)=scalscp*aad(i,1)
1607 aad(i,2)=scalscp*aad(i,2)
1608 bad(i,1)=scalscp*bad(i,1)
1609 bad(i,2)=scalscp*bad(i,2)
1613 print *,'indpdb=',indpdb,' pdbref=',pdbref
1615 ! Read sequence if not taken from the pdb file.
1616 if (iscode.gt.0) then
1617 read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
1619 read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
1621 print *,nres_molec(:),nres
1622 ! Convert sequence to numeric code
1623 do i=1,nres_molec(1)
1626 itype(i,mnum)=rescode(i,sequence(i),iscode,mnum)
1628 do i=nres_molec(1)+1,nres_molec(1)+nres_molec(2)
1631 write (iout,*),i,sequence(i)
1632 itype(i,mnum)=rescode(i,sequence(i),iscode,mnum)
1635 do i=nres_molec(1)+nres_molec(2)+1,nres_molec(1)+nres_molec(2)+nres_molec(5)
1638 itype(i,mnum)=rescode(i,sequence(i),iscode,mnum)
1641 print '(20i4)',(itype(i,molnum(i)),i=1,nres)
1645 if (itype(i,molnum(i)).eq.ntyp1_molec(molnum(i))&
1646 .or. itype(i+1,molnum(i+1)).eq.ntyp1_molec(molnum(i+1))) then
1648 if (itype(i,molnum(i)).eq.ntyp1_molec(molnum(i))) then
1652 else if (iabs(itype(i+1,1)).ne.20) then
1654 else if (iabs(itype(i,1)).ne.20) then
1661 write (iout,*) "ITEL"
1663 write (iout,*) i,itype(i,molnum(i)),itel(i)
1666 print *,'Call Read_Bridge.'
1670 print *,'NNT=',NNT,' NCT=',NCT
1671 if (itype(1,molnum(1)).eq.ntyp1_molec(molnum(1))) nnt=2
1672 if (itype(nres,molnum(nres)).eq.ntyp1_molec(molnum(nres))) nct=nct-1
1673 if (nstart.lt.nnt) nstart=nnt
1674 if (nend.gt.nct .or. nend.eq.0) nend=nct
1675 write (iout,*) "nstart",nstart," nend",nend
1678 ! read(inp,'(a)') pdbfile
1679 ! write (iout,'(2a)') 'PDB data will be read from file ',pdbfile
1680 ! open(ipdbin,file=pdbfile,status='old',err=33)
1682 ! 33 write (iout,'(a)') 'Error opening PDB file.'
1685 ! print *,'Begin reading pdb data'
1687 ! print *,'Finished reading pdb data'
1688 ! write (iout,'(a,i3,a,i3)')'nsup=',nsup,' nstart_sup=',nstart_sup
1690 ! itype_pdb(i)=itype(i)
1693 ! write (iout,'(a,i3)') 'nsup=',nsup
1695 ! if (nsup.le.(nct-nnt+1)) then
1696 ! do i=0,nct-nnt+1-nsup
1697 ! if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),nsup)) then
1702 ! write (iout,'(a)')
1703 ! & 'Error - sequences to be superposed do not match.'
1706 ! do i=0,nsup-(nct-nnt+1)
1707 ! if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),nct-nnt+1))
1709 ! nstart_sup=nstart_sup+i
1714 ! write (iout,'(a)')
1715 ! & 'Error - sequences to be superposed do not match.'
1718 ! write (iout,*) 'nsup=',nsup,' nstart_sup=',nstart_sup,
1719 ! & ' nstart_seq=',nstart_seq
1721 write(iout,*)"przed ini_int_tab"
1723 write(iout,*)"po ini_int_tab"
1724 write(iout,*)"przed setup var"
1726 write(iout,*)"po setup var"
1727 if (ns.gt.0.and.dyn_ss) then
1731 forcon(i-nss)=forcon(i)
1738 ! call hpb_partition
1740 dyn_ss_mask(iss(i))=.true.
1744 write (iout,*) "molread: REFSTR",refstr
1746 if (.not.pdbref) then
1747 call read_angles(inp,*38)
1749 38 write (iout,'(a)') 'Error reading reference structure.'
1751 call mp_stopall(Error_Msg)
1753 stop 'Error reading reference structure'
1762 cref(j,i,kkk)=c(j,i)
1766 call contact(.true.,ncont_ref,icont_ref)
1769 end subroutine molread
1770 !-----------------------------------------------------------------------------
1771 subroutine openunits
1773 ! include 'DIMENSIONS'
1774 use control_data, only: from_cx,from_bx,from_cart
1778 character(len=3) :: liczba
1779 ! include "COMMON.MPI"
1781 ! include 'COMMON.IOUNITS'
1782 ! include 'COMMON.CONTROL'
1783 integer :: lenpre,lenpot !,ilen
1785 character(len=16) :: cformat,cprint
1786 ! character(len=16) ucase
1787 integer :: lenint,lenout
1788 call getenv('INPUT',prefix)
1789 call getenv('OUTPUT',prefout)
1790 call getenv('INTIN',prefintin)
1791 call getenv('COORD',cformat)
1792 call getenv('PRINTCOOR',cprint)
1793 call getenv('SCRATCHDIR',scratchdir)
1796 if (index(ucase(cformat),'CX').gt.0) then
1802 lenout=ilen(prefout)
1803 lenint=ilen(prefintin)
1804 ! Get the names and open the input files
1805 open (inp,file=prefix(:ilen(prefix))//'.inp',status='old')
1807 write (liczba,'(bz,i3.3)') me
1808 outname=prefout(:lenout)//'_clust.out_'//liczba
1810 outname=prefout(:lenout)//'_clust.out'
1813 intinname=prefintin(:lenint)//'.bx'
1814 else if (from_cx) then
1815 intinname=prefintin(:lenint)//'.cx'
1817 intinname=prefintin(:lenint)//'.int'
1819 rmsname=prefintin(:lenint)//'.rms'
1820 open (jplot,file=prefout(:ilen(prefout))//'.tex',&
1822 open (jrms,file=rmsname,status='unknown')
1823 open(iout,file=outname,status='unknown')
1824 ! Get parameter filenames and open the parameter files.
1825 call getenv('BONDPAR',bondname)
1826 open (ibond,file=bondname,status='old')
1827 call getenv('THETPAR',thetname)
1828 open (ithep,file=thetname,status='old')
1829 call getenv('ROTPAR',rotname)
1830 open (irotam,file=rotname,status='old')
1831 call getenv('TORPAR',torname)
1832 open (itorp,file=torname,status='old')
1833 call getenv('TORDPAR',tordname)
1834 open (itordp,file=tordname,status='old')
1835 call getenv('FOURIER',fouriername)
1836 open (ifourier,file=fouriername,status='old')
1837 call getenv('ELEPAR',elename)
1838 open (ielep,file=elename,status='old')
1839 call getenv('SIDEPAR',sidename)
1840 open (isidep,file=sidename,status='old')
1841 call getenv('SIDEP',sidepname)
1842 open (isidep1,file=sidepname,status="old")
1843 call getenv('SCCORPAR',sccorname)
1844 open (isccor,file=sccorname,status="old")
1845 call getenv('BONDPAR_NUCL',bondname_nucl)
1846 open (ibond_nucl,file=bondname_nucl,status='old')
1847 call getenv('THETPAR_NUCL',thetname_nucl)
1848 open (ithep_nucl,file=thetname_nucl,status='old')
1849 call getenv('ROTPAR_NUCL',rotname_nucl)
1850 open (irotam_nucl,file=rotname_nucl,status='old')
1851 call getenv('TORPAR_NUCL',torname_nucl)
1852 open (itorp_nucl,file=torname_nucl,status='old')
1853 call getenv('TORDPAR_NUCL',tordname_nucl)
1854 open (itordp_nucl,file=tordname_nucl,status='old')
1855 call getenv('SIDEPAR_NUCL',sidename_nucl)
1856 open (isidep_nucl,file=sidename_nucl,status='old')
1857 call getenv('SIDEPAR_SCBASE',sidename_scbase)
1858 open (isidep_scbase,file=sidename_scbase,status='old')
1859 call getenv('PEPPAR_PEPBASE',pepname_pepbase)
1860 open (isidep_pepbase,file=pepname_pepbase,status='old')
1861 call getenv('SCPAR_PHOSPH',pepname_scpho)
1862 open (isidep_scpho,file=pepname_scpho,status='old')
1863 call getenv('PEPPAR_PHOSPH',pepname_peppho)
1864 open (isidep_peppho,file=pepname_peppho,status='old')
1867 call getenv('LIPTRANPAR',liptranname)
1868 open (iliptranpar,file=liptranname,status='old')
1869 call getenv('TUBEPAR',tubename)
1870 open (itube,file=tubename,status='old')
1871 call getenv('IONPAR',ionname)
1872 open (iion,file=ionname,status='old')
1876 ! 8/9/01 In the newest version SCp interaction constants are read from a file
1877 ! Use -DOLDSCP to use hard-coded constants instead.
1879 call getenv('SCPPAR',scpname)
1880 open (iscpp,file=scpname,status='old')
1881 call getenv('SCPPAR_NUCL',scpname_nucl)
1882 open (iscpp_nucl,file=scpname_nucl,status='old')
1886 end subroutine openunits
1887 !-----------------------------------------------------------------------------
1889 !-----------------------------------------------------------------------------
1890 subroutine pdboutC(etot,rmsd,tytul)
1892 use energy_data, only: ihpb,jhpb,itype,molnum
1893 use energy, only:boxshift,to_box
1894 ! implicit real*8 (a-h,o-z)
1895 ! include 'DIMENSIONS'
1896 ! include 'COMMON.CONTROL'
1897 ! include 'COMMON.CHAIN'
1898 ! include 'COMMON.INTERACT'
1899 ! include 'COMMON.NAMES'
1900 ! include 'COMMON.IOUNITS'
1901 ! include 'COMMON.HEADER'
1902 ! include 'COMMON.SBRIDGE'
1903 ! include 'COMMON.TEMPFAC'
1904 character(len=50) :: tytul
1905 character(len=1) :: chainid(10)=(/'A','B','C','D','E','F',&
1907 integer :: ica(nres),k,iti1
1908 real(kind=8) :: etot,rmsd
1909 integer :: iatom,ichain,ires,i,j,iti,mnum,mnum1,boxxxshift(3)
1910 real(kind=8) :: boxxxx(3),cbeg(3),Rdist(3)
1911 write (ipdb,'(3a,1pe15.5,a,0pf7.2)') 'REMARK ',tytul(:20),&
1912 ' ENERGY ',etot,' RMS ',rmsd
1916 ! boxxxshift(1)=int(c(1,nnt)/boxxsize)
1917 ! if (c(1,nnt).lt.0) boxxxshift(1)=boxxxshift(1)-1
1918 ! boxxxshift(2)=int(c(2,nnt)/boxzsize)
1919 ! if (c(2,nnt).lt.0) boxxxshift(2)=boxxxshift(2)-1
1920 ! boxxxshift(3)=int(c(3,nnt)/boxzsize)
1921 ! if (c(3,nnt).lt.0) boxxxshift(3)=boxxxshift(3)-1
1923 ! if (boxxxshift(i).gt.2) boxxxshift(i)=2
1924 ! if (boxxxshift(i).lt.-2) boxxxshift(i)=-2
1930 call to_box(cbeg(1),cbeg(2),cbeg(3))
1939 iti1=itype(i+1,mnum)
1941 if ((iti.eq.ntyp1_molec(mnum)).and.(iti1.eq.ntyp1_molec(mnum1))) cycle
1943 if (iti.eq.ntyp1_molec(mnum)) then
1946 write (ipdb,'(a)') 'TER'
1956 call to_box(c(1,i),c(2,i),c(3,i))
1959 Rdist(k) = boxshift(c(k,i) - cbeg(k),boxxxx(k))
1960 c(k,i)=cbeg(k)+Rdist(k)
1962 if ((itype(i,molnum(i)).ne.10).and.(molnum(i).ne.5)) then
1963 call to_box(c(1,i+nres),c(2,i+nres),c(3,i+nres))
1965 Rdist(k) = boxshift(c(k,i+nres) - cbeg(k),boxxxx(k))
1966 c(k,i+nres)=cbeg(k)+Rdist(k)
1969 write (ipdb,10) iatom,restyp(iti,mnum),chainid(ichain),&
1970 ires,(c(j,i),j=1,3),1.0d0,tempfac(1,i)
1971 if ((iti.ne.10).and.(mnum.ne.5)) then
1973 write (ipdb,20) iatom,restyp(iti,mnum),chainid(ichain),&
1974 ires,(c(j,nres+i),j=1,3),1.0d0,tempfac(2,i)
1978 write (ipdb,'(a)') 'TER'
1982 if ((itype(i,mnum).eq.ntyp1).or.(mnum.eq.5)) cycle
1983 if (itype(i,mnum).eq.10 .and. itype(i+1,mnum1).ne.ntyp1_molec(mnum1)) then
1984 write (ipdb,30) ica(i),ica(i+1)
1985 else if (itype(i,mnum).ne.10 .and. itype(i+1,mnum1).ne.ntyp1_molec(mnum1)) then
1986 write (ipdb,30) ica(i),ica(i+1),ica(i)+1
1987 else if (itype(i,mnum).ne.10 .and. itype(i+1,mnum1).eq.ntyp1_molec(mnum1)) then
1988 write (ipdb,30) ica(i),ica(i)+1
1991 if ((itype(nct,molnum(nct)).ne.10).and.(molnum(nct).ne.5)) then
1992 write (ipdb,30) ica(nct),ica(nct)+1
1995 write (ipdb,30) ica(ihpb(i)-nres)+1,ica(jhpb(i)-nres)+1
1997 write (ipdb,'(a6)') 'ENDMDL'
1998 10 FORMAT ('ATOM',I7,' CA ',A3,1X,A1,I4,4X,3F8.3,2f6.2)
1999 20 FORMAT ('ATOM',I7,' CB ',A3,1X,A1,I4,4X,3F8.3,2f6.2)
2000 30 FORMAT ('CONECT',8I5)
2002 end subroutine pdboutC
2003 !-----------------------------------------------------------------------------
2004 subroutine cartout(igr,i,etot,free,rmsd,plik)
2005 ! implicit real*8 (a-h,o-z)
2006 ! include 'DIMENSIONS'
2007 ! include 'sizesclu.dat'
2008 ! include 'COMMON.IOUNITS'
2009 ! include 'COMMON.CHAIN'
2010 ! include 'COMMON.VAR'
2011 ! include 'COMMON.LOCAL'
2012 ! include 'COMMON.INTERACT'
2013 ! include 'COMMON.NAMES'
2014 ! include 'COMMON.GEO'
2015 ! include 'COMMON.CLUSTER'
2016 integer :: igr,i,j,k
2017 real(kind=8) :: etot,free,rmsd
2018 character(len=80) :: plik
2019 open (igeom,file=plik,position='append')
2020 write (igeom,'(2e15.5,f10.5,$)') etot,free,rmsd
2021 write (igeom,'(i4,$)') &
2022 nss_all(i),(ihpb_all(j,i),jhpb_all(j,i),j=1,nss_all(i))
2023 write (igeom,'(i10)') iscore(i)
2024 write (igeom,'(8f10.5)') &
2025 ((allcart(k,j,i),k=1,3),j=1,nres),&
2026 ((allcart(k,j+nres,i),k=1,3),j=nnt,nct)
2028 end subroutine cartout
2029 !------------------------------------------------------------------------------
2030 ! subroutine alloc_clust_arrays(n_conf)
2035 ! allocate(diss(maxdist)) !(maxdist)
2036 !el allocate(energy(0:maxconf),totfree(0:maxconf)) !(0:maxconf)
2037 ! allocatable :: enetb !(max_ene,maxstr_proc)
2038 !el allocate(entfac(maxconf)) !(maxconf)
2039 ! allocatable :: totfree_gr !(maxgr)
2040 !el allocate(rcutoff(max_cut+1)) !(max_cut+1)
2042 ! allocatable :: licz,iass !(maxgr)
2043 ! allocatable :: nconf !(maxgr,maxingr)
2044 ! allocatable :: iass_tot !(maxgr,max_cut)
2045 ! allocatable :: list_conf !(maxconf)
2047 !el allocatable :: allcart !(3,maxres2,maxstr_proc)
2048 !el allocate(rmstb(maxconf)) !(maxconf)
2049 !el allocate(mult(nres)) !(maxres)
2050 !el allocatable :: nss_all !(maxstr_proc)
2051 !el allocatable :: ihpb_all,jhpb_all !(maxss,maxstr_proc)
2052 ! allocate(icc(n_conf),iscore(n_conf)) !(maxconf)
2055 ! allocatable :: tempfac !(2,maxres)
2058 !el allocate(beta_h(maxT)) !(maxT)
2060 ! end subroutine alloc_clust_arrays
2061 !-----------------------------------------------------------------------------
2062 !-----------------------------------------------------------------------------