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)
567 if (from_cart .and. .not. from_bx .and. .not. from_cx) then
569 read (intin,*,end=13,err=11) energy(icon),totfree(icon),&
571 nss_all(icon),(ihpb_all(ii,icon),jhpb_all(i,icon),&
572 i=1,nss_all(icon)),iscore(icon)
574 read (intin,*,end=13,err=11) energy(icon),rmstb(icon),&
575 nss_all(icon),(ihpb_all(ii,icon),jhpb_all(i,icon),&
576 i=1,nss_all(icon)),iscore(icon)
578 read (intin,'(8f10.5)',end=13,err=10) &
579 ((allcart(j,i,icon),j=1,3),i=1,nres),&
580 ((allcart(j,i+nres,icon),j=1,3),i=nnt,nct)
581 print *,icon,energy(icon),nss_all(icon),rmstb(icon)
583 read(intin,'(a80)',end=13,err=12) lineh
584 read(lineh(:5),*,err=8) ic
586 read(lineh(6:),*,err=8) energy(icon)
588 read(lineh(6:),*,err=8) energy(icon)
592 print *,'error, assuming e=1d10',lineh
596 !old read(lineh(18:),*,end=13,err=11) nss_all(icon)
597 ii = index(lineh(15:)," ")+15
598 read(lineh(ii:),*,end=13,err=11) nss_all(icon)
599 IF (NSS_all(icon).LT.9) THEN
600 read (lineh(20:),*,end=102) &
601 (IHPB_all(I,icon),JHPB_all(I,icon),I=1,NSS_all(icon)),&
604 read (lineh(20:),*,end=102) &
605 (IHPB_all(I,icon),JHPB_all(I,icon),I=1,8)
606 read (intin,*) (IHPB_all(I,icon),JHPB_all(I,icon),&
607 I=9,NSS_all(icon)),iscore(icon)
612 PRINT *,'IC:',IC,' ENERGY:',ENERGY(ICON)
613 call read_angles(intin,*13)
615 phiall(i,icon)=phi(i)
616 thetall(i,icon)=theta(i)
617 alphall(i,icon)=alph(i)
618 omall(i,icon)=omeg(i)
624 ! CALCULATE DISTANCES
626 10 print *,'something wrong with angles'
628 11 print *,'something wrong with NSS',nss
630 12 print *,'something wrong with header'
637 open (icbase,file=bprotfiles,status="unknown",&
638 form="unformatted",access="direct",recl=lenrec)
639 ! Read conformations from binary DA files (one per batch) and write them to
640 ! a binary DA scratchfile.
644 write (liczba,'(bz,i3.3)') me
645 IF (ME.EQ.MASTER) THEN
646 ! Only the master reads the database; it'll send it to the other procs
654 open (intin,file=intinname,status="old",form="unformatted",&
655 access="direct",recl=lenrec_in)
657 else if (from_cx) then
658 #if (defined(AIX) && !defined(JUBL))
659 call xdrfopen_(ixdrf,intinname, "r", iret)
661 call xdrfopen(ixdrf,intinname, "r", iret)
664 write (iout,*) "xdrfopen: iret",iret
666 write (iout,*) "Error: coordinate file ",&
667 intinname(:ilen(intinname))," does not exist."
670 call MPI_ABORT(MPI_COMM_WORLD,IERROR,ERRCODE)
675 write (iout,*) "Error: coordinate format not specified"
678 call MPI_ABORT(MPI_COMM_WORLD,IERROR,ERRCODE)
686 write (iout,*) "Opening file ",intinname(:ilen(intinname))
687 write (iout,*) "lenrec",lenrec_in
691 ! write (iout,*) "maxconf",maxconf
695 !el if (i.gt.maxconf) then
696 !el write (iout,*) "Error: too many conformations ",&
697 !el "(",maxconf,") maximum."
699 !el call MPI_Abort(MPI_COMM_WORLD,errcode,ierror)
703 ! write (iout,*) "i",i
706 read(intin,err=101,end=101) &
707 ((csingle(l,k),l=1,3),k=1,nres),&
708 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
709 nss,(ihpb(k),jhpb(k),k=1,nss),&
711 entfac(jj+1),rmstb(jj+1),iscor
718 #if (defined(AIX) && !defined(JUBL))
719 call xdrf3dfcoord_(ixdrf, csingle, itmp, prec, iret)
720 if (iret.eq.0) goto 101
721 call xdrfint_(ixdrf, nss, iret)
722 if (iret.eq.0) goto 101
724 call xdrfint_(ixdrf, ihpb(j), iret)
725 if (iret.eq.0) goto 101
726 call xdrfint_(ixdrf, jhpb(j), iret)
727 if (iret.eq.0) goto 101
729 call xdrffloat_(ixdrf,reini,iret)
730 if (iret.eq.0) goto 101
731 call xdrffloat_(ixdrf,refree,iret)
732 if (iret.eq.0) goto 101
733 call xdrffloat_(ixdrf,rmsdev,iret)
734 if (iret.eq.0) goto 101
735 call xdrfint_(ixdrf,iscor,iret)
736 if (iret.eq.0) goto 101
738 ! write (iout,*) "calling xdrf3dfcoord"
739 call xdrf3dfcoord(ixdrf, csingle, itmp, prec, iret)
740 ! write (iout,*) "iret",iret
742 if (iret.eq.0) goto 101
743 call xdrfint(ixdrf, nss, iret)
744 ! write (iout,*) "iret",iret
745 ! write (iout,*) "nss",nss
747 if (iret.eq.0) goto 101
749 call xdrfint(ixdrf, ihpb(k), iret)
750 if (iret.eq.0) goto 101
751 call xdrfint(ixdrf, jhpb(k), iret)
752 if (iret.eq.0) goto 101
754 call xdrffloat(ixdrf,reini,iret)
755 if (iret.eq.0) goto 101
756 call xdrffloat(ixdrf,refree,iret)
757 if (iret.eq.0) goto 101
758 call xdrffloat(ixdrf,rmsdev,iret)
759 if (iret.eq.0) goto 101
760 call xdrfint(ixdrf,iscor,iret)
761 if (iret.eq.0) goto 101
773 c(l,nres+k)=csingle(l,nres+k-nnt+1)
778 write (iout,'(5hREAD ,i5,3f15.4,i10)') &
779 jj+1,energy(jj+1),entfac(jj+1),&
781 write (iout,*) "Conformation",jjj+1,jj+1
782 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
783 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
786 call add_new_cconf(jjj,jj,jj_old,icount,Next)
789 write (iout,*) i-1," conformations read from DA file ",&
790 intinname(:ilen(intinname))
791 write (iout,*) jj," conformations read so far"
795 #if (defined(AIX) && !defined(JUBL))
796 call xdrfclose_(ixdrf, iret)
798 call xdrfclose(ixdrf, iret)
803 write (iout,*) "jj_old",jj_old," jj",jj
805 call write_and_send_cconf(icount,jj_old,jj,Next)
806 call MPI_Send(0,1,MPI_INTEGER,Next,570,&
807 MPI_COMM_WORLD,IERROR)
810 call write_and_send_cconf(icount,jj_old,jj,Next)
812 t_acq = tcpu() - t_acq
814 write (iout,*) "Processor",me,&
815 " time for conformation read/send",t_acq
817 ! A worker gets the confs from the master and sends them to its neighbor
819 call receive_and_pass_cconf(icount,jj_old,jj,&
821 t_acq = tcpu() - t_acq
828 write(iout,*)"A total of",ncon," conformations read."
830 allocate(enetb(1:max_ene,ncon)) !(max_ene,maxstr_proc)
832 ! Check if everyone has the same number of conformations
833 call MPI_Allgather(ncon,1,MPI_INTEGER,&
834 ntot_all(0),1,MPI_INTEGER,MPI_Comm_World,IERROR)
838 if (ncon.ne.ntot_all(i)) then
839 write (iout,*) "Number of conformations at processor",i,&
840 " differs from that at processor",me,&
848 write (iout,*) "Number of conformations read by processors"
851 write (iout,'(8i10)') i,ntot_all(i)
853 write (iout,*) "Calculation terminated."
859 1111 write(iout,*) "Error opening coordinate file ",&
860 intinname(:ilen(intinname))
863 end subroutine read_coords
864 !------------------------------------------------------------------------------
865 subroutine add_new_cconf(jjj,jj,jj_old,icount,Next)
867 use geometry_data, only: vbld,rad2deg,theta,phi,alph,omeg,deg2rad
868 use energy_data, only: itel,itype,dsc,max_ene,molnum
869 use control_data, only: symetr
870 use geometry, only: int_from_cart1
872 ! include "DIMENSIONS"
873 ! include "sizesclu.dat"
874 ! include "COMMON.CLUSTER"
875 ! include "COMMON.CONTROL"
876 ! include "COMMON.CHAIN"
877 ! include "COMMON.INTERACT"
878 ! include "COMMON.LOCAL"
879 ! include "COMMON.IOUNITS"
880 ! include "COMMON.NAMES"
881 ! include "COMMON.VAR"
882 ! include "COMMON.SBRIDGE"
883 ! include "COMMON.GEO"
884 integer :: i,j,jj,jjj,jj_old,icount,k,kk,l,ii,ib,&
885 nn,nn1,inan,Next,itj,chalen,mnum
886 real(kind=8) :: etot,energia(0:max_ene)
888 chalen=int((nct-nnt+2)/symetr)
889 call int_from_cart1(.false.)
892 write (iout,*) "Check atom",j
894 if (itype(j,mnum).eq.ntyp1_molec(mnum)) cycle
895 if (itype(j+1,molnum(j+1)).eq.ntyp1_molec(molnum(j+1))) cycle
897 if ((vbld(j).lt.2.0d0 .or. vbld(j).gt.5.0d0).and.(mnum.ne.5)) then
899 if (itel(j).ne.0 .and. itel(j-1).ne.0) then
900 write (iout,*) "Conformation",jjj,jj+1
901 write (iout,*) "Bad CA-CA bond length",j," ",vbld(j),itel(j),&
903 write (iout,*) "The Cartesian geometry is:"
904 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
905 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
906 write (iout,*) "The internal geometry is:"
907 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
908 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
909 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
910 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
911 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
912 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
914 "This conformation WILL NOT be added to the database."
924 if (itype(j,1).ne.10 .and. (itype(j,mnum).ne.ntyp1_molec(mnum)) &
925 .and. (vbld(nres+j)-dsc(iabs(itj))) &
927 write (iout,*) "Conformation",jjj,jj+1
928 write (iout,*) "Bad CA-SC bond length",j," ",vbld(nres+j)
929 write (iout,*) "The Cartesian geometry is:"
930 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
931 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
932 write (iout,*) "The internal geometry is:"
933 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
934 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
935 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
936 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
937 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
938 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
940 "This conformation WILL NOT be added to the database."
945 if (theta(j).le.0.0d0) then
947 "Zero theta angle(s) in conformation",jjj,jj+1
948 write (iout,*) "The Cartesian geometry is:"
949 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
950 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
951 write (iout,*) "The internal geometry is:"
952 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
953 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
954 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
955 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
956 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
957 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
959 "This conformation WILL NOT be added to the database."
962 if (theta(j).gt.179.97*deg2rad) theta(j)=179.97*deg2rad
966 write (iout,*) "Conformation",jjj,jj
967 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
968 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
969 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
970 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
971 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
972 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
973 write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct)
974 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
975 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
976 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
977 write (iout,'(e15.5,16i5)') entfac(icount+1)
978 ! & iscore(icount+1,0)
981 call store_cconf_from_file(jj,icount)
982 if (icount.eq.maxstr_proc) then
984 write (iout,* ) "jj_old",jj_old," jj",jj
986 call write_and_send_cconf(icount,jj_old,jj,Next)
991 end subroutine add_new_cconf
992 !------------------------------------------------------------------------------
993 subroutine store_cconf_from_file(jj,icount)
995 use energy_data, only: ihpb,jhpb
997 ! include "DIMENSIONS"
998 ! include "sizesclu.dat"
999 ! include "COMMON.CLUSTER"
1000 ! include "COMMON.CHAIN"
1001 ! include "COMMON.SBRIDGE"
1002 ! include "COMMON.INTERACT"
1003 ! include "COMMON.IOUNITS"
1004 ! include "COMMON.VAR"
1005 integer :: i,j,jj,icount
1006 ! Store the conformation that has been read in
1009 allcart(j,i,icount)=c(j,i)
1014 ihpb_all(i,icount)=ihpb(i)
1015 jhpb_all(i,icount)=jhpb(i)
1018 end subroutine store_cconf_from_file
1019 !------------------------------------------------------------------------------
1020 subroutine write_and_send_cconf(icount,jj_old,jj,Next)
1023 ! include "DIMENSIONS"
1024 ! include "sizesclu.dat"
1029 ! include "COMMON.MPI"
1031 ! include "COMMON.CHAIN"
1032 ! include "COMMON.SBRIDGE"
1033 ! include "COMMON.INTERACT"
1034 ! include "COMMON.IOUNITS"
1035 ! include "COMMON.CLUSTER"
1036 ! include "COMMON.VAR"
1037 integer :: icount,jj_old,jj,Next
1038 ! Write the structures to a scratch file
1040 ! Master sends the portion of conformations that have been read in to the neighbor
1042 write (iout,*) "Processor",me," entered WRITE_AND_SEND_CONF"
1045 call MPI_Send(icount,1,MPI_INTEGER,Next,570,MPI_COMM_WORLD,IERROR)
1046 call MPI_Send(nss_all(1),icount,MPI_INTEGER,&
1047 Next,571,MPI_COMM_WORLD,IERROR)
1048 call MPI_Send(ihpb_all(1,1),icount,MPI_INTEGER,&
1049 Next,572,MPI_COMM_WORLD,IERROR)
1050 call MPI_Send(jhpb_all(1,1),icount,MPI_INTEGER,&
1051 Next,573,MPI_COMM_WORLD,IERROR)
1052 call MPI_Send(rmstb(jj_old),icount,MPI_DOUBLE_PRECISION,&
1053 Next,577,MPI_COMM_WORLD,IERROR)
1054 call MPI_Send(entfac(jj_old),icount,MPI_DOUBLE_PRECISION,&
1055 Next,579,MPI_COMM_WORLD,IERROR)
1056 call MPI_Send(allcart(1,1,1),3*icount*2*nres,&
1057 MPI_REAL,Next,580,MPI_COMM_WORLD,IERROR)
1059 call dawrite_ccoords(jj_old,jj,icbase)
1061 end subroutine write_and_send_cconf
1062 !------------------------------------------------------------------------------
1064 subroutine receive_and_pass_cconf(icount,jj_old,jj,Previous,Next)
1068 ! include "DIMENSIONS"
1069 ! include "sizesclu.dat"
1071 integer :: IERROR !,STATUS(MPI_STATUS_SIZE)
1072 ! include "COMMON.MPI"
1073 ! include "COMMON.CHAIN"
1074 ! include "COMMON.SBRIDGE"
1075 ! include "COMMON.INTERACT"
1076 ! include "COMMON.IOUNITS"
1077 ! include "COMMON.VAR"
1078 ! include "COMMON.GEO"
1079 ! include "COMMON.CLUSTER"
1080 integer :: i,j,k,icount,jj_old,jj,Previous,Next
1083 write (iout,*) "Processor",me," entered RECEIVE_AND_PASS_CONF"
1086 do while (icount.gt.0)
1087 call MPI_Recv(icount,1,MPI_INTEGER,Previous,570,MPI_COMM_WORLD,&
1089 call MPI_Send(icount,1,MPI_INTEGER,Next,570,MPI_COMM_WORLD,&
1092 write (iout,*) "Processor",me," icount",icount
1094 if (icount.eq.0) return
1095 call MPI_Recv(nss_all(1),icount,MPI_INTEGER,&
1096 Previous,571,MPI_COMM_WORLD,STATUS,IERROR)
1097 call MPI_Send(nss_all(1),icount,MPI_INTEGER,&
1098 Next,571,MPI_COMM_WORLD,IERROR)
1099 call MPI_Recv(ihpb_all(1,1),icount,MPI_INTEGER,&
1100 Previous,572,MPI_COMM_WORLD,STATUS,IERROR)
1101 call MPI_Send(ihpb_all(1,1),icount,MPI_INTEGER,&
1102 Next,572,MPI_COMM_WORLD,IERROR)
1103 call MPI_Recv(jhpb_all(1,1),icount,MPI_INTEGER,&
1104 Previous,573,MPI_COMM_WORLD,STATUS,IERROR)
1105 call MPI_Send(jhpb_all(1,1),icount,MPI_INTEGER,&
1106 Next,573,MPI_COMM_WORLD,IERROR)
1107 call MPI_Recv(rmstb(jj_old),icount,MPI_DOUBLE_PRECISION,&
1108 Previous,577,MPI_COMM_WORLD,STATUS,IERROR)
1109 call MPI_Send(rmstb(jj_old),icount,MPI_DOUBLE_PRECISION,&
1110 Next,577,MPI_COMM_WORLD,IERROR)
1111 call MPI_Recv(entfac(jj_old),icount,MPI_DOUBLE_PRECISION,&
1112 Previous,579,MPI_COMM_WORLD,STATUS,IERROR)
1113 call MPI_Send(entfac(jj_old),icount,MPI_DOUBLE_PRECISION,&
1114 Next,579,MPI_COMM_WORLD,IERROR)
1115 call MPI_Recv(allcart(1,1,1),3*icount*2*nres,&
1116 MPI_REAL,Previous,580,MPI_COMM_WORLD,STATUS,IERROR)
1117 call MPI_Send(allcart(1,1,1),3*icount*2*nres,&
1118 MPI_REAL,Next,580,MPI_COMM_WORLD,IERROR)
1120 call dawrite_ccoords(jj_old,jj,icbase)
1123 write (iout,*) "Processor",me," received",icount," conformations"
1125 write (iout,'(8f10.4)') (allcart(l,k,i),l=1,3,k=1,nres)
1126 write (iout,'(8f10.4)')((allcart(l,k,i+nres),l=1,3,k=nnt,nct)
1127 write (iout,'(e15.5,16i5)') entfac(i)
1132 end subroutine receive_and_pass_cconf
1134 !------------------------------------------------------------------------------
1135 subroutine daread_ccoords(istart_conf,iend_conf)
1138 ! include "DIMENSIONS"
1139 ! include "sizesclu.dat"
1143 ! include "COMMON.MPI"
1145 ! include "COMMON.CHAIN"
1146 ! include "COMMON.CLUSTER"
1147 ! include "COMMON.IOUNITS"
1148 ! include "COMMON.INTERACT"
1149 ! include "COMMON.VAR"
1150 ! include "COMMON.SBRIDGE"
1151 ! include "COMMON.GEO"
1152 integer :: istart_conf,iend_conf
1153 integer :: i,j,ij,ii,iii
1155 character(len=16) :: form,acc
1156 character(len=32) :: nam
1158 ! Read conformations off a DA scratchfile.
1161 write (iout,*) "DAREAD_COORDS"
1162 write (iout,*) "istart_conf",istart_conf," iend_conf",iend_conf
1163 inquire(unit=icbase,name=nam,recl=len,form=form,access=acc)
1164 write (iout,*) "len=",len," form=",form," acc=",acc
1165 write (iout,*) "nam=",nam
1168 do ii=istart_conf,iend_conf
1169 ij = ii - istart_conf + 1
1172 write (iout,*) "Reading binary file, record",iii," ii",ii
1175 read(icbase,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres),&
1176 ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),&
1177 nss_all(ij),(ihpb_all(i,ij),jhpb_all(i,ij),i=1,nss),&
1178 entfac(ii),rmstb(ii)
1180 write (iout,*) ii,iii,ij,entfac(ii)
1181 write (iout,'(8f10.5)') ((allcart(j,i,ij),j=1,3),i=1,nres)
1182 write (iout,'(8f10.4)') ((allcart(j,i,ij),j=1,3),&
1183 i=nnt+nres,nct+nres)
1184 write (iout,'(2e15.5)') entfac(ij)
1185 write (iout,'(16i5)') nss_all(ij),(ihpb_all(i,ij),&
1186 jhpb_all(i,ij),i=1,nss)
1191 end subroutine daread_ccoords
1192 !------------------------------------------------------------------------------
1193 subroutine dawrite_ccoords(istart_conf,iend_conf,unit_out)
1196 ! include "DIMENSIONS"
1197 ! include "sizesclu.dat"
1201 ! include "COMMON.MPI"
1203 ! include "COMMON.CHAIN"
1204 ! include "COMMON.INTERACT"
1205 ! include "COMMON.IOUNITS"
1206 ! include "COMMON.VAR"
1207 ! include "COMMON.SBRIDGE"
1208 ! include "COMMON.GEO"
1209 ! include "COMMON.CLUSTER"
1210 integer :: istart_conf,iend_conf
1211 integer :: i,j,ii,ij,iii,unit_out
1213 character(len=16) :: form,acc
1214 character(len=32) :: nam
1216 ! Write conformations to a DA scratchfile.
1219 write (iout,*) "DAWRITE_COORDS"
1220 write (iout,*) "istart_conf",istart_conf," iend_conf",iend_conf
1221 write (iout,*) "lenrec",lenrec
1222 inquire(unit=unit_out,name=nam,recl=len,form=form,access=acc)
1223 write (iout,*) "len=",len," form=",form," acc=",acc
1224 write (iout,*) "nam=",nam
1227 do ii=istart_conf,iend_conf
1229 ij = ii - istart_conf + 1
1231 write (iout,*) "Writing binary file, record",iii," ii",ii
1234 write(unit_out,rec=iii) ((allcart(j,i,ij),j=1,3),i=1,nres),&
1235 ((allcart(j,i,ij),j=1,3),i=nnt+nres,nct+nres),&
1236 nss_all(ij),(ihpb_all(i,ij),jhpb_all(i,ij),i=1,nss_all(ij)),&
1237 entfac(ii),rmstb(ii)
1239 write (iout,'(8f10.5)') ((allcart(j,i,ij),j=1,3),i=1,nres)
1240 write (iout,'(8f10.4)') ((allcart(j,i,ij),j=1,3),i=nnt+nres,&
1242 write (iout,'(2e15.5)') entfac(ij)
1243 write (iout,'(16i5)') nss_all(ij),(ihpb(i,ij),jhpb(i,ij),i=1,&
1249 end subroutine dawrite_ccoords
1250 !-----------------------------------------------------------------------------
1252 !-----------------------------------------------------------------------------
1253 subroutine read_control
1255 ! Read molecular data
1257 use energy_data, only: rescale_mode,distchainmax,ipot !,temp0
1258 use control_data, only: titel,outpdb,outmol2,refstr,pdbref,&
1259 iscode,symetr,punch_dist,print_dist,nstart,nend,&
1260 caonly,iopt,efree,lprint_cart,lprint_int,rlamb_ele,&
1261 r_cut_ele,nclust,tor_mode,scelemode
1263 ! include 'DIMENSIONS'
1264 ! include 'sizesclu.dat'
1265 ! include 'COMMON.IOUNITS'
1266 ! include 'COMMON.TIME1'
1267 ! include 'COMMON.SBRIDGE'
1268 ! include 'COMMON.CONTROL'
1269 ! include 'COMMON.CLUSTER'
1270 ! include 'COMMON.CHAIN'
1271 ! include 'COMMON.HEADER'
1272 ! include 'COMMON.FFIELD'
1273 ! include 'COMMON.FREE'
1274 ! include 'COMMON.INTERACT'
1275 character(len=320) :: controlcard !,ucase
1277 ! include 'COMMON.INFO'
1281 read (INP,'(a80)') titel
1282 call card_concat(controlcard,.true.)
1284 call readi(controlcard,'NRES',nres_molec(1),0)
1285 call readi(controlcard,"NRES_NUCL",nres_molec(2),0)
1286 call readi(controlcard,"NRES_CAT",nres_molec(5),0)
1289 nres=nres_molec(i)+nres
1292 ! call alloc_clust_arrays
1293 allocate(rcutoff(max_cut+1)) !(max_cut+1)
1294 allocate(beta_h(maxT)) !(maxT)
1295 allocate(mult(nres)) !(maxres)
1298 call readi(controlcard,'RESCALE',rescale_mode,2)
1299 call reada(controlcard,'DISTCHAINMAX',distchainmax,50.0d0)
1300 write (iout,*) "DISTCHAINMAX",distchainmax
1301 call readi(controlcard,'PDBOUT',outpdb,0)
1302 call readi(controlcard,'MOL2OUT',outmol2,0)
1303 refstr=(index(controlcard,'REFSTR').gt.0)
1304 write (iout,*) "REFSTR",refstr
1305 pdbref=(index(controlcard,'PDBREF').gt.0)
1306 iscode=index(controlcard,'ONE_LETTER')
1307 tree=(index(controlcard,'MAKE_TREE').gt.0)
1308 min_var=(index(controlcard,'MINVAR').gt.0)
1309 plot_tree=(index(controlcard,'PLOT_TREE').gt.0)
1310 punch_dist=(index(controlcard,'PUNCH_DIST').gt.0)
1311 call readi(controlcard,'NCUT',ncut,1)
1312 call readi(controlcard,'SYM',symetr,1)
1313 write (iout,*) 'sym', symetr
1314 call readi(controlcard,'NSTART',nstart,0)
1315 call reada(controlcard,'BOXX',boxxsize,100.0d0)
1316 call reada(controlcard,'BOXY',boxysize,100.0d0)
1317 call reada(controlcard,'BOXZ',boxzsize,100.0d0)
1318 call readi(controlcard,'NCLUST',nclust,5)
1319 ! ions=index(controlcard,"IONS").gt.0
1320 call readi(controlcard,"SCELEMODE",scelemode,0)
1321 print *,"SCELE",scelemode
1322 call readi(controlcard,'TORMODE',tor_mode,0)
1323 !C if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0) then
1324 write(iout,*) "torsional and valence angle mode",tor_mode
1325 call reada(controlcard,"R_CUT_ELE",r_cut_ele,15.0d0)
1326 call reada(controlcard,"LAMBDA_ELE",rlamb_ele,0.3d0)
1327 write(iout,*) "R_CUT_ELE=",r_cut_ele
1328 call readi(controlcard,'NEND',nend,0)
1329 call reada(controlcard,'ECUT',ecut,10.0d0)
1330 call reada(controlcard,'PROB',prob_limit,0.99d0)
1331 write (iout,*) "Probability limit",prob_limit
1332 lgrp=(index(controlcard,'LGRP').gt.0)
1333 caonly=(index(controlcard,'CA_ONLY').gt.0)
1334 print_dist=(index(controlcard,'PRINT_DIST').gt.0)
1335 call multreada(controlcard,'CUTOFF',rcutoff,ncut,-1.0d0)
1336 call readi(controlcard,'IOPT',iopt,2)
1337 lside = index(controlcard,"SIDE").gt.0
1338 efree = index(controlcard,"EFREE").gt.0
1339 call readi(controlcard,'NTEMP',nT,1)
1340 write (iout,*) "nT",nT
1341 !el call reada(controlcard,'TEMP0',temp0,300.0d0) !el
1342 call multreada(controlcard,'TEMPER',beta_h,nT,300.0d0)
1343 write (iout,*) "nT",nT
1344 write (iout,*) 'beta_h',(beta_h(i),i=1,nT)
1346 beta_h(i)=1.0d0/(1.987D-3*beta_h(i))
1348 write (iout,*) 'beta_h',(beta_h(i),i=1,nT)
1349 lprint_cart=index(controlcard,"PRINT_CART") .gt.0
1350 lprint_int=index(controlcard,"PRINT_INT") .gt.0
1353 end subroutine read_control
1354 !-----------------------------------------------------------------------------
1357 ! Read molecular data.
1359 use geometry_data, only: nsup,cref,nres0,nstart_sup,nstart_seq,dc
1360 use energy_data!, only: wsc,wscp,welec,wbond,wstrain,wtor,wtor_d,&
1361 ! wang,wscloc,wcorr,wcorr4,wcorr5,wcorr6,wsccor,wel_loc,&
1362 ! wturn3,wturn4,wturn6,wvdwpp,weights
1363 use control_data, only: titel,nstart,nend,pdbref,refstr,iscode,&
1365 use geometry, only: chainbuild,alloc_geo_arrays
1366 use energy, only: alloc_ener_arrays
1367 use io_base, only: rescode
1368 use control, only: setup_var,init_int_table
1369 use conform_compar, only: contact
1371 ! include 'DIMENSIONS'
1372 ! include 'COMMON.IOUNITS'
1373 ! include 'COMMON.GEO'
1374 ! include 'COMMON.VAR'
1375 ! include 'COMMON.INTERACT'
1376 ! include 'COMMON.LOCAL'
1377 ! include 'COMMON.NAMES'
1378 ! include 'COMMON.CHAIN'
1379 ! include 'COMMON.FFIELD'
1380 ! include 'COMMON.SBRIDGE'
1381 ! include 'COMMON.HEADER'
1382 ! include 'COMMON.CONTROL'
1383 ! include 'COMMON.CONTACTS'
1384 ! include 'COMMON.TIME1'
1386 ! include 'COMMON.INFO'
1388 character(len=4) :: sequence(nres) !(maxres)
1389 character(len=800) :: weightcard
1391 real(kind=8) :: x(6*nres) !(maxvar)
1392 integer :: itype_pdb(nres) !(maxres)
1394 integer :: i,j,kkk,mnum
1398 !el allocate(weights(n_ene))
1399 allocate(weights(max_ene))
1400 call alloc_geo_arrays
1401 call alloc_ener_arrays
1402 !-----------------------------
1403 allocate(c(3,2*nres+2)) !(3,maxres2+2) maxres2=2*maxres
1404 allocate(dc(3,0:2*nres+2)) !(3,0:maxres2)
1405 allocate(itype(nres+2,5)) !(maxres)
1406 allocate(molnum(nres+1))
1407 allocate(itel(nres+2))
1421 !--------------------------
1422 ! Read weights of the subsequent energy terms.
1423 call card_concat(weightcard,.true.)
1424 call reada(weightcard,'WSC',wsc,1.0d0)
1425 call reada(weightcard,'WLONG',wsc,wsc)
1426 call reada(weightcard,'WSCP',wscp,1.0d0)
1427 call reada(weightcard,'WELEC',welec,1.0D0)
1428 call reada(weightcard,'WVDWPP',wvdwpp,welec)
1429 call reada(weightcard,'WEL_LOC',wel_loc,1.0D0)
1430 call reada(weightcard,'WCORR4',wcorr4,0.0D0)
1431 call reada(weightcard,'WCORR5',wcorr5,0.0D0)
1432 call reada(weightcard,'WCORR6',wcorr6,0.0D0)
1433 call reada(weightcard,'WTURN3',wturn3,1.0D0)
1434 call reada(weightcard,'WTURN4',wturn4,1.0D0)
1435 call reada(weightcard,'WTURN6',wturn6,1.0D0)
1436 call reada(weightcard,'WSTRAIN',wstrain,1.0D0)
1437 call reada(weightcard,'WSCCOR',wsccor,1.0D0)
1438 call reada(weightcard,'WBOND',wbond,1.0D0)
1439 call reada(weightcard,'WTOR',wtor,1.0D0)
1440 call reada(weightcard,'WTORD',wtor_d,1.0D0)
1441 call reada(weightcard,'WANG',wang,1.0D0)
1442 call reada(weightcard,'WSCLOC',wscloc,1.0D0)
1443 call reada(weightcard,'SCAL14',scal14,0.4D0)
1444 call reada(weightcard,'SCALSCP',scalscp,1.0d0)
1445 call reada(weightcard,'CUTOFF',cutoff_corr,7.0d0)
1446 call reada(weightcard,'DELT_CORR',delt_corr,0.5d0)
1447 call reada(weightcard,'WCATCAT',wcatcat,0.0d0)
1448 call reada(weightcard,'WCATPROT',wcatprot,0.0d0)
1449 call reada(weightcard,'TEMP0',temp0,300.0d0) !!! el
1450 if (index(weightcard,'SOFT').gt.0) ipot=6
1451 call reada(weightcard,'WVDWPP_NUCL',wvdwpp_nucl,0.0D0)
1452 call reada(weightcard,'WELPP',welpp,0.0d0)
1453 call reada(weightcard,'WVDWPSB',wvdwpsb,0.0d0)
1454 call reada(weightcard,'WELPSB',welpsb,0.0D0)
1455 call reada(weightcard,'WVDWSB',wvdwsb,0.0d0)
1456 call reada(weightcard,'WELSB',welsb,0.0D0)
1457 call reada(weightcard,'WBOND_NUCL',wbond_nucl,0.0D0)
1458 call reada(weightcard,'WANG_NUCL',wang_nucl,0.0D0)
1459 call reada(weightcard,'WSBLOC',wsbloc,0.0D0)
1460 call reada(weightcard,'WTOR_NUCL',wtor_nucl,0.0D0)
1461 ! print *,"KUR..",wtor_nucl
1462 call reada(weightcard,'WTORD_NUCL',wtor_d_nucl,0.0D0)
1463 call reada(weightcard,'WCORR_NUCL',wcorr_nucl,0.0D0)
1464 call reada(weightcard,'WCORR3_NUC',wcorr3_nucl,0.0D0)
1465 ! 12/1/95 Added weight for the multi-body term WCORR
1466 call reada(weightcard,'WCORRH',wcorr,1.0D0)
1467 call reada(weightcard,'WSCBASE',wscbase,0.0D0)
1468 call reada(weightcard,'WPEPBASE',wpepbase,1.0d0)
1469 call reada(weightcard,'WSCPHO',wscpho,0.0d0)
1470 call reada(weightcard,'WPEPPHO',wpeppho,0.0d0)
1473 if (wcorr4.gt.0.0d0) wcorr=wcorr4
1492 !el weights(19)=wsccor !!!!!!!!!!!!!!!!
1494 weights(26)=wvdwpp_nucl
1500 weights(32)=wbond_nucl
1501 weights(33)=wang_nucl
1503 weights(35)=wtor_nucl
1504 weights(36)=wtor_d_nucl
1505 weights(37)=wcorr_nucl
1506 weights(38)=wcorr3_nucl
1508 weights(42)=wcatprot
1514 write (iout,10) wsc,wscp,welec,wvdwpp,wbond,wang,wscloc,wtor,&
1515 wtor_d,wstrain,wel_loc,wcorr,wcorr5,wcorr6,wturn3,&
1516 wturn4,wturn6,wsccor
1517 10 format (/'Energy-term weights (unscaled):'// &
1518 'WSCC= ',f10.6,' (SC-SC)'/ &
1519 'WSCP= ',f10.6,' (SC-p)'/ &
1520 'WELEC= ',f10.6,' (p-p electr)'/ &
1521 'WVDWPP= ',f10.6,' (p-p VDW)'/ &
1522 'WBOND= ',f10.6,' (stretching)'/ &
1523 'WANG= ',f10.6,' (bending)'/ &
1524 'WSCLOC= ',f10.6,' (SC local)'/ &
1525 'WTOR= ',f10.6,' (torsional)'/ &
1526 'WTORD= ',f10.6,' (double torsional)'/ &
1527 'WSTRAIN=',f10.6,' (SS bridges & dist. cnstr.)'/ &
1528 'WEL_LOC=',f10.6,' (multi-body 3-rd order)'/ &
1529 'WCORR4= ',f10.6,' (multi-body 4th order)'/ &
1530 'WCORR5= ',f10.6,' (multi-body 5th order)'/ &
1531 'WCORR6= ',f10.6,' (multi-body 6th order)'/ &
1532 'WTURN3= ',f10.6,' (turns, 3rd order)'/ &
1533 'WTURN4= ',f10.6,' (turns, 4th order)'/ &
1534 'WTURN6= ',f10.6,' (turns, 6th order)'/ &
1535 'WSCCOR= ',f10.6,' (SC-backbone torsinal correalations)')
1537 if (wcorr4.gt.0.0d0) then
1538 write (iout,'(/2a/)') 'Local-electrostatic type correlation ',&
1539 'between contact pairs of peptide groups'
1540 write (iout,'(2(a,f5.3/))') &
1541 'Cutoff on 4-6th order correlation terms: ',cutoff_corr,&
1542 'Range of quenching the correlation terms:',2*delt_corr
1543 else if (wcorr.gt.0.0d0) then
1544 write (iout,'(/2a/)') 'Hydrogen-bonding correlation ',&
1545 'between contact pairs of peptide groups'
1547 write (iout,'(a,f8.3)') &
1548 'Scaling factor of 1,4 SC-p interactions:',scal14
1549 write (iout,'(a,f8.3)') &
1550 'General scaling factor of SC-p interactions:',scalscp
1551 r0_corr=cutoff_corr-delt_corr
1553 aad(i,1)=scalscp*aad(i,1)
1554 aad(i,2)=scalscp*aad(i,2)
1555 bad(i,1)=scalscp*bad(i,1)
1556 bad(i,2)=scalscp*bad(i,2)
1560 print *,'indpdb=',indpdb,' pdbref=',pdbref
1562 ! Read sequence if not taken from the pdb file.
1563 if (iscode.gt.0) then
1564 read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
1566 read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
1568 ! Convert sequence to numeric code
1569 do i=1,nres_molec(1)
1572 itype(i,mnum)=rescode(i,sequence(i),iscode,mnum)
1574 do i=nres_molec(1)+1,nres_molec(1)+nres_molec(2)
1577 itype(i,mnum)=rescode(i,sequence(i),iscode,mnum)
1579 do i=nres_molec(1)+nres_molec(2)+1,nres_molec(1)+nres_molec(2)+nres_molec(5)
1582 itype(i,mnum)=rescode(i,sequence(i),iscode,mnum)
1585 print '(20i4)',(itype(i,molnum(i)),i=1,nres)
1589 if (itype(i,molnum(i)).eq.ntyp1_molec(molnum(i))&
1590 .or. itype(i+1,molnum(i+1)).eq.ntyp1_molec(molnum(i+1))) then
1592 if (itype(i,molnum(i)).eq.ntyp1_molec(molnum(i))) then
1596 else if (iabs(itype(i+1,1)).ne.20) then
1598 else if (iabs(itype(i,1)).ne.20) then
1605 write (iout,*) "ITEL"
1607 write (iout,*) i,itype(i,molnum(i)),itel(i)
1610 print *,'Call Read_Bridge.'
1614 print *,'NNT=',NNT,' NCT=',NCT
1615 if (itype(1,molnum(1)).eq.ntyp1_molec(molnum(1))) nnt=2
1616 if (itype(nres,molnum(nres)).eq.ntyp1_molec(molnum(nres))) nct=nct-1
1617 if (nstart.lt.nnt) nstart=nnt
1618 if (nend.gt.nct .or. nend.eq.0) nend=nct
1619 write (iout,*) "nstart",nstart," nend",nend
1622 ! read(inp,'(a)') pdbfile
1623 ! write (iout,'(2a)') 'PDB data will be read from file ',pdbfile
1624 ! open(ipdbin,file=pdbfile,status='old',err=33)
1626 ! 33 write (iout,'(a)') 'Error opening PDB file.'
1629 ! print *,'Begin reading pdb data'
1631 ! print *,'Finished reading pdb data'
1632 ! write (iout,'(a,i3,a,i3)')'nsup=',nsup,' nstart_sup=',nstart_sup
1634 ! itype_pdb(i)=itype(i)
1637 ! write (iout,'(a,i3)') 'nsup=',nsup
1639 ! if (nsup.le.(nct-nnt+1)) then
1640 ! do i=0,nct-nnt+1-nsup
1641 ! if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),nsup)) then
1646 ! write (iout,'(a)')
1647 ! & 'Error - sequences to be superposed do not match.'
1650 ! do i=0,nsup-(nct-nnt+1)
1651 ! if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),nct-nnt+1))
1653 ! nstart_sup=nstart_sup+i
1658 ! write (iout,'(a)')
1659 ! & 'Error - sequences to be superposed do not match.'
1662 ! write (iout,*) 'nsup=',nsup,' nstart_sup=',nstart_sup,
1663 ! & ' nstart_seq=',nstart_seq
1665 write(iout,*)"przed ini_int_tab"
1667 write(iout,*)"po ini_int_tab"
1668 write(iout,*)"przed setup var"
1670 write(iout,*)"po setup var"
1671 write (iout,*) "molread: REFSTR",refstr
1673 if (.not.pdbref) then
1674 call read_angles(inp,*38)
1676 38 write (iout,'(a)') 'Error reading reference structure.'
1678 call mp_stopall(Error_Msg)
1680 stop 'Error reading reference structure'
1689 cref(j,i,kkk)=c(j,i)
1693 call contact(.true.,ncont_ref,icont_ref)
1696 end subroutine molread
1697 !-----------------------------------------------------------------------------
1698 subroutine openunits
1700 ! include 'DIMENSIONS'
1701 use control_data, only: from_cx,from_bx,from_cart
1705 character(len=3) :: liczba
1706 ! include "COMMON.MPI"
1708 ! include 'COMMON.IOUNITS'
1709 ! include 'COMMON.CONTROL'
1710 integer :: lenpre,lenpot !,ilen
1712 character(len=16) :: cformat,cprint
1713 ! character(len=16) ucase
1714 integer :: lenint,lenout
1715 call getenv('INPUT',prefix)
1716 call getenv('OUTPUT',prefout)
1717 call getenv('INTIN',prefintin)
1718 call getenv('COORD',cformat)
1719 call getenv('PRINTCOOR',cprint)
1720 call getenv('SCRATCHDIR',scratchdir)
1723 if (index(ucase(cformat),'CX').gt.0) then
1729 lenout=ilen(prefout)
1730 lenint=ilen(prefintin)
1731 ! Get the names and open the input files
1732 open (inp,file=prefix(:ilen(prefix))//'.inp',status='old')
1734 write (liczba,'(bz,i3.3)') me
1735 outname=prefout(:lenout)//'_clust.out_'//liczba
1737 outname=prefout(:lenout)//'_clust.out'
1740 intinname=prefintin(:lenint)//'.bx'
1741 else if (from_cx) then
1742 intinname=prefintin(:lenint)//'.cx'
1744 intinname=prefintin(:lenint)//'.int'
1746 rmsname=prefintin(:lenint)//'.rms'
1747 open (jplot,file=prefout(:ilen(prefout))//'.tex',&
1749 open (jrms,file=rmsname,status='unknown')
1750 open(iout,file=outname,status='unknown')
1751 ! Get parameter filenames and open the parameter files.
1752 call getenv('BONDPAR',bondname)
1753 open (ibond,file=bondname,status='old')
1754 call getenv('THETPAR',thetname)
1755 open (ithep,file=thetname,status='old')
1756 call getenv('ROTPAR',rotname)
1757 open (irotam,file=rotname,status='old')
1758 call getenv('TORPAR',torname)
1759 open (itorp,file=torname,status='old')
1760 call getenv('TORDPAR',tordname)
1761 open (itordp,file=tordname,status='old')
1762 call getenv('FOURIER',fouriername)
1763 open (ifourier,file=fouriername,status='old')
1764 call getenv('ELEPAR',elename)
1765 open (ielep,file=elename,status='old')
1766 call getenv('SIDEPAR',sidename)
1767 open (isidep,file=sidename,status='old')
1768 call getenv('SIDEP',sidepname)
1769 open (isidep1,file=sidepname,status="old")
1770 call getenv('SCCORPAR',sccorname)
1771 open (isccor,file=sccorname,status="old")
1772 call getenv('BONDPAR_NUCL',bondname_nucl)
1773 open (ibond_nucl,file=bondname_nucl,status='old')
1774 call getenv('THETPAR_NUCL',thetname_nucl)
1775 open (ithep_nucl,file=thetname_nucl,status='old')
1776 call getenv('ROTPAR_NUCL',rotname_nucl)
1777 open (irotam_nucl,file=rotname_nucl,status='old')
1778 call getenv('TORPAR_NUCL',torname_nucl)
1779 open (itorp_nucl,file=torname_nucl,status='old')
1780 call getenv('TORDPAR_NUCL',tordname_nucl)
1781 open (itordp_nucl,file=tordname_nucl,status='old')
1782 call getenv('SIDEPAR_NUCL',sidename_nucl)
1783 open (isidep_nucl,file=sidename_nucl,status='old')
1784 call getenv('SIDEPAR_SCBASE',sidename_scbase)
1785 open (isidep_scbase,file=sidename_scbase,status='old')
1786 call getenv('PEPPAR_PEPBASE',pepname_pepbase)
1787 open (isidep_pepbase,file=pepname_pepbase,status='old')
1788 call getenv('SCPAR_PHOSPH',pepname_scpho)
1789 open (isidep_scpho,file=pepname_scpho,status='old')
1790 call getenv('PEPPAR_PHOSPH',pepname_peppho)
1791 open (isidep_peppho,file=pepname_peppho,status='old')
1794 call getenv('LIPTRANPAR',liptranname)
1795 open (iliptranpar,file=liptranname,status='old')
1796 call getenv('TUBEPAR',tubename)
1797 open (itube,file=tubename,status='old')
1798 call getenv('IONPAR',ionname)
1799 open (iion,file=ionname,status='old')
1803 ! 8/9/01 In the newest version SCp interaction constants are read from a file
1804 ! Use -DOLDSCP to use hard-coded constants instead.
1806 call getenv('SCPPAR',scpname)
1807 open (iscpp,file=scpname,status='old')
1808 call getenv('SCPPAR_NUCL',scpname_nucl)
1809 open (iscpp_nucl,file=scpname_nucl,status='old')
1813 end subroutine openunits
1814 !-----------------------------------------------------------------------------
1816 !-----------------------------------------------------------------------------
1817 subroutine pdboutC(etot,rmsd,tytul)
1819 use energy_data, only: ihpb,jhpb,itype,molnum
1820 ! implicit real*8 (a-h,o-z)
1821 ! include 'DIMENSIONS'
1822 ! include 'COMMON.CONTROL'
1823 ! include 'COMMON.CHAIN'
1824 ! include 'COMMON.INTERACT'
1825 ! include 'COMMON.NAMES'
1826 ! include 'COMMON.IOUNITS'
1827 ! include 'COMMON.HEADER'
1828 ! include 'COMMON.SBRIDGE'
1829 ! include 'COMMON.TEMPFAC'
1830 character(len=50) :: tytul
1831 character(len=1) :: chainid(10)=(/'A','B','C','D','E','F',&
1833 integer :: ica(nres)
1834 real(kind=8) :: etot,rmsd
1835 integer :: iatom,ichain,ires,i,j,iti,mnum,mnum1,boxxxshift(3)
1836 real(kind=8) :: boxxxx(3)
1837 write (ipdb,'(3a,1pe15.5,a,0pf7.2)') 'REMARK ',tytul(:20),&
1838 ' ENERGY ',etot,' RMS ',rmsd
1842 boxxxshift(1)=int(c(1,nnt)/boxxsize)
1843 ! if (c(1,nnt).lt.0) boxxxshift(1)=boxxxshift(1)-1
1844 boxxxshift(2)=int(c(2,nnt)/boxzsize)
1845 ! if (c(2,nnt).lt.0) boxxxshift(2)=boxxxshift(2)-1
1846 boxxxshift(3)=int(c(3,nnt)/boxzsize)
1847 ! if (c(3,nnt).lt.0) boxxxshift(3)=boxxxshift(3)-1
1856 if (iti.eq.ntyp1_molec(mnum)) then
1859 write (ipdb,'(a)') 'TER'
1866 if ((c(j,i).gt.0).and.(c(j,nnt).lt.0)) then
1867 c(j,i)=dmod(c(j,i),boxxxx(j))+(boxxxshift(j)-1)*boxxxx(j)
1868 else if ((c(j,i).lt.0).and.(c(j,nnt).gt.0)) then
1869 c(j,i)=dmod(c(j,i),boxxxx(j))+(boxxxshift(j)+1)*boxxxx(j)
1871 c(j,i)=dmod(c(j,i),boxxxx(j))+(boxxxshift(j))*boxxxx(j)
1875 write (ipdb,10) iatom,restyp(iti,mnum),chainid(ichain),&
1876 ires,(c(j,i),j=1,3),1.0d0,tempfac(1,i)
1877 if ((iti.ne.10).and.(mnum.ne.5)) then
1879 write (ipdb,20) iatom,restyp(iti,mnum),chainid(ichain),&
1880 ires,(c(j,nres+i),j=1,3),1.0d0,tempfac(2,i)
1884 write (ipdb,'(a)') 'TER'
1888 if ((itype(i,mnum).eq.ntyp1).or.(mnum.eq.5)) cycle
1889 if (itype(i,mnum).eq.10 .and. itype(i+1,mnum1).ne.ntyp1_molec(mnum1)) then
1890 write (ipdb,30) ica(i),ica(i+1)
1891 else if (itype(i,mnum).ne.10 .and. itype(i+1,mnum1).ne.ntyp1_molec(mnum1)) then
1892 write (ipdb,30) ica(i),ica(i+1),ica(i)+1
1893 else if (itype(i,mnum).ne.10 .and. itype(i+1,mnum1).eq.ntyp1_molec(mnum1)) then
1894 write (ipdb,30) ica(i),ica(i)+1
1897 if ((itype(nct,molnum(nct)).ne.10).and.(molnum(nct).ne.5)) then
1898 write (ipdb,30) ica(nct),ica(nct)+1
1901 write (ipdb,30) ica(ihpb(i)-nres)+1,ica(jhpb(i)-nres)+1
1903 write (ipdb,'(a6)') 'ENDMDL'
1904 10 FORMAT ('ATOM',I7,' CA ',A3,1X,A1,I4,4X,3F8.3,2f6.2)
1905 20 FORMAT ('ATOM',I7,' CB ',A3,1X,A1,I4,4X,3F8.3,2f6.2)
1906 30 FORMAT ('CONECT',8I5)
1908 end subroutine pdboutC
1909 !-----------------------------------------------------------------------------
1910 subroutine cartout(igr,i,etot,free,rmsd,plik)
1911 ! implicit real*8 (a-h,o-z)
1912 ! include 'DIMENSIONS'
1913 ! include 'sizesclu.dat'
1914 ! include 'COMMON.IOUNITS'
1915 ! include 'COMMON.CHAIN'
1916 ! include 'COMMON.VAR'
1917 ! include 'COMMON.LOCAL'
1918 ! include 'COMMON.INTERACT'
1919 ! include 'COMMON.NAMES'
1920 ! include 'COMMON.GEO'
1921 ! include 'COMMON.CLUSTER'
1922 integer :: igr,i,j,k
1923 real(kind=8) :: etot,free,rmsd
1924 character(len=80) :: plik
1925 open (igeom,file=plik,position='append')
1926 write (igeom,'(2e15.5,f10.5,$)') etot,free,rmsd
1927 write (igeom,'(i4,$)') &
1928 nss_all(i),(ihpb_all(j,i),jhpb_all(j,i),j=1,nss_all(i))
1929 write (igeom,'(i10)') iscore(i)
1930 write (igeom,'(8f10.5)') &
1931 ((allcart(k,j,i),k=1,3),j=1,nres),&
1932 ((allcart(k,j+nres,i),k=1,3),j=nnt,nct)
1934 end subroutine cartout
1935 !------------------------------------------------------------------------------
1936 ! subroutine alloc_clust_arrays(n_conf)
1941 ! allocate(diss(maxdist)) !(maxdist)
1942 !el allocate(energy(0:maxconf),totfree(0:maxconf)) !(0:maxconf)
1943 ! allocatable :: enetb !(max_ene,maxstr_proc)
1944 !el allocate(entfac(maxconf)) !(maxconf)
1945 ! allocatable :: totfree_gr !(maxgr)
1946 !el allocate(rcutoff(max_cut+1)) !(max_cut+1)
1948 ! allocatable :: licz,iass !(maxgr)
1949 ! allocatable :: nconf !(maxgr,maxingr)
1950 ! allocatable :: iass_tot !(maxgr,max_cut)
1951 ! allocatable :: list_conf !(maxconf)
1953 !el allocatable :: allcart !(3,maxres2,maxstr_proc)
1954 !el allocate(rmstb(maxconf)) !(maxconf)
1955 !el allocate(mult(nres)) !(maxres)
1956 !el allocatable :: nss_all !(maxstr_proc)
1957 !el allocatable :: ihpb_all,jhpb_all !(maxss,maxstr_proc)
1958 ! allocate(icc(n_conf),iscore(n_conf)) !(maxconf)
1961 ! allocatable :: tempfac !(2,maxres)
1964 !el allocate(beta_h(maxT)) !(maxT)
1966 ! end subroutine alloc_clust_arrays
1967 !-----------------------------------------------------------------------------
1968 !-----------------------------------------------------------------------------