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 include "COMMON.SBRIDGE"
22 include "COMMON.SHIELD"
23 character*800 controlcard
24 integer i,j,k,ii,n_ene_found
25 integer ind,itype1,itype2,itypf,itypsc,itypp
32 call card_concat(controlcard,.true.)
33 call readi(controlcard,"N_ENE",n_ene,max_ene)
34 if (n_ene.gt.max_ene) then
35 write (iout,*) "Error: parameter out of range: N_ENE",n_ene,
39 call readi(controlcard,"NPARMSET",nparmset,1)
40 separate_parset = index(controlcard,"SEPARATE_PARSET").gt.0
41 call readi(controlcard,"IPARMPRINT",iparmprint,1)
42 write (iout,*) "PARMPRINT",iparmprint
43 if (nparmset.gt.max_parm) then
44 write (iout,*) "Error: parameter out of range: NPARMSET",
48 call readi(controlcard,"MAXIT",maxit,5000)
49 call reada(controlcard,"FIMIN",fimin,1.0d-3)
50 call readi(controlcard,"ENSEMBLES",ensembles,0)
51 hamil_rep=index(controlcard,"HAMIL_REP").gt.0
52 write (iout,*) "Number of energy parameter sets",nparmset
53 call multreadi(controlcard,"ISAMPL",isampl,nparmset,1)
55 if (isampl(i).eq.0) then
56 write (iout,*) "ERROR: isampl is 0 for parmset",i
61 write (iout,*) "MaxSlice",MaxSlice
62 call readi(controlcard,"NSLICE",nslice,1)
64 if (nslice.gt.MaxSlice) then
65 write (iout,*) "Error: parameter out of range: NSLICE",nslice,
69 write (iout,*) "Frequency of storing conformations",
70 & (isampl(i),i=1,nparmset)
71 write (iout,*) "Maxit",maxit," Fimin",fimin
72 call readi(controlcard,"NQ",nQ,1)
74 write (iout,*) "Error: parameter out of range: NQ",nq,
79 energy_dec=(index(controlcard,'ENERGY_DEC').gt.0)
80 rmsrgymap = (index(controlcard,'RMSRGYMAP').gt.0)
81 if (index(controlcard,"CLASSIFY").gt.0) indpdb=1
82 call reada(controlcard,"DELTA",delta,1.0d-2)
83 call readi(controlcard,"EINICHECK",einicheck,2)
84 call reada(controlcard,"DELTRMS",deltrms,5.0d-2)
85 call reada(controlcard,"DELTRGY",deltrgy,5.0d-2)
86 call readi(controlcard,"RESCALE",rescale_mode,1)
87 check_conf=index(controlcard,"NO_CHECK_CONF").eq.0
88 call readi(controlcard,'TORMODE',tor_mode,0)
89 C if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0) then
90 write(iout,*) "torsional and valence angle mode",tor_mode
91 call reada(controlcard,'DISTCHAINMAX',distchainmax,50.0d0)
92 call reada(controlcard,'BOXX',boxxsize,100.0d0)
93 call reada(controlcard,'BOXY',boxysize,100.0d0)
94 call reada(controlcard,'BOXZ',boxzsize,100.0d0)
95 c Cutoff range for interactions
96 call reada(controlcard,"R_CUT",r_cut,15.0d0)
97 call reada(controlcard,"LAMBDA",rlamb,0.3d0)
98 call reada(controlcard,"LIPTHICK",lipthick,0.0d0)
99 call reada(controlcard,"LIPAQBUF",lipbufthick,0.0d0)
100 if (lipthick.gt.0.0d0) then
101 bordliptop=(boxzsize+lipthick)/2.0
102 bordlipbot=bordliptop-lipthick
104 if ((bordliptop.gt.boxzsize).or.(bordlipbot.lt.0.0))
105 & write(iout,*) "WARNING WRONG SIZE OF LIPIDIC PHASE"
106 buflipbot=bordlipbot+lipbufthick
107 bufliptop=bordliptop-lipbufthick
108 if ((lipbufthick*2.0d0).gt.lipthick)
109 &write(iout,*) "WARNING WRONG SIZE OF LIP AQ BUF"
111 write(iout,*) "bordliptop=",bordliptop
112 write(iout,*) "bordlipbot=",bordlipbot
113 write(iout,*) "bufliptop=",bufliptop
114 write(iout,*) "buflipbot=",buflipbot
115 call readi(controlcard,'SYM',symetr,1)
116 write (iout,*) "DISTCHAINMAX",distchainmax
117 write (iout,*) "delta",delta
118 write (iout,*) "einicheck",einicheck
119 write (iout,*) "rescale_mode",rescale_mode
121 bxfile=index(controlcard,"BXFILE").gt.0
122 cxfile=index(controlcard,"CXFILE").gt.0
123 if (nslice .eq. 1 .and. .not.bxfile .and. .not.cxfile)
125 histfile=index(controlcard,"HISTFILE").gt.0
126 histout=index(controlcard,"HISTOUT").gt.0
127 entfile=index(controlcard,"ENTFILE").gt.0
128 zscfile=index(controlcard,"ZSCFILE").gt.0
129 with_dihed_constr = index(controlcard,"WITH_DIHED_CONSTR").gt.0
130 write (iout,*) "with_dihed_constr ",with_dihed_constr
131 with_theta_constr = index(controlcard,"WITH_THETA_CONSTR").gt.0
132 write (iout,*) "with_theta_constr ",with_theta_constr
133 call readi(controlcard,'SHIELD',shield_mode,0)
134 write(iout,*) "shield_mode",shield_mode
136 call readi(controlcard,'TORMODE',tor_mode,0)
137 write(iout,*) "torsional and valence angle mode",tor_mode
138 if (shield_mode.gt.0) then
140 C VSolvSphere the volume of solving sphere
142 C rpp(1,1) is the energy r0 for peptide group contact and will be used for it
143 C there will be no distinction between proline peptide group and normal peptide
144 C group in case of shielding parameters
145 VSolvSphere=4.0/3.0*pi*rpp(1,1)**3
146 VSolvSphere_div=VSolvSphere-4.0/3.0*pi*(rpp(1,1)/2.0)**3
147 write (iout,*) VSolvSphere,VSolvSphere_div
148 C long axis of side chain
150 C long_r_sidechain(i)=vbldsc0(1,i)
151 C short_r_sidechain(i)=sigma0(i)
156 call readi(controlcard,'CONSTR_DIST',constr_dist,0)
157 dyn_ss=index(controlcard,"DYN_SS").gt.0
158 adaptive = index(controlcard,"ADAPTIVE").gt.0
161 c------------------------------------------------------------------------------
162 subroutine read_efree(*)
164 C Read molecular data
168 include 'DIMENSIONS.ZSCOPT'
169 include 'DIMENSIONS.COMPAR'
170 include 'DIMENSIONS.FREE'
171 include 'COMMON.IOUNITS'
172 include 'COMMON.TIME1'
173 include 'COMMON.SBRIDGE'
174 include 'COMMON.CONTROL'
175 include 'COMMON.CHAIN'
176 include 'COMMON.HEADER'
178 include 'COMMON.FREE'
179 character*320 controlcard,ucase
180 integer iparm,ib,i,j,npars
192 call card_concat(controlcard,.true.)
193 call readi(controlcard,'NT',nT_h(iparm),1)
194 write (iout,*) "iparm",iparm," nt",nT_h(iparm)
196 if (nT_h(iparm).gt.MaxT_h) then
197 write (iout,*) "Error: parameter out of range: NT",nT_h(iparm),
201 replica(iparm)=index(controlcard,"REPLICA").gt.0
202 umbrella(iparm)=index(controlcard,"UMBRELLA").gt.0
203 read_iset(iparm)=index(controlcard,"READ_ISET").gt.0
204 write (iout,*) "nQ",nQ," nT",nT_h(iparm)," replica ",
205 & replica(iparm)," umbrella ",umbrella(iparm),
206 & " read_iset",read_iset(iparm)
209 call card_concat(controlcard,.true.)
210 call readi(controlcard,'NR',nR(ib,iparm),1)
211 if (umbrella(iparm)) then
214 nRR(ib,iparm)=nR(ib,iparm)
216 if (nR(ib,iparm).gt.MaxR) then
217 write (iout,*) "Error: parameter out of range: NR",
221 call reada(controlcard,'TEMP',beta_h(ib,iparm),298.0d0)
222 beta_h(ib,iparm)=1.0d0/(beta_h(ib,iparm)*1.987D-3)
223 call multreada(controlcard,'FI',f(1,ib,iparm),nR(ib,iparm),
226 call card_concat(controlcard,.true.)
227 call multreada(controlcard,'KH',KH(1,i,ib,iparm),nQ,
229 call multreada(controlcard,'Q0',Q0(1,i,ib,iparm),nQ,
234 write (iout,*) "ib",ib," beta_h",
235 & 1.0d0/(0.001987*beta_h(ib,iparm))
236 write (iout,*) "nR",nR(ib,iparm)
237 write (iout,*) "fi",(f(i,ib,iparm),i=1,nR(ib,iparm))
239 write (iout,*) "i",i," Kh",(Kh(j,i,ib,iparm),j=1,nQ),
240 & "q0",(q0(j,i,ib,iparm),j=1,nQ)
252 nR(ib,iparm)=nR(ib,1)
253 if (umbrella(iparm)) then
256 nRR(ib,iparm)=nR(ib,1)
258 beta_h(ib,iparm)=beta_h(ib,1)
260 f(i,ib,iparm)=f(i,ib,1)
262 KH(j,i,ib,iparm)=KH(j,i,ib,1)
263 Q0(j,i,ib,iparm)=Q0(j,i,ib,1)
266 replica(iparm)=replica(1)
267 umbrella(iparm)=umbrella(1)
268 read_iset(iparm)=read_iset(1)
276 c-----------------------------------------------------------------------------
277 subroutine read_protein_data(*)
280 include "DIMENSIONS.ZSCOPT"
281 include "DIMENSIONS.FREE"
284 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
287 include "COMMON.CHAIN"
288 include "COMMON.IOUNITS"
289 include "COMMON.PROT"
290 include "COMMON.PROTFILES"
291 include "COMMON.NAMES"
292 include "COMMON.FREE"
293 include "COMMON.OBCINKA"
295 character*16000 controlcard
296 integer i,ii,ib,iR,iparm,ilen,iroof,nthr,npars
306 C Read names of files with conformation data.
307 if (replica(iparm)) then
313 do ii=1,nRR(ib,iparm)
314 write (iout,*) "Parameter set",iparm," temperature",ib,
317 call card_concat(controlcard,.true.)
318 write (iout,*) controlcard(:ilen(controlcard))
319 call readi(controlcard,"NFILE_BIN",nfile_bin(ii,ib,iparm),0)
320 call readi(controlcard,"NFILE_ASC",nfile_asc(ii,ib,iparm),0)
321 call readi(controlcard,"NFILE_CX",nfile_cx(ii,ib,iparm),0)
322 call readi(controlcard,"REC_START",rec_start(ii,ib,iparm),1)
323 call readi(controlcard,"REC_END",rec_end(ii,ib,iparm),
324 & maxstr*isampl(iparm)+rec_start(ii,ib,iparm)-1)
325 call reada(controlcard,"TIME_START",
326 & time_start_collect(ii,ib,iparm),0.0d0)
327 call reada(controlcard,"TIME_END",time_end_collect(ii,ib,iparm),
329 write (iout,*) "rec_start",rec_start(ii,ib,iparm),
330 & " rec_end",rec_end(ii,ib,iparm)
331 write (iout,*) "time_start",time_start_collect(ii,ib,iparm),
332 & " time_end",time_end_collect(ii,ib,iparm)
334 if (replica(iparm)) then
335 call readi(controlcard,"TOTRAJ",totraj(ii,iparm),1)
336 write (iout,*) "Number of trajectories",totraj(ii,iparm)
339 if (nfile_bin(ii,ib,iparm).lt.2
340 & .and. nfile_asc(ii,ib,iparm).eq.0
341 & .and. nfile_cx(ii,ib,iparm).eq.0) then
342 write (iout,*) "Error - no action specified!"
345 if (nfile_bin(ii,ib,iparm).gt.0) then
346 call card_concat(controlcard,.false.)
347 call split_string(controlcard,protfiles(1,1,ii,ib,iparm),
348 & maxfile_prot,nfile_bin(ii,ib,iparm))
350 write(iout,*)"nfile_bin",nfile_bin(ii,ib,iparm)
351 write(iout,*) (protfiles(i,1,ii,ib,iparm),
352 & i=1,nfile_bin(ii,ib,iparm))
355 if (nfile_asc(ii,ib,iparm).gt.0) then
356 call card_concat(controlcard,.false.)
357 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
358 & maxfile_prot,nfile_asc(ii,ib,iparm))
360 write(iout,*) "nfile_asc(ii,ib,iparm)",nfile_asc(ii,ib,iparm)
361 write(iout,*) (protfiles(i,2,ii,ib,iparm),
362 & i=1,nfile_asc(ii,ib,iparm))
364 else if (nfile_cx(ii,ib,iparm).gt.0) then
365 call card_concat(controlcard,.false.)
366 call split_string(controlcard,protfiles(1,2,ii,ib,iparm),
367 & maxfile_prot,nfile_cx(ii,ib,iparm))
369 write(iout,*) "nfile_cx(ii,ib,iparm)",nfile_cx(ii,ib,iparm)
370 write(iout,*) (protfiles(i,2,ii,ib,iparm),
371 & i=1,nfile_cx(ii,ib,iparm))
382 c-------------------------------------------------------------------------------
383 subroutine opentmp(islice,iunit,bprotfile_temp)
386 include "DIMENSIONS.ZSCOPT"
387 include "DIMENSIONS.FREE"
390 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
393 include "COMMON.IOUNITS"
394 include "COMMON.PROTFILES"
395 include "COMMON.PROT"
396 include "COMMON.FREE"
397 character*64 bprotfile_temp
398 character*3 liczba,liczba2
405 write (liczba1,'(bz,i2.2)') islice
406 write (liczba,'(bz,i3.3)') me
408 c write (iout,*) "separate_parset ",separate_parset,
410 if (separate_parset) then
411 write (liczba2,'(bz,i3.3)') myparm
412 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
413 & prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
414 open (iunit,file=bprotfile_temp,status="unknown",
415 & form="unformatted",access="direct",recl=lenrec)
417 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
418 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
419 open (iunit,file=bprotfile_temp,status="unknown",
420 & form="unformatted",access="direct",recl=lenrec)
423 bprotfile_temp = scratchdir(:ilen(scratchdir))//
424 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
425 open (iunit,file=bprotfile_temp,status="unknown",
426 & form="unformatted",access="direct",recl=lenrec)
428 c write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
433 c-------------------------------------------------------------------------------
434 subroutine read_database(*)
437 include "DIMENSIONS.ZSCOPT"
438 include "DIMENSIONS.FREE"
441 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
444 include "COMMON.CHAIN"
445 include "COMMON.IOUNITS"
446 include "COMMON.PROTFILES"
447 include "COMMON.NAMES"
450 include "COMMON.ENEPS"
451 include "COMMON.PROT"
452 include "COMMON.INTERACT"
453 include "COMMON.FREE"
454 include "COMMON.SBRIDGE"
455 include "COMMON.OBCINKA"
456 real*4 csingle(3,maxres2)
457 character*64 nazwa,bprotfile_temp
460 integer i,j,ii,jj(maxslice),k,kk(maxslice),l,
461 & ll(maxslice),mm(maxslice),if
462 integer nrec,nlines,iscor,iunit,islice
463 double precision energ
466 double precision rmsdev,energia(0:max_ene),efree,eini,temp
467 double precision prop(maxQ)
468 integer ntot_all(maxslice,0:maxprocs-1), maxslice_buff
469 integer iparm,ib,iib,ir,nprop,nthr,npars
470 double precision etot,time
474 lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
475 lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
477 write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,
487 write (iout,*) "nparmset",nparmset
495 if (replica(iparm)) then
502 do iR=1,nRR(ib,iparm)
504 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
510 IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
511 c Read conformations from binary DA files (one per batch) and write them to
512 c a binary DA scratchfile.
513 write (liczba,'(bz,i3.3)') me
514 do if=1,nfile_bin(iR,ib,iparm)
515 nazwa=protfiles(if,1,iR,ib,iparm)
516 & (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
517 open (ientin,file=nazwa,status="old",form="unformatted",
518 & access="direct",recl=lenrec2,err=1111)
521 call opentmp(islice,ientout,bprotfile_temp)
522 call bxread(nazwa,ii,jj(islice),kk(islice),ll(islice),
523 & mm(islice),iR,ib,iparm)
530 IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
531 c Read conformations from multiple ASCII int files and write them to a binary
533 do if=1,nfile_asc(iR,ib,iparm)
534 nazwa=protfiles(if,2,iR,ib,iparm)
535 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
536 open(unit=ientin,file=nazwa,status='old',err=1111)
537 write(iout,*) "reading ",nazwa(:ilen(nazwa))
539 call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
542 IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
543 c Read conformations from cx files and write them to a binary
545 do if=1,nfile_cx(iR,ib,iparm)
546 nazwa=protfiles(if,2,iR,ib,iparm)
547 & (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
548 write(iout,*) "reading ",nazwa(:ilen(nazwa))
550 print *,"Calling cxread"
551 call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,
554 write (iout,*) "exit cxread"
560 stot(islice)=stot(islice)+jj(islice)
565 write (iout,*) "IPARM",iparm
568 if (nslice.eq.1) then
570 write (liczba,'(bz,i3.3)') me
571 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
572 & prefix(:ilen(prefix))//liczba//".xbin.tmp"
574 bprotfile_temp = scratchdir(:ilen(scratchdir))//
575 & "/"//prefix(:ilen(prefix))//".xbin.tmp"
577 write(iout,*) mm(1)," conformations read",ll(1),
578 & " conformations written to ",
579 & bprotfile_temp(:ilen(bprotfile_temp))
582 write (liczba1,'(bz,i2.2)') islice
584 write (liczba,'(bz,i3.3)') me
585 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"//
586 & prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
588 bprotfile_temp = scratchdir(:ilen(scratchdir))//
589 & "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
591 write(iout,*) mm(islice)," conformations read",ll(islice),
592 & " conformations written to ",
593 & bprotfile_temp(:ilen(bprotfile_temp))
598 c Check if everyone has the same number of conformations
600 c call MPI_ALLgather(MPI_IN_PLACE,stot(1),MPI_DATATYPE_NULL,
601 c & ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
603 maxslice_buff=maxslice
605 call MPI_Allgather(stot(1),maxslice_buff,MPI_INTEGER,
606 & ntot_all(1,0),maxslice,MPI_INTEGER,MPI_Comm_World,IERROR)
611 if (stot(islice).ne.ntot_all(islice,i)) then
612 write (iout,*) "Number of conformations at processor",i,
613 & " differs from that at processor",me,
614 & stot(islice),ntot_all(islice,i)," slice",islice
622 write (iout,*) "Numbers of conformations read by processors"
625 write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
627 write (iout,*) "Calculation terminated."
632 ntot(islice)=stot(islice)
636 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
640 c------------------------------------------------------------------------------
641 subroutine card_concat(card,to_upper)
643 include 'DIMENSIONS.ZSCOPT'
644 include "COMMON.IOUNITS"
646 character*80 karta,ucase
650 read (inp,'(a)') karta
651 if (to_upper) karta=ucase(karta)
653 do while (karta(80:80).eq.'&')
654 card=card(:ilen(card)+1)//karta(:79)
655 read (inp,'(a)') karta
656 if (to_upper) karta=ucase(karta)
658 card=card(:ilen(card)+1)//karta
661 c------------------------------------------------------------------------------
662 subroutine readi(rekord,lancuch,wartosc,default)
664 character*(*) rekord,lancuch
665 integer wartosc,default
668 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
673 iread=iread+ilen(lancuch)+1
674 read (rekord(iread:),*) wartosc
677 c----------------------------------------------------------------------------
678 subroutine reada(rekord,lancuch,wartosc,default)
680 character*(*) rekord,lancuch
682 double precision wartosc,default
685 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
690 iread=iread+ilen(lancuch)+1
691 read (rekord(iread:),*) wartosc
694 c----------------------------------------------------------------------------
695 subroutine multreadi(rekord,lancuch,tablica,dim,default)
698 integer tablica(dim),default
699 character*(*) rekord,lancuch
706 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
707 if (iread.eq.0) return
708 iread=iread+ilen(lancuch)+1
709 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
712 c----------------------------------------------------------------------------
713 subroutine multreada(rekord,lancuch,tablica,dim,default)
716 double precision tablica(dim),default
717 character*(*) rekord,lancuch
724 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
725 if (iread.eq.0) return
726 iread=iread+ilen(lancuch)+1
727 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
730 c----------------------------------------------------------------------------
731 subroutine reads(rekord,lancuch,wartosc,default)
733 character*(*) rekord,lancuch,wartosc,default
735 integer ilen,lenlan,lenrec,iread,ireade
741 iread=index(rekord,lancuch(:lenlan)//"=")
742 c print *,"rekord",rekord," lancuch",lancuch
743 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
749 c print *,"iread",iread
750 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
751 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
753 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
755 c print *,"iread",iread
756 if (iread.gt.lenrec) then
761 c print *,"ireade",ireade
762 do while (ireade.lt.lenrec .and.
763 & .not.iblnk(rekord(ireade:ireade)))
766 wartosc=rekord(iread:ireade)
769 c----------------------------------------------------------------------------
770 subroutine multreads(rekord,lancuch,tablica,dim,default)
773 character*(*) rekord,lancuch,tablica(dim),default
775 integer ilen,lenlan,lenrec,iread,ireade
784 iread=index(rekord,lancuch(:lenlan)//"=")
785 c print *,"rekord",rekord," lancuch",lancuch
786 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
787 if (iread.eq.0) return
790 c print *,"iread",iread
791 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
792 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
794 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
796 c print *,"iread",iread
797 if (iread.gt.lenrec) return
799 c print *,"ireade",ireade
800 do while (ireade.lt.lenrec .and.
801 & .not.iblnk(rekord(ireade:ireade)))
804 tablica(i)=rekord(iread:ireade)
808 c----------------------------------------------------------------------------
809 subroutine split_string(rekord,tablica,dim,nsub)
811 integer dim,nsub,i,ii,ll,kk
812 character*(*) tablica(dim)
823 C Find the start of term name
825 do while (ii.le.ll .and. rekord(ii:ii).eq." ")
828 C Parse the name into TABLICA(i) until blank found
829 do while (ii.le.ll .and. rekord(ii:ii).ne." ")
831 tablica(i)(kk:kk)=rekord(ii:ii)
834 if (kk.gt.0) nsub=nsub+1
839 c--------------------------------------------------------------------------------
840 integer function iroof(n,m)
842 if (ii*m .lt. n) ii=ii+1
846 c-------------------------------------------------------------------------------
847 subroutine read_dist_constr
848 implicit real*8 (a-h,o-z)
850 include 'COMMON.CONTROL'
851 include 'COMMON.CHAIN'
852 include 'COMMON.IOUNITS'
853 include 'COMMON.SBRIDGE'
854 integer ifrag_(2,100),ipair_(2,100)
855 double precision wfrag_(100),wpair_(100)
856 character*500 controlcard
859 write (iout,*) "Calling read_dist_constr"
860 c write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
862 call card_concat(controlcard,.true.)
863 call readi(controlcard,"NFRAG",nfrag_,0)
864 call readi(controlcard,"NPAIR",npair_,0)
865 call readi(controlcard,"NDIST",ndist_,0)
866 call reada(controlcard,'DIST_CUT',dist_cut,5.0d0)
867 call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0)
868 call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0)
869 call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0)
870 call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0)
871 normalize = index(controlcard,"NORMALIZE").gt.0
872 c write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_
873 c write (iout,*) "IFRAG"
875 c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
877 c write (iout,*) "IPAIR"
879 c write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i)
883 if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup
884 if (ifrag_(2,i).gt.nstart_sup+nsup-1)
885 & ifrag_(2,i)=nstart_sup+nsup-1
886 c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
888 if (wfrag_(i).gt.0.0d0) then
889 do j=ifrag_(1,i),ifrag_(2,i)-1
891 c write (iout,*) "j",j," k",k
893 if (constr_dist.eq.1) then
898 forcon(nhpb)=wfrag_(i)
899 else if (constr_dist.eq.2) then
900 if (ddjk.le.dist_cut) then
905 forcon(nhpb)=wfrag_(i)
912 forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2)
914 write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",
915 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
921 if (wpair_(i).gt.0.0d0) then
929 do j=ifrag_(1,ii),ifrag_(2,ii)
930 do k=ifrag_(1,jj),ifrag_(2,jj)
934 forcon(nhpb)=wpair_(i)
936 write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
937 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
942 write (iout,*) "Distance restraints as read from input"
944 if (constr_dist.eq.11) then
945 read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(i),dhpb1(i),
946 & ibecarb(i),forcon(nhpb+1),fordepth(nhpb+1)
947 c fordepth(nhpb+1)=fordepth(nhpb+1)/forcon(nhpb+1)
948 if (forcon(nhpb+1).le.0.0d0.or.fordepth(nhpb+1).le.0.0d0)cycle
950 write (iout,'(a,4i5,2f8.2,2f10.5)') "+dist.restr ",
951 & nhpb,ihpb(nhpb),jhpb(nhpb),ibecarb(nhpb),dhpb(nhpb),
952 & dhpb1(nhpb),forcon(nhpb),fordepth(nhpb)
953 if (ibecarb(i).gt.0) then
959 read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(i),dhpb1(i),
960 & ibecarb(i),forcon(nhpb+1)
961 if (forcon(nhpb+1).gt.0.0d0) then
963 if (ibecarb(i).gt.0) then
967 if (dhpb(nhpb).eq.0.0d0)
968 & dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
970 write (iout,'(a,4i5,f8.2,f10.1)') "+dist.restr ",
971 & nhpb,ihpb(nhpb),jhpb(nhpb),ibecarb(i),dhpb(nhpb),forcon(nhpb)
973 C read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1)
974 C if (forcon(nhpb+1).gt.0.0d0) then
976 C dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
978 if (constr_dist.eq.11 .and. normalize) then
979 fordepthmax=fordepth(1)
981 if (fordepth(i).gt.fordepthmax) fordepthmax=fordepth(i)
984 fordepth(i)=fordepth(i)/fordepthmax
986 write (iout,'(/a/4a5,2a8,2a10)')
987 & "Normalized Lorenzian-like distance restraints",
988 & " Nr"," res1"," res2"," beta"," d1"," d2"," k"," V"
990 write (iout,'(4i5,2f8.2,2f10.5)')i,ihpb(i),jhpb(i),ibecarb(i),
991 & dhpb(i),dhpb1(i),forcon(i),fordepth(i)