1 subroutine read_general_data(*)
4 include "DIMENSIONS.ZSCOPT"
8 integer ierror,kolor,klucz
10 include "COMMON.WEIGHTS"
11 include "COMMON.NAMES"
12 include "COMMON.VMCPAR"
13 include "COMMON.TORSION"
14 include "COMMON.INTERACT"
15 include "COMMON.ENERGIES"
16 include "COMMON.MINPAR"
17 include "COMMON.IOUNITS"
18 include "COMMON.TIME1"
19 include "COMMON.PROTFILES"
20 include "COMMON.CHAIN"
21 include "COMMON.CLASSES"
22 include "COMMON.THERMAL"
23 include "COMMON.FFIELD"
24 include "COMMON.OPTIM"
25 include "COMMON.CONTROL"
26 include "COMMON.SCCOR"
27 include "COMMON.SPLITELE"
28 character*800 controlcard
29 integer i,j,k,l,ii,n_ene_found,ist1,iet1,ist2,iet2,ls,le
30 integer ind,itype1,itype2,itypf,itypsc,itypp,itypt1,itypt2
31 integer ilen,lenpot,lenpre
33 character*4 liczba,liczba1
40 double precision xchuj,weitemp,weitemp_low,weitemp_up
45 c write (iout,*) "Enter read_general_data"
53 C Read first record: seed and number of energy components
54 call card_concat(controlcard,.true.)
55 c write (iout,*) "card1",controlcard
57 call readi(controlcard,"ISEED",iseed,0)
58 if (iseed.eq.0) stop "Seed is zero!"
59 c print *,me," iseed",iseed
60 call readi(controlcard,"NPARMSET",nparmset,1)
61 c print *,me," nparmset",nparmset
63 c Split processor pool if multiple parameter sets are treated
64 if (nparmset.eq.1) then
65 WHAM_COMM = MPI_COMM_WORLD
66 c print *,me," opening ",outname," opened"
67 open(iout,file=outname,status='unknown')
69 c print *,me," outname ",outname," opened"
71 if (nprocs.lt.nparmset) then
73 & "*** Cannot split parameter sets for fewer processors than sets",
75 call MPI_Finalize(ierror)
78 c write (iout,*) "nparmset",nparmset
79 nprocs = nprocs/nparmset
81 klucz = mod(me,nprocs)
82 c write (*,*) "My old rank",me," kolor",kolor," klucz",klucz
83 call MPI_Comm_split(MPI_COMM_WORLD,kolor,klucz,WHAM_COMM,ierror)
84 call MPI_Comm_size(WHAM_COMM,nprocs,ierror)
85 call MPI_Comm_rank(WHAM_COMM,me,ierror)
86 c write (*,*) "My new rank",me," comm size",nprocs
87 c write (*,*) "MPI_COMM_WORLD",MPI_COMM_WORLD,
88 c & " WHAM_COMM",WHAM_COMM
90 c write (*,*) "My parameter set is",myparm
91 write(liczba,'(bz,i4.4)') me
92 write(liczba1,'(bz,i4.4)') myparm
93 outname=prefix(:lenpre)//'.out_par'//liczba1//'_'//
94 & pot(:lenpot)//liczba
95 open(iout,file=outname,status='unknown')
98 c print *,me," checkpoint 1"
99 energy_dec=(index(controlcard,'ENERGY_DEC').gt.0)
100 call readi(controlcard,'TORMODE',tor_mode,0)
101 c print *,me," tor_mode",tor_mode
102 call readi(controlcard,'SHIELD',shield_mode,0)
103 call readi(controlcard,"N_ENE",n_ene,max_ene)
104 c print *,"iseed",iseed," n_ene",n_ene
105 call readi(controlcard,"NPARMSET",nparmset,1)
106 geom_and_wham_weights =
107 & index(controlcard,"GEOM_AND_WHAM_WEIGHTS").gt.0
108 c write (iout,*) "GEOM_AND_WHAM_WEIGHTS",geom_and_wham_weights
109 if (iseed.gt.0) iseed=-iseed
111 out_newe=index(controlcard,"OUT_NEWE").gt.0
112 unres_pdb = index(controlcard,"UNRES_PDB").gt.0
113 c write (iout,*) "UNRES_PDB ",unres_pdb
114 c Energy calculation settings section
115 call readi(controlcard,'TORMODE',tor_mode,0)
116 C if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0) then
117 call reada(controlcard,'BOXX',boxxsize,100.0d0)
118 call reada(controlcard,'BOXY',boxysize,100.0d0)
119 call reada(controlcard,'BOXZ',boxzsize,100.0d0)
120 c print *,me," checkpoint 2"
121 c Cutoff range for interactions
122 call reada(controlcard,"R_CUT",r_cut,15.0d0)
123 call reada(controlcard,"LAMBDA",rlamb,0.3d0)
124 call reada(controlcard,"LIPTHICK",lipthick,0.0d0)
125 call reada(controlcard,"LIPAQBUF",lipbufthick,0.0d0)
126 if (lipthick.gt.0.0d0) then
127 bordliptop=(boxzsize+lipthick)/2.0
128 bordlipbot=bordliptop-lipthick
130 if ((bordliptop.gt.boxzsize).or.(bordlipbot.lt.0.0))
131 & write(iout,*) "WARNING WRONG SIZE OF LIPIDIC PHASE"
132 buflipbot=bordlipbot+lipbufthick
133 bufliptop=bordliptop-lipbufthick
134 if ((lipbufthick*2.0d0).gt.lipthick)
135 &write(iout,*) "WARNING WRONG SIZE OF LIP AQ BUF"
137 c write(iout,*) "bordliptop=",bordliptop
138 c write(iout,*) "bordlipbot=",bordlipbot
139 c write(iout,*) "bufliptop=",bufliptop
140 c write(iout,*) "buflipbot=",buflipbot
141 call readi(controlcard,'SYM',symetr,1)
142 c print *,me," checkpoint 3"
143 c write (iout,*) "n_ene",n_ene
145 c Energy-term-weights section
147 C Read third record: weights
151 c print *,me," checkpoint 4"
152 C Read fourth record: masks
153 call card_concat(controlcard,.true.)
154 c write (iout,*) "card2",controlcard
156 key = "MASK_"//wname(i)(2:ilen(wname(i)))
157 call readi(controlcard,key(:ilen(key)),imask(i),0)
159 c print *,me," checkpoint 5"
160 C Read fifth record: lower limits of weights
161 call card_concat(controlcard,.true.)
162 c write (iout,*) "card3",controlcard
164 key = "WLOW_"//wname(i)(2:ilen(wname(i)))
165 call reada(controlcard,key(:ilen(key)),ww_low(i),ww(i))
167 C Read sixth record: upper limits of weights
168 call card_concat(controlcard,.true.)
169 c write (iout,*) "card4",controlcard
171 key = "WUP_"//wname(i)(2:ilen(wname(i)))
172 call reada(controlcard,key(:ilen(key)),ww_up(i),ww(i))
174 C Read seventh record: VMC parameters
175 call card_concat(controlcard,.true.)
176 c write (iout,*) "card5",controlcard
177 call readi(controlcard,"MODE",mode,3)
178 call readi(controlcard,"NSCANCYCLE",nscancycle,3)
179 call readi(controlcard,"MAXSTEP_SCAN",maxstep_scan,50)
180 call reada(controlcard,"RSTEP_SCAN",step_scan,1.0d-1)
181 call readi(controlcard,"READ_STAT",read_stat,3)
182 call readi(controlcard,"RESCALE_MODE",rescale_mode,2)
183 c init_ene = index(controlcard,"INIT_ENE").gt.0 .and. read_stat.gt.1
185 call readi(controlcard,"NMCM",nmcm,0)
186 call readi(controlcard,"MAXSCAN",maxscan,0)
187 call readi(controlcard,"MAXMIN",maxmin,100)
188 call readi(controlcard,"MAXFUN",maxfun,100)
189 call reada(controlcard,"TOLF",tolf,1.0d-6)
190 call reada(controlcard,"RTOLF",rtolf,1.0d-6)
192 if (index(controlcard,"OUT_MINIM").gt.0) out_minim=iout
194 if (index(controlcard,"PRINT_INI").gt.0) print_ini=1
196 if (index(controlcard,"PRINT_FIN").gt.0) print_fin=1
198 if (index(controlcard,"PRINT_STAT").gt.0) print_stat=1
199 call reada(controlcard,"RSTIME",rstime,9.0d2)
200 call reads(controlcard,"MINIMIZER",minimizer,"SUMSL")
201 call readi(controlcard,"OPT_MODE",opt_mode,0)
202 mod_other_params=index(controlcard,"OPTIMIZE_OTHER").gt.0
203 if (read_stat.eq.0 .and. mod_other_params) then
204 write (iout,*) "Error: only optimization of energy-term",
205 & " weights can be performed with READ_STAT=",read_stat
209 if (index(controlcard,"RESTART").gt.0) then
214 c print *,me," checkpoint 6"
217 c-------------------------------------------------------------------------------------------------
218 subroutine read_pmf_data(*)
221 include "DIMENSIONS.ZSCOPT"
225 integer ierror,kolor,klucz
227 include "COMMON.WEIGHTS"
228 include "COMMON.NAMES"
229 include "COMMON.VMCPAR"
230 include "COMMON.TORSION"
231 include "COMMON.INTERACT"
232 include "COMMON.ENERGIES"
233 include "COMMON.MINPAR"
234 include "COMMON.IOUNITS"
235 include "COMMON.TIME1"
236 include "COMMON.PROTFILES"
237 include "COMMON.CHAIN"
238 include "COMMON.CLASSES"
239 include "COMMON.THERMAL"
240 include "COMMON.FFIELD"
241 include "COMMON.OPTIM"
242 include "COMMON.CONTROL"
243 include "COMMON.SCCOR"
244 include "COMMON.SPLITELE"
245 character*800 controlcard
246 integer i,j,k,l,ii,n_ene_found,ist1,iet1,ist2,iet2,ls,le
247 integer ind,itype1,itype2,itypf,itypsc,itypp,itypt1,itypt2
248 integer ilen,lenpot,lenpre
250 character*4 liczba,liczba1
257 double precision xchuj,weitemp,weitemp_low,weitemp_up
262 c A 2/11/18 Read PMF parameters
263 call card_concat(controlcard,.true.)
264 torsion_pmf=index(controlcard,"TORSION_PMF").gt.0
265 turn_pmf=index(controlcard,"TURN_PMF").gt.0
266 eello_pmf=index(controlcard,"EELLO_PMF").gt.0
267 write (iout,*) "TORSION_PMF", TORSION_PMF," TURN_PMF ",TURN_PMF,
268 & " EELLO_PMF",EELLO_PMF
269 call reada(controlcard,"WELLO_PMF",wello_pmf,1.0d0)
270 call reada(controlcard,"WELLO_PMF_LOW",wello_pmf_low,0.0d0)
271 call reada(controlcard,"WELLO_PMF_UP",wello_pmf_up,5.0d0)
272 call reada(controlcard,"WTURN_PMF",wturn_pmf,1.0d0)
273 call reada(controlcard,"WTURN_PMF_LOW",wturn_pmf_low,0.0d0)
274 call reada(controlcard,"WTURN_PMF_UP",wturn_pmf_up,5.0d0)
275 call reada(controlcard,"WPMF",wpmf,1.0d0)
276 write (iout,*) "nloctyp",nloctyp
277 call multreada(controlcard,"E0",e0(0,-nloctyp+1),
278 & (ntyp+1)*(2*nloctyp-1),0.0d0)
279 call multreada(controlcard,"E0_LOW",e0_low(0,-nloctyp+1),
280 & (ntyp+1)*(2*nloctyp-1),
282 call multreada(controlcard,"E0_UP",e0_up(0,-nloctyp+1),
283 & (ntyp+1)*(2*nloctyp-1),
285 write (iout,*) "WTURN_PMF",WTURN_PMF,WTURN_PMF_LOW,WTURN_PMF_UP
286 write (iout,*) "WELLO_PMF",WELLO_PMF,WELLO_PMF_LOW,WELLO_PMF_UP
290 write (iout,"(2i5,3f15.5)") i,j,e0(i,j),e0_low(i,j),e0_up(i,j)
296 c-----------------------------------------------------------------------
297 subroutine read_optim_parm(*)
300 include "DIMENSIONS.ZSCOPT"
304 integer ierror,kolor,klucz
306 include "COMMON.WEIGHTS"
307 include "COMMON.NAMES"
308 include "COMMON.VMCPAR"
309 include "COMMON.TORSION"
310 include "COMMON.LOCAL"
311 include "COMMON.INTERACT"
312 include "COMMON.ENERGIES"
313 include "COMMON.MINPAR"
314 include "COMMON.IOUNITS"
315 include "COMMON.TIME1"
316 include "COMMON.PROTFILES"
317 include "COMMON.CHAIN"
318 include "COMMON.CLASSES"
319 include "COMMON.THERMAL"
320 include "COMMON.FFIELD"
321 include "COMMON.OPTIM"
322 include "COMMON.CONTROL"
323 include "COMMON.SCCOR"
324 character*800 controlcard
326 integer i,j,k,l,ii,n_ene_found,ist1,iet1,ist2,iet2,ls,le
327 integer ind,ind1,ind2,itype1,itype2,itypf,itypsc,itypp,
328 & itypt1,itypt2,masktemp,iblock,iaux,itypa
329 integer ilen,lenpot,lenpre
331 double precision aux,vb_low,vb_up,vb_rel
332 character*4 liczba,liczba1
339 double precision xchuj,weitemp,weitemp_low,weitemp_up
341 character*3 typf,typa
344 integer ind_shield /25/
348 write (iout,*) "MOD_OTHER_PARAMS ",mod_other_params
352 mask_tor(0,itypt1,itypt2,iblock)=0
353 weitor(0,itypt1,itypt2,iblock)=1.0d0
354 weitor_low(0,itypt1,itypt2,iblock)=1.0d0
355 weitor_up(0,itypt1,itypt2,iblock)=1.0d0
361 mask_tor(l,itypt1,itypt2,iblock)=0
362 weitor(l,itypt1,itypt2,iblock)=1.0
363 weitor_low(l,itypt1,itypt2,iblock)=1.0
364 weitor_up(l,itypt1,itypt2,iblock)=1.0
370 write (iout,*) "ntyp",ntyp
373 write (iout,*) "itypt1",itypt1," itypt2",itypt2,
374 & "weitor",weitor(0,itypt1,itypt2,1),eitor(0,itypt1,itypt2,2)
378 if (mod_other_params) then
380 c &"Internal parameters of UNRES energy components will be optimized"
381 call card_concat(controlcard,.true.)
382 write (iout,*) "mod_side ",controlcard
384 mod_side = (index(controlcard,"MOD_SIDE").gt.0)
388 call card_concat(controlcard,.true.)
389 do while ( index(controlcard,"END").eq.0 )
390 call split_string(controlcard,item(1),4,nitem)
391 if (nitem.eq.1 .or. item(2)(:1).eq."*") then
392 nsingle_sc=nsingle_sc+1
393 ityp_ssc(nsingle_sc)=rescode(1,item(1),0,.false.)
395 epss_low(nsingle_sc)=0.02d0
397 read (item(3),*) epss_low(nsingle_sc)
400 epss_up(nsingle_sc)=0.0d0
402 read (item(4),*) epss_up(nsingle_sc)
406 ityp_psc(1,npair_sc)=rescode(1,item(1),0,.false.)
407 ityp_psc(2,npair_sc)=rescode(1,item(2),0,.false.)
409 epsp_low(npair_sc)=0.02d0
411 read (item(3),*) epsp_low(npair_sc)
414 epsp_up(npair_sc)=0.0d0
416 read (item(4),*) epsp_up(npair_sc)
419 call card_concat(controlcard,.true.)
421 if (nsingle_sc+npair_sc.eq.0) mod_side=.false.
422 call card_concat(controlcard,.true.)
425 & index(controlcard,"MOD_SIDE_OTHER").gt.0
426 write (iout,*) "mod_side_other ",controlcard,mod_side_other
427 if (mod_side_other) then
428 mod_side_other=.false.
434 call card_concat(controlcard,.true.)
435 c write (iout,*) "mod_side_oter ",controlcard
436 do while ( index(controlcard,"END").eq.0 )
437 call reads(controlcard,"RESKIND",reskind," ")
438 itypsc=rescode(1,reskind,0,.false.)
439 if (itypsc.lt.1 .or. itypsc.gt.20) then
440 write (iout,*) "Error in SC optimization data;",
441 & " SC type must not be dummy, type is" ,restyp," ",itypsc
442 write (*,*) "Error in SC optimization data;",
443 & " SC type must not be dummy, type is" ,restyp," ",itypsc
446 call readi(controlcard,"MASK_SIGMA0",mask_sigma(1,itypsc),0)
447 call readi(controlcard,"MASK_SIGMAII",mask_sigma(2,itypsc),0)
448 call readi(controlcard,"MASK_CHIP",mask_sigma(3,itypsc),0)
449 call readi(controlcard,"MASK_ALP",mask_sigma(4,itypsc),0)
450 call reada(controlcard,"SIGMA0_LOW",sigma_low(1,itypsc),0.d0)
451 call reada(controlcard,"SIGMA0_UP",sigma_up(1,itypsc),0.0d0)
452 call reada(controlcard,"SIGMAII_LOW",sigma_low(2,itypsc),
454 call reada(controlcard,"SIGMAII_UP",sigma_up(2,itypsc),0.0d0)
455 call reada(controlcard,"CHIP_LOW",sigma_low(3,itypsc),0.d0)
456 call reada(controlcard,"CHIP_UP",sigma_up(3,itypsc),0.0d0)
457 call reada(controlcard,"ALP_LOW",sigma_low(4,itypsc),0.d0)
458 call reada(controlcard,"ALP_UP",sigma_up(4,itypsc),0.0d0)
460 if (sigma_low(k,itypsc).eq.0.0d0 .and.
461 & sigma_up(k,itypsc).eq.0.0d0) mask_sigma(k,itypsc)=0
464 mod_side_other=mod_side_other.or.mask_sigma(k,itypsc).gt.0
466 write (iout,'(a4,i3,4(i2,2f8.3))') reskind,itypsc,
467 & (mask_sigma(k,itypsc),
468 & sigma_low(k,itypsc),sigma_up(k,itypsc),k=1,4)
469 call card_concat(controlcard,.true.)
470 c write (iout,*) "mod_side_oter ",controlcard
472 write (iout,*) "Optimization flags of one-body SC parameters"
474 write (iout,'(a4,i3,4(i2,2f8.3))') restyp(i),i,
475 & (mask_sigma(k,i),sigma_low(k,i),sigma_up(k,i),k=1,4)
477 call card_concat(controlcard,.true.)
479 c write (iout,*) "mod_side_other ",mod_side_other
480 c write (iout,*) "mod_tor ",controlcard
482 mod_tor = index(controlcard,"MOD_TOR").gt.0
485 do i=-ntortyp,ntortyp
486 do j=-ntortyp,ntortyp
487 mask_tor(0,i,j,iblock)=0
488 weitor(0,i,j,iblock)=1.0d0
489 weitor_low(0,i,j,iblock)=0.0d0
490 weitor_up(0,i,j,iblock)=2.0d0
494 call card_concat(controlcard,.true.)
495 write (iout,*) controlcard
496 do while ( index(controlcard,"END").eq.0 )
497 call split_string(controlcard,item(1),7,nitem)
499 write (*,*) "Error in torsional optimization data;",
500 & " must specify both residues and type"
508 write (iout,*) "item3 ",item(3)," item4 ",item(4),
511 if (nitem.gt.2) read(item(3),*) masktemp
512 if (nitem.gt.3) read(item(4),*) weitemp
513 if (nitem.gt.4) read(item(5),*) weitemp_low
514 if (nitem.gt.5) read(item(6),*) weitemp_up
515 if (nitem.gt.6) read(item(7),*) iblock
516 write (iout,*) controlcard
517 write (iout,*) item(1)," ",item(2),weitemp,
518 & weitemp_low,weitemp_up
519 if (index(item(1),"*").gt.0) then
523 ist1=itortyp(rescode(1,item(1),0,.false.))
526 if (index(item(2),"*").gt.0) then
530 ist2=itortyp(rescode(1,item(2),0,.false.))
533 c write (iout,*) "ist1",ist1," iet1",iet1," ist2",ist2,
538 c write (iout,*) "itypt1",itypt1," itypt2",itypt2,
539 c & weitemp,weitemp_low,weitemp_up
540 mask_tor(0,itypt1,itypt2,iblock)=masktemp
541 weitor(0,itypt1,itypt2,iblock)=weitemp
542 weitor_low(0,itypt1,itypt2,iblock)=weitemp_low
543 weitor_up(0,itypt1,itypt2,iblock)=weitemp_up
544 c write (iout,*) "itypt1",itypt1," itypt2",itypt2,
545 c & mask_tor(0,itypt1,itypt2,iblock),
546 c & weitor(0,itypt1,itypt2,iblock),
547 c & weitor_low(0,itypt1,itypt2,iblock),
548 c & weitor_up(0,itypt1,itypt2,iblock)
551 call card_concat(controlcard,.true.)
552 write (iout,*) controlcard
555 if (tor_mode.gt.1) then
556 write (iout,*) "TORMODE is",tor_mode,
557 & " torsional constants are computed from energy map expansion."
561 write (iout,*) "Initialized torsional parameters:"
563 do itypt1=-ntortyp,ntortyp
564 do itypt2=-ntortyp,ntortyp
565 write (iout,'(3i5,3f10.5)') itypt1,itypt2,
566 & mask_tor(0,itypt1,itypt2,iblock),
567 & weitor(0,itypt1,itypt2,iblock),
568 & weitor_low(0,itypt1,itypt2,iblock),
569 & weitor_up(0,itypt1,itypt2,iblock)
574 if (tor_mode.eq.1) then
575 c Exchange indices because the residue order in new torsionals is reversed
577 do itypt1=-ntortyp,ntortyp
578 do itypt2=itypt1+1,ntortyp
579 iaux=mask_tor(0,itypt1,itypt2,iblock)
580 mask_tor(0,itypt1,itypt2,iblock)=
581 & mask_tor(0,itypt2,itypt1,iblock)
582 mask_tor(0,itypt2,itypt1,iblock)=iaux
583 aux=weitor(0,itypt1,itypt2,iblock)
584 weitor(0,itypt1,itypt2,iblock)=
585 & weitor(0,itypt2,itypt1,iblock)
586 weitor(0,itypt2,itypt1,iblock)=aux
587 aux=weitor_low(0,itypt1,itypt2,iblock)
588 weitor_low(0,itypt1,itypt2,iblock)=
589 & weitor_low(0,itypt2,itypt1,iblock)
590 weitor_low(0,itypt2,itypt1,iblock)=aux
591 aux=weitor_up(0,itypt1,itypt2,iblock)
592 weitor_up(0,itypt1,itypt2,iblock)=
593 & weitor_up(0,itypt2,itypt1,iblock)
594 weitor_up(0,itypt2,itypt1,iblock)=aux
599 call card_concat(controlcard,.true.)
601 write (iout,*) "mod_sccor ",controlcard
603 mod_sccor = index(controlcard,"MOD_SCCOR").gt.0
605 call card_concat(controlcard,.true.)
608 do i=-nsccortyp,nsccortyp
609 do j=-nsccortyp,nsccortyp
610 mask_tor(l,i,j,iblock)=0
611 weitor(l,i,j,iblock)=1.0d0
612 weitor_low(l,i,j,iblock)=0.0d0
613 weitor_up(l,i,j,iblock)=2.0d0
618 do while ( index(controlcard,"END").eq.0 )
619 call split_string(controlcard,item(1),7,nitem)
621 write (*,*) "Error in torsional optimization data;",
622 & " must specify both residues and type"
628 if (nitem.gt.3) read(item(4),*) masktemp
629 if (nitem.gt.4) read(item(5),*) weitemp
630 if (nitem.gt.5) read(item(6),*) weitemp_low
631 if (nitem.gt.6) read(item(7),*) weitemp_up
632 if (index(item(1),"*").gt.0) then
636 ist1=isccortyp(rescode(1,item(1),0,.false.))
639 if (index(item(2),"*").gt.0) then
643 ist2=isccortyp(rescode(1,item(2),0,.false.))
646 if (index(item(3),"*").gt.0) then
656 mask_tor(l,itypt1,itypt2,1)=masktemp
657 weitor(l,itypt1,itypt2,1)=weitemp
658 weitor_low(l,itypt1,itypt2,1)=weitemp_low
659 weitor_up(l,itypt1,itypt2,1)=weitemp_up
663 call card_concat(controlcard,.true.)
665 call card_concat(controlcard,.true.)
668 write (iout,*) "Optimizable sidechain-torsional parameters:"
669 do itypt1=-nsccortyp,nsccortyp
670 do itypt2=-nsccortyp,nsccortyp
672 if (mask_tor(l,itypt1,itypt2,1).gt.0)
673 & write (iout,'(4i5,3f10.5)') itypt1,itypt2,l,
674 & mask_tor(l,itypt1,itypt2,1),weitor(l,itypt1,itypt2,1),
675 & weitor_low(l,itypt1,itypt2,1),weitor_up(l,itypt1,itypt2,1)
680 mod_ang=tor_mode.gt.0 .and. index(controlcard,"MOD_ANGLE").gt.0
681 write (iout,*) "mod_angle ",controlcard
682 write (iout,*) "mod_ang",mod_ang
689 call card_concat(controlcard,.true.)
690 do while (index(controlcard,'END').eq.0)
691 write (iout,'(a)') "angle: ",controlcard
692 call reads(controlcard,"TYPE",typa," ")
693 itypa=rescode(1,typa,0,.false.)
694 c write (iout,*) "TYPA ",typa," itypa",itypa
695 if (itypa.eq.0 .or. itypa.gt.ntyp) then
696 write (*,*) "Error specifying angle parms"
701 write (iout,*) "bend type",itypa
702 call reada(controlcard,"VB_LOW",vb_low,-1.0d5)
703 call reada(controlcard,"VB_UP",vb_up,1.0d5)
704 call reada(controlcard,"VB_REL",vb_rel,0.0d0)
705 write (iout,*) "VB_LOW",vb_low," VB_UP",vb_up," VB_REL",vb_rel
706 do i=1,nbend_kcc_TB(itypa)
707 if (vb_rel.gt.0) then
708 write (iout,*) "v1bend_chyb=",v1bend_chyb(i,itypa)
709 v1bend_low(i,itypa)=v1bend_chyb(i,itypa)*(1.0d0
710 & -dsign(vb_rel,v1bend_chyb(i,itypa)))
711 v1bend_up(i,itypa)=v1bend_chyb(i,itypa)*(1.0d0
712 & +dsign(vb_rel,v1bend_chyb(i,itypa)))
714 v1bend_low(i,itypa)=vb_low
715 v1bend_up(i,itypa)=vb_up
718 call card_concat(controlcard,.true.)
721 call card_concat(controlcard,.true.)
723 write (iout,*) "Boundaries of angle potential coefficients"
724 write (iout,'(3a)') "Index"," Low bound"," Up bound"
726 if (mask_ang(i).eq.0) cycle
727 write (iout,'(2a)') 'Type ',restyp(iloctyp(i))
728 do j=1,nbend_kcc_TB(i)
729 write (iout,'(i5,2f15.1)') j,v1bend_low(j,i),v1bend_up(j,i)
735 write (iout,*) "mod_fourier ",controlcard
737 mod_fourier(nloctyp)=index(controlcard,"MOD_FOURIER")
739 if (mod_fourier(nloctyp).gt.0) then
740 mod_fourier(nloctyp)=0
754 mask_eenew(ii,j,k,i)=0
762 call card_concat(controlcard,.true.)
763 do while ( index(controlcard,"END").eq.0 )
764 c write(iout,*) controlcard
765 call reads(controlcard,"TYPF",typf," ")
766 itypf=rescode(1,typf,0,.false.)
767 c write (iout,*) "TYPF ",typf," itypf",itypf
768 if (itypf.eq.0 .or. itypf.gt.ntyp) then
769 write (*,*) "Error specifying fourier parms"
772 itypf=itype2loc(itypf)
773 write (iout,*) "local type",itypf
775 if (index(controlcard,"B1_LOW").gt.0) then
777 call readi(controlcard,"IND",ind,0)
778 call readi(controlcard,"COEFF",ii,0)
779 if (ind.eq.0 .or. ii.eq.0 .or. ii.gt.3) then
781 & "Error specifying B1, components undefined",ind,ii
784 mask_bnew1(ii,ind,itypf)=1
785 call reada(controlcard,"B1_LOW",bnew1_low(ii,ind,itypf),
787 call reada(controlcard,"B1_UP",bnew1_up(ii,ind,itypf),
789 mod_fourier(itypf)=mod_fourier(itypf)
790 & +mask_bnew1(ii,ind,itypf)
792 else if (index(controlcard,"B2_LOW").gt.0) then
794 mask_bnew2(ii,ind,itypf)=1
795 call readi(controlcard,"IND",ind,0)
796 call readi(controlcard,"COEFF",ii,0)
797 if (ind.eq.0 .or. ii.eq.0 .or. ii.gt.3) then
799 & "Error specifying B2, components undefined",ind,ii
802 call reada(controlcard,"B2_LOW",bnew2_low(ii,ind,itypf),
804 call reada(controlcard,"B2_UP",bnew2_up(ii,ind,itypf),
806 mod_fourier(itypf)=mod_fourier(itypf)
807 & +mask_bnew2(ii,ind,itypf)
809 else if (index(controlcard,"C_LOW").gt.0) then
811 call readi(controlcard,"IND",ind,0)
812 call readi(controlcard,"COEFF",ii,0)
813 if (ind.eq.0 .or. ii.eq.0 .or. ii.gt.3) then
815 & "Error specifying C, components undefined",ind,ii
818 mask_ccnew(ii,ind,itypf)=1
819 call reada(controlcard,"C_LOW",ccnew_low(ii,ind,itypf),.1d0)
820 call reada(controlcard,"C_UP",ccnew_up(ii,ind,itypf),0.0d0)
821 mod_fourier(itypf)=mod_fourier(itypf)
822 & +mask_ccnew(ii,ind,itypf)
824 else if (index(controlcard,"D_LOW").gt.0) then
826 call readi(controlcard,"IND",ind,0)
827 call readi(controlcard,"COEFF",ii,0)
828 if (ind.eq.0 .or. ii.eq.0 .or. ii.gt.3) then
830 & "Error specifying D, components undefined",ind,ii
833 mask_ddnew(ii,ind,itypf)=1
834 call reada(controlcard,"D_LOW",ddnew_low(ii,ind,itypf),
836 call reada(controlcard,"D_UP",ddnew_up(ii,ind,itypf),
838 mod_fourier(itypf)=mod_fourier(itypf)
839 & +mask_ddnew(ii,ind,itypf)
841 else if (index(controlcard,"E0_LOW").gt.0) then
843 call readi(controlcard,"COEFF",ii,0)
844 if (ii.eq.0 .or. ii.gt.3) then
846 & "Error specifying E0, components undefined",ii
849 mask_e0new(ii,itypf)=1
850 call reada(controlcard,"E0_LOW",e0new_low(ii,itypf),
852 call reada(controlcard,"E0_UP",e0new_up(ii,itypf),
854 mod_fourier(itypf)=mod_fourier(itypf)
855 & +mask_e0new(ii,itypf)
857 else if (index(controlcard,"E_LOW").gt.0) then
859 call readi(controlcard,"IND1",ind1,0)
860 call readi(controlcard,"IND2",ind2,0)
861 call readi(controlcard,"COEFF",ii,0)
862 if (ind1.eq.0 .or. ind2.eq.0 .or. ii.eq.0 .or. ii.gt.2) then
864 & "Error specifying E, components undefined",ind1,ind2,ii
867 mask_eenew(ii,ind1,ind2,itypf)=1
868 call reada(controlcard,"E_LOW",
869 & eenew_low(ii,ind1,ind2,itypf),0.1d0)
870 call reada(controlcard,"E_UP",
871 & eenew_up(ii,ind1,ind2,itypf),0.0d0)
872 mod_fourier(itypf)=mod_fourier(itypf)
873 & +mask_eenew(ii,ind1,ind2,itypf)
875 call card_concat(controlcard,.true.)
877 call card_concat(controlcard,.true.)
878 write (iout,*) "mod_fouriertor card ",controlcard
879 mod_fouriertor(nloctyp)=index(controlcard,"MOD_FOURIERTOR")
880 write (iout,*) "mod_fouriertor value",mod_fouriertor(nloctyp)
881 write (2,*) "SPLIT_FOURIERTOR",SPLIT_FOURIERTOR,
882 & " tor_mode",tor_mode
883 IF (SPLIT_FOURIERTOR .and. tor_mode.eq.2
884 & .and. mod_fouriertor(nloctyp).gt.0) THEN
889 mask_bnew1tor(ii,j,i)=0
890 mask_bnew2tor(ii,j,i)=0
891 mask_ccnewtor(ii,j,i)=0
892 mask_ddnewtor(ii,j,i)=0
898 mask_eenewtor(ii,j,k,i)=0
903 mask_e0newtor(ii,i)=0
906 call card_concat(controlcard,.true.)
907 do while ( index(controlcard,"END").eq.0 )
908 c write(iout,*) controlcard
909 call reads(controlcard,"TYPF",typf," ")
910 itypf=rescode(1,typf,0,.false.)
911 c write (iout,*) "TYPF ",typf," itypf",itypf
912 if (itypf.eq.0 .or. itypf.gt.ntyp) then
913 write (*,*) "Error specifying fourier parms"
916 itypf=itype2loc(itypf)
917 write (iout,*) "local type",itypf
919 if (index(controlcard,"B1_LOW").gt.0) then
921 call readi(controlcard,"IND",ind,0)
922 call readi(controlcard,"COEFF",ii,0)
923 if (ind.eq.0 .or. ii.eq.0 .or. ii.gt.3) then
925 & "Error specifying B1, components undefined",ind,ii
928 mask_bnew1tor(ii,ind,itypf)=1
929 call reada(controlcard,"B1_LOW",bnew1tor_low(ii,ind,itypf),
931 call reada(controlcard,"B1_UP",bnew1tor_up(ii,ind,itypf),
933 mod_fouriertor(itypf)=mod_fouriertor(itypf)
934 & +mask_bnew1tor(ii,ind,itypf)
936 else if (index(controlcard,"B2_LOW").gt.0) then
938 mask_bnew2tor(ii,ind,itypf)=1
939 call readi(controlcard,"IND",ind,0)
940 call readi(controlcard,"COEFF",ii,0)
941 if (ind.eq.0 .or. ii.eq.0 .or. ii.gt.3) then
943 & "Error specifying B2, components undefined",ind,ii
946 call reada(controlcard,"B2_LOW",bnew2tor_low(ii,ind,itypf),
948 call reada(controlcard,"B2_UP",bnew2tor_up(ii,ind,itypf),
950 mod_fouriertor(itypf)=mod_fouriertor(itypf)
951 & +mask_bnew2tor(ii,ind,itypf)
953 else if (index(controlcard,"C_LOW").gt.0) then
955 call readi(controlcard,"IND",ind,0)
956 call readi(controlcard,"COEFF",ii,0)
957 if (ind.eq.0 .or. ii.eq.0 .or. ii.gt.3) then
959 & "Error specifying C, components undefined",ind,ii
962 mask_ccnewtor(ii,ind,itypf)=1
963 call reada(controlcard,"C_LOW",ccnewtor_low(ii,ind,itypf),
965 call reada(controlcard,"C_UP",ccnewtor_up(ii,ind,itypf),
967 mod_fouriertor(itypf)=mod_fouriertor(itypf)
968 & +mask_ccnewtor(ii,ind,itypf)
970 else if (index(controlcard,"D_LOW").gt.0) then
972 call readi(controlcard,"IND",ind,0)
973 call readi(controlcard,"COEFF",ii,0)
974 if (ind.eq.0 .or. ii.eq.0 .or. ii.gt.3) then
976 & "Error specifying D, components undefined",ind,ii
979 mask_ddnewtor(ii,ind,itypf)=1
980 call reada(controlcard,"D_LOW",ddnewtor_low(ii,ind,itypf),
982 call reada(controlcard,"D_UP",ddnewtor_up(ii,ind,itypf),
984 mod_fouriertor(itypf)=mod_fouriertor(itypf)
985 & +mask_ddnewtor(ii,ind,itypf)
987 else if (index(controlcard,"E0_LOW").gt.0) then
989 call readi(controlcard,"COEFF",ii,0)
990 if (ii.eq.0 .or. ii.gt.3) then
992 & "Error specifying E0, components undefined",ii
995 mask_e0newtor(ii,itypf)=1
996 call reada(controlcard,"E0_LOW",e0newtor_low(ii,itypf),
998 call reada(controlcard,"E0_UP",e0newtor_up(ii,itypf),
1000 mod_fouriertor(itypf)=mod_fouriertor(itypf)
1001 & +mask_e0newtor(ii,itypf)
1003 else if (index(controlcard,"E_LOW").gt.0) then
1005 call readi(controlcard,"IND1",ind1,0)
1006 call readi(controlcard,"IND2",ind2,0)
1007 call readi(controlcard,"COEFF",ii,0)
1008 if (ind1.eq.0 .or. ind2.eq.0 .or. ii.eq.0 .or. ii.gt.2) then
1010 & "Error specifying E, components undefined",ind1,ind2,ii
1013 mask_eenewtor(ii,ind1,ind2,itypf)=1
1014 call reada(controlcard,"E_LOW",
1015 & eenewtor_low(ii,ind1,ind2,itypf),0.1d0)
1016 call reada(controlcard,"E_UP",
1017 & eenewtor_up(ii,ind1,ind2,itypf),0.0d0)
1018 mod_fouriertor(itypf)=mod_fouriertor(itypf)
1019 & +mask_eenewtor(ii,ind1,ind2,itypf)
1021 call card_concat(controlcard,.true.)
1023 call card_concat(controlcard,.true.)
1024 write (2,*) "End read FOURIERTOR ",controlcard
1026 ENDIF ! SPLIT_FOURIERTOR
1029 do itypf=0,nloctyp-1
1030 write (iout,*) "itypf",itypf," mod_fourier",
1031 & mod_fourier(itypf)
1032 mod_fourier(nloctyp)=mod_fourier(nloctyp)
1033 & +mod_fourier(itypf)
1035 write (iout,*) "mod_fourier all",mod_fourier(nloctyp)
1036 do itypf=0,nloctyp-1
1037 write (iout,*) "itypf",itypf," mod_fouriertor",
1038 & mod_fouriertor(itypf)
1039 mod_fouriertor(nloctyp)=mod_fouriertor(nloctyp)
1040 & +mod_fouriertor(itypf)
1042 write (iout,*) "mod_fouriertor all",mod_fouriertor(nloctyp)
1044 if (mod_fourier(nloctyp).gt.0) then
1045 call card_concat(controlcard,.true.)
1046 do while ( index(controlcard,"END").eq.0 )
1047 call readi(controlcard,"TYPF",typf," ")
1048 itypf=rescode(1,typf,0,.false.)
1049 if (itypf.eq.0 .or. itypf.gt.ntyp) then
1050 write (*,*) "Error specifying fourier parms"
1053 itypf=itype2loc(itypf)
1054 call readi(controlcard,"IND",ind,0)
1055 call reada(controlcard,"B_LOW",b_low(ind,itypf),0.1d0)
1056 call reada(controlcard,"B_UP",b_up(ind,itypf),0.0d0)
1057 mask_fourier(ind,itypf)=1
1058 mod_fourier(itypf)=mod_fourier(itypf)
1059 & +mask_fourier(ind,itypf)
1060 mod_fourier(nloctyp)=mod_fourier(nloctyp)
1061 & +mask_fourier(ind,itypf)
1062 call card_concat(controlcard,.true.)
1064 call card_concat(controlcard,.true.)
1066 do itypf=0,nloctyp-1
1067 write (iout,*) "itypf",itypf," mod_fourier",mod_fourier(itypf)
1068 mod_fourier(nloctyp)=mod_fourier(nloctyp)+mod_fourier(itypf)
1070 write (iout,*) "mod_fourier all",mod_fourier(nloctyp)
1079 mod_elec=index(controlcard,"MOD_ELEC").gt.0
1082 call card_concat(controlcard,.true.)
1083 do while ( index(controlcard,"END").eq.0 )
1084 c write (iout,*) "controlcard ",controlcard
1085 call readi(controlcard,"TYPE1",itype1,0)
1086 call readi(controlcard,"TYPE2",itype2,0)
1087 c write (iout,*) "itype1",itype1," itype2",itype2
1088 if (itype1.eq.0 .or. itype1.gt.2 .or. itype2.eq.0
1089 & .or. itype2.gt.2) then
1090 write (iout,*) "Error specifying elec parms"
1093 if (index(controlcard,"EPP").gt.0) then
1095 mask_elec(itype1,itype2,1)=1
1096 mask_elec(itype2,itype1,1)=1
1097 call reada(controlcard,"EPP_LOW",epp_low(itype1,itype2),
1099 epp_low(itype2,itype1)=epp_low(itype1,itype2)
1100 call reada(controlcard,"EPP_UP",epp_up(itype1,itype2),
1102 epp_up(itype2,itype1)=epp_up(itype1,itype2)
1104 if (index(controlcard,"RPP").gt.0) then
1106 mask_elec(itype1,itype2,2)=1
1107 mask_elec(itype2,itype1,2)=1
1108 call reada(controlcard,"RPP_LOW",rpp_low(itype1,itype2),
1110 rpp_low(itype2,itype1)=rpp_low(itype1,itype2)
1111 call reada(controlcard,"RPP_UP",rpp_up(itype1,itype2),
1113 rpp_up(itype2,itype1)=rpp_up(itype1,itype2)
1115 if (index(controlcard,"ELPP6").gt.0) then
1117 mask_elec(itype1,itype2,3)=1
1118 mask_elec(itype2,itype1,3)=1
1119 call reada(controlcard,"ELPP6_LOW",
1120 & elpp6_low(itype1,itype2),0.1d0)
1121 elpp6_low(itype2,itype1)=elpp6_low(itype1,itype2)
1122 call reada(controlcard,"ELPP6_UP",
1123 & elpp6_up(itype1,itype2),0.0d0)
1124 elpp6_up(itype2,itype1)=elpp6_up(itype1,itype2)
1126 if (index(controlcard,"ELPP3").gt.0) then
1128 mask_elec(itype1,itype2,4)=1
1129 mask_elec(itype2,itype1,4)=1
1130 call reada(controlcard,"ELPP3_LOW",
1131 & elpp3_low(itype1,itype2),0.1d0)
1132 elpp3_low(itype2,itype1)=elpp3_low(itype1,itype2)
1133 call reada(controlcard,"ELPP3_UP",
1134 & elpp3_up(itype1,itype2),0.0d0)
1135 elpp3_up(itype2,itype1)=elpp3_up(itype1,itype2)
1137 call card_concat(controlcard,.true.)
1139 call card_concat(controlcard,.true.)
1148 mod_scp=index(controlcard,"MOD_SCP").gt.0
1151 call card_concat(controlcard,.true.)
1152 do while (index(controlcard,"END").eq.0)
1153 if (index(controlcard,"EPSCP").gt.0) then
1155 call readi(controlcard,"ITYPSC",itypsc,0)
1156 call readi(controlcard,"ITYPP",itypp,0)
1157 if (itypsc.gt.20 .or. itypp.eq.0 .or. itypp.gt.2) then
1158 write (iout,*) "Error specifying scp parms"
1161 mask_scp(itypsc,itypp,1)=1
1162 call reada(controlcard,"EPSCP_LOW",
1163 & epscp_low(itypsc,itypp),0.1d0)
1164 call reada(controlcard,"EPSCP_UP",
1165 & epscp_up(itypsc,itypp),0.0d0)
1167 if (index(controlcard,"RSCP").gt.0) then
1169 call readi(controlcard,"ITYPSC",itypsc,0)
1170 call readi(controlcard,"ITYPP",itypp,0)
1171 mask_scp(itypsc,itypp,2)=1
1172 call readi(controlcard,"ITYPSC",itypsc,0)
1173 call readi(controlcard,"ITYPP",itypp,0)
1174 if (itypsc.gt.20 .or. itypp.eq.0 .or. itypp.gt.2) then
1175 write (iout,*) "Error specifying scp parms"
1178 call reada(controlcard,"RSCP_LOW",
1179 & rscp_low(itypsc,itypp),0.1d0)
1180 c write(iout,*)itypsc,itypp,rscp_low(itypsc,itypp)
1181 call reada(controlcard,"RSCP_UP",
1182 & rscp_up(itypsc,itypp),0.0d0)
1183 c write(iout,*)itypsc,itypp,rscp_up(itypsc,itypp)
1185 call card_concat(controlcard,.true.)
1188 write (iout,*) "END ",controlcard
1191 & mod_fourier(nloctyp).gt.0 .or. mod_elec .or. mod_scp
1192 & .or. mod_sccor .or. mod_ang .or. imask(ind_shield).eq.1
1193 if (read_stat.lt.2. .and. mod_other_params) then
1194 write (iout,*)"Error - only weights and sidechain parameters",
1195 & " can be optimized with READ_STAT=",read_stat
1199 init_ene = init_ene .or. read_stat.eq.2
1201 write (iout,*) "End read: MOD_OTHER_PARAMS ",mod_other_params
1202 c write (*,*) "MOD_SIDE ",mod_side," MOD_FOURIER",
1203 c & mod_fourier(nloctyp),
1204 c & " MOD_ELEC ",mod_elec," MOD_SCP ",mod_scp
1205 init_ene=init_ene .or. mod_other_params
1206 c write (iout,*) "init_ene",init_ene
1211 c-----------------------------------------------------------------------------
1212 subroutine print_general_data(*)
1214 include "DIMENSIONS"
1215 include "DIMENSIONS.ZSCOPT"
1218 include "COMMON.MPI"
1219 integer ierror,kolor,klucz
1221 include "COMMON.WEIGHTS"
1222 include "COMMON.NAMES"
1223 include "COMMON.VMCPAR"
1224 include "COMMON.TORSION"
1225 include "COMMON.LOCAL"
1226 include "COMMON.INTERACT"
1227 include "COMMON.ENERGIES"
1228 include "COMMON.MINPAR"
1229 include "COMMON.IOUNITS"
1230 include "COMMON.TIME1"
1231 include "COMMON.PROTFILES"
1232 include "COMMON.CHAIN"
1233 include "COMMON.CLASSES"
1234 include "COMMON.THERMAL"
1235 include "COMMON.FFIELD"
1236 include "COMMON.OPTIM"
1237 include "COMMON.CONTROL"
1238 include "COMMON.SCCOR"
1239 character*800 controlcard
1240 integer i,j,k,l,ii,n_ene_found,itypt,jtypt
1241 integer ind,itype1,itype2,itypf,itypsc,itypp
1242 integer ilen,lenpot,lenpre
1244 character*4 liczba,liczba1
1246 write (iout,*) "Single-point target function evaluation"
1247 else if (mode.eq.2) then
1248 write (iout,*) "Grid search of the parameter space"
1249 else if (mode.eq.3) then
1250 write (iout,*) "Minimization of the target function"
1252 write (iout,*) "Wrong MODE",mode,", should be 1, 2, or 3"
1256 write (iout,*) "RESCALE_MODE",rescale_mode
1257 c mod_other_params=index(controlcard,"OPTIMIZE_OTHER").gt.0
1258 c if (read_stat.eq.0 .and. mod_other_params) then
1259 c write (iout,*) "Error: only optimization of energy-term",
1260 c & " weights can be performed with READ_STAT=",read_stat
1265 write (iout,*) "Parameters of the SUMSL procedure:"
1266 write (iout,'(a,t15,i5)') "MAXMIN",maxmin
1267 write (iout,'(a,t15,i5)') "MAXFUN",maxfun
1268 write (iout,'(a,t15,e15.7)') "TOLF",tolf
1269 write (iout,'(a,t15,e15.7)') "RTOLF",rtolf
1271 if (mod_other_params) then
1273 &"Internal parameters of UNRES energy components will be optimized"
1274 c call card_concat(controlcard,.true.)
1275 c mod_side = (index(controlcard,"MOD_SIDE").gt.0)
1277 write (iout,*) "SC epsilons to be optimized:"
1278 write (iout,*) "Single: eps(i,j)=eps(i,j)+(x(i)+x(j))/2"
1280 write (iout,*) restyp(ityp_ssc(i)),epss_low(i),epss_up(i)
1282 write (iout,*) "Pairs: eps(i,j)=eps(i,j)+x(i,j):"
1284 write (iout,*) restyp(ityp_psc(1,i)),
1285 & restyp(ityp_psc(2,i)),epsp_low(i),epsp_up(i)
1289 write (iout,*)"Torsional weights/coefficients to be optimized"
1290 write(iout,'(a)') "TYP1 TYP2 WEIGHT",
1291 & " LOWER BOUND UPPER BOUND"
1292 do itypt=-nsccortyp,nsccortyp
1293 do jtypt=-nsccortyp,nsccortyp
1295 if (mask_tor(l,itypt,jtypt,1).gt.0) then
1296 write(iout,'(3a4,3f10.5)')l,restyp(itypt),restyp(jtypt),
1297 & weitor(l,itypt,jtypt,1),weitor_low(l,itypt,jtypt,1),
1298 & weitor_up(l,itypt,jtypt,1)
1304 c 7/8/17 AL: Optimization of the bending parameters
1306 write (iout,*) "Boundaries of angle potential coefficients"
1307 write (iout,'(3a)') "Index"," Low bound"," Up bound"
1309 if (mask_ang(i).eq.0) cycle
1310 write (iout,'(2a)') 'Type ',restyp(iloctyp(i))
1311 do j=1,nbend_kcc_TB(i)
1312 write (iout,'(i5,2f15.3)') j,v1bend_low(j,i),v1bend_up(j,i)
1316 if (mod_fourier(nloctyp).gt.0) then
1317 write (iout,*) "Fourier coefficients to be optimized"
1318 do itypf=0,nloctyp-1
1319 if (mod_fourier(itypf).gt.0) then
1320 write (iout,'(3a,i2)') "Type ",restyp(iloctyp(itypf)),
1321 & " number of coeffs",mod_fourier(itypf)
1322 write(iout,'(a)') "NAME COEFF LOWER BOUND UPPER BOUND"
1326 if (mask_bnew1(k,j,itypf).gt.0)
1327 & write (iout,'(2hB1,2i2,f10.5)') k,j,bnew1(k,j,itypf),
1328 & bnew1_low(k,j,itypf),bnew1_up(k,j,itypf)
1333 if (mask_bnew2(k,j,itypf).gt.0)
1334 & write (iout,'(2hB2,2i2,3f11.5)') k,j,bnew2(k,j,itypf),
1335 & bnew2_low(k,j,itypf),bnew2_up(k,j,itypf)
1340 if (mask_ccnew(k,j,itypf).gt.0)
1341 & write (iout,'(2hCC,2i2,3f11.5)') k,j,ccnew(k,j,itypf),
1342 & ccnew_low(k,j,itypf),ccnew_up(k,j,itypf)
1347 if (mask_ddnew(k,j,itypf).gt.0)
1348 & write (iout,'(2hDD,2i2,3f11.5)') k,j,ddnew(k,j,itypf),
1349 & ddnew_low(k,j,itypf),ddnew_up(k,j,itypf)
1353 if (mask_e0new(j,itypf).gt.0)
1354 & write (iout,'(2hE0,i2,3f11.5)') j,e0new(j,itypf),
1355 & e0new_low(j,itypf),e0new_up(j,itypf)
1360 if (mask_eenew(l,k,j,itypf).gt.0)
1361 & write (iout,'(2hEE,3i2,3f11.5)')
1362 & l,k,j,eenew(l,k,j,itypf),eenew_low(l,k,j,itypf),
1363 & eenew_up(l,k,j,itypf)
1369 if (mask_fourier(i,itypf).gt.0) then
1370 write (iout,'(1hB,i2,3f11.5)')
1371 & i,b(i,itypf),b_low(i,itypf),b_up(i,itypf)
1380 write (iout,*) "Electrostatic parameters to be optimized"
1383 if (mask_elec(itype1,itype2,1).gt.0)
1384 & write (iout,'(2i3," EPP",f8.3," EPP_LOW",f8.3,
1386 & itype1,itype2,epp(itype1,itype2),epp_low(itype1,itype2),
1387 & epp_up(itype1,itype2)
1388 if (mask_elec(itype1,itype2,2).gt.0)
1389 & write (iout,'(2i3," RPP",f8.3," RPP_LOW",f8.3,
1391 & itype1,itype2,rpp(itype1,itype2),rpp_low(itype1,itype2),
1392 & rpp_up(itype1,itype2)
1393 if (mask_elec(itype1,itype2,3).gt.0)
1394 & write (iout,'(2i3," ELPP6",f8.3," ELPP6_LOW",f8.3,
1395 & " ELPP6_UP",f8.3)')
1396 & itype1,itype2,elpp6(itype1,itype2),
1397 & elpp6_low(itype1,itype2),elpp6_up(itype1,itype2)
1398 if (mask_elec(itype1,itype2,4).gt.0)
1399 & write (iout,'(2i3," ELPP3",f8.3," ELPP3_LOW",f8.3,
1400 & " ELPP3_UP",f8.3)')
1401 & itype1,itype2,elpp3(itype1,itype2),
1402 & elpp3_low(itype1,itype2),elpp3_up(itype1,itype2)
1409 write (iout,*) "SCP parameters to be optimized:"
1410 if (mask_scp(0,1,1)+mask_scp(0,2,1)+mask_scp(0,1,2)+
1411 & mask_scp(0,2,2) .gt. 0) then
1413 & "SCP parameters are assumed to depend on peptide group type only"
1415 if (mask_scp(0,j,1).gt.0)
1416 & write (iout,'(i3," EPSCP",f8.3," EPSCP_LOW",f8.3,
1417 & " EPSCP_UP",f8.3)') j,eps_scp(1,j),epscp_low(0,j),
1419 if (mask_scp(0,j,2).gt.0)
1420 & write (iout,'(i3," RSCP",f8.3," RSCP_LOW",f8.3,
1421 & " RSCP_UP",f8.3)') j,rscp(1,j),rscp_low(0,j),
1427 if (mask_scp(i,j,1).gt.0)
1428 & write (iout,'(2i3," EPSCP",f8.3," EPSCP_LOW",f8.3,
1429 & " EPSCP_UP",f8.3)') i,j,eps_scp(i,j),epscp_low(i,j),
1431 if (mask_scp(i,j,2).gt.0)
1432 & write (iout,'(2i3," RSCP",f8.3," RSCP_LOW",f8.3,
1433 & " RSCP_UP",f8.3)') i,j,rscp(i,j),rscp_low(i,j),
1440 write (iout,*) "MOD_OTHER_PARAMS ",mod_other_params
1441 write (iout,*) "MOD_SIDE ",mod_side," MOD_FOURIER",
1442 & mod_fourier(nloctyp),
1443 & " MOD_ELEC ",mod_elec," MOD_SCP ",mod_scp," mod_ang",mod_ang
1444 init_ene=init_ene .or. mod_other_params
1445 write (iout,*) "init_ene",init_ene
1450 c-----------------------------------------------------------------------------
1451 subroutine read_protein_data(*)
1453 include "DIMENSIONS"
1454 include "DIMENSIONS.ZSCOPT"
1457 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
1458 include "COMMON.MPI"
1460 include "COMMON.CONTROL"
1461 include "COMMON.CHAIN"
1462 include "COMMON.CLASSES"
1463 include "COMMON.COMPAR"
1464 include "COMMON.ENERGIES"
1465 include "COMMON.IOUNITS"
1466 include "COMMON.PROTFILES"
1467 include "COMMON.PROTNAME"
1468 include "COMMON.VMCPAR"
1469 include "COMMON.OPTIM"
1470 include "COMMON.WEIGHTS"
1471 include "COMMON.NAMES"
1472 include "COMMON.ALLPROT"
1473 include "COMMON.THERMAL"
1474 include "COMMON.TIME1"
1475 include "COMMON.WEIGHTDER"
1476 character*64 nazwa,key
1477 character*16000 controlcard
1478 character*16000 all_protfiles
1479 integer i,j,k,kk,l,iprot,ii,if,ib,ibatch,nn,nn1,iww,maskcheck,
1480 & il,inat,ilevel,weightread,jfrag,jclass,nfragm,iparm
1481 integer nrec,nlines,iscor
1482 double precision energ,wangnorm(maxT),wq(maxT),wrms(maxT),
1483 & wrgy(maxT),wsign(maxT),wangnat(maxT),wqnat(maxT),wrmsnat(maxT),
1484 & wrgynat(maxT),wcorangnorm(maxT),wcorq(maxT),
1485 & wcorrms(maxT),wcorrgy(maxT),wcorsign(maxT),wcorangnat(maxT),
1486 & wcorqnat(maxT),wcorrmsnat(maxT),wcorrgynat(maxT),
1487 & angnormlow(maxT),qlow(maxT),rmslow(maxT),
1488 & rgylow(maxT),signlow(maxT),angnatlow(maxT),
1489 & qnatlow(maxT),rmsnatlow(maxT),rgynatlow(maxT),
1490 & angcorlow(maxT),qcorlow(maxT),rmscorlow(maxT),rgycorlow(maxT),
1491 & signcorlow(maxT),angcornatlow(maxT),qcornatlow(maxT),
1492 & rmscornatlow(maxT),rgycornatlow(maxT),signcornatlow(maxT),
1493 & angnormup(maxT),qup(maxT),rmsup(maxT),rgyup(maxT),signup(maxT),
1494 & angnatup(maxT),qnatup(maxT),rmsnatup(maxT),rgynatup(maxT),
1495 & angcorup(maxT),qcorup(maxT),rmscorup(maxT),rgycorup(maxT),
1497 & angcornatup(maxT),qcornatup(maxT),rmscornatup(maxT),
1498 & rgycornatup(maxT),signcornatup(MaxT)
1501 double precision ebia(maxprot),rjunk
1503 character*64 zeros /
1504 &'0000000000000000000000000000000000000000000000000000000000000000'
1507 c print *,"Processor",me," calls read_protein_data"
1509 C Read seventh record: general data of proteins used for calibration
1510 call card_concat(controlcard,.true.)
1511 c write(2, *)controlcard
1512 call readi(controlcard,"NPROT",nprot,0)
1513 pdbref=index(controlcard,"PDBREF").gt.0
1514 print_refstr=index(controlcard,"PRINT_REFSTR").gt.0
1515 if (nprot.eq.0) then
1516 write(iout,*) "Error: at least one protein must be specified."
1522 read (inp,'(a)') protname(iprot)
1524 write (iout,*) "Reading data of protein",iprot," named ",
1526 call card_concat(controlcard,.true.)
1527 call reada(controlcard,"ENECUT_MIN",enecut_min(iprot),15.0d0)
1528 call reada(controlcard,"ENECUT_MAX",enecut_max(iprot),100.0d0)
1529 call reada(controlcard,"ENECUT",enecut(iprot),enecut_min(iprot))
1530 if (enecut(iprot).lt.enecut_min(iprot))
1531 & enecut(iprot)=enecut_min(iprot)
1532 if (enecut_max(iprot).le.enecut_min(iprot))
1533 & enecut_max(iprot)=2*enecut_min(iprot)
1534 write (iout,'(3(a,f9.1))') "ENECUT",enecut(iprot)," ENECUT_MIN",
1535 & enecut_min(iprot)," ENECUT_MAX",enecut_max(iprot)
1536 call readi(controlcard,"HEFAC",hefac(iprot),50)
1537 call readi(controlcard,"HTFAC",htfac(iprot),50)
1538 call readi(controlcard,"HELOW",hemax_low(iprot),20)
1539 call readi(controlcard,"HTLOW",htmax_low(iprot),20)
1540 write (iout,*) "iprot",iprot,
1541 & " hefac",hefac(iprot)," helow",hemax_low(iprot),
1542 & " htfac",htfac(iprot)," htlow",htmax_low(iprot)
1543 c 7/27/2013 AL Read maximum likelihood data
1544 call card_concat(controlcard,.true.)
1545 call readi(controlcard,"NBETA_L",nbeta(1,iprot),0)
1546 write (iout,*) "NBETA_L",nbeta(1,iprot)
1547 caonly(iprot)=index(controlcard,"CAONLY").gt.0
1548 sconly(iprot)=index(controlcard,"SCONLY").gt.0
1549 rmscomp(iprot)=index(controlcard,"RMSCOMP").gt.0
1550 anglecomp(iprot)=index(controlcard,"ANGLECOMP").gt.0
1551 sidecomp(iprot)=index(controlcard,"SIDECOMP").gt.0
1552 call reada(controlcard,"SIGMA",sigma2(iprot),4.0d0)
1553 call reada(controlcard,"SIGMAANG",sigmaang2(iprot),4.0d0)
1554 call reada(controlcard,"SIGMASIDE",sigmaside2(iprot),4.0d0)
1555 write (iout,*) "RMSCOMP",rmscomp(iprot)," SIGMA",sigma2(iprot),
1556 & " CAONLY ",caonly(iprot)," SCONLY",sconly(iprot)
1557 write (iout,*) "ANGLECOMP",anglecomp(iprot),
1558 & " SIGMAANG",sigmaang2(iprot)
1559 write (iout,*) "SIDECOMP",sidecomp(iprot),
1560 & " SIGMASIDE",sigmaside2(iprot)
1561 do ib=1,nbeta(1,iprot)
1562 read(inp,*)betaT(ib,1,iprot),weilik(ib,iprot),
1564 write (iout,*) i,betaT(ib,1,iprot),weilik(ib,iprot),
1567 c 7/27/2013 AL Read heat-capacity data
1568 call card_concat(controlcard,.true.)
1569 call readi(controlcard,"NBETA_CV",nbeta(2,iprot),0)
1570 call reada(controlcard,"WCV",wcv(iprot),1.0d0)
1571 call reada(controlcard,"BASE",heatbase(iprot),0.0d0)
1572 write (iout,*) "NBETA_CV",nbeta(2,iprot)," WCV",wcv(iprot)
1573 do ib=1,nbeta(2,iprot)
1574 read(inp,*) betaT(ib,2,iprot),target_cv(ib,iprot),
1576 weiCv(ib,iprot)=weiCv(ib,iprot)*wcv(iprot)
1577 write (iout,*) betaT(ib,2,iprot),target_cv(ib,iprot),
1580 write (iout,*) "Experimental heat capacity curve"
1581 do ib=1,nbeta(2,iprot)
1582 write (iout,'(f6.2,2f10.5)') betaT(ib,2,iprot),
1583 & target_cv(ib,iprot),weiCv(ib,iprot)
1585 write (iout,'(a,f10.5)') "Baseline",heatbase(iprot)
1587 c Conformation-dependent averages
1588 call card_concat(controlcard,.true.)
1589 call readi(controlcard,"NATLIKE",natlike(iprot),0)
1590 do i=1,natlike(iprot)
1591 call card_concat(controlcard,.true.)
1592 call reada(controlcard,"WNAT",wnat(i,iprot),1.0d0)
1593 call readi(controlcard,"NUMNAT",numnat(i,iprot),1)
1594 call readi(controlcard,"NATDIM",natdim(i,iprot),1)
1595 do ib=1,nbeta(i+2,iprot)
1596 read(inp,*)betaT(ib,i+2,iprot),(weinat(k,ib,i,iprot),
1597 & nuexp(k,ib,i,iprot),k=1,natdim(i,iprot))
1600 do i=1,natlike(iprot)+2
1601 do j=1,nbeta(i,iprot)
1602 betaT(j,i,iprot)=1.0d0/(Rgas*betaT(j,i,iprot))
1603 write (2,*) "R i",i," j",j," beta",betaT(j,i,iprot)
1608 C Read names of files with the data for protein IPROT
1609 call card_concat(controlcard,.false.)
1611 if (iparm.eq.myparm) then
1612 call split_string(controlcard,protfiles(1,iprot),
1613 & maxfile_prot,nfile_prot(iprot))
1615 write(iout,*)"iprot",iprot," nfile",nfile_prot(iprot)
1617 & (protfiles(i,iprot),i=1,nfile_prot(iprot))
1623 c Read molecular information of the protein
1624 call molread_zs(iprot)
1625 c 3/31/04 AL Read the reference structures of protein iprot
1626 c print *,"Calling read_ref_structure"
1627 call read_ref_structure(iprot,*11)
1629 write (iout,*) "Protein",iprot," left READ_REF_STRUCTURE"
1635 c-------------------------------------------------------------------------------
1636 subroutine read_database(*)
1638 include "DIMENSIONS"
1639 include "DIMENSIONS.ZSCOPT"
1642 integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
1643 include "COMMON.MPI"
1645 include "COMMON.CHAIN"
1646 include "COMMON.INTERACT"
1647 include "COMMON.CLASSES"
1648 include "COMMON.ENERGIES"
1649 include "COMMON.IOUNITS"
1650 include "COMMON.PROTFILES"
1651 include "COMMON.PROTNAME"
1652 include "COMMON.VMCPAR"
1653 include "COMMON.NAMES"
1654 include "COMMON.ALLPROT"
1655 include "COMMON.WEIGHTS"
1656 include "COMMON.WEIGHTDER"
1657 include "COMMON.VAR"
1658 include "COMMON.SBRIDGE"
1659 include "COMMON.GEO"
1660 include "COMMON.COMPAR"
1661 include "COMMON.OPTIM"
1663 character*16000 controlcard
1664 character*16000 all_protfiles
1665 character*4 liczba,liczba1
1666 integer i,j,jj,jjj,jj_old,icount,k,kk,l,iprot,ii,if,ib,ibatch,
1668 integer ixdrf,iret,itmp
1669 integer nrec,nlines,iscor
1670 double precision energ,t_acq,tcpu
1673 double precision rjunk
1674 integer ntot_all(0:maxprot,0:maxprocs-1)
1676 double precision energia(0:max_ene),etot
1677 real*4 csingle(3,maxres2),reini,refree,rmsdev,prec
1678 integer Previous,Next
1679 c print *,"Processor",me," calls read_protein_data"
1681 if (me.eq.master) then
1682 Previous=MPI_PROC_NULL
1686 if (me.eq.nprocs-1) then
1692 c Set the scratchfile names
1693 write (liczba,'(bz,i4.4)') me
1694 write (liczba1,'(bz,i4.4)') myparm
1696 c 1/27/05 AL Change stored coordinates to single precision and don't store
1697 c energy components in the binary databases.
1698 lenrec(iprot)=12*(nres_zs(iprot)+nct_zs(iprot)-nnt_zs(iprot)+1)
1699 & +4*(2*nss_zs(1,iprot)+1)+8*natlike(iprot)*maxdimnat+16
1700 c 4/13/04 AL Add space for similarity measure
1701 lenrec_ene(iprot) = (2*nntyp+3*n_ene+2)*8
1702 & +8*natlike(iprot)*maxdimnat
1705 write (iout,*) "Protein i"," lenrec",lenrec(iprot)
1706 write (iout,*) "lenrec_ene",lenrec_ene(iprot)
1709 bprotfiles(iprot) = scratchdir(:ilen(scratchdir))//
1710 & "/"//protname(iprot)(:ilen(protname(iprot)))//
1711 & "_par"//liczba1//"_"//liczba//".xbin"
1712 benefiles(iprot) = scratchdir(:ilen(scratchdir))//
1713 & "/"//protname(iprot)(:ilen(protname(iprot)))//
1714 & "_par"//liczba1//"_"//liczba//".enebin"
1715 c write (iout,*) "scratchfile ",
1716 c & bprotfiles(iprot)(:ilen(bprotfiles(iprot)))
1725 call restore_molinfo(iprot)
1726 open (ientout,file=bprotfiles(iprot),status="unknown",
1727 & form="unformatted",access="direct",recl=lenrec(iprot))
1728 c Change AL 12/30/2017
1729 if (.not.mod_other_params)
1730 & open (istat,file=benefiles(iprot),status="unknown",
1731 & form="unformatted",access="direct",recl=lenrec_ene(iprot))
1732 c Read conformations from binary DA files (one per batch) and write them to
1733 c a binary DA scratchfile.
1736 write (liczba,'(bz,i4.4)') me
1738 IF (ME.EQ.MASTER) THEN
1739 c Only the master reads the database; it'll send it to the other procs
1745 do if=1,nfile_prot(iprot)
1746 nazwa=protfiles(if,iprot)
1747 & (:ilen(protfiles(if,iprot)))//".cx"
1749 write (iout,*) "Opening file ",nazwa(:ilen(nazwa))
1751 #if (defined(AIX) && !defined(JUBL))
1752 call xdrfopen_(ixdrf,nazwa, "r", iret)
1754 call xdrfopen(ixdrf,nazwa, "r", iret)
1756 if (iret.eq.0) goto 1111
1760 #if (defined(AIX) && !defined(JUBL))
1761 call xdrf3dfcoord_(ixdrf, csingle, itmp, prec, iret)
1762 if (iret.eq.0) goto 102
1763 call xdrfint_(ixdrf, nss, iret)
1764 if (iret.eq.0) goto 102
1766 call xdrfint_(ixdrf, ihpb(j), iret)
1767 if (iret.eq.0) goto 102
1768 call xdrfint_(ixdrf, jhpb(j), iret)
1769 if (iret.eq.0) goto 102
1771 call xdrffloat_(ixdrf,reini,iret)
1772 if (iret.eq.0) goto 102
1773 call xdrffloat_(ixdrf,refree,iret)
1774 if (iret.eq.0) goto 102
1775 call xdrfint_(ixdrf,natlik,iret)
1776 if (iret.eq.0) goto 102
1778 call xdrfint(ixdrf,natliktemp(j),iret)
1779 if (iret.eq.0) goto 102
1780 do k=1,natliktemp(j)
1781 call xdrffloat(ixdrf,nutemp(k,j),iret)
1782 if (iret.eq.0) goto 102
1786 c write (0,*) "me",me," iprot",iprot," i",i
1787 call xdrf3dfcoord(ixdrf, csingle, itmp, prec, iret)
1788 if (iret.eq.0) goto 102
1789 call xdrfint(ixdrf, nss, iret)
1790 if (iret.eq.0) goto 102
1792 call xdrfint(ixdrf, ihpb(k), iret)
1793 if (iret.eq.0) goto 102
1794 call xdrfint(ixdrf, jhpb(k), iret)
1795 if (iret.eq.0) goto 102
1797 call xdrffloat(ixdrf,reini,iret)
1798 if (iret.eq.0) goto 102
1799 call xdrffloat(ixdrf,refree,iret)
1800 if (iret.eq.0) goto 102
1802 call xdrfint(ixdrf,natlik,iret)
1803 if (iret.eq.0) goto 102
1805 call xdrfint(ixdrf,natliktemp(j),iret)
1806 if (iret.eq.0) goto 102
1807 do k=1,natliktemp(j)
1808 call xdrffloat(ixdrf,nutemp(k,j),iret)
1809 if (iret.eq.0) goto 102
1814 call xdrffloat(ixdrf,rmsdev,iret)
1815 if (iret.eq.0) goto 102
1816 c write (2,*) "rmsdev",rmsdev
1817 call xdrfint(ixdrf,iscor,iret)
1818 if (iret.eq.0) goto 102
1819 c write (2,*) "iscor",iscor
1822 eini(jj+1,iprot)=reini
1823 entfac(jj+1,iprot)=refree
1831 c(l,nres+k)=csingle(l,nres+k-nnt+1)
1835 write (iout,'(5hREAD ,i5,2f15.4)')
1836 & jj+1,eini(jj+1,iprot),entfac(jj+1,iprot)
1837 write (iout,*) "natlik",natlik
1839 write (iout,*) "natliketemp",natliktemp(l)
1840 write(iout,*) (nutemp(k,l),k=1,natliktemp(l))
1842 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1843 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1846 call add_new_cconf(jjj,jj,jj_old,icount,iprot,
1850 write (iout,*) "Protein ",protname(iprot),
1851 & i-1," conformations read from DA file ",
1852 & nazwa(:ilen(nazwa))
1853 write (iout,*) jj," conformations read so far"
1854 #if (defined(AIX) && !defined(JUBL))
1855 call xdrfclose_(ixdrf,nazwa,iret)
1857 call xdrfclose(ixdrf,nazwa,iret)
1859 c print *,"file ",nazwa(:ilen(nazwa))," closed"
1863 write (iout,*) "jj_old",jj_old," jj",jj
1865 call write_and_send_cconf(icount,jj_old,jj,iprot,Next)
1866 if (icount.gt.0) call MPI_Send(0,1,MPI_INTEGER,Next,570,
1870 call write_and_send_cconf(icount,jj_old,jj,iprot,Next)
1872 t_acq = tcpu() - t_acq
1873 write (iout,*) "Processor",me," protein",iprot,
1874 & " batch",ibatch," time for conformation read/send",t_acq
1877 c A worker gets the confs from the master and sends them to its neighbor
1879 call receive_and_pass_cconf(icount,jj_old,jj,iprot,
1881 t_acq = tcpu() - t_acq
1882 write (iout,*) "Processor",me," protein",iprot,
1884 & " time for conformation receive/send",t_acq
1888 write (iout,*) "Protein",
1889 & protname(iprot)(:ilen(protname(iprot))),", ",ntot(iprot),
1890 & " conformatons read ",ntot(iprot)," conformations written to ",
1891 & bprotfiles(iprot)(:ilen(bprotfiles(iprot)))
1892 ntot(0) = ntot(0)+ntot(iprot)
1897 write(iout,*)"A total of",ntot(0)," conformations read."
1899 c Check if everyone has the same number of conformations
1900 call MPI_Allgather(ntot(0),maxprot+1,MPI_INTEGER,
1901 & ntot_all(0,0),maxprot+1,MPI_INTEGER,MPI_Comm_World,IERROR)
1906 if (ntot(j).ne.ntot_all(j,i)) then
1907 write (iout,*) "Number of conformations at processor",i,
1908 & " for protein",j," differs from that at processor",me,
1909 & ntot(j),ntot_all(j,i)
1916 c----------- Temporary; reading probs from external file
1917 open (88,file='1LE1_wham_last_2.rms',status='old')
1919 read (88,*) ii,weirms(i,1)
1922 write (iout,*) "i",i," weirms",weirms(i,1)
1925 call MPI_Bcast(weirms(1,1), ntot(1), MPI_Double_Precision,
1926 & Master, MPI_COMM_WORLD, IERROR)
1927 c----------- end temportary stuff
1931 write (iout,*) "Number of conformations read by processors"
1932 write (iout,'(10x,7a10)') (protname(i),i=0,nprot)
1935 write (iout,'(8i10)') i,(ntot_all(j,i),j=0,nprot)
1937 write (iout,*) "Calculation terminated."
1943 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
1946 call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
1950 c------------------------------------------------------------------------------
1951 subroutine add_new_cconf(jjj,jj,jj_old,icount,iprot,Next)
1953 include "DIMENSIONS"
1954 include "DIMENSIONS.ZSCOPT"
1955 include "COMMON.CHAIN"
1956 include "COMMON.INTERACT"
1957 include "COMMON.LOCAL"
1958 include "COMMON.CLASSES"
1959 include "COMMON.ENERGIES"
1960 include "COMMON.IOUNITS"
1961 include "COMMON.PROTFILES"
1962 include "COMMON.PROTNAME"
1963 include "COMMON.VMCPAR"
1964 include "COMMON.WEIGHTS"
1965 include "COMMON.NAMES"
1966 include "COMMON.ALLPROT"
1967 include "COMMON.WEIGHTDER"
1968 include "COMMON.VAR"
1969 include "COMMON.SBRIDGE"
1970 include "COMMON.GEO"
1971 include "COMMON.COMPAR"
1975 integer i,j,jj,jjj,jj_old,icount,k,kk,l,iprot,ii,ib,ibatch,
1976 & nn,nn1,inan,Next,itj
1977 double precision etot,energia(0:max_ene)
1979 c if (entfac(jj+1,iprot).gt.-10.0d0
1980 c & .or. entfac(jj+1,iprot).lt.-150.0d0) then
1981 c write (iout,*) "Entropy factor out of range for conformation",
1982 c & jjj,entfac(jj+1,iprot),", conformation skipped."
1985 call int_from_cart1(.false.)
1987 if (itype(j-1).eq.ntyp1 .or. itype(j).eq.ntyp1) cycle
1988 if (vbld(j).lt.2.0d0 .or. vbld(j).gt.6.5d0) then
1989 write (iout,*) "nnt",nnt," nct",nct
1990 write (iout,*) "Bad CA-CA bond length",j," ",vbld(j)
1991 write (iout,*) "Bad CA-SC bond length",ii," ",vbld(nres+j)
1992 write (iout,*) "The Cartesian geometry is:"
1993 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1994 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1995 write (iout,*) "The internal geometry is:"
1996 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
1997 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
1998 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
1999 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
2000 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
2001 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
2003 & "This conformation WILL NOT be added to the database."
2009 if (itype(j).ne.10 .and. itype(j).ne.ntyp1
2010 & .and. (vbld(nres+j)-dsc(itj)).gt.2.0d0) then
2011 write (iout,*) "nnt",nnt," nct",nct
2012 write (iout,*) "Bad CA-SC bond length",ii," ",vbld(nres+j)
2013 write (iout,*) "The Cartesian geometry is:"
2014 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
2015 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
2016 write (iout,*) "The internal geometry is:"
2017 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
2018 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
2019 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
2020 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
2021 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
2022 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
2024 & "This conformation WILL NOT be added to the database."
2029 if (itype(j).eq.ntyp1 .or. itype(j-1).eq.ntyp1
2030 & .or. itype(j-1).eq.ntyp1) cycle
2031 if (theta(j).le.0.0d0) then
2033 & "Zero theta angle(s) in conformation",ii
2034 write (iout,*) "The Cartesian geometry is:"
2035 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
2036 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
2037 write (iout,*) "The internal geometry is:"
2038 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
2039 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
2040 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
2041 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
2042 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
2043 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
2045 & "This conformation WILL NOT be added to the database."
2048 if (theta(j).gt.179.97*deg2rad) theta(j)=179.97*deg2rad
2050 if (.not. init_ene) then
2053 etot=etot+ww(j)*enetb(icount+1,j,iprot)
2059 if (isnan(etot).ne.0) inan=1
2061 if (isnan(etot)) inan=1
2065 idumm=proc_proc(etot,inan)
2067 call proc_proc(etot,inan)
2074 write (iout,*) "NaNs detected in some of the energy",
2075 & " components for protein",iprot," batch",ibatch,
2076 & " conformation",jjj
2077 write (iout,*) "etot",etot
2078 write (iout,*) "The Cartesian geometry is:"
2079 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
2080 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
2081 write (iout,*) "The internal geometry is:"
2082 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
2083 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
2084 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
2085 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
2086 write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct)
2087 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
2088 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
2089 write (iout,*) "The components of the energy are:"
2092 energia(k)=enetb(jj+1,k,iprot)
2094 call enerprint(energia(0))
2096 & "This conformation WILL NOT be added to the database."
2101 write (iout,'(e15.5,16i5)') entfac(icount+1,iprot),
2102 & iscore(icount+1,0,iprot),ibatch
2103 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
2104 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
2105 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
2106 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
2107 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
2108 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
2109 write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct)
2110 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
2111 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
2112 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
2113 write (iout,'(5e15.5)') (enetb(icount+1,j,iprot),j=1,n_ene)
2114 c write (iout,'(2e15.5)') (eneps(1,j,icount+1,iprot),
2115 c & eneps(2,j,icount+1,iprot),j=1,nntyp)
2117 c write (iout,*) "First nu in readrtms"
2120 list_conf(jj,iprot)=jj
2121 call store_cconf_from_file(jj,icount,iprot)
2122 do j=1,natlike(iprot)
2124 if (k.eq.numnat(j,iprot)) then
2125 do l=1,natdim(j,iprot)
2126 nu(l,k,j,iprot)=nutemp(l,k)
2132 c write (iout,*) "End First nu in readrtms"
2134 if (icount.eq.maxstr_proc) then
2136 write (iout,* ) "jj_old",jj_old," jj",jj
2137 write (iout,*) "Calling write_and_send_cconf"
2140 call write_and_send_cconf(icount,jj_old,jj,iprot,Next)
2143 write (iout,*) "Exited write_and_send_cconf"
2151 c------------------------------------------------------------------------------
2152 subroutine store_cconf_from_file(jj,icount,iprot)
2154 include "DIMENSIONS"
2155 include "DIMENSIONS.ZSCOPT"
2156 include "COMMON.CHAIN"
2157 include "COMMON.SBRIDGE"
2158 include "COMMON.INTERACT"
2159 include "COMMON.IOUNITS"
2160 include "COMMON.CLASSES"
2161 include "COMMON.ALLPROT"
2162 include "COMMON.VAR"
2163 integer i,j,jj,icount,ibatch,iprot
2164 c Store the conformation that has been read in
2167 c_zs(j,i,icount,iprot)=c(j,i)
2170 nss_zs(icount,iprot)=nss
2172 ihpb_zs(i,icount,iprot)=ihpb(i)
2173 jhpb_zs(i,icount,iprot)=jhpb(i)
2177 c------------------------------------------------------------------------------
2178 subroutine write_and_send_cconf(icount,jj_old,jj,iprot,Next)
2180 include "DIMENSIONS"
2181 include "DIMENSIONS.ZSCOPT"
2185 include "COMMON.MPI"
2187 include "COMMON.WEIGHTS"
2188 include "COMMON.CHAIN"
2189 include "COMMON.SBRIDGE"
2190 include "COMMON.INTERACT"
2191 include "COMMON.IOUNITS"
2192 include "COMMON.CLASSES"
2193 include "COMMON.VAR"
2194 include "COMMON.ALLPROT"
2195 include "COMMON.ENERGIES"
2196 include "COMMON.WEIGHTDER"
2197 include "COMMON.OPTIM"
2198 include "COMMON.COMPAR"
2199 integer icount,jj_old,jj,Next,iprot
2200 c Write the structures to a scratch file
2202 c Master sends the portion of conformations that have been read in to the neighbor
2204 write (iout,*) "Processor",me," entered WRITE_AND_SEND_CONF"
2207 call MPI_Send(icount,1,MPI_INTEGER,Next,570,WHAM_COMM,IERROR)
2209 write (iout,*) "Processor",me," Next",next," sent icount=",icount
2210 write (iout,*) "Processor",me," jj_old",jj_old," jj",jj
2213 if (icount.gt.0) then
2214 call MPI_Send(nss_zs(1,iprot),icount,MPI_INTEGER,
2215 & Next,571,WHAM_COMM,IERROR)
2216 call MPI_Send(ihpb_zs(1,1,iprot),icount,MPI_INTEGER,
2217 & Next,572,WHAM_COMM,IERROR)
2218 if (.not. init_ene) call MPI_Send(enetb(jj_old,1,iprot),
2220 & MPI_DOUBLE_PRECISION,Next,576,WHAM_COMM,IERROR)
2221 call MPI_Send(nu(1,1,jj_old,iprot),
2222 & maxdimnat*natlike(iprot)*icount,
2223 & MPI_DOUBLE_PRECISION,
2224 & Next,577,WHAM_COMM,IERROR)
2225 call MPI_Send(eini(jj_old,iprot),icount,
2226 & MPI_DOUBLE_PRECISION,Next,578,WHAM_COMM,IERROR)
2227 call MPI_Send(entfac(jj_old,iprot),icount,MPI_DOUBLE_PRECISION,
2228 & Next,579,WHAM_COMM,IERROR)
2229 call MPI_Send(c_zs(1,1,1,iprot),3*icount*maxres2,
2230 & MPI_REAL,Next,580,WHAM_COMM,IERROR)
2231 if (.not. init_ene) call MPI_Send(eneps(1,1,1,iprot),
2233 & MPI_DOUBLE_PRECISION,Next,582,WHAM_COMM,IERROR)
2237 call dawrite_ccoords(iprot,jj_old,jj,ientout)
2238 c Change AL 20/12/2017
2239 if (.not. mod_other_params)
2240 &call dawrite_ene(iprot,jj_old,jj,istat)
2243 c------------------------------------------------------------------------------
2245 subroutine receive_and_pass_cconf(icount,jj_old,jj,iprot,Previous,
2248 include "DIMENSIONS"
2249 include "DIMENSIONS.ZSCOPT"
2251 integer IERROR,STATUS(MPI_STATUS_SIZE)
2252 include "COMMON.MPI"
2253 include "COMMON.CHAIN"
2254 include "COMMON.SBRIDGE"
2255 include "COMMON.INTERACT"
2256 include "COMMON.IOUNITS"
2257 include "COMMON.CLASSES"
2258 include "COMMON.ALLPROT"
2259 include "COMMON.ENERGIES"
2260 include "COMMON.VAR"
2261 include "COMMON.GEO"
2262 include "COMMON.WEIGHTS"
2263 include "COMMON.WEIGHTDER"
2264 include "COMMON.COMPAR"
2265 include "COMMON.OPTIM"
2266 integer i,j,k,icount,jj_old,jj,iprot,Previous,Next
2269 write (iout,*) "Processor",me," entered RECEIVE_AND_PASS_CONF"
2272 do while (icount.gt.0)
2273 c call MPI_Probe(Previous,570,WHAM_COMM,STATUS,IERROR)
2274 call MPI_Recv(icount,1,MPI_INTEGER,Previous,570,WHAM_COMM,
2277 write (iout,*)"Processor",me," previous",previous," icount",icount
2280 call MPI_Send(icount,1,MPI_INTEGER,Next,570,WHAM_COMM,
2283 write (iout,*) "Processor",me," icount",icount
2286 if (icount.eq.0) return
2287 call MPI_Recv(nss_zs(1,iprot),icount,MPI_INTEGER,
2288 & Previous,571,WHAM_COMM,STATUS,IERROR)
2289 call MPI_Send(nss_zs(1,iprot),icount,MPI_INTEGER,
2290 & Next,571,WHAM_COMM,IERROR)
2291 call MPI_Recv(ihpb_zs(1,1,iprot),icount,MPI_INTEGER,
2292 & Previous,572,WHAM_COMM,STATUS,IERROR)
2293 call MPI_Send(ihpb_zs(1,1,iprot),icount,MPI_INTEGER,
2294 & Next,572,WHAM_COMM,IERROR)
2295 if (.not. init_ene) then
2296 call MPI_Recv(enetb(jj_old,1,iprot),maxstr*n_ene,
2297 & MPI_DOUBLE_PRECISION,Previous,576,WHAM_COMM,STATUS,IERROR)
2298 call MPI_Send(enetb(jj_old,1,iprot),maxstr*n_ene,
2299 & MPI_DOUBLE_PRECISION,Next,576,WHAM_COMM,STATUS,IERROR)
2301 call MPI_Recv(nu(1,1,jj_old,iprot),
2302 & maxdimnat*natlike(iprot)*icount,
2303 & MPI_DOUBLE_PRECISION,
2304 & Previous,577,WHAM_COMM,STATUS,IERROR)
2305 call MPI_Send(nu(1,1,jj_old,iprot),
2306 & maxdimnat*natlike(iprot)*icount,
2307 & MPI_DOUBLE_PRECISION,
2308 & Next,577,WHAM_COMM,IERROR)
2309 call MPI_Recv(eini(jj_old,iprot),icount,
2310 & MPI_DOUBLE_PRECISION,Previous,578,WHAM_COMM,STATUS,IERROR)
2311 call MPI_Send(eini(jj_old,iprot),icount,
2312 & MPI_DOUBLE_PRECISION,Next,578,WHAM_COMM,IERROR)
2313 call MPI_Recv(entfac(jj_old,iprot),icount,MPI_DOUBLE_PRECISION,
2314 & Previous,579,WHAM_COMM,STATUS,IERROR)
2315 call MPI_Send(entfac(jj_old,iprot),icount,MPI_DOUBLE_PRECISION,
2316 & Next,579,WHAM_COMM,IERROR)
2317 call MPI_Recv(c_zs(1,1,1,iprot),3*icount*maxres2,
2318 & MPI_REAL,Previous,580,WHAM_COMM,STATUS,IERROR)
2319 call MPI_Send(c_zs(1,1,1,iprot),3*icount*maxres2,
2320 & MPI_REAL,Next,580,WHAM_COMM,IERROR)
2321 if (.not. init_ene) then
2322 call MPI_Recv(eneps(1,1,1,iprot),2*icount*nntyp,
2323 & MPI_DOUBLE_PRECISION,Previous,582,WHAM_COMM,STATUS,IERROR)
2324 call MPI_Send(eneps(1,1,1,iprot),2*icount*nntyp,
2325 & MPI_DOUBLE_PRECISION,Next,582,WHAM_COMM,IERROR)
2329 list_conf(i,iprot)=i
2331 call dawrite_ccoords(iprot,jj_old,jj,ientout)
2332 c Change AL 12/20/2017
2333 if (.not. mod_other_params)
2334 &call dawrite_ene(iprot,jj_old,jj,istat)
2337 write (iout,*) "Protein",iprot
2338 write (iout,*) "Processor",me," received",icount," conformations"
2340 write (iout,'(8f10.4)') ((c_zs(l,k,i,iprot),l=1,3,k=1,nres)
2341 write (iout,'(8f10.4)')((c_zs(l,k,i+nres,iprot),l=1,3,k=nnt,nct)
2342 write (iout,'(16i5)') nss_zs(i,iprot),(ihpb_zs(k,i,iprot),
2343 & jhpb_zs(k,i,iprot),k=1,nss_zs(i,iprot))
2344 write (iout,'(5e15.5)') (enetb(i,j,iprot),j=1,n_ene)
2345 write (iout,'(2e15.5)') (eneps(1,j,i,iprot),
2346 & eneps(2,j,i,iprot),j=1,nntyp)
2347 write (iout,'(e15.5,16i5)') entfac(i),iscore(i,0,iprot),
2355 c------------------------------------------------------------------------------
2356 subroutine read_thermal
2358 include "DIMENSIONS"
2359 include "DIMENSIONS.ZSCOPT"
2360 include "COMMON.CHAIN"
2361 include "COMMON.SBRIDGE"
2362 include "COMMON.INTERACT"
2363 include "COMMON.IOUNITS"
2364 include "COMMON.CLASSES"
2365 include "COMMON.VAR"
2366 include "COMMON.THERMAL"
2367 character*800 controlcard
2368 call card_concat(controlcard,.true.)
2369 call readi(controlcard,"NGRIDT",NGridT,200)
2370 call reada(controlcard,"DELTAT",deltaT,5.0d0)
2371 call reada(controlcard,"T0",GridT0,2.0d2)
2372 write (iout,*) "Parameters of thermal curves"
2373 write (iout,*) "NGRIDT",NGridT," DELTAT",deltaT," T0",GridT0
2376 c------------------------------------------------------------------------------
2377 subroutine daread_ene(iprot,istart_conf,iend_conf)
2379 include "DIMENSIONS"
2380 include "DIMENSIONS.ZSCOPT"
2383 include "COMMON.MPI"
2385 include "COMMON.CHAIN"
2386 include "COMMON.CLASSES"
2387 include "COMMON.ENERGIES"
2388 include "COMMON.IOUNITS"
2389 include "COMMON.PROTFILES"
2390 include "COMMON.ALLPROT"
2391 include "COMMON.WEIGHTDER"
2392 include "COMMON.COMPAR"
2393 include "COMMON.VMCPAR"
2394 integer iprot,istart_conf,iend_conf
2395 integer i,ii,iii,j,k
2397 write (iout,*) "Calling DAREAD_ENE"
2399 c write (iout,*) "Reading: n_ene",n_ene
2401 c write (iout,*) "DAREAD_ENE",istart_conf,iend_conf
2403 c Read conformations off a DA scratchfile.
2405 do ii=istart_conf,iend_conf
2406 iii=list_conf(ii,iprot)
2407 i = ii - istart_conf + 1
2408 read(ientin,rec=iii) (enetb(i,j,iprot),j=1,n_ene),
2409 & (enetb_orig(i,j,iprot),j=1,n_ene),
2410 & (enetb_oorig(i,j,iprot),j=1,n_ene),
2411 & eini(ii,iprot),entfac(ii,iprot),
2412 & (eneps(1,j,i,iprot),eneps(2,j,i,iprot),j=1,nntyp),
2413 & ((nu(k,j,i,iprot),k=1,maxdimnat),j=1,natlike(iprot))
2415 write (iout,'(3i5,3e15.4,i5,i10)') ii,iii,i,eini(ii,iprot),
2417 write (iout,'(20(1pe12.4)') (enetb(i,j,iprot),j=1,n_ene)
2418 write (iout,'(18(1pe12.4))')
2419 & ((nu(k,j,i,iprot)k=1,maxdimnat),j=1,natlike(iprot))
2426 c------------------------------------------------------------------------------
2427 subroutine dawrite_ene(iprot,istart_conf,iend_conf,unit_out)
2429 include "DIMENSIONS"
2430 include "DIMENSIONS.ZSCOPT"
2433 include "COMMON.MPI"
2435 include "COMMON.CHAIN"
2436 include "COMMON.CLASSES"
2437 include "COMMON.ENERGIES"
2438 include "COMMON.IOUNITS"
2439 include "COMMON.PROTFILES"
2440 include "COMMON.ALLPROT"
2441 include "COMMON.WEIGHTDER"
2442 include "COMMON.VMCPAR"
2443 include "COMMON.COMPAR"
2444 integer iprot,istart_conf,iend_conf,unit_out
2445 integer i,ii,iii,j,k
2446 c write (iout,*) "Writing: n_ene",n_ene
2448 c write (iout,*) "DAWRITE_ENE",istart_conf,iend_conf
2450 c Write conformations to a DA scratchfile.
2452 do ii=istart_conf,iend_conf
2453 iii=list_conf(ii,iprot)
2454 i = ii - istart_conf + 1
2455 write(unit_out,rec=iii) (enetb(i,j,iprot),j=1,n_ene),
2456 & (enetb_orig(i,j,iprot),j=1,n_ene),
2457 & (enetb_oorig(i,j,iprot),j=1,n_ene),
2458 & eini(ii,iprot),entfac(ii,iprot),
2459 & (eneps(1,j,i,iprot),eneps(2,j,i,iprot),j=1,nntyp),
2460 & ((nu(k,j,i,iprot),k=1,maxdimnat),j=1,natlike(iprot))
2462 write (iout,'(3i5,3e15.4,i5,i10)') ii,iii,i,eini(ii,iprot),
2464 write (iout,'(20(1pe12.4)') (enetb(i,j,iprot),j=1,n_ene)
2465 write (iout,'(18(1pe12.4))')
2466 & ((nu(kj,i,iprot),k=1,maxdimnat),j=1,natlike(iprot))
2473 c------------------------------------------------------------------------------
2474 subroutine daread_ccoords(iprot,istart_conf,iend_conf)
2476 include "DIMENSIONS"
2477 include "DIMENSIONS.ZSCOPT"
2480 include "COMMON.MPI"
2482 include "COMMON.CHAIN"
2483 include "COMMON.CLASSES"
2484 include "COMMON.ENERGIES"
2485 include "COMMON.IOUNITS"
2486 include "COMMON.PROTFILES"
2487 include "COMMON.ALLPROT"
2488 include "COMMON.INTERACT"
2489 include "COMMON.VAR"
2490 include "COMMON.SBRIDGE"
2491 include "COMMON.GEO"
2492 include "COMMON.COMPAR"
2493 include "COMMON.VMCPAR"
2494 include "COMMON.WEIGHTDER"
2495 integer iprot,istart_conf,iend_conf
2496 integer i,j,k,ij,ii,iii
2498 real*4 csingle(3,maxres2)
2499 character*16 form,acc
2502 c Read conformations off a DA scratchfile.
2505 write (iout,*) "DAREAD_COORDS"
2506 write (iout,*) "istart_conf",istart_conf," iend_conf",iend_conf
2507 write (iout,*) "lenrec",lenrec(iprot)
2508 inquire(unit=ientin,name=nam,recl=len,form=form,access=acc)
2509 write (iout,*) "len=",len," form=",form," acc=",acc
2510 write (iout,*) "nam=",nam
2513 do ii=istart_conf,iend_conf
2514 iii=list_conf(ii,iprot)
2515 ij = ii - istart_conf + 1
2517 write (iout,*) "Reading binary file, record",iii," ii",ii
2520 read(icbase,rec=iii) ((csingle(j,i),j=1,3),i=1,nres),
2521 & ((csingle(j,i),j=1,3),i=nnt+nres,nct+nres),
2522 & nss,(ihpb(i),jhpb(i),i=1,nss),eini(ii,iprot),
2524 & ((nu(k,i,ij,iprot),k=1,maxdimnat),i=1,natlike(iprot))
2531 write (iout,*) "iprot",iprot," ii",ii
2532 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
2533 write (iout,'(8f10.4)') ((c(j,i),j=1,3),i=nnt+nres,nct+nres)
2534 write (iout,'(2e15.5)') eini(ii,iprot),entfac(ii,iprot)
2535 write (iout,'(16i5)') nss,(ihpb(i),jhpb(i),i=1,nss)
2538 call store_ccoords(iprot,ii-istart_conf+1)
2542 c------------------------------------------------------------------------------
2543 subroutine dawrite_ccoords(iprot,istart_conf,iend_conf,unit_out)
2545 include "DIMENSIONS"
2546 include "DIMENSIONS.ZSCOPT"
2549 include "COMMON.MPI"
2551 include "COMMON.CHAIN"
2552 include "COMMON.INTERACT"
2553 include "COMMON.CLASSES"
2554 include "COMMON.ENERGIES"
2555 include "COMMON.IOUNITS"
2556 include "COMMON.PROTFILES"
2557 include "COMMON.ALLPROT"
2558 include "COMMON.VAR"
2559 include "COMMON.SBRIDGE"
2560 include "COMMON.GEO"
2561 include "COMMON.COMPAR"
2562 include "COMMON.WEIGHTDER"
2563 include "COMMON.VMCPAR"
2564 real*4 csingle(3,maxres2)
2565 integer iprot,istart_conf,iend_conf
2566 integer i,j,k,ii,ij,iii,unit_out
2568 character*16 form,acc
2571 c Write conformations to a DA scratchfile.
2574 write (iout,*) "DAWRITE_COORDS"
2575 write (iout,*) "istart_conf",istart_conf," iend_conf",iend_conf
2576 write (iout,*) "lenrec",lenrec(iprot)
2577 inquire(unit=unit_out,name=nam,recl=len,form=form,access=acc)
2578 write (iout,*) "len=",len," form=",form," acc=",acc
2579 write (iout,*) "nam=",nam
2582 do ii=istart_conf,iend_conf
2583 iii=list_conf(ii,iprot)
2584 ij = ii - istart_conf + 1
2585 call restore_ccoords(iprot,ii-istart_conf+1)
2587 write (iout,*) "Writing binary file, record",iii," ii",ii
2595 write(unit_out,rec=iii) ((csingle(j,i),j=1,3),i=1,nres),
2596 & ((csingle(j,i),j=1,3),i=nnt+nres,nct+nres),
2597 & nss,(ihpb(i),jhpb(i),i=1,nss),eini(ii,iprot),
2599 & ((nu(k,i,ij,iprot),k=1,maxdimnat),i=1,natlike(iprot))
2601 write (iout,*) "kbatch",kbatch(ii,iprot)," iscore",
2602 & iscore(ii,0,iprot)
2603 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
2604 write (iout,'(8f10.4)') ((c(j,i),j=1,3),i=nnt+nres,nct+nres)
2605 write (iout,'(2e15.5)') eini(ii,iprot),entfac(ii,iprot)
2606 write (iout,'(16i5)') nss,(ihpb(i),jhpb(i),i=1,nss)
2612 c------------------------------------------------------------------------------
2613 subroutine store_ccoords(iprot,ii)
2615 include "DIMENSIONS"
2616 include "DIMENSIONS.ZSCOPT"
2617 include "COMMON.VAR"
2618 include "COMMON.CHAIN"
2619 include "COMMON.ALLPROT"
2620 include "COMMON.IOUNITS"
2621 include "COMMON.GEO"
2622 include "COMMON.SBRIDGE"
2623 double precision thetnorm
2624 integer iprot,i,j,ii
2625 do i=1,nres_zs(iprot)
2627 c_zs(j,i,ii,iprot)=c(j,i)
2630 do i=nnt_zs(iprot),nct_zs(iprot)
2632 c_zs(j,i+nres,ii,iprot)=c(j,i+nres)
2635 c 5/7/02 AL: store sbridge info
2636 nss_zs(ii,iprot)=nss
2638 ihpb_zs(i,ii,iprot)=ihpb(i)
2639 jhpb_zs(i,ii,iprot)=jhpb(i)
2643 c------------------------------------------------------------------------------
2644 subroutine restore_ccoords(iprot,ii)
2646 include "DIMENSIONS"
2647 include "DIMENSIONS.ZSCOPT"
2648 include "COMMON.INTERACT"
2649 include "COMMON.VAR"
2650 include "COMMON.ALLPROT"
2651 include "COMMON.SBRIDGE"
2652 include "COMMON.CHAIN"
2653 include "COMMON.IOUNITS"
2654 integer iprot,i,j,ii
2655 do i=1,nres_zs(iprot)
2657 c(j,i)=c_zs(j,i,ii,iprot)
2660 do i=nnt_zs(iprot),nct_zs(iprot)
2662 c(j,i+nres)=c_zs(j,i+nres,ii,iprot)
2665 c 5/7/02 AL: restore sbridge info
2666 nss=nss_zs(ii,iprot)
2668 ihpb(i)=ihpb_zs(i,ii,iprot)+nres
2669 jhpb(i)=jhpb_zs(i,ii,iprot)+nres
2674 write (iout,*) "restore_ccoords",ii
2675 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
2676 write (iout,'(8f10.4)') ((c(j,i),j=1,3),i=nnt+nres,nct+nres)
2677 write (iout,'(16i5)') nss,(ihpb(i),jhpb(i),i=1,nss)
2682 c------------------------------------------------------------------------------
2683 subroutine card_concat(card,to_upper)
2685 include 'DIMENSIONS.ZSCOPT'
2686 include "COMMON.IOUNITS"
2688 character*80 karta,ucase
2692 read (inp,'(a)') karta
2693 if (to_upper) karta=ucase(karta)
2695 do while (karta(80:80).eq.'&')
2696 card=card(:ilen(card)+1)//karta(:79)
2697 read (inp,'(a)') karta
2698 if (to_upper) karta=ucase(karta)
2700 card=card(:ilen(card)+1)//karta
2703 c------------------------------------------------------------------------------
2704 subroutine readi(rekord,lancuch,wartosc,default)
2706 character*(*) rekord,lancuch
2707 integer wartosc,default
2710 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
2711 if (iread.eq.0) then
2715 iread=iread+ilen(lancuch)+1
2716 read (rekord(iread:),*) wartosc
2719 c----------------------------------------------------------------------------
2720 subroutine reada(rekord,lancuch,wartosc,default)
2722 character*(*) rekord,lancuch
2724 double precision wartosc,default
2727 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
2728 if (iread.eq.0) then
2732 iread=iread+ilen(lancuch)+1
2733 read (rekord(iread:),*) wartosc
2736 c----------------------------------------------------------------------------
2737 subroutine multreadi(rekord,lancuch,tablica,dim,default)
2740 integer tablica(dim),default
2741 character*(*) rekord,lancuch
2748 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
2749 if (iread.eq.0) return
2750 iread=iread+ilen(lancuch)+1
2751 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
2754 c----------------------------------------------------------------------------
2755 subroutine multreada(rekord,lancuch,tablica,dim,default)
2758 double precision tablica(dim),default
2759 character*(*) rekord,lancuch
2766 iread=index(rekord,lancuch(:ilen(lancuch))//"=")
2767 if (iread.eq.0) return
2768 iread=iread+ilen(lancuch)+1
2769 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim)
2772 c----------------------------------------------------------------------------
2773 subroutine reads(rekord,lancuch,wartosc,default)
2775 character*(*) rekord,lancuch,wartosc,default
2777 integer ilen,lenlan,lenrec,iread,ireade
2781 lenlan=ilen(lancuch)
2783 iread=index(rekord,lancuch(:lenlan)//"=")
2784 c print *,"rekord",rekord," lancuch",lancuch
2785 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
2786 if (iread.eq.0) then
2790 iread=iread+lenlan+1
2791 c print *,"iread",iread
2792 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
2793 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
2795 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
2797 c print *,"iread",iread
2798 if (iread.gt.lenrec) then
2803 c print *,"ireade",ireade
2804 do while (ireade.lt.lenrec .and.
2805 & .not.iblnk(rekord(ireade:ireade)))
2808 wartosc=rekord(iread:ireade)
2811 c----------------------------------------------------------------------------
2812 subroutine multreads(rekord,lancuch,tablica,dim,default)
2815 character*(*) rekord,lancuch,tablica(dim),default
2817 integer ilen,lenlan,lenrec,iread,ireade
2824 lenlan=ilen(lancuch)
2826 iread=index(rekord,lancuch(:lenlan)//"=")
2827 c print *,"rekord",rekord," lancuch",lancuch
2828 c print *,"iread",iread," lenlan",lenlan," lenrec",lenrec
2829 if (iread.eq.0) return
2830 iread=iread+lenlan+1
2832 c print *,"iread",iread
2833 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
2834 do while (iread.le.lenrec .and. iblnk(rekord(iread:iread)))
2836 c print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread))
2838 c print *,"iread",iread
2839 if (iread.gt.lenrec) return
2841 c print *,"ireade",ireade
2842 do while (ireade.lt.lenrec .and.
2843 & .not.iblnk(rekord(ireade:ireade)))
2846 tablica(i)=rekord(iread:ireade)
2850 c----------------------------------------------------------------------------
2851 subroutine split_string(rekord,tablica,dim,nsub)
2853 integer dim,nsub,i,ii,ll,kk
2854 character*(*) tablica(dim)
2855 character*(*) rekord
2865 C Find the start of term name
2867 do while (ii.le.ll .and. rekord(ii:ii).eq." ")
2870 C Parse the name into TABLICA(i) until blank found
2871 do while (ii.le.ll .and. rekord(ii:ii).ne." ")
2873 tablica(i)(kk:kk)=rekord(ii:ii)
2876 if (kk.gt.0) nsub=nsub+1
2877 if (ii.gt.ll) return
2881 c-------------------------------------------------------------------------------
2882 integer function iroof(n,m)
2884 if (ii*m .lt. n) ii=ii+1