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 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 call readi(controlcard,'CONSTR_HOMOL',constr_homology,0)
97 write (iout,*) "with_homology_constr ",with_dihed_constr,
98 & " CONSTR_HOMOLOGY",constr_homology
99 refstr = index(controlcard,'REFSTR').gt.0
100 pdbref = index(controlcard,'PDBREF').gt.0
101 dyn_ss=(index(controlcard,'DYN_SS').gt.0)
102 C /06/28/2013 Adasko: dyn_ss is keyword allowing to break and create bond
103 C disulfide bond. Note that in conterary to dynamics this in
104 C CONTROLCARD. The bond are read in molread_zs.F
108 c------------------------------------------------------------------------------
109 subroutine read_efree(*)
111 C Read molecular data
115 include 'DIMENSIONS.ZSCOPT'
116 include 'DIMENSIONS.COMPAR'
117 include 'DIMENSIONS.FREE'
118 include 'COMMON.IOUNITS'
119 include 'COMMON.TIME1'
120 include 'COMMON.SBRIDGE'
121 include 'COMMON.CONTROL'
122 include 'COMMON.CHAIN'
123 include 'COMMON.HEADER'
125 include 'COMMON.FREE'
126 character*320 controlcard,ucase
127 integer iparm,ib,i,j,npars
139 call card_concat(controlcard,.true.)
140 call readi(controlcard,'NT',nT_h(iparm),1)
141 write (iout,*) "iparm",iparm," nt",nT_h(iparm)
143 if (nT_h(iparm).gt.MaxT_h) then
144 write (iout,*) "Error: parameter out of range: NT",nT_h(iparm),
148 replica(iparm)=index(controlcard,"REPLICA").gt.0
149 umbrella(iparm)=index(controlcard,"UMBRELLA").gt.0
150 read_iset(iparm)=index(controlcard,"READ_ISET").gt.0
151 write (iout,*) "nQ",nQ," nT",nT_h(iparm)," replica ",
152 & replica(iparm)," umbrella ",umbrella(iparm),
153 & " read_iset",read_iset(iparm)
156 call card_concat(controlcard,.true.)
157 call readi(controlcard,'NR',nR(ib,iparm),1)
158 if (umbrella(iparm)) then
161 nRR(ib,iparm)=nR(ib,iparm)
163 if (nR(ib,iparm).gt.MaxR) then
164 write (iout,*) "Error: parameter out of range: NR",
168 call reada(controlcard,'TEMP',beta_h(ib,iparm),298.0d0)
169 beta_h(ib,iparm)=1.0d0/(beta_h(ib,iparm)*1.987D-3)
170 call multreada(controlcard,'FI',f(1,ib,iparm),nR(ib,iparm),
173 call card_concat(controlcard,.true.)
174 call multreada(controlcard,'KH',KH(1,i,ib,iparm),nQ,
176 call multreada(controlcard,'Q0',Q0(1,i,ib,iparm),nQ,
181 write (iout,*) "ib",ib," beta_h",
182 & 1.0d0/(0.001987*beta_h(ib,iparm))
183 write (iout,*) "nR",nR(ib,iparm)
184 write (iout,*) "fi",(f(i,ib,iparm),i=1,nR(ib,iparm))
186 write (iout,*) "i",i," Kh",(Kh(j,i,ib,iparm),j=1,nQ),
187 & "q0",(q0(j,i,ib,iparm),j=1,nQ)
199 nR(ib,iparm)=nR(ib,1)
200 if (umbrella(iparm)) then
203 nRR(ib,iparm)=nR(ib,1)
205 beta_h(ib,iparm)=beta_h(ib,1)
207 f(i,ib,iparm)=f(i,ib,1)
209 KH(j,i,ib,iparm)=KH(j,i,ib,1)
210 Q0(j,i,ib,iparm)=Q0(j,i,ib,1)
213 replica(iparm)=replica(1)
214 umbrella(iparm)=umbrella(1)
215 read_iset(iparm)=read_iset(1)
223 c-----------------------------------------------------------------------------
224 subroutine read_protein_data(*)
227 include "DIMENSIONS.ZSCOPT"
228 include "DIMENSIONS.FREE"
231 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
234 include "COMMON.CHAIN"
235 include "COMMON.IOUNITS"
236 include "COMMON.PROT"
237 include "COMMON.PROTFILES"
238 include "COMMON.NAMES"
239 include "COMMON.FREE"
240 include "COMMON.OBCINKA"
242 character*16000 controlcard
243 integer i,ii,ib,iR,iparm,ilen,iroof,nthr,npars
253 C Read names of files with conformation data.
254 if (replica(iparm)) then
260 do ii=1,nRR(ib,iparm)
261 write (iout,*) "Parameter set",iparm," temperature",ib,
264 call card_concat(controlcard,.true.)
265 write (iout,*) controlcard(:ilen(controlcard))
266 call readi(controlcard,"NFILE_BIN",nfile_bin(ii,ib,iparm),0)
267 call readi(controlcard,"NFILE_ASC",nfile_asc(ii,ib,iparm),0)
268 call readi(controlcard,"NFILE_CX",nfile_cx(ii,ib,iparm),0)
269 call readi(controlcard,"REC_START",rec_start(ii,ib,iparm),1)
270 call readi(controlcard,"REC_END",rec_end(ii,ib,iparm),
271 & maxstr*isampl(iparm)+rec_start(ii,ib,iparm)-1)
272 call reada(controlcard,"TIME_START",
273 & time_start_collect(ii,ib,iparm),0.0d0)
274 call reada(controlcard,"TIME_END",time_end_collect(ii,ib,iparm),
276 write (iout,*) "rec_start",rec_start(ii,ib,iparm),
277 & " rec_end",rec_end(ii,ib,iparm)
278 write (iout,*) "time_start",time_start_collect(ii,ib,iparm),
279 & " time_end",time_end_collect(ii,ib,iparm)
281 if (replica(iparm)) then
282 call readi(controlcard,"TOTRAJ",totraj(ii,iparm),1)
283 write (iout,*) "Number of trajectories",totraj(ii,iparm)
286 if (nfile_bin(ii,ib,iparm).lt.2
287 & .and. nfile_asc(ii,ib,iparm).eq.0
288 & .and. nfile_cx(ii,ib,iparm).eq.0) then
289 write (iout,*) "Error - no action specified!"
292 if (nfile_bin(ii,ib,iparm).gt.0) then
293 call card_concat(controlcard,.false.)
294 call split_string(controlcard,protfiles(1,1,ii,ib,iparm),
295 & maxfile_prot,nfile_bin(ii,ib,iparm))
297 write(iout,*)"nfile_bin",nfile_bin(ii,ib,iparm)
298 write(iout,*) (protfiles(i,1,ii,ib,iparm),
299 & i=1,nfile_bin(ii,ib,iparm))
302 if (nfile_asc(ii,ib,iparm).gt.0) then
303 call card_concat(controlcard,.false.)
304 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
305 & maxfile_prot,nfile_asc(ii,ib,iparm))
307 write(iout,*) "nfile_asc(ii,ib,iparm)",nfile_asc(ii,ib,iparm)
308 write(iout,*) (protfiles(i,2,ii,ib,iparm),
309 & i=1,nfile_asc(ii,ib,iparm))
311 else if (nfile_cx(ii,ib,iparm).gt.0) then
312 call card_concat(controlcard,.false.)
313 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
314 & maxfile_prot,nfile_cx(ii,ib,iparm))
316 write(iout,*) "nfile_cx(ii,ib,iparm)",nfile_cx(ii,ib,iparm)
317 write(iout,*) (protfiles(i,2,ii,ib,iparm),
318 & i=1,nfile_cx(ii,ib,iparm))
329 c-------------------------------------------------------------------------------
330 subroutine opentmp(islice,iunit,bprotfile_temp)
333 include "DIMENSIONS.ZSCOPT"
334 include "DIMENSIONS.FREE"
337 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
340 include "COMMON.IOUNITS"
341 include "COMMON.PROTFILES"
342 include "COMMON.PROT"
343 include "COMMON.FREE"
344 character*64 bprotfile_temp
345 character*3 liczba,liczba2
352 write (liczba1,'(bz,i2.2)') islice
353 write (liczba,'(bz,i3.3)') me
355 c write (iout,*) "separate_parset ",separate_parset,
357 if (separate_parset) then
358 write (liczba2,'(bz,i3.3)') myparm
359 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
360 & prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
361 open (iunit,file=bprotfile_temp,status="unknown",
362 & form="unformatted",access="direct",recl=lenrec)
364 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
365 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
366 open (iunit,file=bprotfile_temp,status="unknown",
367 & form="unformatted",access="direct",recl=lenrec)
370 bprotfile_temp = scratchdir(:ilen(scratchdir))//
371 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
372 open (iunit,file=bprotfile_temp,status="unknown",
373 & form="unformatted",access="direct",recl=lenrec)
375 c write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
380 c-------------------------------------------------------------------------------
381 subroutine read_database(*)
384 include "DIMENSIONS.ZSCOPT"
385 include "DIMENSIONS.FREE"
388 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
391 include "COMMON.CHAIN"
392 include "COMMON.IOUNITS"
393 include "COMMON.PROTFILES"
394 include "COMMON.NAMES"
397 include "COMMON.ENEPS"
398 include "COMMON.PROT"
399 include "COMMON.INTERACT"
400 include "COMMON.FREE"
401 include "COMMON.SBRIDGE"
402 include "COMMON.OBCINKA"
403 real*4 csingle(3,maxres2)
404 character*64 nazwa,bprotfile_temp
407 integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
408 & ll(maxslice),mm(maxslice),if
409 integer nrec,nlines,iscor,iunit,islice
410 double precision energ
413 double precision rmsdev,energia(0:max_ene),efree,eini,temp
414 double precision prop(maxQ)
415 integer ntot_all(maxslice,0:maxprocs-1)
416 integer iparm,ib,iib,ir,nprop,nthr,npars
417 double precision etot,time
421 lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
422 lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
424 write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
434 write (iout,*) "nparmset",nparmset
442 if (replica(iparm)) then
449 do iR=1,nRR(ib,iparm)
451 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
457 IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
458 c Read conformations from binary DA files (one per batch) and write them to
459 c a binary DA scratchfile.
460 write (liczba,'(bz,i3.3)') me
461 do if=1,nfile_bin(iR,ib,iparm)
462 nazwa=protfiles(if,1,iR,ib,iparm)
463 & (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
464 open (ientin,file=nazwa,status="old",form="unformatted",
465 & access="direct",recl=lenrec2,err=1111)
468 call opentmp(islice,ientout,bprotfile_temp)
469 call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
470 & mm(islice),iR,ib,iparm)
477 IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
478 c Read conformations from multiple ASCII int files and write them to a binary
480 do if=1,nfile_asc(iR,ib,iparm)
481 nazwa=protfiles(if,2,iR,ib,iparm)
482 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
483 open(unit=ientin,file=nazwa,status='old',err=1111)
484 write(iout,*) "reading ",nazwa(:ilen(nazwa))
486 call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
489 IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
490 c Read conformations from cx files and write them to a binary
492 do if=1,nfile_cx(iR,ib,iparm)
493 nazwa=protfiles(if,2,iR,ib,iparm)
494 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
495 write(iout,*) "reading ",nazwa(:ilen(nazwa))
497 print *,"Calling cxread"
498 call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
501 write (iout,*) "exit cxread"
507 stot(islice)=stot(islice)+jj(islice)
512 write (iout,*) "IPARM",iparm
515 if (nslice.eq.1) then
517 write (liczba,'(bz,i3.3)') me
518 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
519 & prefix(:ilen(prefix))//liczba//".xbin.tmp"
521 bprotfile_temp = scratchdir(:ilen(scratchdir))//
522 & "/"//prefix(:ilen(prefix))//".xbin.tmp"
524 write(iout,*) mm(1)," conformations read",ll(1),
525 & " conformations written to ",
526 & bprotfile_temp(:ilen(bprotfile_temp))
529 write (liczba1,'(bz,i2.2)') islice
531 write (liczba,'(bz,i3.3)') me
532 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
533 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
535 bprotfile_temp = scratchdir(:ilen(scratchdir))//
536 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
538 write(iout,*) mm(islice)," conformations read",ll(islice),
539 & " conformations written to ",
540 & bprotfile_temp(:ilen(bprotfile_temp))
545 c Check if everyone has the same number of conformations
546 call MPI_Allgather(stot(1),maxslice,MPI_INTEGER,
547 & ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
552 if (stot(islice).ne.ntot_all(islice,i)) then
553 write (iout,*) "Number of conformations at processor",i,
554 & " differs from that at processor",me,
555 & stot(islice),ntot_all(islice,i)," slice",islice
563 write (iout,*) "Numbers of conformations read by processors"
566 write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
568 write (iout,*) "Calculation terminated."
573 ntot(islice)=stot(islice)
577 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
581 c------------------------------------------------------------------------------
582 subroutine card_concat(card,to_upper)
584 include 'DIMENSIONS.ZSCOPT'
585 include "COMMON.IOUNITS"
587 character*80 karta,ucase
591 read (inp,'(a)') karta
592 if (to_upper) karta=ucase(karta)
594 do while (karta(80:80).eq.'&')
595 card=card(:ilen(card)+1)//karta(:79)
596 read (inp,'(a)') karta
597 if (to_upper) karta=ucase(karta)
599 card=card(:ilen(card)+1)//karta
602 c------------------------------------------------------------------------------
603 subroutine readi(rekord,lancuch,wartosc,default)
605 character*(*) rekord,lancuch
606 integer wartosc,default
609 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
614 iread=iread+ilen(lancuch)+1
615 read (rekord(iread:),*) wartosc
618 c----------------------------------------------------------------------------
619 subroutine reada(rekord,lancuch,wartosc,default)
621 character*(*) rekord,lancuch
623 double precision wartosc,default
626 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
631 iread=iread+ilen(lancuch)+1
632 read (rekord(iread:),*) wartosc
635 c----------------------------------------------------------------------------
636 subroutine multreadi(rekord,lancuch,tablica,dim,default)
639 integer tablica(dim),default
640 character*(*) rekord,lancuch
647 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
648 if (iread.eq.0) return
649 iread=iread+ilen(lancuch)+1
650 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
653 c----------------------------------------------------------------------------
654 subroutine multreada(rekord,lancuch,tablica,dim,default)
657 double precision tablica(dim),default
658 character*(*) rekord,lancuch
665 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
666 if (iread.eq.0) return
667 iread=iread+ilen(lancuch)+1
668 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
671 c----------------------------------------------------------------------------
672 subroutine reads(rekord,lancuch,wartosc,default)
674 character*(*) rekord,lancuch,wartosc,default
676 integer ilen,lenlan,lenrec,iread,ireade
682 iread=index(rekord,lancuch(:lenlan)//"=")
683 c print *,"rekord",rekord," lancuch",lancuch
684 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
690 c print *,"iread",iread
691 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
692 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
694 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
696 c print *,"iread",iread
697 if (iread.gt.lenrec) then
702 c print *,"ireade",ireade
703 do while (ireade.lt.lenrec .and.
704 & .not.iblnk(rekord(ireade:ireade)))
707 wartosc=rekord(iread:ireade)
710 c----------------------------------------------------------------------------
711 subroutine multreads(rekord,lancuch,tablica,dim,default)
714 character*(*) rekord,lancuch,tablica(dim),default
716 integer ilen,lenlan,lenrec,iread,ireade
725 iread=index(rekord,lancuch(:lenlan)//"=")
726 c print *,"rekord",rekord," lancuch",lancuch
727 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
728 if (iread.eq.0) return
731 c print *,"iread",iread
732 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
733 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
735 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
737 c print *,"iread",iread
738 if (iread.gt.lenrec) return
740 c print *,"ireade",ireade
741 do while (ireade.lt.lenrec .and.
742 & .not.iblnk(rekord(ireade:ireade)))
745 tablica(i)=rekord(iread:ireade)
749 c----------------------------------------------------------------------------
750 subroutine split_string(rekord,tablica,dim,nsub)
752 integer dim,nsub,i,ii,ll,kk
753 character*(*) tablica(dim)
764 C Find the start of term name
766 do while (ii.le.ll .and. rekord(ii:ii).eq." ")
769 C Parse the name into TABLICA(i) until blank found
770 do while (ii.le.ll .and. rekord(ii:ii).ne." ")
772 tablica(i)(kk:kk)=rekord(ii:ii)
775 if (kk.gt.0) nsub=nsub+1
780 c--------------------------------------------------------------------------------
781 integer function iroof(n,m)
783 if (ii*m .lt. n) ii=ii+1