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 character*800 controlcard
21 integer i,j,k,ii,n_ene_found
22 integer ind,itype1,itype2,itypf,itypsc,itypp
29 call card_concat(controlcard,.true.)
30 call readi(controlcard,"N_ENE",n_ene,max_ene)
31 if (n_ene.gt.max_ene) then
32 write (iout,*) "Error: parameter out of range: N_ENE",n_ene,
36 call readi(controlcard,"NPARMSET",nparmset,1)
37 separate_parset = index(controlcard,"SEPARATE_PARSET").gt.0
38 call readi(controlcard,"IPARMPRINT",iparmprint,1)
39 write (iout,*) "PARMPRINT",iparmprint
40 if (nparmset.gt.max_parm) then
41 write (iout,*) "Error: parameter out of range: NPARMSET",
45 energy_dec=index(controlcard,"ENERGY_DEC").gt.0
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 readi(controlcard,"NGRIDT",NGridT,400)
75 call reada(controlcard,"STARTGRIDT",StartGridT,200.0d0)
76 call reada(controlcard,"DELTA_T",Delta_T,1.0d0)
77 call reada(controlcard,"DELTRGY",deltrgy,5.0d-2)
78 call readi(controlcard,"RESCALE",rescale_mode,1)
79 check_conf=index(controlcard,"NO_CHECK_CONF").eq.0
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 call readi(controlcard,'CONSTR_DIST',constr_dist,0)
94 write (iout,*) "with_dihed_constr ",with_dihed_constr,
95 & " CONSTR_DIST",constr_dist
96 refstr = index(controlcard,'REFSTR').gt.0
97 pdbref = index(controlcard,'PDBREF').gt.0
101 c------------------------------------------------------------------------------
102 subroutine read_efree(*)
104 C Read molecular data
108 include 'DIMENSIONS.ZSCOPT'
109 include 'DIMENSIONS.COMPAR'
110 include 'DIMENSIONS.FREE'
111 include 'COMMON.IOUNITS'
112 include 'COMMON.TIME1'
113 include 'COMMON.SBRIDGE'
114 include 'COMMON.CONTROL'
115 include 'COMMON.CHAIN'
116 include 'COMMON.HEADER'
118 include 'COMMON.FREE'
119 character*320 controlcard,ucase
120 integer iparm,ib,i,j,npars
132 call card_concat(controlcard,.true.)
133 call readi(controlcard,'NT',nT_h(iparm),1)
134 write (iout,*) "iparm",iparm," nt",nT_h(iparm)
136 if (nT_h(iparm).gt.MaxT_h) then
137 write (iout,*) "Error: parameter out of range: NT",nT_h(iparm),
141 replica(iparm)=index(controlcard,"REPLICA").gt.0
142 umbrella(iparm)=index(controlcard,"UMBRELLA").gt.0
143 read_iset(iparm)=index(controlcard,"READ_ISET").gt.0
144 write (iout,*) "nQ",nQ," nT",nT_h(iparm)," replica ",
145 & replica(iparm)," umbrella ",umbrella(iparm),
146 & " read_iset",read_iset(iparm)
149 call card_concat(controlcard,.true.)
150 call readi(controlcard,'NR',nR(ib,iparm),1)
151 if (umbrella(iparm)) then
154 nRR(ib,iparm)=nR(ib,iparm)
156 if (nR(ib,iparm).gt.MaxR) then
157 write (iout,*) "Error: parameter out of range: NR",
161 call reada(controlcard,'TEMP',beta_h(ib,iparm),298.0d0)
162 beta_h(ib,iparm)=1.0d0/(beta_h(ib,iparm)*1.987D-3)
163 call multreada(controlcard,'FI',f(1,ib,iparm),nR(ib,iparm),
166 call card_concat(controlcard,.true.)
167 call multreada(controlcard,'KH',KH(1,i,ib,iparm),nQ,
169 call multreada(controlcard,'Q0',Q0(1,i,ib,iparm),nQ,
174 write (iout,*) "ib",ib," beta_h",
175 & 1.0d0/(0.001987*beta_h(ib,iparm))
176 write (iout,*) "nR",nR(ib,iparm)
177 write (iout,*) "fi",(f(i,ib,iparm),i=1,nR(ib,iparm))
179 write (iout,*) "i",i," Kh",(Kh(j,i,ib,iparm),j=1,nQ),
180 & "q0",(q0(j,i,ib,iparm),j=1,nQ)
192 nR(ib,iparm)=nR(ib,1)
193 if (umbrella(iparm)) then
196 nRR(ib,iparm)=nR(ib,1)
198 beta_h(ib,iparm)=beta_h(ib,1)
200 f(i,ib,iparm)=f(i,ib,1)
202 KH(j,i,ib,iparm)=KH(j,i,ib,1)
203 Q0(j,i,ib,iparm)=Q0(j,i,ib,1)
206 replica(iparm)=replica(1)
207 umbrella(iparm)=umbrella(1)
208 read_iset(iparm)=read_iset(1)
216 c-----------------------------------------------------------------------------
217 subroutine read_protein_data(*)
220 include "DIMENSIONS.ZSCOPT"
221 include "DIMENSIONS.FREE"
224 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
227 include "COMMON.CHAIN"
228 include "COMMON.IOUNITS"
229 include "COMMON.PROT"
230 include "COMMON.PROTFILES"
231 include "COMMON.NAMES"
232 include "COMMON.FREE"
233 include "COMMON.OBCINKA"
235 character*16000 controlcard
236 integer i,ii,ib,iR,iparm,ilen,iroof,nthr,npars
246 C Read names of files with conformation data.
247 if (replica(iparm)) then
253 do ii=1,nRR(ib,iparm)
254 write (iout,*) "Parameter set",iparm," temperature",ib,
257 call card_concat(controlcard,.true.)
258 write (iout,*) controlcard(:ilen(controlcard))
259 call readi(controlcard,"NFILE_BIN",nfile_bin(ii,ib,iparm),0)
260 call readi(controlcard,"NFILE_ASC",nfile_asc(ii,ib,iparm),0)
261 call readi(controlcard,"NFILE_CX",nfile_cx(ii,ib,iparm),0)
262 call readi(controlcard,"REC_START",rec_start(ii,ib,iparm),1)
263 call readi(controlcard,"REC_END",rec_end(ii,ib,iparm),
264 & maxstr*isampl(iparm)+rec_start(ii,ib,iparm)-1)
265 call reada(controlcard,"TIME_START",
266 & time_start_collect(ii,ib,iparm),0.0d0)
267 call reada(controlcard,"TIME_END",time_end_collect(ii,ib,iparm),
269 write (iout,*) "rec_start",rec_start(ii,ib,iparm),
270 & " rec_end",rec_end(ii,ib,iparm)
271 write (iout,*) "time_start",time_start_collect(ii,ib,iparm),
272 & " time_end",time_end_collect(ii,ib,iparm)
274 if (replica(iparm)) then
275 call readi(controlcard,"TOTRAJ",totraj(ii,iparm),1)
276 write (iout,*) "Number of trajectories",totraj(ii,iparm)
279 if (nfile_bin(ii,ib,iparm).lt.2
280 & .and. nfile_asc(ii,ib,iparm).eq.0
281 & .and. nfile_cx(ii,ib,iparm).eq.0) then
282 write (iout,*) "Error - no action specified!"
285 if (nfile_bin(ii,ib,iparm).gt.0) then
286 call card_concat(controlcard,.false.)
287 call split_string(controlcard,protfiles(1,1,ii,ib,iparm),
288 & maxfile_prot,nfile_bin(ii,ib,iparm))
290 write(iout,*)"nfile_bin",nfile_bin(ii,ib,iparm)
291 write(iout,*) (protfiles(i,1,ii,ib,iparm),
292 & i=1,nfile_bin(ii,ib,iparm))
295 if (nfile_asc(ii,ib,iparm).gt.0) then
296 call card_concat(controlcard,.false.)
297 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
298 & maxfile_prot,nfile_asc(ii,ib,iparm))
300 write(iout,*) "nfile_asc(ii,ib,iparm)",nfile_asc(ii,ib,iparm)
301 write(iout,*) (protfiles(i,2,ii,ib,iparm),
302 & i=1,nfile_asc(ii,ib,iparm))
304 else if (nfile_cx(ii,ib,iparm).gt.0) then
305 call card_concat(controlcard,.false.)
306 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
307 & maxfile_prot,nfile_cx(ii,ib,iparm))
309 write(iout,*) "nfile_cx(ii,ib,iparm)",nfile_cx(ii,ib,iparm)
310 write(iout,*) (protfiles(i,2,ii,ib,iparm),
311 & i=1,nfile_cx(ii,ib,iparm))
322 c-------------------------------------------------------------------------------
323 subroutine opentmp(islice,iunit,bprotfile_temp)
326 include "DIMENSIONS.ZSCOPT"
327 include "DIMENSIONS.FREE"
330 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
333 include "COMMON.IOUNITS"
334 include "COMMON.PROTFILES"
335 include "COMMON.PROT"
336 include "COMMON.FREE"
337 character*64 bprotfile_temp
338 character*3 liczba,liczba2
345 write (liczba1,'(bz,i2.2)') islice
346 write (liczba,'(bz,i3.3)') me
348 c write (iout,*) "separate_parset ",separate_parset,
350 if (separate_parset) then
351 write (liczba2,'(bz,i3.3)') myparm
352 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
353 & prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
354 open (iunit,file=bprotfile_temp,status="unknown",
355 & form="unformatted",access="direct",recl=lenrec)
357 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
358 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
359 open (iunit,file=bprotfile_temp,status="unknown",
360 & form="unformatted",access="direct",recl=lenrec)
363 bprotfile_temp = scratchdir(:ilen(scratchdir))//
364 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
365 open (iunit,file=bprotfile_temp,status="unknown",
366 & form="unformatted",access="direct",recl=lenrec)
368 c write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
373 c-------------------------------------------------------------------------------
374 subroutine read_database(*)
377 include "DIMENSIONS.ZSCOPT"
378 include "DIMENSIONS.FREE"
381 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
384 include "COMMON.CHAIN"
385 include "COMMON.IOUNITS"
386 include "COMMON.PROTFILES"
387 include "COMMON.NAMES"
390 include "COMMON.ENEPS"
391 include "COMMON.PROT"
392 include "COMMON.INTERACT"
393 include "COMMON.FREE"
394 include "COMMON.SBRIDGE"
395 include "COMMON.OBCINKA"
396 real*4 csingle(3,maxres2)
397 character*64 nazwa,bprotfile_temp
400 integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
401 & ll(maxslice),mm(maxslice),if
402 integer nrec,nlines,iscor,iunit,islice
403 double precision energ
406 double precision rmsdev,energia(0:max_ene),efree,eini,temp
407 double precision prop(maxQ)
408 integer ntot_all(maxslice,0:maxprocs-1)
409 integer iparm,ib,iib,ir,nprop,nthr,npars
410 double precision etot,time
414 lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
415 lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
417 write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
427 write (iout,*) "nparmset",nparmset
435 if (replica(iparm)) then
442 do iR=1,nRR(ib,iparm)
444 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
450 IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
451 c Read conformations from binary DA files (one per batch) and write them to
452 c a binary DA scratchfile.
453 write (liczba,'(bz,i3.3)') me
454 do if=1,nfile_bin(iR,ib,iparm)
455 nazwa=protfiles(if,1,iR,ib,iparm)
456 & (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
457 open (ientin,file=nazwa,status="old",form="unformatted",
458 & access="direct",recl=lenrec2,err=1111)
461 call opentmp(islice,ientout,bprotfile_temp)
462 call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
463 & mm(islice),iR,ib,iparm)
470 IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
471 c Read conformations from multiple ASCII int files and write them to a binary
473 do if=1,nfile_asc(iR,ib,iparm)
474 nazwa=protfiles(if,2,iR,ib,iparm)
475 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
476 open(unit=ientin,file=nazwa,status='old',err=1111)
477 write(iout,*) "reading ",nazwa(:ilen(nazwa))
479 call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
482 IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
483 c Read conformations from cx files and write them to a binary
485 do if=1,nfile_cx(iR,ib,iparm)
486 nazwa=protfiles(if,2,iR,ib,iparm)
487 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
488 write(iout,*) "reading ",nazwa(:ilen(nazwa))
490 print *,"Calling cxread"
491 call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
494 write (iout,*) "exit cxread"
500 stot(islice)=stot(islice)+jj(islice)
505 write (iout,*) "IPARM",iparm
508 if (nslice.eq.1) then
510 write (liczba,'(bz,i3.3)') me
511 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
512 & prefix(:ilen(prefix))//liczba//".xbin.tmp"
514 bprotfile_temp = scratchdir(:ilen(scratchdir))//
515 & "/"//prefix(:ilen(prefix))//".xbin.tmp"
517 write(iout,*) mm(1)," conformations read",ll(1),
518 & " conformations written to ",
519 & bprotfile_temp(:ilen(bprotfile_temp))
522 write (liczba1,'(bz,i2.2)') islice
524 write (liczba,'(bz,i3.3)') me
525 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
526 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
528 bprotfile_temp = scratchdir(:ilen(scratchdir))//
529 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
531 write(iout,*) mm(islice)," conformations read",ll(islice),
532 & " conformations written to ",
533 & bprotfile_temp(:ilen(bprotfile_temp))
538 c Check if everyone has the same number of conformations
539 call MPI_Allgather(stot(1),maxslice,MPI_INTEGER,
540 & ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
545 if (stot(islice).ne.ntot_all(islice,i)) then
546 write (iout,*) "Number of conformations at processor",i,
547 & " differs from that at processor",me,
548 & stot(islice),ntot_all(islice,i)," slice",islice
556 write (iout,*) "Numbers of conformations read by processors"
559 write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
561 write (iout,*) "Calculation terminated."
566 ntot(islice)=stot(islice)
570 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
574 c------------------------------------------------------------------------------
575 subroutine card_concat(card,to_upper)
577 include 'DIMENSIONS.ZSCOPT'
578 include "COMMON.IOUNITS"
580 character*80 karta,ucase
584 read (inp,'(a)') karta
585 if (to_upper) karta=ucase(karta)
587 do while (karta(80:80).eq.'&')
588 card=card(:ilen(card)+1)//karta(:79)
589 read (inp,'(a)') karta
590 if (to_upper) karta=ucase(karta)
592 card=card(:ilen(card)+1)//karta
595 c------------------------------------------------------------------------------
596 subroutine readi(rekord,lancuch,wartosc,default)
598 character*(*) rekord,lancuch
599 integer wartosc,default
602 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
607 iread=iread+ilen(lancuch)+1
608 read (rekord(iread:),*) wartosc
611 c----------------------------------------------------------------------------
612 subroutine reada(rekord,lancuch,wartosc,default)
614 character*(*) rekord,lancuch
616 double precision wartosc,default
619 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
624 iread=iread+ilen(lancuch)+1
625 read (rekord(iread:),*) wartosc
628 c----------------------------------------------------------------------------
629 subroutine multreadi(rekord,lancuch,tablica,dim,default)
632 integer tablica(dim),default
633 character*(*) rekord,lancuch
640 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
641 if (iread.eq.0) return
642 iread=iread+ilen(lancuch)+1
643 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
646 c----------------------------------------------------------------------------
647 subroutine multreada(rekord,lancuch,tablica,dim,default)
650 double precision tablica(dim),default
651 character*(*) rekord,lancuch
658 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
659 if (iread.eq.0) return
660 iread=iread+ilen(lancuch)+1
661 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
664 c----------------------------------------------------------------------------
665 subroutine reads(rekord,lancuch,wartosc,default)
667 character*(*) rekord,lancuch,wartosc,default
669 integer ilen,lenlan,lenrec,iread,ireade
675 iread=index(rekord,lancuch(:lenlan)//"=")
676 c print *,"rekord",rekord," lancuch",lancuch
677 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
683 c print *,"iread",iread
684 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
685 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
687 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
689 c print *,"iread",iread
690 if (iread.gt.lenrec) then
695 c print *,"ireade",ireade
696 do while (ireade.lt.lenrec .and.
697 & .not.iblnk(rekord(ireade:ireade)))
700 wartosc=rekord(iread:ireade)
703 c----------------------------------------------------------------------------
704 subroutine multreads(rekord,lancuch,tablica,dim,default)
707 character*(*) rekord,lancuch,tablica(dim),default
709 integer ilen,lenlan,lenrec,iread,ireade
718 iread=index(rekord,lancuch(:lenlan)//"=")
719 c print *,"rekord",rekord," lancuch",lancuch
720 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
721 if (iread.eq.0) return
724 c print *,"iread",iread
725 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
726 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
728 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
730 c print *,"iread",iread
731 if (iread.gt.lenrec) return
733 c print *,"ireade",ireade
734 do while (ireade.lt.lenrec .and.
735 & .not.iblnk(rekord(ireade:ireade)))
738 tablica(i)=rekord(iread:ireade)
742 c----------------------------------------------------------------------------
743 subroutine split_string(rekord,tablica,dim,nsub)
745 integer dim,nsub,i,ii,ll,kk
746 character*(*) tablica(dim)
757 C Find the start of term name
759 do while (ii.le.ll .and. rekord(ii:ii).eq." ")
762 C Parse the name into TABLICA(i) until blank found
763 do while (ii.le.ll .and. rekord(ii:ii).ne." ")
765 tablica(i)(kk:kk)=rekord(ii:ii)
768 if (kk.gt.0) nsub=nsub+1
773 c--------------------------------------------------------------------------------
774 integer function iroof(n,m)
776 if (ii*m .lt. n) ii=ii+1