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.SPLITELE"
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 reada(controlcard,'BOXX',boxxsize,100.0d0)
79 call reada(controlcard,'BOXY',boxysize,100.0d0)
80 call reada(controlcard,'BOXZ',boxzsize,100.0d0)
81 c Cutoff range for interactions
82 call reada(controlcard,"R_CUT",r_cut,15.0d0)
83 call reada(controlcard,"LAMBDA",rlamb,0.3d0)
84 call readi(controlcard,'SYM',symetr,1)
85 write (iout,*) "DISTCHAINMAX",distchainmax
86 write (iout,*) "delta",delta
87 write (iout,*) "einicheck",einicheck
88 write (iout,*) "rescale_mode",rescale_mode
90 bxfile=index(controlcard,"BXFILE").gt.0
91 cxfile=index(controlcard,"CXFILE").gt.0
92 if (nslice .eq. 1 .and. .not.bxfile .and. .not.cxfile)
94 histfile=index(controlcard,"HISTFILE").gt.0
95 histout=index(controlcard,"HISTOUT").gt.0
96 entfile=index(controlcard,"ENTFILE").gt.0
97 zscfile=index(controlcard,"ZSCFILE").gt.0
98 with_dihed_constr = index(controlcard,"WITH_DIHED_CONSTR").gt.0
99 write (iout,*) "with_dihed_constr ",with_dihed_constr
100 call readi(controlcard,'CONSTR_DIST',constr_dist,0)
103 c------------------------------------------------------------------------------
104 subroutine read_efree(*)
106 C Read molecular data
110 include 'DIMENSIONS.ZSCOPT'
111 include 'DIMENSIONS.COMPAR'
112 include 'DIMENSIONS.FREE'
113 include 'COMMON.IOUNITS'
114 include 'COMMON.TIME1'
115 include 'COMMON.SBRIDGE'
116 include 'COMMON.CONTROL'
117 include 'COMMON.CHAIN'
118 include 'COMMON.HEADER'
120 include 'COMMON.FREE'
121 character*320 controlcard,ucase
122 integer iparm,ib,i,j,npars
134 call card_concat(controlcard,.true.)
135 call readi(controlcard,'NT',nT_h(iparm),1)
136 write (iout,*) "iparm",iparm," nt",nT_h(iparm)
138 if (nT_h(iparm).gt.MaxT_h) then
139 write (iout,*) "Error: parameter out of range: NT",nT_h(iparm),
143 replica(iparm)=index(controlcard,"REPLICA").gt.0
144 umbrella(iparm)=index(controlcard,"UMBRELLA").gt.0
145 read_iset(iparm)=index(controlcard,"READ_ISET").gt.0
146 write (iout,*) "nQ",nQ," nT",nT_h(iparm)," replica ",
147 & replica(iparm)," umbrella ",umbrella(iparm),
148 & " read_iset",read_iset(iparm)
151 call card_concat(controlcard,.true.)
152 call readi(controlcard,'NR',nR(ib,iparm),1)
153 if (umbrella(iparm)) then
156 nRR(ib,iparm)=nR(ib,iparm)
158 if (nR(ib,iparm).gt.MaxR) then
159 write (iout,*) "Error: parameter out of range: NR",
163 call reada(controlcard,'TEMP',beta_h(ib,iparm),298.0d0)
164 beta_h(ib,iparm)=1.0d0/(beta_h(ib,iparm)*1.987D-3)
165 call multreada(controlcard,'FI',f(1,ib,iparm),nR(ib,iparm),
168 call card_concat(controlcard,.true.)
169 call multreada(controlcard,'KH',KH(1,i,ib,iparm),nQ,
171 call multreada(controlcard,'Q0',Q0(1,i,ib,iparm),nQ,
176 write (iout,*) "ib",ib," beta_h",
177 & 1.0d0/(0.001987*beta_h(ib,iparm))
178 write (iout,*) "nR",nR(ib,iparm)
179 write (iout,*) "fi",(f(i,ib,iparm),i=1,nR(ib,iparm))
181 write (iout,*) "i",i," Kh",(Kh(j,i,ib,iparm),j=1,nQ),
182 & "q0",(q0(j,i,ib,iparm),j=1,nQ)
194 nR(ib,iparm)=nR(ib,1)
195 if (umbrella(iparm)) then
198 nRR(ib,iparm)=nR(ib,1)
200 beta_h(ib,iparm)=beta_h(ib,1)
202 f(i,ib,iparm)=f(i,ib,1)
204 KH(j,i,ib,iparm)=KH(j,i,ib,1)
205 Q0(j,i,ib,iparm)=Q0(j,i,ib,1)
208 replica(iparm)=replica(1)
209 umbrella(iparm)=umbrella(1)
210 read_iset(iparm)=read_iset(1)
218 c-----------------------------------------------------------------------------
219 subroutine read_protein_data(*)
222 include "DIMENSIONS.ZSCOPT"
223 include "DIMENSIONS.FREE"
226 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
229 include "COMMON.CHAIN"
230 include "COMMON.IOUNITS"
231 include "COMMON.PROT"
232 include "COMMON.PROTFILES"
233 include "COMMON.NAMES"
234 include "COMMON.FREE"
235 include "COMMON.OBCINKA"
237 character*16000 controlcard
238 integer i,ii,ib,iR,iparm,ilen,iroof,nthr,npars
248 C Read names of files with conformation data.
249 if (replica(iparm)) then
255 do ii=1,nRR(ib,iparm)
256 write (iout,*) "Parameter set",iparm," temperature",ib,
259 call card_concat(controlcard,.true.)
260 write (iout,*) controlcard(:ilen(controlcard))
261 call readi(controlcard,"NFILE_BIN",nfile_bin(ii,ib,iparm),0)
262 call readi(controlcard,"NFILE_ASC",nfile_asc(ii,ib,iparm),0)
263 call readi(controlcard,"NFILE_CX",nfile_cx(ii,ib,iparm),0)
264 call readi(controlcard,"REC_START",rec_start(ii,ib,iparm),1)
265 call readi(controlcard,"REC_END",rec_end(ii,ib,iparm),
266 & maxstr*isampl(iparm)+rec_start(ii,ib,iparm)-1)
267 call reada(controlcard,"TIME_START",
268 & time_start_collect(ii,ib,iparm),0.0d0)
269 call reada(controlcard,"TIME_END",time_end_collect(ii,ib,iparm),
271 write (iout,*) "rec_start",rec_start(ii,ib,iparm),
272 & " rec_end",rec_end(ii,ib,iparm)
273 write (iout,*) "time_start",time_start_collect(ii,ib,iparm),
274 & " time_end",time_end_collect(ii,ib,iparm)
276 if (replica(iparm)) then
277 call readi(controlcard,"TOTRAJ",totraj(ii,iparm),1)
278 write (iout,*) "Number of trajectories",totraj(ii,iparm)
281 if (nfile_bin(ii,ib,iparm).lt.2
282 & .and. nfile_asc(ii,ib,iparm).eq.0
283 & .and. nfile_cx(ii,ib,iparm).eq.0) then
284 write (iout,*) "Error - no action specified!"
287 if (nfile_bin(ii,ib,iparm).gt.0) then
288 call card_concat(controlcard,.false.)
289 call split_string(controlcard,protfiles(1,1,ii,ib,iparm),
290 & maxfile_prot,nfile_bin(ii,ib,iparm))
292 write(iout,*)"nfile_bin",nfile_bin(ii,ib,iparm)
293 write(iout,*) (protfiles(i,1,ii,ib,iparm),
294 & i=1,nfile_bin(ii,ib,iparm))
297 if (nfile_asc(ii,ib,iparm).gt.0) then
298 call card_concat(controlcard,.false.)
299 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
300 & maxfile_prot,nfile_asc(ii,ib,iparm))
302 write(iout,*) "nfile_asc(ii,ib,iparm)",nfile_asc(ii,ib,iparm)
303 write(iout,*) (protfiles(i,2,ii,ib,iparm),
304 & i=1,nfile_asc(ii,ib,iparm))
306 else if (nfile_cx(ii,ib,iparm).gt.0) then
307 call card_concat(controlcard,.false.)
308 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
309 & maxfile_prot,nfile_cx(ii,ib,iparm))
311 write(iout,*) "nfile_cx(ii,ib,iparm)",nfile_cx(ii,ib,iparm)
312 write(iout,*) (protfiles(i,2,ii,ib,iparm),
313 & i=1,nfile_cx(ii,ib,iparm))
324 c-------------------------------------------------------------------------------
325 subroutine opentmp(islice,iunit,bprotfile_temp)
328 include "DIMENSIONS.ZSCOPT"
329 include "DIMENSIONS.FREE"
332 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
335 include "COMMON.IOUNITS"
336 include "COMMON.PROTFILES"
337 include "COMMON.PROT"
338 include "COMMON.FREE"
339 character*64 bprotfile_temp
340 character*3 liczba,liczba2
347 write (liczba1,'(bz,i2.2)') islice
348 write (liczba,'(bz,i3.3)') me
350 c write (iout,*) "separate_parset ",separate_parset,
352 if (separate_parset) then
353 write (liczba2,'(bz,i3.3)') myparm
354 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
355 & prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
356 open (iunit,file=bprotfile_temp,status="unknown",
357 & form="unformatted",access="direct",recl=lenrec)
359 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
360 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
361 open (iunit,file=bprotfile_temp,status="unknown",
362 & form="unformatted",access="direct",recl=lenrec)
365 bprotfile_temp = scratchdir(:ilen(scratchdir))//
366 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
367 open (iunit,file=bprotfile_temp,status="unknown",
368 & form="unformatted",access="direct",recl=lenrec)
370 c write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
375 c-------------------------------------------------------------------------------
376 subroutine read_database(*)
379 include "DIMENSIONS.ZSCOPT"
380 include "DIMENSIONS.FREE"
383 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
386 include "COMMON.CHAIN"
387 include "COMMON.IOUNITS"
388 include "COMMON.PROTFILES"
389 include "COMMON.NAMES"
392 include "COMMON.ENEPS"
393 include "COMMON.PROT"
394 include "COMMON.INTERACT"
395 include "COMMON.FREE"
396 include "COMMON.SBRIDGE"
397 include "COMMON.OBCINKA"
398 real*4 csingle(3,maxres2)
399 character*64 nazwa,bprotfile_temp
402 integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
403 & ll(maxslice),mm(maxslice),if
404 integer nrec,nlines,iscor,iunit,islice
405 double precision energ
408 double precision rmsdev,energia(0:max_ene),efree,eini,temp
409 double precision prop(maxQ)
410 integer ntot_all(maxslice,0:maxprocs-1)
411 integer iparm,ib,iib,ir,nprop,nthr,npars
412 double precision etot,time
416 lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
417 lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
419 write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
429 write (iout,*) "nparmset",nparmset
437 if (replica(iparm)) then
444 do iR=1,nRR(ib,iparm)
446 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
452 IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
453 c Read conformations from binary DA files (one per batch) and write them to
454 c a binary DA scratchfile.
455 write (liczba,'(bz,i3.3)') me
456 do if=1,nfile_bin(iR,ib,iparm)
457 nazwa=protfiles(if,1,iR,ib,iparm)
458 & (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
459 open (ientin,file=nazwa,status="old",form="unformatted",
460 & access="direct",recl=lenrec2,err=1111)
463 call opentmp(islice,ientout,bprotfile_temp)
464 call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
465 & mm(islice),iR,ib,iparm)
472 IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
473 c Read conformations from multiple ASCII int files and write them to a binary
475 do if=1,nfile_asc(iR,ib,iparm)
476 nazwa=protfiles(if,2,iR,ib,iparm)
477 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
478 open(unit=ientin,file=nazwa,status='old',err=1111)
479 write(iout,*) "reading ",nazwa(:ilen(nazwa))
481 call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
484 IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
485 c Read conformations from cx files and write them to a binary
487 do if=1,nfile_cx(iR,ib,iparm)
488 nazwa=protfiles(if,2,iR,ib,iparm)
489 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
490 write(iout,*) "reading ",nazwa(:ilen(nazwa))
492 print *,"Calling cxread"
493 call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
496 write (iout,*) "exit cxread"
502 stot(islice)=stot(islice)+jj(islice)
507 write (iout,*) "IPARM",iparm
510 if (nslice.eq.1) then
512 write (liczba,'(bz,i3.3)') me
513 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
514 & prefix(:ilen(prefix))//liczba//".xbin.tmp"
516 bprotfile_temp = scratchdir(:ilen(scratchdir))//
517 & "/"//prefix(:ilen(prefix))//".xbin.tmp"
519 write(iout,*) mm(1)," conformations read",ll(1),
520 & " conformations written to ",
521 & bprotfile_temp(:ilen(bprotfile_temp))
524 write (liczba1,'(bz,i2.2)') islice
526 write (liczba,'(bz,i3.3)') me
527 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
528 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
530 bprotfile_temp = scratchdir(:ilen(scratchdir))//
531 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
533 write(iout,*) mm(islice)," conformations read",ll(islice),
534 & " conformations written to ",
535 & bprotfile_temp(:ilen(bprotfile_temp))
540 c Check if everyone has the same number of conformations
541 call MPI_Allgather(stot(1),maxslice,MPI_INTEGER,
542 & ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
547 if (stot(islice).ne.ntot_all(islice,i)) then
548 write (iout,*) "Number of conformations at processor",i,
549 & " differs from that at processor",me,
550 & stot(islice),ntot_all(islice,i)," slice",islice
558 write (iout,*) "Numbers of conformations read by processors"
561 write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
563 write (iout,*) "Calculation terminated."
568 ntot(islice)=stot(islice)
572 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
576 c------------------------------------------------------------------------------
577 subroutine card_concat(card,to_upper)
579 include 'DIMENSIONS.ZSCOPT'
580 include "COMMON.IOUNITS"
582 character*80 karta,ucase
586 read (inp,'(a)') karta
587 if (to_upper) karta=ucase(karta)
589 do while (karta(80:80).eq.'&')
590 card=card(:ilen(card)+1)//karta(:79)
591 read (inp,'(a)') karta
592 if (to_upper) karta=ucase(karta)
594 card=card(:ilen(card)+1)//karta
597 c------------------------------------------------------------------------------
598 subroutine readi(rekord,lancuch,wartosc,default)
600 character*(*) rekord,lancuch
601 integer wartosc,default
604 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
609 iread=iread+ilen(lancuch)+1
610 read (rekord(iread:),*) wartosc
613 c----------------------------------------------------------------------------
614 subroutine reada(rekord,lancuch,wartosc,default)
616 character*(*) rekord,lancuch
618 double precision wartosc,default
621 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
626 iread=iread+ilen(lancuch)+1
627 read (rekord(iread:),*) wartosc
630 c----------------------------------------------------------------------------
631 subroutine multreadi(rekord,lancuch,tablica,dim,default)
634 integer tablica(dim),default
635 character*(*) rekord,lancuch
642 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
643 if (iread.eq.0) return
644 iread=iread+ilen(lancuch)+1
645 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
648 c----------------------------------------------------------------------------
649 subroutine multreada(rekord,lancuch,tablica,dim,default)
652 double precision tablica(dim),default
653 character*(*) rekord,lancuch
660 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
661 if (iread.eq.0) return
662 iread=iread+ilen(lancuch)+1
663 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
666 c----------------------------------------------------------------------------
667 subroutine reads(rekord,lancuch,wartosc,default)
669 character*(*) rekord,lancuch,wartosc,default
671 integer ilen,lenlan,lenrec,iread,ireade
677 iread=index(rekord,lancuch(:lenlan)//"=")
678 c print *,"rekord",rekord," lancuch",lancuch
679 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
685 c print *,"iread",iread
686 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
687 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
689 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
691 c print *,"iread",iread
692 if (iread.gt.lenrec) then
697 c print *,"ireade",ireade
698 do while (ireade.lt.lenrec .and.
699 & .not.iblnk(rekord(ireade:ireade)))
702 wartosc=rekord(iread:ireade)
705 c----------------------------------------------------------------------------
706 subroutine multreads(rekord,lancuch,tablica,dim,default)
709 character*(*) rekord,lancuch,tablica(dim),default
711 integer ilen,lenlan,lenrec,iread,ireade
720 iread=index(rekord,lancuch(:lenlan)//"=")
721 c print *,"rekord",rekord," lancuch",lancuch
722 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
723 if (iread.eq.0) return
726 c print *,"iread",iread
727 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
728 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
730 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
732 c print *,"iread",iread
733 if (iread.gt.lenrec) return
735 c print *,"ireade",ireade
736 do while (ireade.lt.lenrec .and.
737 & .not.iblnk(rekord(ireade:ireade)))
740 tablica(i)=rekord(iread:ireade)
744 c----------------------------------------------------------------------------
745 subroutine split_string(rekord,tablica,dim,nsub)
747 integer dim,nsub,i,ii,ll,kk
748 character*(*) tablica(dim)
759 C Find the start of term name
761 do while (ii.le.ll .and. rekord(ii:ii).eq." ")
764 C Parse the name into TABLICA(i) until blank found
765 do while (ii.le.ll .and. rekord(ii:ii).ne." ")
767 tablica(i)(kk:kk)=rekord(ii:ii)
770 if (kk.gt.0) nsub=nsub+1
775 c--------------------------------------------------------------------------------
776 integer function iroof(n,m)
778 if (ii*m .lt. n) ii=ii+1