1 subroutine read_general_data(*)
4 include "DIMENSIONS.ZSCOPT"
5 include "DIMENSIONS.FREE"
6 include "COMMON.TORSION"
7 include "COMMON.INTERACT"
8 include "COMMON.IOUNITS"
11 include "COMMON.PROTFILES"
12 include "COMMON.CHAIN"
13 include "COMMON.NAMES"
14 include "COMMON.FFIELD"
15 include "COMMON.ENEPS"
16 include "COMMON.WEIGHTS"
18 include "COMMON.CONTROL"
19 include "COMMON.ENERGIES"
20 include "COMMON.SBRIDGE"
21 character*800 controlcard
22 integer i,j,k,ii,n_ene_found
23 integer ind,itype1,itype2,itypf,itypsc,itypp
30 call card_concat(controlcard,.true.)
31 call readi(controlcard,"N_ENE",n_ene,max_ene)
32 if (n_ene.gt.max_ene) then
33 write (iout,*) "Error: parameter out of range: N_ENE",n_ene,
37 call readi(controlcard,"NPARMSET",nparmset,1)
38 separate_parset = index(controlcard,"SEPARATE_PARSET").gt.0
39 call readi(controlcard,"IPARMPRINT",iparmprint,1)
40 write (iout,*) "PARMPRINT",iparmprint
41 if (nparmset.gt.max_parm) then
42 write (iout,*) "Error: parameter out of range: NPARMSET",
46 call readi(controlcard,"MAXIT",maxit,5000)
47 call reada(controlcard,"FIMIN",fimin,1.0d-3)
48 call readi(controlcard,"ENSEMBLES",ensembles,0)
49 hamil_rep=index(controlcard,"HAMIL_REP").gt.0
50 write (iout,*) "Number of energy parameter sets",nparmset
51 call multreadi(controlcard,"ISAMPL",isampl,nparmset,1)
52 write (iout,*) "MaxSlice",MaxSlice
53 call readi(controlcard,"NSLICE",nslice,1)
55 if (nslice.gt.MaxSlice) then
56 write (iout,*) "Error: parameter out of range: NSLICE",nslice,
60 write (iout,*) "Frequency of storing conformations",
61 & (isampl(i),i=1,nparmset)
62 write (iout,*) "Maxit",maxit," Fimin",fimin
63 call readi(controlcard,"NQ",nQ,1)
65 write (iout,*) "Error: parameter out of range: NQ",nq,
70 if (index(controlcard,"CLASSIFY").gt.0) indpdb=1
71 call reada(controlcard,"DELTA",delta,1.0d-2)
72 call readi(controlcard,"EINICHECK",einicheck,2)
73 call reada(controlcard,"DELTRMS",deltrms,5.0d-2)
74 call reada(controlcard,"DELTRGY",deltrgy,5.0d-2)
75 call readi(controlcard,"RESCALE",rescale_mode,1)
76 check_conf=index(controlcard,"NO_CHECK_CONF").eq.0
77 call reada(controlcard,'DISTCHAINMAX',distchainmax,50.0d0)
78 call readi(controlcard,'SYM',symetr,1)
79 write (iout,*) "DISTCHAINMAX",distchainmax
80 write (iout,*) "delta",delta
81 write (iout,*) "einicheck",einicheck
82 write (iout,*) "rescale_mode",rescale_mode
84 bxfile=index(controlcard,"BXFILE").gt.0
85 cxfile=index(controlcard,"CXFILE").gt.0
86 if (nslice .eq. 1 .and. .not.bxfile .and. .not.cxfile)
88 histfile=index(controlcard,"HISTFILE").gt.0
89 histout=index(controlcard,"HISTOUT").gt.0
90 entfile=index(controlcard,"ENTFILE").gt.0
91 zscfile=index(controlcard,"ZSCFILE").gt.0
92 with_dihed_constr = index(controlcard,"WITH_DIHED_CONSTR").gt.0
93 write (iout,*) "with_dihed_constr ",with_dihed_constr
94 with_theta_constr = index(controlcard,"WITH_THETA_CONSTR").gt.0
95 write (iout,*) "with_theta_constr ",with_theta_constr
96 call readi(controlcard,'CONSTR_DIST',constr_dist,0)
97 dyn_ss=index(controlcard,"DYN_SS").gt.0
100 c------------------------------------------------------------------------------
101 subroutine read_efree(*)
103 C Read molecular data
107 include 'DIMENSIONS.ZSCOPT'
108 include 'DIMENSIONS.COMPAR'
109 include 'DIMENSIONS.FREE'
110 include 'COMMON.IOUNITS'
111 include 'COMMON.TIME1'
112 include 'COMMON.SBRIDGE'
113 include 'COMMON.CONTROL'
114 include 'COMMON.CHAIN'
115 include 'COMMON.HEADER'
117 include 'COMMON.FREE'
118 character*320 controlcard,ucase
119 integer iparm,ib,i,j,npars
131 call card_concat(controlcard,.true.)
132 call readi(controlcard,'NT',nT_h(iparm),1)
133 write (iout,*) "iparm",iparm," nt",nT_h(iparm)
135 if (nT_h(iparm).gt.MaxT_h) then
136 write (iout,*) "Error: parameter out of range: NT",nT_h(iparm),
140 replica(iparm)=index(controlcard,"REPLICA").gt.0
141 umbrella(iparm)=index(controlcard,"UMBRELLA").gt.0
142 read_iset(iparm)=index(controlcard,"READ_ISET").gt.0
143 write (iout,*) "nQ",nQ," nT",nT_h(iparm)," replica ",
144 & replica(iparm)," umbrella ",umbrella(iparm),
145 & " read_iset",read_iset(iparm)
148 call card_concat(controlcard,.true.)
149 call readi(controlcard,'NR',nR(ib,iparm),1)
150 if (umbrella(iparm)) then
153 nRR(ib,iparm)=nR(ib,iparm)
155 if (nR(ib,iparm).gt.MaxR) then
156 write (iout,*) "Error: parameter out of range: NR",
160 call reada(controlcard,'TEMP',beta_h(ib,iparm),298.0d0)
161 beta_h(ib,iparm)=1.0d0/(beta_h(ib,iparm)*1.987D-3)
162 call multreada(controlcard,'FI',f(1,ib,iparm),nR(ib,iparm),
165 call card_concat(controlcard,.true.)
166 call multreada(controlcard,'KH',KH(1,i,ib,iparm),nQ,
168 call multreada(controlcard,'Q0',Q0(1,i,ib,iparm),nQ,
173 write (iout,*) "ib",ib," beta_h",
174 & 1.0d0/(0.001987*beta_h(ib,iparm))
175 write (iout,*) "nR",nR(ib,iparm)
176 write (iout,*) "fi",(f(i,ib,iparm),i=1,nR(ib,iparm))
178 write (iout,*) "i",i," Kh",(Kh(j,i,ib,iparm),j=1,nQ),
179 & "q0",(q0(j,i,ib,iparm),j=1,nQ)
191 nR(ib,iparm)=nR(ib,1)
192 if (umbrella(iparm)) then
195 nRR(ib,iparm)=nR(ib,1)
197 beta_h(ib,iparm)=beta_h(ib,1)
199 f(i,ib,iparm)=f(i,ib,1)
201 KH(j,i,ib,iparm)=KH(j,i,ib,1)
202 Q0(j,i,ib,iparm)=Q0(j,i,ib,1)
205 replica(iparm)=replica(1)
206 umbrella(iparm)=umbrella(1)
207 read_iset(iparm)=read_iset(1)
215 c-----------------------------------------------------------------------------
216 subroutine read_protein_data(*)
219 include "DIMENSIONS.ZSCOPT"
220 include "DIMENSIONS.FREE"
223 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
226 include "COMMON.CHAIN"
227 include "COMMON.IOUNITS"
228 include "COMMON.PROT"
229 include "COMMON.PROTFILES"
230 include "COMMON.NAMES"
231 include "COMMON.FREE"
232 include "COMMON.OBCINKA"
234 character*16000 controlcard
235 integer i,ii,ib,iR,iparm,ilen,iroof,nthr,npars
245 C Read names of files with conformation data.
246 if (replica(iparm)) then
252 do ii=1,nRR(ib,iparm)
253 write (iout,*) "Parameter set",iparm," temperature",ib,
256 call card_concat(controlcard,.true.)
257 write (iout,*) controlcard(:ilen(controlcard))
258 call readi(controlcard,"NFILE_BIN",nfile_bin(ii,ib,iparm),0)
259 call readi(controlcard,"NFILE_ASC",nfile_asc(ii,ib,iparm),0)
260 call readi(controlcard,"NFILE_CX",nfile_cx(ii,ib,iparm),0)
261 call readi(controlcard,"REC_START",rec_start(ii,ib,iparm),1)
262 call readi(controlcard,"REC_END",rec_end(ii,ib,iparm),
263 & maxstr*isampl(iparm)+rec_start(ii,ib,iparm)-1)
264 call reada(controlcard,"TIME_START",
265 & time_start_collect(ii,ib,iparm),0.0d0)
266 call reada(controlcard,"TIME_END",time_end_collect(ii,ib,iparm),
268 write (iout,*) "rec_start",rec_start(ii,ib,iparm),
269 & " rec_end",rec_end(ii,ib,iparm)
270 write (iout,*) "time_start",time_start_collect(ii,ib,iparm),
271 & " time_end",time_end_collect(ii,ib,iparm)
273 if (replica(iparm)) then
274 call readi(controlcard,"TOTRAJ",totraj(ii,iparm),1)
275 write (iout,*) "Number of trajectories",totraj(ii,iparm)
278 if (nfile_bin(ii,ib,iparm).lt.2
279 & .and. nfile_asc(ii,ib,iparm).eq.0
280 & .and. nfile_cx(ii,ib,iparm).eq.0) then
281 write (iout,*) "Error - no action specified!"
284 if (nfile_bin(ii,ib,iparm).gt.0) then
285 call card_concat(controlcard,.false.)
286 call split_string(controlcard,protfiles(1,1,ii,ib,iparm),
287 & maxfile_prot,nfile_bin(ii,ib,iparm))
289 write(iout,*)"nfile_bin",nfile_bin(ii,ib,iparm)
290 write(iout,*) (protfiles(i,1,ii,ib,iparm),
291 & i=1,nfile_bin(ii,ib,iparm))
294 if (nfile_asc(ii,ib,iparm).gt.0) then
295 call card_concat(controlcard,.false.)
296 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
297 & maxfile_prot,nfile_asc(ii,ib,iparm))
299 write(iout,*) "nfile_asc(ii,ib,iparm)",nfile_asc(ii,ib,iparm)
300 write(iout,*) (protfiles(i,2,ii,ib,iparm),
301 & i=1,nfile_asc(ii,ib,iparm))
303 else if (nfile_cx(ii,ib,iparm).gt.0) then
304 call card_concat(controlcard,.false.)
305 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
306 & maxfile_prot,nfile_cx(ii,ib,iparm))
308 write(iout,*) "nfile_cx(ii,ib,iparm)",nfile_cx(ii,ib,iparm)
309 write(iout,*) (protfiles(i,2,ii,ib,iparm),
310 & i=1,nfile_cx(ii,ib,iparm))
321 c-------------------------------------------------------------------------------
322 subroutine opentmp(islice,iunit,bprotfile_temp)
325 include "DIMENSIONS.ZSCOPT"
326 include "DIMENSIONS.FREE"
329 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
332 include "COMMON.IOUNITS"
333 include "COMMON.PROTFILES"
334 include "COMMON.PROT"
335 include "COMMON.FREE"
336 character*64 bprotfile_temp
337 character*3 liczba,liczba2
344 write (liczba1,'(bz,i2.2)') islice
345 write (liczba,'(bz,i3.3)') me
347 c write (iout,*) "separate_parset ",separate_parset,
349 if (separate_parset) then
350 write (liczba2,'(bz,i3.3)') myparm
351 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
352 & prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
353 open (iunit,file=bprotfile_temp,status="unknown",
354 & form="unformatted",access="direct",recl=lenrec)
356 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
357 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
358 open (iunit,file=bprotfile_temp,status="unknown",
359 & form="unformatted",access="direct",recl=lenrec)
362 bprotfile_temp = scratchdir(:ilen(scratchdir))//
363 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
364 open (iunit,file=bprotfile_temp,status="unknown",
365 & form="unformatted",access="direct",recl=lenrec)
367 c write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
372 c-------------------------------------------------------------------------------
373 subroutine read_database(*)
376 include "DIMENSIONS.ZSCOPT"
377 include "DIMENSIONS.FREE"
380 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
383 include "COMMON.CHAIN"
384 include "COMMON.IOUNITS"
385 include "COMMON.PROTFILES"
386 include "COMMON.NAMES"
389 include "COMMON.ENEPS"
390 include "COMMON.PROT"
391 include "COMMON.INTERACT"
392 include "COMMON.FREE"
393 include "COMMON.SBRIDGE"
394 include "COMMON.OBCINKA"
395 real*4 csingle(3,maxres2)
396 character*64 nazwa,bprotfile_temp
399 integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
400 & ll(maxslice),mm(maxslice),if
401 integer nrec,nlines,iscor,iunit,islice
402 double precision energ
405 double precision rmsdev,energia(0:max_ene),efree,eini,temp
406 double precision prop(maxQ)
407 integer ntot_all(maxslice,0:maxprocs-1), maxslice_buff
408 integer iparm,ib,iib,ir,nprop,nthr,npars
409 double precision etot,time
413 lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
414 lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
416 write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
426 write (iout,*) "nparmset",nparmset
434 if (replica(iparm)) then
441 do iR=1,nRR(ib,iparm)
443 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
449 IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
450 c Read conformations from binary DA files (one per batch) and write them to
451 c a binary DA scratchfile.
452 write (liczba,'(bz,i3.3)') me
453 do if=1,nfile_bin(iR,ib,iparm)
454 nazwa=protfiles(if,1,iR,ib,iparm)
455 & (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
456 open (ientin,file=nazwa,status="old",form="unformatted",
457 & access="direct",recl=lenrec2,err=1111)
460 call opentmp(islice,ientout,bprotfile_temp)
461 call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
462 & mm(islice),iR,ib,iparm)
469 IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
470 c Read conformations from multiple ASCII int files and write them to a binary
472 do if=1,nfile_asc(iR,ib,iparm)
473 nazwa=protfiles(if,2,iR,ib,iparm)
474 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
475 open(unit=ientin,file=nazwa,status='old',err=1111)
476 write(iout,*) "reading ",nazwa(:ilen(nazwa))
478 call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
481 IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
482 c Read conformations from cx files and write them to a binary
484 do if=1,nfile_cx(iR,ib,iparm)
485 nazwa=protfiles(if,2,iR,ib,iparm)
486 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
487 write(iout,*) "reading ",nazwa(:ilen(nazwa))
489 print *,"Calling cxread"
490 call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
493 write (iout,*) "exit cxread"
499 stot(islice)=stot(islice)+jj(islice)
504 write (iout,*) "IPARM",iparm
507 if (nslice.eq.1) then
509 write (liczba,'(bz,i3.3)') me
510 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
511 & prefix(:ilen(prefix))//liczba//".xbin.tmp"
513 bprotfile_temp = scratchdir(:ilen(scratchdir))//
514 & "/"//prefix(:ilen(prefix))//".xbin.tmp"
516 write(iout,*) mm(1)," conformations read",ll(1),
517 & " conformations written to ",
518 & bprotfile_temp(:ilen(bprotfile_temp))
521 write (liczba1,'(bz,i2.2)') islice
523 write (liczba,'(bz,i3.3)') me
524 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
525 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
527 bprotfile_temp = scratchdir(:ilen(scratchdir))//
528 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
530 write(iout,*) mm(islice)," conformations read",ll(islice),
531 & " conformations written to ",
532 & bprotfile_temp(:ilen(bprotfile_temp))
537 c Check if everyone has the same number of conformations
539 c call MPI_ALLgather(MPI_IN_PLACE,stot(1),MPI_DATATYPE_NULL,
540 c & ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
542 maxslice_buff=maxslice
544 call MPI_Allgather(stot(1),maxslice_buff,MPI_INTEGER,
545 & ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
550 if (stot(islice).ne.ntot_all(islice,i)) then
551 write (iout,*) "Number of conformations at processor",i,
552 & " differs from that at processor",me,
553 & stot(islice),ntot_all(islice,i)," slice",islice
561 write (iout,*) "Numbers of conformations read by processors"
564 write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
566 write (iout,*) "Calculation terminated."
571 ntot(islice)=stot(islice)
575 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
579 c------------------------------------------------------------------------------
580 subroutine card_concat(card,to_upper)
582 include 'DIMENSIONS.ZSCOPT'
583 include "COMMON.IOUNITS"
585 character*80 karta,ucase
589 read (inp,'(a)') karta
590 if (to_upper) karta=ucase(karta)
592 do while (karta(80:80).eq.'&')
593 card=card(:ilen(card)+1)//karta(:79)
594 read (inp,'(a)') karta
595 if (to_upper) karta=ucase(karta)
597 card=card(:ilen(card)+1)//karta
600 c------------------------------------------------------------------------------
601 subroutine readi(rekord,lancuch,wartosc,default)
603 character*(*) rekord,lancuch
604 integer wartosc,default
607 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
612 iread=iread+ilen(lancuch)+1
613 read (rekord(iread:),*) wartosc
616 c----------------------------------------------------------------------------
617 subroutine reada(rekord,lancuch,wartosc,default)
619 character*(*) rekord,lancuch
621 double precision wartosc,default
624 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
629 iread=iread+ilen(lancuch)+1
630 read (rekord(iread:),*) wartosc
633 c----------------------------------------------------------------------------
634 subroutine multreadi(rekord,lancuch,tablica,dim,default)
637 integer tablica(dim),default
638 character*(*) rekord,lancuch
645 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
646 if (iread.eq.0) return
647 iread=iread+ilen(lancuch)+1
648 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
651 c----------------------------------------------------------------------------
652 subroutine multreada(rekord,lancuch,tablica,dim,default)
655 double precision tablica(dim),default
656 character*(*) rekord,lancuch
663 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
664 if (iread.eq.0) return
665 iread=iread+ilen(lancuch)+1
666 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
669 c----------------------------------------------------------------------------
670 subroutine reads(rekord,lancuch,wartosc,default)
672 character*(*) rekord,lancuch,wartosc,default
674 integer ilen,lenlan,lenrec,iread,ireade
680 iread=index(rekord,lancuch(:lenlan)//"=")
681 c print *,"rekord",rekord," lancuch",lancuch
682 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
688 c print *,"iread",iread
689 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
690 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
692 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
694 c print *,"iread",iread
695 if (iread.gt.lenrec) then
700 c print *,"ireade",ireade
701 do while (ireade.lt.lenrec .and.
702 & .not.iblnk(rekord(ireade:ireade)))
705 wartosc=rekord(iread:ireade)
708 c----------------------------------------------------------------------------
709 subroutine multreads(rekord,lancuch,tablica,dim,default)
712 character*(*) rekord,lancuch,tablica(dim),default
714 integer ilen,lenlan,lenrec,iread,ireade
723 iread=index(rekord,lancuch(:lenlan)//"=")
724 c print *,"rekord",rekord," lancuch",lancuch
725 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
726 if (iread.eq.0) return
729 c print *,"iread",iread
730 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
731 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
733 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
735 c print *,"iread",iread
736 if (iread.gt.lenrec) return
738 c print *,"ireade",ireade
739 do while (ireade.lt.lenrec .and.
740 & .not.iblnk(rekord(ireade:ireade)))
743 tablica(i)=rekord(iread:ireade)
747 c----------------------------------------------------------------------------
748 subroutine split_string(rekord,tablica,dim,nsub)
750 integer dim,nsub,i,ii,ll,kk
751 character*(*) tablica(dim)
762 C Find the start of term name
764 do while (ii.le.ll .and. rekord(ii:ii).eq." ")
767 C Parse the name into TABLICA(i) until blank found
768 do while (ii.le.ll .and. rekord(ii:ii).ne." ")
770 tablica(i)(kk:kk)=rekord(ii:ii)
773 if (kk.gt.0) nsub=nsub+1
778 c--------------------------------------------------------------------------------
779 integer function iroof(n,m)
781 if (ii*m .lt. n) ii=ii+1