2 !-----------------------------------------------------------------------------
7 use energy_data, only:nnt,nct,nss,ihpb,jhpb,iset
8 use geometry_data, only:nres,c
11 ! include "COMMON.MPI"
15 !-----------------------------------------------------------------------------
18 !-----------------------------------------------------------------------------
20 !-----------------------------------------------------------------------------
22 !-------------------------------------------------------------------------------
23 subroutine opentmp(islice,iunit,bprotfile_temp)
25 ! include "DIMENSIONS"
26 ! include "DIMENSIONS.ZSCOPT"
27 ! include "DIMENSIONS.FREE"
28 ! use MPI_data, only:me
31 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
32 ! include "COMMON.MPI"
34 ! include "COMMON.IOUNITS"
35 ! include "COMMON.PROTFILES"
36 ! include "COMMON.PROT"
37 ! include "COMMON.FREE"
38 character(len=64) :: bprotfile_temp
39 character(len=3) :: liczba,liczba2
40 character(len=2) :: liczba1
41 integer :: iunit,islice
45 ! integer :: lenrec,lenrec2
48 ! lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
50 write (liczba1,'(bz,i2.2)') islice
52 write (liczba,'(bz,i3.3)') me
54 ! write (iout,*) "separate_parset ",separate_parset,
56 if (separate_parset) then
57 write (liczba2,'(bz,i3.3)') myparm
58 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
59 prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
60 open (iunit,file=bprotfile_temp,status="unknown",&
61 form="unformatted",access="direct",recl=lenrec)
63 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
64 prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
65 open (iunit,file=bprotfile_temp,status="unknown",&
66 form="unformatted",access="direct",recl=lenrec)
69 bprotfile_temp = scratchdir(:ilen(scratchdir))// &
70 "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
71 open (iunit,file=bprotfile_temp,status="unknown",&
72 form="unformatted",access="direct",recl=lenrec)
74 ! write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
78 end subroutine opentmp
79 !-------------------------------------------------------------------------------
80 subroutine read_database(*)
82 ! use energy_data, only:nct,nnt,nss
84 ! include "DIMENSIONS"
85 ! include "DIMENSIONS.ZSCOPT"
86 ! include "DIMENSIONS.FREE"
87 use MPI_data, only:me,nprocs
90 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
91 ! include "COMMON.MPI"
93 ! include "COMMON.CHAIN"
94 ! include "COMMON.IOUNITS"
95 ! include "COMMON.PROTFILES"
96 ! include "COMMON.NAMES"
97 ! include "COMMON.VAR"
98 ! include "COMMON.GEO"
99 ! include "COMMON.ENEPS"
100 ! include "COMMON.PROT"
101 ! include "COMMON.INTERACT"
102 ! include "COMMON.FREE"
103 ! include "COMMON.SBRIDGE"
104 ! include "COMMON.OBCINKA"
105 real(kind=4) :: csingle(3,nres*2) !(3,maxres2)
106 character(len=64) :: nazwa,bprotfile_temp
107 character(len=3) :: liczba
108 character(len=2) :: liczba1
109 integer :: i,j,ii,jj(nslice),k,kk(nslice),l,&
110 ll(nslice),mm(nslice),if
111 integer :: nrec,nlines,iscor,iunit,islice
112 real(kind=8) :: energ
114 ! external ilen,iroof
115 real(kind=8) :: rmsdev,energia(0:n_ene),efree,eini,temp
116 !el real(kind=8) :: rmsdev,energia(0:max_eneW),efree,eini,temp
117 real(kind=8) :: prop(nQ) !(maxQ)
118 integer :: ntot_all(nslice,0:nprocs-1)!(maxslice,0:maxprocs-1)
119 integer :: iparm,ib,iib,ir,nprop,nthr,npars
120 real(kind=8) :: etot,time
121 integer :: ixdrf,iret
122 logical :: lerr,linit
124 lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
125 lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
127 write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,&
137 write (iout,*) "nparmset",nparmset
145 if (replica(iparm)) then
152 do iR=1,nRR(ib,iparm)
154 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
160 IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
161 ! Read conformations from binary DA files (one per batch) and write them to
162 ! a binary DA scratchfile.
163 write (liczba,'(bz,i3.3)') me
164 do if=1,nfile_bin(iR,ib,iparm)
165 nazwa=protfiles(if,1,iR,ib,iparm) &
166 (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
167 open (ientin,file=nazwa,status="old",form="unformatted",&
168 access="direct",recl=lenrec2,err=1111)
171 call opentmp(islice,ientout,bprotfile_temp)
172 call bxread(nazwa,islice,ii,jj(islice),kk(islice),ll(islice),&
173 mm(islice),iR,ib,iparm)
180 IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
181 ! Read conformations from multiple ASCII int files and write them to a binary
183 do if=1,nfile_asc(iR,ib,iparm)
184 nazwa=protfiles(if,2,iR,ib,iparm) &
185 (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
186 open(unit=ientin,file=nazwa,status='old',err=1111)
187 write(iout,*) "reading ",nazwa(:ilen(nazwa))
189 call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
192 IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
193 ! Read conformations from cx files and write them to a binary
195 do if=1,nfile_cx(iR,ib,iparm)
196 nazwa=protfiles(if,2,iR,ib,iparm) &
197 (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
198 write(iout,*) "reading ",nazwa(:ilen(nazwa))
200 print *,"Calling cxread"
201 call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,&
203 write(iout,*)"after call cxread"
205 write (iout,*) "exit cxread"
209 write(iout,*)"*********************in read database"
213 stot(islice)=stot(islice)+jj(islice)
218 write (iout,*) "IPARM",iparm
221 if (nslice.eq.1) then
223 write (liczba,'(bz,i3.3)') me
224 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
225 prefix(:ilen(prefix))//liczba//".xbin.tmp"
227 bprotfile_temp = scratchdir(:ilen(scratchdir))// &
228 "/"//prefix(:ilen(prefix))//".xbin.tmp"
230 write(iout,*) mm(1)," conformations read",ll(1),&
231 " conformations written to ",&
232 bprotfile_temp(:ilen(bprotfile_temp))
235 write (liczba1,'(bz,i2.2)') islice
237 write (liczba,'(bz,i3.3)') me
238 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
239 prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
241 bprotfile_temp = scratchdir(:ilen(scratchdir))// &
242 "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
244 write(iout,*) mm(islice)," conformations read",ll(islice),&
245 " conformations written to ",&
246 bprotfile_temp(:ilen(bprotfile_temp))
251 ! Check if everyone has the same number of conformations
252 call MPI_Allgather(stot(1),nslice,MPI_INTEGER,&
253 ntot_all(1,0),nslice,MPI_INTEGER,MPI_Comm_World,IERROR)
258 if (stot(islice).ne.ntot_all(islice,i)) then
259 write (iout,*) "Number of conformations at processor",i,&
260 " differs from that at processor",me,&
261 stot(islice),ntot_all(islice,i)," slice",islice
269 write (iout,*) "Numbers of conformations read by processors"
272 write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
274 write (iout,*) "Calculation terminated."
279 ntot(islice)=stot(islice)
281 write(iout,*) "end of read database"
284 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
287 end subroutine read_database
288 !--------------------------------------------------------------------------------
289 integer function iroof(n,m)
292 if (ii*m .lt. n) ii=ii+1
296 !--------------------------------------------------------------------------------
298 !--------------------------------------------------------------------------------
299 subroutine bxread(nazwa,islice,ii,jj,kk,ll,mm,iR,ib,iparm)
301 ! include "DIMENSIONS"
302 ! include "DIMENSIONS.ZSCOPT"
303 ! include "DIMENSIONS.FREE"
304 ! use energy_data, only:nnt,nct,nss,ihpb,jhpbi
305 use MPI_data, only:nprocs
308 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
309 ! include "COMMON.MPI"
311 ! include "COMMON.CHAIN"
312 ! include "COMMON.IOUNITS"
313 ! include "COMMON.PROTFILES"
314 ! include "COMMON.NAMES"
315 ! include "COMMON.VAR"
316 ! include "COMMON.GEO"
317 ! include "COMMON.ENEPS"
318 ! include "COMMON.PROT"
319 ! include "COMMON.INTERACT"
320 ! include "COMMON.FREE"
321 ! include "COMMON.SBRIDGE"
322 real(kind=4) :: csingle(3,nres*2) !(3,maxres2)
323 character(len=64) :: nazwa,bprotfile_temp
324 character(len=3) :: liczba
325 integer :: i,is,ie,j,ii,jj,k,kk,l,ll,mm,if
326 integer :: nrec,nlines,iscor,islice
327 real(kind=8) :: energ
329 ! external ilen,iroof
330 real(kind=8) :: rmsdev,energia(0:n_ene),efree,eini,temp
331 !el real(kind=8) :: rmsdev,energia(0:max_eneW),efree,eini,temp
332 real(kind=8) :: prop(nQ) !(maxQ)
333 integer :: ntot_all(0:nprocs-1)!(0:maxprocs-1)
334 integer :: iparm,ib,iib,ir,nprop,nthr,nrec_slice
335 real(kind=8) :: etot,time
337 nrec_slice=(rec_end(iR,ib,iparm)-rec_start(iR,ib,iparm)+1)/nslice
338 is=rec_start(iR,ib,iparm)+(islice-1)*nrec_slice
339 ie=rec_start(iR,ib,iparm)+islice*nrec_slice-1
340 write (iout,*) "bxread: islice",islice," nslice",nslice,&
341 " nrec_slice",nrec_slice
342 write (iout,*) "is",is," ie",ie,"rec_start",&
343 rec_start(iR,ib,iparm)," rec_end",rec_end(iR,ib,iparm)
345 read(ientin,rec=i+1,err=101) &
346 ((csingle(l,k),l=1,3),k=1,nres),&
347 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
348 nss,(ihpb(k),jhpb(k),k=1,nss),&
349 eini,efree,rmsdev,(prop(j),j=1,nQ),iscor
352 if (mod(kk,isampl(iparm)).eq.0) then
354 write(ientout,rec=jj) &
355 ((csingle(l,k),l=1,3),k=1,nres),&
356 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
357 nss,(ihpb(k),jhpb(k),k=1,nss),&
358 eini,efree,rmsdev,(prop(j),j=1,nQ),iR,ib,iparm
365 call int_from_cart1(.false.)
366 write (iout,*) "Writing conformation, record",jj
367 write (iout,*) "Cartesian coordinates"
368 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
369 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
370 write (iout,*) "Internal coordinates"
371 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
372 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
373 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
374 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
375 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
376 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
377 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
378 write (iout,'(f10.5,i5)') rmsdev,iscor
384 write (iout,*) ii," conformations read from DA file ",&
386 write (iout,*) kk," conformations read so far, slice",islice
387 write (iout,*) jj," conformations stored so far, slice",islice
390 end subroutine bxread
391 !--------------------------------------------------------------------------------
393 !--------------------------------------------------------------------------------
394 subroutine cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,*)
398 use geometry, only:int_from_cart1
399 use geometry_data, only:vbld,rad2deg,theta,phi,alph,omeg
402 ! implicit real*8 (a-h,o-z)
403 ! include 'DIMENSIONS'
404 ! include 'DIMENSIONS.ZSCOPT'
405 ! include 'DIMENSIONS.FREE'
406 integer,parameter :: MaxTraj=2050
407 ! include 'COMMON.CHAIN'
408 ! include 'COMMON.INTERACT'
409 ! include 'COMMON.NAMES'
410 ! include 'COMMON.IOUNITS'
411 ! include 'COMMON.HEADER'
412 ! include 'COMMON.SBRIDGE'
413 ! include 'COMMON.PROTFILES'
414 ! include 'COMMON.OBCINKA'
415 ! include 'COMMON.FREE'
416 ! include 'COMMON.VAR'
417 ! include 'COMMON.GEO'
418 ! include 'COMMON.PROT'
419 character(len=64) :: nazwa,bprotfile_temp
420 real(kind=4) :: rtime,rpotE,ruconst,rt_bath,rprop(nQ) !(2000) !(maxQ)
422 integer :: iret,itmp,itraj,ntraj
423 real(kind=4) :: xoord(3,2*nres+2),prec
424 integer :: nstep(0:MaxTraj-1)
427 integer :: ii,jj(nslice),kk(nslice),ll(nslice),mm(nslice) !(maxslice)
428 integer :: is(nSlice),ie(nSlice),nrec_slice
429 real(kind=8) :: ts(nSlice),te(nSlice),time_slice
430 integer :: iR,ib,iparm,i,j,it,islice,nprop_prev
431 integer :: k,l,iib,islice1,nprop
432 real(kind=8) :: efree,rmsdev
435 ! logical :: conf_check
444 call set_slices(is,ie,ts,te,iR,ib,iparm)
455 #if (defined(AIX) && !defined(JUBL))
456 call xdrfopen_(ixdrf,nazwa, "r", iret)
458 call xdrfopen(ixdrf,nazwa, "r", iret)
460 if (iret.eq.0) return 1
463 call opentmp(islice1,ientout,bprotfile_temp)
467 #if (defined(AIX) && !defined(JUBL))
468 call xdrffloat_(ixdrf, rtime, iret)
469 print *,"rtime",rtime," iret",iret !d
470 call xdrffloat_(ixdrf, rpotE, iret)
471 write (iout,*) "rpotE",rpotE," iret",iret !d
473 call xdrffloat_(ixdrf, ruconst, iret)
474 call xdrffloat_(ixdrf, rt_bath, iret)
475 call xdrfint_(ixdrf, nss, iret)
477 call xdrfint_(ixdrf, ihpb(j), iret)
478 call xdrfint_(ixdrf, jhpb(j), iret)
480 call xdrfint_(ixdrf, nprop, iret)
481 if (umbrella(iparm) .or. read_iset(iparm) .or. hamil_rep) &
482 call xdrfint(ixdrf, iset, iret)
484 call xdrffloat_(ixdrf, rprop(i), iret)
487 call xdrffloat(ixdrf, rtime, iret)
488 call xdrffloat(ixdrf, rpotE, iret)
489 write (iout,*) "rpotE",rpotE," iret",iret !d
491 call xdrffloat(ixdrf, ruconst, iret)
492 call xdrffloat(ixdrf, rt_bath, iret)
493 call xdrfint(ixdrf, nss, iret)
495 call xdrfint(ixdrf, ihpb(j), iret)
496 call xdrfint(ixdrf, jhpb(j), iret)
498 call xdrfint(ixdrf, nprop, iret)
499 write (iout,*) "nprop",nprop !d
500 if (it.gt.0 .and. nprop.ne.nprop_prev) then
501 write (iout,*) "Warning previous nprop",nprop_prev,&
508 if (umbrella(iparm) .or. read_iset(iparm) .or. hamil_rep) &
509 call xdrfint(ixdrf, iset, iret)
511 call xdrffloat(ixdrf, rprop(i), iret)
515 itraj=mod(it,totraj(iR,iparm))
518 write (iout,*) "ii",ii," itraj",itraj," it",it
520 if (iset.eq.0) iset = 1
523 if (itraj.gt.ntraj) ntraj=itraj
524 nstep(itraj)=nstep(itraj)+1
525 ! rprop(2)=dsqrt(rprop(2))
526 ! rprop(3)=dsqrt(rprop(3))
528 write (iout,*) "umbrella ",umbrella
529 write (iout,*) rtime,rpotE,rt_bath,nss,&
530 (ihpb(j),jhpb(j),j=1,nss),(rprop(j),j=1,nprop)
531 write (iout,*) "nprop",nprop," iset",iset," myparm",myparm
536 #if (defined(AIX) && !defined(JUBL))
537 call xdrf3dfcoord_(ixdrf, xoord, itmp, prec, iret)
539 call xdrf3dfcoord(ixdrf, xoord, itmp, prec, iret)
542 write (iout,'(10f8.3)') ((xoord(j,i),j=1,3),i=1,2*nres+2)
546 if (itmp .ne. nres + nct - nnt + 1) then
547 write (iout,*) "Error: inconsistent sizes",itmp,nres+nct-nnt+1
553 write (iout,*) "calling slice" !d
555 islice=slice(nstep(itraj),time,is,ie,ts,te)
556 write (iout,*) "islice",islice !d
566 c(j,i+nres+nnt-1)=xoord(j,i+nres)
570 if (islice.gt.0 .and. islice.le.nslice .and. (.not.separate_parset &
571 .or. iset.eq.myparm)) then
573 kk(islice)=kk(islice)+1
574 mm(islice)=mm(islice)+1
575 if (mod(nstep(itraj),isampl(iparm)).eq.0 .and. &
576 conf_check(ll(islice)+1,1)) then
577 if (replica(iparm)) then
578 rt_bath=1.0d0/(rt_bath*1.987D-3)
580 if (abs(real(beta_h(i,iparm))-rt_bath).lt.1.0e-4) then
586 if (i.gt.nT_h(iparm)) then
587 write (iout,*) "Error - temperature of conformation",&
588 ii,1.0d0/(rt_bath*1.987D-3),&
589 " does not match any of the list"
591 1.0d0/(rt_bath*1.987D-3),&
592 (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
595 ! call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
597 kk(islice)=kk(islice)-1
598 mm(islice)=mm(islice)-1
606 jj(islice)=jj(islice)+1
607 if (umbrella(iparm)) then
608 snk(iset,iib,iparm,islice)=snk(iset,iib,iparm,islice)+1
609 else if (hamil_rep) then
610 snk(1,iib,iparm,islice)=snk(1,iib,iparm,islice)+1
612 snk(iR,iib,iparm,islice)=snk(iR,iib,iparm,islice)+1
614 ll(islice)=ll(islice)+1
616 write (iout,*) "Writing conformation, record",ll(islice)
617 write (iout,*) "ib",ib," iib",iib
618 write (iout,*) "ntraj",ntraj," itraj",itraj,&
619 " nstep",nstep(itraj)
620 write (iout,*) "pote",rpotE," time",rtime
621 ! if (replica(iparm)) then
622 ! write (iout,*) "TEMP",1.0d0/(rt_bath*1.987D-3)
623 ! write (iout,*) "TEMP list"
625 ! & (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
627 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
628 ! write (iout,*) "nres",nres," nnt",nnt," nct",nct," nss",nss
629 ! write (iout,*) "length",nres*4+(nct-nnt+1)*4+4+2*nss*4
632 if (islice.ne.islice1) then
633 ! write (iout,*) "islice",islice," islice1",islice1
635 ! write (iout,*) "Closing file ",
636 ! & bprotfile_temp(:ilen(bprotfile_temp))
637 call opentmp(islice,ientout,bprotfile_temp)
638 ! write (iout,*) "Opening file ",
639 ! & bprotfile_temp(:ilen(bprotfile_temp))
642 if (umbrella(iparm)) then
643 write(ientout,rec=ll(islice)) &
644 ((xoord(l,k),l=1,3),k=1,nres),&
645 ((xoord(l,k),l=1,3),k=nres+1,nres+nct-nnt+1),&
646 nss,(ihpb(k),jhpb(k),k=1,nss),&
647 rpotE+0.0d0,efree,rmsdev,(rprop(i)+0.0d0,i=1,nQ),&
649 else if (hamil_rep) then
650 write(ientout,rec=ll(islice)) &
651 ((xoord(l,k),l=1,3),k=1,nres),&
652 ((xoord(l,k),l=1,3),k=nres+1,nres+nct-nnt+1),&
653 nss,(ihpb(k),jhpb(k),k=1,nss),&
654 rpotE+0.0d0,efree,rmsdev,(rprop(i)+0.0d0,i=1,nQ),&
657 write(ientout,rec=ll(islice)) &
658 ((xoord(l,k),l=1,3),k=1,nres),&
659 ((xoord(l,k),l=1,3),k=nres+1,nres+nct-nnt+1),&
660 nss,(ihpb(k),jhpb(k),k=1,nss),&
661 rpotE+0.0d0,efree,rmsdev,(rprop(i)+0.0d0,i=1,nQ),&
665 call int_from_cart1(.false.)
666 write (iout,*) "Writing conformation, record",ll(islice)
667 write (iout,*) "Cartesian coordinates"
668 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
669 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
670 write (iout,*) "Internal coordinates"
671 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
672 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
673 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
674 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
675 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
676 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
677 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
678 ! write (iout,'(8f10.5)') (rprop(j),j=1,nQ)
679 write (iout,'(16i5)') iscor
689 #if (defined(AIX) && !defined(JUBL))
690 call xdrfclose_(ixdrf, iret)
692 call xdrfclose(ixdrf, iret)
694 write (iout,'(i10," trajectories found in file.")') ntraj+1
695 write (iout,'(a)') "Numbers of steps in trajectories:"
696 write (iout,'(8i10)') (nstep(i),i=0,ntraj)
697 write (iout,*) ii," conformations read from file",&
700 write (iout,*) mm(islice)," conformations read so far, slice",&
702 write (iout,*) ll(islice),&
703 " conformations stored so far, slice",islice
708 end subroutine cxread
709 !--------------------------------------------------------------------------------
711 !--------------------------------------------------------------------------------
712 subroutine xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
716 ! include "DIMENSIONS"
717 ! include "DIMENSIONS.ZSCOPT"
718 ! include "DIMENSIONS.FREE"
719 use MPI_data, only:nprocs
722 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
723 ! include "COMMON.MPI"
725 integer,parameter :: MaxTraj=2050
726 ! include "COMMON.CHAIN"
727 ! include "COMMON.IOUNITS"
728 ! include "COMMON.PROTFILES"
729 ! include "COMMON.NAMES"
730 ! include "COMMON.VAR"
731 ! include "COMMON.GEO"
732 ! include "COMMON.ENEPS"
733 ! include "COMMON.PROT"
734 ! include "COMMON.INTERACT"
735 ! include "COMMON.FREE"
736 ! include "COMMON.SBRIDGE"
737 ! include "COMMON.OBCINKA"
738 real(kind=4) :: csingle(3,nres*2)
739 character(len=64) :: nazwa,bprotfile_temp
740 integer :: i,j,k,l,ii,jj(nslice),kk(nslice),ll(nslice),&
741 mm(nslice) !(maxslice)
742 integer :: iscor,islice,islice1 !el,slice
743 real(kind=8) :: energ
745 ! external ilen,iroof
746 real(kind=8) :: rmsdev,energia(0:n_ene),efree,eini,temp
747 !el real(kind=8) :: rmsdev,energia(0:max_eneW),efree,eini,temp
748 real(kind=8) :: prop(nQ) !(maxQ)
749 integer :: ntot_all(0:nprocs-1)!(0:maxprocs-1)
750 integer :: iparm,ib,iib,ir,nprop,nthr
751 real(kind=8) :: etot,time,ts(nslice),te(nslice)
752 integer :: is(nslice),ie(nslice),itraj,ntraj,it,iset
753 integer :: nstep(0:MaxTraj-1)
756 call set_slices(is,ie,ts,te,iR,ib,iparm)
766 call opentmp(islice1,ientout,bprotfile_temp)
768 if (replica(iparm)) then
769 if (hamil_rep .or. umbrella(iparm)) then
770 read (ientin,*,end=1112,err=1112) time,eini,&
771 etot,temp,nss,(ihpb(j),jhpb(j),j=1,nss),&
772 nprop,(prop(j),j=1,nprop),iset
774 read (ientin,*,end=1112,err=1112) time,eini,&
775 etot,temp,nss,(ihpb(j),jhpb(j),j=1,nss),&
776 nprop,(prop(j),j=1,nprop)
778 temp=1.0d0/(temp*1.987D-3)
779 ! write (iout,*) time,eini,etot,nss,
780 ! & (ihpb(j),jhpb(j),j=1,nss),(prop(j),j=1,nprop)
783 if (beta_h(i,iparm).eq.temp) then
789 if (i.gt.nT_h(iparm)) then
790 write (iout,*) "Error - temperature of conformation",&
791 ii,1.0d0/(temp*1.987D-3),&
792 " does not match any of the list"
794 1.0d0/(temp*1.987D-3),&
795 (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
798 call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
802 read (ientin,*,end=1112,err=1112) time,eini,&
803 etot,nss,(ihpb(j),jhpb(j),j=1,nss),&
804 nprop,(prop(j),j=1,nprop)
807 itraj=mod(it,totraj(iR,iparm))
808 ! write (*,*) "ii",ii," itraj",itraj
811 if (itraj.gt.ntraj) ntraj=itraj
812 nstep(itraj)=nstep(itraj)+1
813 islice=slice(nstep(itraj),time,is,ie,ts,te)
814 read (ientin,'(8f10.5)',end=1112,err=1112) &
815 ((csingle(l,k),l=1,3),k=1,nres),&
816 ((csingle(l,k+nres),l=1,3),k=nnt,nct)
818 if (islice.gt.0 .and. islice.le.nslice) then
820 kk(islice)=kk(islice)+1
821 mm(islice)=mm(islice)+1
822 if (mod(nstep(itraj),isampl(iparm)).eq.0) then
823 jj(islice)=jj(islice)+1
825 snk(iR,iib,iset,islice)=snk(iR,iib,iset,islice)+1
826 else if (umbrella(iparm)) then
827 snk(iset,iib,iparm,islice)=snk(iset,iib,iparm,islice)+1
829 snk(iR,iib,iparm,islice)=snk(iR,iib,iparm,islice)+1
831 ll(islice)=ll(islice)+1
832 ! write (iout,*) ii,kk,jj,ll,eini,(prop(j),j=1,nprop)
834 ! write (iout,*) "Writing conformation, record",ll(islice)
835 ! write (iout,*) "ib",ib," iib",iib
836 if (replica(iparm)) then
837 write (iout,*) "TEMP",1.0d0/(temp*1.987D-3)
838 write (iout,*) "TEMP list"
840 (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
844 ! write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
845 ! write (iout,*) "nres",nres," nnt",nnt," nct",nct," nss",nss
846 ! write (iout,*) "length",nres*4+(nct-nnt+1)*4+4+2*nss*4
848 if (islice.ne.islice1) then
849 ! write (iout,*) "islice",islice," islice1",islice1
851 ! write (iout,*) "Closing file ",
852 ! & bprotfile_temp(:ilen(bprotfile_temp))
853 call opentmp(islice,ientout,bprotfile_temp)
854 ! write (iout,*) "Opening file ",
855 ! & bprotfile_temp(:ilen(bprotfile_temp))
859 write(ientout,rec=ll(islice)) &
860 ((csingle(l,k),l=1,3),k=1,nres),&
861 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
862 nss,(ihpb(k),jhpb(k),k=1,nss),&
863 eini,efree,rmsdev,(prop(i),i=1,nQ),iR,iib,iparm
870 call int_from_cart1(.false.)
871 write (iout,*) "Writing conformation, record",ll(islice)
872 write (iout,*) "Cartesian coordinates"
873 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
874 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
875 write (iout,*) "Internal coordinates"
876 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
877 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
878 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
879 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
880 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
881 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
882 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
883 ! write (iout,'(8f10.5)') (prop(j),j=1,nQ)
884 write (iout,'(16i5)') iscor
892 write (iout,'(i10," trajectories found in file.")') ntraj+1
893 write (iout,'(a)') "Numbers of steps in trajectories:"
894 write (iout,'(8i10)') (nstep(i),i=0,ntraj)
895 write (iout,*) ii," conformations read from file",&
897 write (iout,*) mm(islice)," conformations read so far, slice",&
899 write (iout,*) ll(islice)," conformations stored so far, slice",&
904 !--------------------------------------------------------------------------------
906 !--------------------------------------------------------------------------------
907 subroutine write_dbase(islice,*)
910 use control_data, only:indpdb
912 use conform_compar, only:conf_compar
914 ! include "DIMENSIONS"
915 ! include "DIMENSIONS.ZSCOPT"
916 ! include "DIMENSIONS.FREE"
917 ! include "DIMENSIONS.COMPAR"
918 use geometry, only:int_from_cart1
921 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
922 ! include "COMMON.MPI"
924 ! include "COMMON.CONTROL"
925 ! include "COMMON.CHAIN"
926 ! include "COMMON.IOUNITS"
927 ! include "COMMON.PROTFILES"
928 ! include "COMMON.NAMES"
929 ! include "COMMON.VAR"
930 ! include "COMMON.SBRIDGE"
931 ! include "COMMON.GEO"
932 ! include "COMMON.FFIELD"
933 ! include "COMMON.ENEPS"
934 ! include "COMMON.LOCAL"
935 ! include "COMMON.WEIGHTS"
936 ! include "COMMON.INTERACT"
937 ! include "COMMON.FREE"
938 ! include "COMMON.ENERGIES"
939 ! include "COMMON.COMPAR"
940 ! include "COMMON.PROT"
941 ! include "COMMON.CONTACTS1"
942 character(len=64) :: nazwa
943 character(len=80) :: bxname,cxname
944 character(len=64) :: bprotfile_temp
945 character(len=3) :: liczba,licz
946 character(len=2) :: licz2
947 integer :: i,itj,ii,iii,j,k,l
948 integer :: ixdrf,iret
949 integer :: iscor,islice
950 real(kind=8) :: rmsdev,efree,eini
951 real(kind=4) :: csingle(3,nres*2)
952 real(kind=8) :: energ
954 ! external ilen,iroof
955 integer :: ir,ib,iparm
956 integer :: isecstr(nres)
957 write (licz2,'(bz,i2.2)') islice
958 call opentmp(islice,ientout,bprotfile_temp)
959 write (iout,*) "bprotfile_temp ",bprotfile_temp
961 if (.not.bxfile .and. .not. cxfile .and. indpdb.eq.0 &
962 .and. ensembles.eq.0) then
963 close(ientout,status="delete")
967 write (liczba,'(bz,i3.3)') me
968 if (bxfile .or. cxfile .or. ensembles.gt.0) then
969 if (.not.separate_parset) then
970 bxname = prefix(:ilen(prefix))//liczba//".bx"
972 write (licz,'(bz,i3.3)') myparm
973 bxname = prefix(:ilen(prefix))//liczba//"_par"//licz//".bx"
975 open (ientin,file=bxname,status="unknown",&
976 form="unformatted",access="direct",recl=lenrec1)
979 if (bxfile .or. cxfile .or. ensembles.gt.0) then
980 if (nslice.eq.1) then
981 bxname = prefix(:ilen(prefix))//".bx"
983 bxname = prefix(:ilen(prefix))// &
984 "_slice_"//licz2//".bx"
986 open (ientin,file=bxname,status="unknown",&
987 form="unformatted",access="direct",recl=lenrec1)
988 write (iout,*) "Calculating energies; writing geometry",&
989 " and energy components to ",bxname(:ilen(bxname))
991 #if (defined(AIX) && !defined(JUBL))
992 call xdrfopen_(ixdrf,cxname, "w", iret)
994 call xdrfopen(ixdrf,cxname, "w", iret)
997 write (iout,*) "Error opening cxfile ",cxname(:ilen(cxname))
1002 if (indpdb.gt.0) then
1003 if (nslice.eq.1) then
1005 if (.not.separate_parset) then
1006 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot)) &
1009 write (licz,'(bz,i3.3)') myparm
1010 statname=prefix(:ilen(prefix))//'_par'//licz//'_'// &
1011 pot(:ilen(pot))//liczba//'.stat'
1015 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))//'.stat'
1019 if (.not.separate_parset) then
1020 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))// &
1021 "_slice_"//licz2//liczba//'.stat'
1023 write (licz,'(bz,i3.3)') myparm
1024 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))// &
1025 '_par'//licz//"_slice_"//licz2//liczba//'.stat'
1028 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot)) &
1029 //"_slice_"//licz2//'.stat'
1032 open(istat,file=statname,status="unknown")
1040 read(ientout,rec=i,err=101) &
1041 ((csingle(l,k),l=1,3),k=1,nres),&
1042 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
1043 nss,(ihpb(k),jhpb(k),k=1,nss),&
1044 eini,efree,rmsdev,(q(k,i),k=1,nQ),iR,ib,iparm
1045 ! write (iout,*) iR,ib,iparm,eini,efree
1051 call int_from_cart1(.false.)
1053 ! write (iout,*) "Calling conf_compar",i
1055 anatemp= 1.0d0/(beta_h(ib,iparm)*1.987D-3)
1056 if (indpdb.gt.0) then
1057 call conf_compar(i,.false.,.true.)
1059 ! call elecont(.false.,ncont,icont,nnt,nct)
1060 ! call secondary2(.false.,.false.,ncont,icont,isecstr)
1062 ! write (iout,*) "Exit conf_compar",i
1064 if (bxfile .or.cxfile .or. ensembles.gt.0) write (ientin,rec=i) &
1065 ((csingle(l,k),l=1,3),k=1,nres),&
1066 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
1067 nss,(ihpb(k),jhpb(k),k=1,nss),&
1068 ! & potE(i,iparm),-entfac(i),rms_nat,iscore
1069 potE(i,nparmset),-entfac(i),rms_nat,iscore
1070 ! write (iout,'(2i5,3e15.5)') i,me,potE(i,iparm),-entfac(i)
1072 if (cxfile) call cxwrite(ixdrf,csingle,potE(i,nparmset),&
1073 -entfac(i),rms_nat,iscore)
1076 close(ientout,status="delete")
1078 if (bxfile .or. cxfile .or. ensembles.gt.0) close(ientin)
1080 call MPI_Barrier(WHAM_COMM,IERROR)
1081 if (me.ne.Master .or. .not.bxfile .and. .not. cxfile &
1082 .and. ensembles.eq.0) return
1084 if (bxfile .or. ensembles.gt.0) then
1085 if (nslice.eq.1) then
1086 if (.not.separate_parset) then
1087 bxname = prefix(:ilen(prefix))//".bx"
1089 write (licz,'(bz,i3.3)') myparm
1090 bxname = prefix(:ilen(prefix))//"_par"//licz//".bx"
1093 if (.not.separate_parset) then
1094 bxname = prefix(:ilen(prefix))//"_slice_"//licz2//".bx"
1096 write (licz,'(bz,i3.3)') myparm
1097 bxname = prefix(:ilen(prefix))//"par_"//licz// &
1098 "_slice_"//licz2//".bx"
1101 open (ientout,file=bxname,status="unknown",&
1102 form="unformatted",access="direct",recl=lenrec1)
1103 write (iout,*) "Master is creating binary database ",&
1104 bxname(:ilen(bxname))
1107 if (nslice.eq.1) then
1108 if (.not.separate_parset) then
1109 cxname = prefix(:ilen(prefix))//".cx"
1111 cxname = prefix(:ilen(prefix))//"_par"//licz//".cx"
1114 if (.not.separate_parset) then
1115 cxname = prefix(:ilen(prefix))// &
1116 "_slice_"//licz2//".cx"
1118 cxname = prefix(:ilen(prefix))//"_par"//licz// &
1119 "_slice_"//licz2//".cx"
1122 #if (defined(AIX) && !defined(JUBL))
1123 call xdrfopen_(ixdrf,cxname, "w", iret)
1125 call xdrfopen(ixdrf,cxname, "w", iret)
1128 write (iout,*) "Error opening cxfile ",cxname(:ilen(cxname))
1133 write (liczba,'(bz,i3.3)') j
1134 if (separate_parset) then
1135 write (licz,'(bz,i3.3)') myparm
1136 bxname = prefix(:ilen(prefix))//liczba//"_par"//licz//".bx"
1138 bxname = prefix(:ilen(prefix))//liczba//".bx"
1140 open (ientin,file=bxname,status="unknown",&
1141 form="unformatted",access="direct",recl=lenrec1)
1142 write (iout,*) "Master is reading conformations from ",&
1143 bxname(:ilen(bxname))
1145 ! write (iout,*) "j",j," indstart",indstart(j)," indend",indend(j)
1147 do i=indstart(j),indend(j)
1149 read(ientin,rec=iii,err=101) &
1150 ((csingle(l,k),l=1,3),k=1,nres),&
1151 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
1152 nss,(ihpb(k),jhpb(k),k=1,nss),&
1153 eini,efree,rmsdev,iscor
1154 if (bxfile .or. ensembles.gt.0) then
1155 write (ientout,rec=i) &
1156 ((csingle(l,k),l=1,3),k=1,nres),&
1157 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
1158 nss,(ihpb(k),jhpb(k),k=1,nss),&
1159 eini,efree,rmsdev,iscor
1161 if(cxfile)call cxwrite(ixdrf,csingle,eini,efree,rmsdev,iscor)
1168 call int_from_cart1(.false.)
1169 write (iout,'(2i5,3e15.5)') i,iii,eini,efree
1170 write (iout,*) "The Cartesian geometry is:"
1171 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1172 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1173 write (iout,*) "The internal geometry is:"
1174 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
1175 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
1176 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
1177 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
1178 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
1179 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
1180 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
1181 write (iout,'(f10.5,i5)') rmsdev,iscor
1184 write (iout,*) iii," conformations (from",indstart(j)," to",&
1185 indend(j),") read from ",&
1186 bxname(:ilen(bxname))
1187 close (ientin,status="delete")
1189 if (bxfile .or. cxfile .or. ensembles.gt.0) close (ientout)
1190 #if (defined(AIX) && !defined(JUBL))
1191 if (cxfile) call xdrfclose_(ixdrf,cxname,iret)
1193 if (cxfile) call xdrfclose(ixdrf,cxname,iret)
1197 101 write (iout,*) "Error in scratchfile."
1200 end subroutine write_dbase
1201 !-------------------------------------------------------------------------------
1202 subroutine cxwrite(ixdrf,csingle,eini,efree,rmsdev,iscor)
1204 ! include "DIMENSIONS"
1205 ! include "DIMENSIONS.ZSCOPT"
1206 ! include "DIMENSIONS.FREE"
1207 ! include "DIMENSIONS.COMPAR"
1210 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
1211 ! include "COMMON.MPI"
1213 ! include "COMMON.CONTROL"
1214 ! include "COMMON.CHAIN"
1215 ! include "COMMON.IOUNITS"
1216 ! include "COMMON.PROTFILES"
1217 ! include "COMMON.NAMES"
1218 ! include "COMMON.VAR"
1219 ! include "COMMON.SBRIDGE"
1220 ! include "COMMON.GEO"
1221 ! include "COMMON.FFIELD"
1222 ! include "COMMON.ENEPS"
1223 ! include "COMMON.LOCAL"
1224 ! include "COMMON.WEIGHTS"
1225 ! include "COMMON.INTERACT"
1226 ! include "COMMON.FREE"
1227 ! include "COMMON.ENERGIES"
1228 ! include "COMMON.COMPAR"
1229 ! include "COMMON.PROT"
1230 integer :: i,j,itmp,iscor,iret,ixdrf
1231 real(kind=8) :: rmsdev,efree,eini
1232 real(kind=4) :: csingle(3,nres*2),xoord(3,2*nres+2)
1233 real(kind=4) :: prec
1235 ! write (iout,*) "cxwrite"
1240 xoord(j,i)=csingle(j,i)
1245 xoord(j,nres+i-nnt+1)=csingle(j,i+nres)
1251 ! write (iout,*) "itmp",itmp
1253 #if (defined(AIX) && !defined(JUBL))
1254 call xdrf3dfcoord_(ixdrf, xoord, itmp, prec, iret)
1256 ! write (iout,*) "xdrf3dfcoord"
1258 call xdrfint_(ixdrf, nss, iret)
1260 call xdrfint_(ixdrf, ihpb(j), iret)
1261 call xdrfint_(ixdrf, jhpb(j), iret)
1263 call xdrffloat_(ixdrf,real(eini),iret)
1264 call xdrffloat_(ixdrf,real(efree),iret)
1265 call xdrffloat_(ixdrf,real(rmsdev),iret)
1266 call xdrfint_(ixdrf,iscor,iret)
1268 call xdrf3dfcoord(ixdrf, xoord, itmp, prec, iret)
1270 call xdrfint(ixdrf, nss, iret)
1272 call xdrfint(ixdrf, ihpb(j), iret)
1273 call xdrfint(ixdrf, jhpb(j), iret)
1275 call xdrffloat(ixdrf,real(eini),iret)
1276 call xdrffloat(ixdrf,real(efree),iret)
1277 call xdrffloat(ixdrf,real(rmsdev),iret)
1278 call xdrfint(ixdrf,iscor,iret)
1282 end subroutine cxwrite
1283 !-------------------------------------------------------------------------------
1285 !-------------------------------------------------------------------------------
1286 subroutine set_slices(is,ie,ts,te,iR,ib,iparm)
1288 ! include 'DIMENSIONS'
1289 ! include 'DIMENSIONS.ZSCOPT'
1290 ! include 'DIMENSIONS.FREE'
1291 ! include 'COMMON.IOUNITS'
1292 ! include 'COMMON.PROTFILES'
1293 ! include 'COMMON.OBCINKA'
1294 ! include 'COMMON.PROT'
1295 integer :: islice,iR,ib,iparm
1296 integer :: is(MaxSlice),ie(MaxSlice),nrec_slice
1297 real(kind=8) :: ts(MaxSlice),te(MaxSlice),time_slice
1300 if (time_end_collect(iR,ib,iparm).ge.1.0d10) then
1301 ts(islice)=time_start_collect(iR,ib,iparm)
1302 te(islice)=time_end_collect(iR,ib,iparm)
1303 nrec_slice=(rec_end(iR,ib,iparm)- &
1304 rec_start(iR,ib,iparm)+1)/nslice
1305 is(islice)=rec_start(iR,ib,iparm)+(islice-1)*nrec_slice
1306 ie(islice)=rec_start(iR,ib,iparm)+islice*nrec_slice-1
1308 time_slice=(time_end_collect(iR,ib,iparm) &
1309 -time_start_collect(iR,ib,iparm))/nslice
1310 ts(islice)=time_start_collect(iR,ib,iparm)+(islice-1)* &
1312 te(islice)=time_start_collect(iR,ib,iparm)+islice*time_slice
1313 is(islice)=rec_start(iR,ib,iparm)
1314 ie(islice)=rec_end(iR,ib,iparm)
1318 write (iout,*) "nrec_slice",nrec_slice," time_slice",time_slice
1319 write (iout,*) "is",(is(islice),islice=1,nslice)
1320 write (iout,*) "ie",(ie(islice),islice=1,nslice)
1321 write (iout,*) "rec_start",&
1322 rec_start(iR,ib,iparm)," rec_end",rec_end(iR,ib,iparm)
1323 write (iout,*) "ts",(ts(islice),islice=1,nslice)
1324 write (iout,*) "te",(te(islice),islice=1,nslice)
1325 write (iout,*) "time_start",&
1326 time_start_collect(iR,ib,iparm)," time_end",&
1327 time_end_collect(iR,ib,iparm)
1331 end subroutine set_slices
1332 !-----------------------------------------------------------------------------
1333 integer function slice(irecord,time,is,ie,ts,te)
1335 ! include 'DIMENSIONS'
1336 ! include 'DIMENSIONS.ZSCOPT'
1337 ! include 'DIMENSIONS.FREE'
1338 ! include 'COMMON.IOUNITS'
1339 ! include 'COMMON.PROTFILES'
1340 ! include 'COMMON.OBCINKA'
1341 ! include 'COMMON.PROT'
1342 integer :: is(MaxSlice),ie(MaxSlice),nrec_slice
1343 real(kind=8) :: ts(MaxSlice),te(MaxSlice),time_slice
1344 integer :: i,ii,irecord
1345 real(kind=8) :: time
1347 ! write (iout,*) "within slice nslice",nslice
1349 if (irecord.lt.is(1) .or. time.lt.ts(1)) then
1353 do while (ii.le.nslice .and. &
1354 (irecord.lt.is(ii) .or. irecord.gt.ie(ii) .or. &
1355 time.lt.ts(ii) .or. time.gt.te(ii)) )
1356 ! write (iout,*) "ii",ii,time,ts(ii)
1361 ! write (iout,*) "end: ii",ii
1366 !-----------------------------------------------------------------------------
1368 !-----------------------------------------------------------------------------
1369 logical function conf_check(ii,iprint)
1371 use names, only:ntyp1
1373 use energy_data, only:itype,dsc
1374 use geometry, only:int_from_cart1
1376 ! include "DIMENSIONS"
1377 ! include "DIMENSIONS.ZSCOPT"
1378 ! include "DIMENSIONS.FREE"
1382 ! include "COMMON.MPI"
1384 ! include "COMMON.CHAIN"
1385 ! include "COMMON.IOUNITS"
1386 ! include "COMMON.PROTFILES"
1387 ! include "COMMON.NAMES"
1388 ! include "COMMON.VAR"
1389 ! include "COMMON.SBRIDGE"
1390 ! include "COMMON.GEO"
1391 ! include "COMMON.FFIELD"
1392 ! include "COMMON.ENEPS"
1393 ! include "COMMON.LOCAL"
1394 ! include "COMMON.WEIGHTS"
1395 ! include "COMMON.INTERACT"
1396 ! include "COMMON.FREE"
1397 ! include "COMMON.ENERGIES"
1398 ! include "COMMON.CONTROL"
1399 ! include "COMMON.TORCNSTR"
1403 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
1405 integer :: j,k,l,ii,itj,iprint
1406 if (.not. check_conf) then
1410 call int_from_cart1(.false.)
1412 if (itype(j-1).ne.ntyp1 .and. itype(j).ne.ntyp1 .and. &
1413 (vbld(j).lt.2.0d0 .or. vbld(j).gt.5.0d0)) then
1415 write (iout,*) "Bad CA-CA bond length",j," ",vbld(j),&
1416 " for conformation",ii
1417 if (iprint.gt.1) then
1418 write (iout,*) "The Cartesian geometry is:"
1419 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1420 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1421 write (iout,*) "The internal geometry is:"
1422 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
1423 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
1424 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
1425 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
1426 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
1427 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
1429 if (iprint.gt.0) write (iout,*) &
1430 "This conformation WILL NOT be added to the database."
1437 if (itype(j).ne.10 .and.itype(j).ne.ntyp1 .and. &
1438 (vbld(nres+j)-dsc(iabs(itj))).gt.2.0d0) then
1440 write (iout,*) "Bad CA-SC bond length",j," ",vbld(nres+j),&
1441 " for conformation",ii
1442 if (iprint.gt.1) then
1443 write (iout,*) "The Cartesian geometry is:"
1444 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1445 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1446 write (iout,*) "The internal geometry is:"
1447 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
1448 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
1449 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
1450 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
1451 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
1452 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
1454 if (iprint.gt.0) write (iout,*) &
1455 "This conformation WILL NOT be added to the database."
1461 if (theta(j).le.0.0d0) then
1463 write (iout,*) "Zero theta angle(s) in conformation",ii
1464 if (iprint.gt.1) then
1465 write (iout,*) "The Cartesian geometry is:"
1466 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1467 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1468 write (iout,*) "The internal geometry is:"
1469 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
1470 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
1471 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
1472 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
1473 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
1474 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
1476 if (iprint.gt.0) write (iout,*) &
1477 "This conformation WILL NOT be added to the database."
1481 if (theta(j).gt.179.97*deg2rad) theta(j)=179.97*deg2rad
1484 ! write (iout,*) "conf_check passed",ii
1486 end function conf_check
1487 !-----------------------------------------------------------------------------
1488 end module io_database