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 call readi(controlcard,'CONSTR_DIST',constr_dist,0)
95 dyn_ss=index(controlcard,"DYN_SS").gt.0
98 c------------------------------------------------------------------------------
99 subroutine read_efree(*)
101 C Read molecular data
105 include 'DIMENSIONS.ZSCOPT'
106 include 'DIMENSIONS.COMPAR'
107 include 'DIMENSIONS.FREE'
108 include 'COMMON.IOUNITS'
109 include 'COMMON.TIME1'
110 include 'COMMON.SBRIDGE'
111 include 'COMMON.CONTROL'
112 include 'COMMON.CHAIN'
113 include 'COMMON.HEADER'
115 include 'COMMON.FREE'
116 character*320 controlcard,ucase
117 integer iparm,ib,i,j,npars
129 call card_concat(controlcard,.true.)
130 call readi(controlcard,'NT',nT_h(iparm),1)
131 write (iout,*) "iparm",iparm," nt",nT_h(iparm)
133 if (nT_h(iparm).gt.MaxT_h) then
134 write (iout,*) "Error: parameter out of range: NT",nT_h(iparm),
138 replica(iparm)=index(controlcard,"REPLICA").gt.0
139 umbrella(iparm)=index(controlcard,"UMBRELLA").gt.0
140 read_iset(iparm)=index(controlcard,"READ_ISET").gt.0
141 write (iout,*) "nQ",nQ," nT",nT_h(iparm)," replica ",
142 & replica(iparm)," umbrella ",umbrella(iparm),
143 & " read_iset",read_iset(iparm)
146 call card_concat(controlcard,.true.)
147 call readi(controlcard,'NR',nR(ib,iparm),1)
148 if (umbrella(iparm)) then
151 nRR(ib,iparm)=nR(ib,iparm)
153 if (nR(ib,iparm).gt.MaxR) then
154 write (iout,*) "Error: parameter out of range: NR",
158 call reada(controlcard,'TEMP',beta_h(ib,iparm),298.0d0)
159 beta_h(ib,iparm)=1.0d0/(beta_h(ib,iparm)*1.987D-3)
160 call multreada(controlcard,'FI',f(1,ib,iparm),nR(ib,iparm),
163 call card_concat(controlcard,.true.)
164 call multreada(controlcard,'KH',KH(1,i,ib,iparm),nQ,
166 call multreada(controlcard,'Q0',Q0(1,i,ib,iparm),nQ,
171 write (iout,*) "ib",ib," beta_h",
172 & 1.0d0/(0.001987*beta_h(ib,iparm))
173 write (iout,*) "nR",nR(ib,iparm)
174 write (iout,*) "fi",(f(i,ib,iparm),i=1,nR(ib,iparm))
176 write (iout,*) "i",i," Kh",(Kh(j,i,ib,iparm),j=1,nQ),
177 & "q0",(q0(j,i,ib,iparm),j=1,nQ)
189 nR(ib,iparm)=nR(ib,1)
190 if (umbrella(iparm)) then
193 nRR(ib,iparm)=nR(ib,1)
195 beta_h(ib,iparm)=beta_h(ib,1)
197 f(i,ib,iparm)=f(i,ib,1)
199 KH(j,i,ib,iparm)=KH(j,i,ib,1)
200 Q0(j,i,ib,iparm)=Q0(j,i,ib,1)
203 replica(iparm)=replica(1)
204 umbrella(iparm)=umbrella(1)
205 read_iset(iparm)=read_iset(1)
213 c-----------------------------------------------------------------------------
214 subroutine read_protein_data(*)
217 include "DIMENSIONS.ZSCOPT"
218 include "DIMENSIONS.FREE"
221 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
224 include "COMMON.CHAIN"
225 include "COMMON.IOUNITS"
226 include "COMMON.PROT"
227 include "COMMON.PROTFILES"
228 include "COMMON.NAMES"
229 include "COMMON.FREE"
230 include "COMMON.OBCINKA"
232 character*16000 controlcard
233 integer i,ii,ib,iR,iparm,ilen,iroof,nthr,npars
243 C Read names of files with conformation data.
244 if (replica(iparm)) then
250 do ii=1,nRR(ib,iparm)
251 write (iout,*) "Parameter set",iparm," temperature",ib,
254 call card_concat(controlcard,.true.)
255 write (iout,*) controlcard(:ilen(controlcard))
256 call readi(controlcard,"NFILE_BIN",nfile_bin(ii,ib,iparm),0)
257 call readi(controlcard,"NFILE_ASC",nfile_asc(ii,ib,iparm),0)
258 call readi(controlcard,"NFILE_CX",nfile_cx(ii,ib,iparm),0)
259 call readi(controlcard,"REC_START",rec_start(ii,ib,iparm),1)
260 call readi(controlcard,"REC_END",rec_end(ii,ib,iparm),
261 & maxstr*isampl(iparm)+rec_start(ii,ib,iparm)-1)
262 call reada(controlcard,"TIME_START",
263 & time_start_collect(ii,ib,iparm),0.0d0)
264 call reada(controlcard,"TIME_END",time_end_collect(ii,ib,iparm),
266 write (iout,*) "rec_start",rec_start(ii,ib,iparm),
267 & " rec_end",rec_end(ii,ib,iparm)
268 write (iout,*) "time_start",time_start_collect(ii,ib,iparm),
269 & " time_end",time_end_collect(ii,ib,iparm)
271 if (replica(iparm)) then
272 call readi(controlcard,"TOTRAJ",totraj(ii,iparm),1)
273 write (iout,*) "Number of trajectories",totraj(ii,iparm)
276 if (nfile_bin(ii,ib,iparm).lt.2
277 & .and. nfile_asc(ii,ib,iparm).eq.0
278 & .and. nfile_cx(ii,ib,iparm).eq.0) then
279 write (iout,*) "Error - no action specified!"
282 if (nfile_bin(ii,ib,iparm).gt.0) then
283 call card_concat(controlcard,.false.)
284 call split_string(controlcard,protfiles(1,1,ii,ib,iparm),
285 & maxfile_prot,nfile_bin(ii,ib,iparm))
287 write(iout,*)"nfile_bin",nfile_bin(ii,ib,iparm)
288 write(iout,*) (protfiles(i,1,ii,ib,iparm),
289 & i=1,nfile_bin(ii,ib,iparm))
292 if (nfile_asc(ii,ib,iparm).gt.0) then
293 call card_concat(controlcard,.false.)
294 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
295 & maxfile_prot,nfile_asc(ii,ib,iparm))
297 write(iout,*) "nfile_asc(ii,ib,iparm)",nfile_asc(ii,ib,iparm)
298 write(iout,*) (protfiles(i,2,ii,ib,iparm),
299 & i=1,nfile_asc(ii,ib,iparm))
301 else if (nfile_cx(ii,ib,iparm).gt.0) then
302 call card_concat(controlcard,.false.)
303 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
304 & maxfile_prot,nfile_cx(ii,ib,iparm))
306 write(iout,*) "nfile_cx(ii,ib,iparm)",nfile_cx(ii,ib,iparm)
307 write(iout,*) (protfiles(i,2,ii,ib,iparm),
308 & i=1,nfile_cx(ii,ib,iparm))
319 c-------------------------------------------------------------------------------
320 subroutine opentmp(islice,iunit,bprotfile_temp)
323 include "DIMENSIONS.ZSCOPT"
324 include "DIMENSIONS.FREE"
327 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
330 include "COMMON.IOUNITS"
331 include "COMMON.PROTFILES"
332 include "COMMON.PROT"
333 include "COMMON.FREE"
334 character*64 bprotfile_temp
335 character*3 liczba,liczba2
342 write (liczba1,'(bz,i2.2)') islice
343 write (liczba,'(bz,i3.3)') me
345 c write (iout,*) "separate_parset ",separate_parset,
347 if (separate_parset) then
348 write (liczba2,'(bz,i3.3)') myparm
349 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
350 & prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
351 open (iunit,file=bprotfile_temp,status="unknown",
352 & form="unformatted",access="direct",recl=lenrec)
354 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
355 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
356 open (iunit,file=bprotfile_temp,status="unknown",
357 & form="unformatted",access="direct",recl=lenrec)
360 bprotfile_temp = scratchdir(:ilen(scratchdir))//
361 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
362 open (iunit,file=bprotfile_temp,status="unknown",
363 & form="unformatted",access="direct",recl=lenrec)
365 c write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
370 c-------------------------------------------------------------------------------
371 subroutine read_database(*)
374 include "DIMENSIONS.ZSCOPT"
375 include "DIMENSIONS.FREE"
378 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
381 include "COMMON.CHAIN"
382 include "COMMON.IOUNITS"
383 include "COMMON.PROTFILES"
384 include "COMMON.NAMES"
387 include "COMMON.ENEPS"
388 include "COMMON.PROT"
389 include "COMMON.INTERACT"
390 include "COMMON.FREE"
391 include "COMMON.SBRIDGE"
392 include "COMMON.OBCINKA"
393 real*4 csingle(3,maxres2)
394 character*64 nazwa,bprotfile_temp
397 integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
398 & ll(maxslice),mm(maxslice),if
399 integer nrec,nlines,iscor,iunit,islice
400 double precision energ
403 double precision rmsdev,energia(0:max_ene),efree,eini,temp
404 double precision prop(maxQ)
405 integer ntot_all(maxslice,0:maxprocs-1), maxslice_buff
406 integer iparm,ib,iib,ir,nprop,nthr,npars
407 double precision etot,time
411 lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
412 lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
414 write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
424 write (iout,*) "nparmset",nparmset
432 if (replica(iparm)) then
439 do iR=1,nRR(ib,iparm)
441 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
447 IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
448 c Read conformations from binary DA files (one per batch) and write them to
449 c a binary DA scratchfile.
450 write (liczba,'(bz,i3.3)') me
451 do if=1,nfile_bin(iR,ib,iparm)
452 nazwa=protfiles(if,1,iR,ib,iparm)
453 & (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
454 open (ientin,file=nazwa,status="old",form="unformatted",
455 & access="direct",recl=lenrec2,err=1111)
458 call opentmp(islice,ientout,bprotfile_temp)
459 call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
460 & mm(islice),iR,ib,iparm)
467 IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
468 c Read conformations from multiple ASCII int files and write them to a binary
470 do if=1,nfile_asc(iR,ib,iparm)
471 nazwa=protfiles(if,2,iR,ib,iparm)
472 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
473 open(unit=ientin,file=nazwa,status='old',err=1111)
474 write(iout,*) "reading ",nazwa(:ilen(nazwa))
476 call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
479 IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
480 c Read conformations from cx files and write them to a binary
482 do if=1,nfile_cx(iR,ib,iparm)
483 nazwa=protfiles(if,2,iR,ib,iparm)
484 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
485 write(iout,*) "reading ",nazwa(:ilen(nazwa))
487 print *,"Calling cxread"
488 call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
491 write (iout,*) "exit cxread"
497 stot(islice)=stot(islice)+jj(islice)
502 write (iout,*) "IPARM",iparm
505 if (nslice.eq.1) then
507 write (liczba,'(bz,i3.3)') me
508 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
509 & prefix(:ilen(prefix))//liczba//".xbin.tmp"
511 bprotfile_temp = scratchdir(:ilen(scratchdir))//
512 & "/"//prefix(:ilen(prefix))//".xbin.tmp"
514 write(iout,*) mm(1)," conformations read",ll(1),
515 & " conformations written to ",
516 & bprotfile_temp(:ilen(bprotfile_temp))
519 write (liczba1,'(bz,i2.2)') islice
521 write (liczba,'(bz,i3.3)') me
522 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
523 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
525 bprotfile_temp = scratchdir(:ilen(scratchdir))//
526 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
528 write(iout,*) mm(islice)," conformations read",ll(islice),
529 & " conformations written to ",
530 & bprotfile_temp(:ilen(bprotfile_temp))
535 c Check if everyone has the same number of conformations
537 c call MPI_ALLgather(MPI_IN_PLACE,stot(1),MPI_DATATYPE_NULL,
538 c & ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
540 maxslice_buff=maxslice
542 call MPI_Allgather(stot(1),maxslice_buff,MPI_INTEGER,
543 & ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
548 if (stot(islice).ne.ntot_all(islice,i)) then
549 write (iout,*) "Number of conformations at processor",i,
550 & " differs from that at processor",me,
551 & stot(islice),ntot_all(islice,i)," slice",islice
559 write (iout,*) "Numbers of conformations read by processors"
562 write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
564 write (iout,*) "Calculation terminated."
569 ntot(islice)=stot(islice)
573 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
577 c------------------------------------------------------------------------------
578 subroutine card_concat(card,to_upper)
580 include 'DIMENSIONS.ZSCOPT'
581 include "COMMON.IOUNITS"
583 character*80 karta,ucase
587 read (inp,'(a)') karta
588 if (to_upper) karta=ucase(karta)
590 do while (karta(80:80).eq.'&')
591 card=card(:ilen(card)+1)//karta(:79)
592 read (inp,'(a)') karta
593 if (to_upper) karta=ucase(karta)
595 card=card(:ilen(card)+1)//karta
598 c------------------------------------------------------------------------------
599 subroutine readi(rekord,lancuch,wartosc,default)
601 character*(*) rekord,lancuch
602 integer wartosc,default
605 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
610 iread=iread+ilen(lancuch)+1
611 read (rekord(iread:),*) wartosc
614 c----------------------------------------------------------------------------
615 subroutine reada(rekord,lancuch,wartosc,default)
617 character*(*) rekord,lancuch
619 double precision wartosc,default
622 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
627 iread=iread+ilen(lancuch)+1
628 read (rekord(iread:),*) wartosc
631 c----------------------------------------------------------------------------
632 subroutine multreadi(rekord,lancuch,tablica,dim,default)
635 integer tablica(dim),default
636 character*(*) rekord,lancuch
643 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
644 if (iread.eq.0) return
645 iread=iread+ilen(lancuch)+1
646 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
649 c----------------------------------------------------------------------------
650 subroutine multreada(rekord,lancuch,tablica,dim,default)
653 double precision tablica(dim),default
654 character*(*) rekord,lancuch
661 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
662 if (iread.eq.0) return
663 iread=iread+ilen(lancuch)+1
664 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
667 c----------------------------------------------------------------------------
668 subroutine reads(rekord,lancuch,wartosc,default)
670 character*(*) rekord,lancuch,wartosc,default
672 integer ilen,lenlan,lenrec,iread,ireade
678 iread=index(rekord,lancuch(:lenlan)//"=")
679 c print *,"rekord",rekord," lancuch",lancuch
680 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
686 c print *,"iread",iread
687 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
688 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
690 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
692 c print *,"iread",iread
693 if (iread.gt.lenrec) then
698 c print *,"ireade",ireade
699 do while (ireade.lt.lenrec .and.
700 & .not.iblnk(rekord(ireade:ireade)))
703 wartosc=rekord(iread:ireade)
706 c----------------------------------------------------------------------------
707 subroutine multreads(rekord,lancuch,tablica,dim,default)
710 character*(*) rekord,lancuch,tablica(dim),default
712 integer ilen,lenlan,lenrec,iread,ireade
721 iread=index(rekord,lancuch(:lenlan)//"=")
722 c print *,"rekord",rekord," lancuch",lancuch
723 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
724 if (iread.eq.0) return
727 c print *,"iread",iread
728 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
729 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
731 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
733 c print *,"iread",iread
734 if (iread.gt.lenrec) return
736 c print *,"ireade",ireade
737 do while (ireade.lt.lenrec .and.
738 & .not.iblnk(rekord(ireade:ireade)))
741 tablica(i)=rekord(iread:ireade)
745 c----------------------------------------------------------------------------
746 subroutine split_string(rekord,tablica,dim,nsub)
748 integer dim,nsub,i,ii,ll,kk
749 character*(*) tablica(dim)
760 C Find the start of term name
762 do while (ii.le.ll .and. rekord(ii:ii).eq." ")
765 C Parse the name into TABLICA(i) until blank found
766 do while (ii.le.ll .and. rekord(ii:ii).ne." ")
768 tablica(i)(kk:kk)=rekord(ii:ii)
771 if (kk.gt.0) nsub=nsub+1
776 c--------------------------------------------------------------------------------
777 integer function iroof(n,m)
779 if (ii*m .lt. n) ii=ii+1