2 !-----------------------------------------------------------------------------
7 use energy_data, only:nnt,nct,nss,ihpb,jhpb
9 use geometry_data, only:nres,c
12 ! include "COMMON.MPI"
16 !-----------------------------------------------------------------------------
19 !-----------------------------------------------------------------------------
21 !-----------------------------------------------------------------------------
23 !-------------------------------------------------------------------------------
24 subroutine opentmp(islice,iunit,bprotfile_temp)
26 ! include "DIMENSIONS"
27 ! include "DIMENSIONS.ZSCOPT"
28 ! include "DIMENSIONS.FREE"
29 ! use MPI_data, only:me
32 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
33 ! include "COMMON.MPI"
35 ! include "COMMON.IOUNITS"
36 ! include "COMMON.PROTFILES"
37 ! include "COMMON.PROT"
38 ! include "COMMON.FREE"
39 character(len=64) :: bprotfile_temp
40 character(len=3) :: liczba,liczba2
41 character(len=2) :: liczba1
42 integer :: iunit,islice
46 ! integer :: lenrec,lenrec2
49 ! lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
51 write (liczba1,'(bz,i2.2)') islice
53 write (liczba,'(bz,i3.3)') me
55 ! write (iout,*) "separate_parset ",separate_parset,
57 if (separate_parset) then
58 write (liczba2,'(bz,i3.3)') myparm
59 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
60 prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
61 open (iunit,file=bprotfile_temp,status="unknown",&
62 form="unformatted",access="direct",recl=lenrec)
64 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
65 prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
66 open (iunit,file=bprotfile_temp,status="unknown",&
67 form="unformatted",access="direct",recl=lenrec)
70 bprotfile_temp = scratchdir(:ilen(scratchdir))// &
71 "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
72 open (iunit,file=bprotfile_temp,status="unknown",&
73 form="unformatted",access="direct",recl=lenrec)
75 ! write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
79 end subroutine opentmp
80 !-------------------------------------------------------------------------------
81 subroutine read_database(*)
83 ! use energy_data, only:nct,nnt,nss
85 ! include "DIMENSIONS"
86 ! include "DIMENSIONS.ZSCOPT"
87 ! include "DIMENSIONS.FREE"
88 use MPI_data, only:me,nprocs
91 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
92 ! include "COMMON.MPI"
94 ! include "COMMON.CHAIN"
95 ! include "COMMON.IOUNITS"
96 ! include "COMMON.PROTFILES"
97 ! include "COMMON.NAMES"
98 ! include "COMMON.VAR"
99 ! include "COMMON.GEO"
100 ! include "COMMON.ENEPS"
101 ! include "COMMON.PROT"
102 ! include "COMMON.INTERACT"
103 ! include "COMMON.FREE"
104 ! include "COMMON.SBRIDGE"
105 ! include "COMMON.OBCINKA"
106 real(kind=4) :: csingle(3,nres*2) !(3,maxres2)
107 character(len=64) :: nazwa,bprotfile_temp
108 character(len=3) :: liczba
109 character(len=2) :: liczba1
110 integer :: i,j,ii,jj(nslice),k,kk(nslice),l,&
111 ll(nslice),mm(nslice),if
112 integer :: nrec,nlines,iscor,iunit,islice
113 real(kind=8) :: energ
115 ! external ilen,iroof
116 real(kind=8) :: rmsdev,energia(0:n_ene),efree,eini,temp
117 !el real(kind=8) :: rmsdev,energia(0:max_eneW),efree,eini,temp
118 real(kind=8) :: prop(nQ) !(maxQ)
119 integer :: ntot_all(nslice,0:nprocs-1)!(maxslice,0:maxprocs-1)
120 integer :: iparm,ib,iib,ir,nprop,nthr,npars
121 real(kind=8) :: etot,time
122 integer :: ixdrf,iret
123 logical :: lerr,linit
125 lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
126 lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
128 write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,&
138 write (iout,*) "nparmset",nparmset
146 if (replica(iparm)) then
153 do iR=1,nRR(ib,iparm)
155 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
161 IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
162 ! Read conformations from binary DA files (one per batch) and write them to
163 ! a binary DA scratchfile.
164 write (liczba,'(bz,i3.3)') me
165 do if=1,nfile_bin(iR,ib,iparm)
166 nazwa=protfiles(if,1,iR,ib,iparm) &
167 (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
168 open (ientin,file=nazwa,status="old",form="unformatted",&
169 access="direct",recl=lenrec2,err=1111)
172 call opentmp(islice,ientout,bprotfile_temp)
173 call bxread(nazwa,islice,ii,jj(islice),kk(islice),ll(islice),&
174 mm(islice),iR,ib,iparm)
181 IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
182 ! Read conformations from multiple ASCII int files and write them to a binary
184 do if=1,nfile_asc(iR,ib,iparm)
185 nazwa=protfiles(if,2,iR,ib,iparm) &
186 (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
187 open(unit=ientin,file=nazwa,status='old',err=1111)
188 write(iout,*) "reading ",nazwa(:ilen(nazwa))
190 call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
193 IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
194 ! Read conformations from cx files and write them to a binary
196 do if=1,nfile_cx(iR,ib,iparm)
197 nazwa=protfiles(if,2,iR,ib,iparm) &
198 (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
199 write(iout,*) "reading ",nazwa(:ilen(nazwa))
201 print *,"Calling cxread"
202 call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,&
204 print *,"after call cxread"
205 write(iout,*)"after call cxread"
207 write (iout,*) "exit cxread"
211 write(iout,*)"*********************in read database"
215 stot(islice)=stot(islice)+jj(islice)
220 write (iout,*) "IPARM",iparm
223 if (nslice.eq.1) then
225 write (liczba,'(bz,i3.3)') me
226 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
227 prefix(:ilen(prefix))//liczba//".xbin.tmp"
229 bprotfile_temp = scratchdir(:ilen(scratchdir))// &
230 "/"//prefix(:ilen(prefix))//".xbin.tmp"
232 write(iout,*) mm(1)," conformations read",ll(1),&
233 " conformations written to ",&
234 bprotfile_temp(:ilen(bprotfile_temp))
237 write (liczba1,'(bz,i2.2)') islice
239 write (liczba,'(bz,i3.3)') me
240 bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
241 prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
243 bprotfile_temp = scratchdir(:ilen(scratchdir))// &
244 "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
246 write(iout,*) mm(islice)," conformations read",ll(islice),&
247 " conformations written to ",&
248 bprotfile_temp(:ilen(bprotfile_temp))
253 ! Check if everyone has the same number of conformations
254 call MPI_Allgather(stot(1),nslice,MPI_INTEGER,&
255 ntot_all(1,0),nslice,MPI_INTEGER,MPI_Comm_World,IERROR)
260 if (stot(islice).ne.ntot_all(islice,i)) then
261 write (iout,*) "Number of conformations at processor",i,&
262 " differs from that at processor",me,&
263 stot(islice),ntot_all(islice,i)," slice",islice
271 write (iout,*) "Numbers of conformations read by processors"
274 write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
276 write (iout,*) "Calculation terminated."
281 ntot(islice)=stot(islice)
283 write(iout,*) "end of read database"
286 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
289 end subroutine read_database
290 !--------------------------------------------------------------------------------
291 integer function iroof(n,m)
294 if (ii*m .lt. n) ii=ii+1
298 !--------------------------------------------------------------------------------
300 !--------------------------------------------------------------------------------
301 subroutine bxread(nazwa,islice,ii,jj,kk,ll,mm,iR,ib,iparm)
303 ! include "DIMENSIONS"
304 ! include "DIMENSIONS.ZSCOPT"
305 ! include "DIMENSIONS.FREE"
306 ! use energy_data, only:nnt,nct,nss,ihpb,jhpbi
307 use MPI_data, only:nprocs
310 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
311 ! include "COMMON.MPI"
313 ! include "COMMON.CHAIN"
314 ! include "COMMON.IOUNITS"
315 ! include "COMMON.PROTFILES"
316 ! include "COMMON.NAMES"
317 ! include "COMMON.VAR"
318 ! include "COMMON.GEO"
319 ! include "COMMON.ENEPS"
320 ! include "COMMON.PROT"
321 ! include "COMMON.INTERACT"
322 ! include "COMMON.FREE"
323 ! include "COMMON.SBRIDGE"
324 real(kind=4) :: csingle(3,nres*2) !(3,maxres2)
325 character(len=64) :: nazwa,bprotfile_temp
326 character(len=3) :: liczba
327 integer :: i,is,ie,j,ii,jj,k,kk,l,ll,mm,if
328 integer :: nrec,nlines,iscor,islice
329 real(kind=8) :: energ
331 ! external ilen,iroof
332 real(kind=8) :: rmsdev,energia(0:n_ene),efree,eini,temp
333 !el real(kind=8) :: rmsdev,energia(0:max_eneW),efree,eini,temp
334 real(kind=8) :: prop(nQ) !(maxQ)
335 integer :: ntot_all(0:nprocs-1)!(0:maxprocs-1)
336 integer :: iparm,ib,iib,ir,nprop,nthr,nrec_slice
337 real(kind=8) :: etot,time
339 nrec_slice=(rec_end(iR,ib,iparm)-rec_start(iR,ib,iparm)+1)/nslice
340 is=rec_start(iR,ib,iparm)+(islice-1)*nrec_slice
341 ie=rec_start(iR,ib,iparm)+islice*nrec_slice-1
342 write (iout,*) "bxread: islice",islice," nslice",nslice,&
343 " nrec_slice",nrec_slice
344 write (iout,*) "is",is," ie",ie,"rec_start",&
345 rec_start(iR,ib,iparm)," rec_end",rec_end(iR,ib,iparm)
347 read(ientin,rec=i+1,err=101) &
348 ((csingle(l,k),l=1,3),k=1,nres),&
349 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
350 nss,(ihpb(k),jhpb(k),k=1,nss),&
351 eini,efree,rmsdev,(prop(j),j=1,nQ),iscor
354 if (mod(kk,isampl(iparm)).eq.0) then
356 write(ientout,rec=jj) &
357 ((csingle(l,k),l=1,3),k=1,nres),&
358 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
359 nss,(ihpb(k),jhpb(k),k=1,nss),&
360 eini,efree,rmsdev,(prop(j),j=1,nQ),iR,ib,iparm
367 call int_from_cart1(.false.)
368 write (iout,*) "Writing conformation, record",jj
369 write (iout,*) "Cartesian coordinates"
370 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
371 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
372 write (iout,*) "Internal coordinates"
373 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
374 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
375 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
376 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
377 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
378 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
379 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
380 write (iout,'(f10.5,i5)') rmsdev,iscor
386 write (iout,*) ii," conformations read from DA file ",&
388 write (iout,*) kk," conformations read so far, slice",islice
389 write (iout,*) jj," conformations stored so far, slice",islice
392 end subroutine bxread
393 !--------------------------------------------------------------------------------
395 !--------------------------------------------------------------------------------
396 subroutine cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,*)
400 use geometry, only:int_from_cart1
401 use geometry_data, only:vbld,rad2deg,theta,phi,alph,omeg
404 ! implicit real*8 (a-h,o-z)
405 ! include 'DIMENSIONS'
406 ! include 'DIMENSIONS.ZSCOPT'
407 ! include 'DIMENSIONS.FREE'
408 integer,parameter :: MaxTraj=2050
409 ! include 'COMMON.CHAIN'
410 ! include 'COMMON.INTERACT'
411 ! include 'COMMON.NAMES'
412 ! include 'COMMON.IOUNITS'
413 ! include 'COMMON.HEADER'
414 ! include 'COMMON.SBRIDGE'
415 ! include 'COMMON.PROTFILES'
416 ! include 'COMMON.OBCINKA'
417 ! include 'COMMON.FREE'
418 ! include 'COMMON.VAR'
419 ! include 'COMMON.GEO'
420 ! include 'COMMON.PROT'
421 character(len=64) :: nazwa,bprotfile_temp
422 real(kind=4) :: rtime,rpotE,ruconst,rt_bath,rprop(nQ) !(2000) !(maxQ)
424 integer :: iret,itmp,itraj,ntraj
425 real(kind=4) :: xoord(3,2*nres+2),prec
426 integer :: nstep(0:MaxTraj-1)
429 integer :: ii,jj(nslice),kk(nslice),ll(nslice),mm(nslice) !(maxslice)
430 integer :: is(nSlice),ie(nSlice),nrec_slice
431 real(kind=8) :: ts(nSlice),te(nSlice),time_slice
432 integer :: iR,ib,iparm,i,j,it,islice,nprop_prev
433 integer :: k,l,iib,islice1,nprop
434 real(kind=8) :: efree,rmsdev
437 ! logical :: conf_check
446 call set_slices(is,ie,ts,te,iR,ib,iparm)
457 #if (defined(AIX) && !defined(JUBL))
458 call xdrfopen_(ixdrf,nazwa, "r", iret)
460 call xdrfopen(ixdrf,nazwa, "r", iret)
462 if (iret.eq.0) return 1
465 call opentmp(islice1,ientout,bprotfile_temp)
469 #if (defined(AIX) && !defined(JUBL))
470 call xdrffloat_(ixdrf, rtime, iret)
471 print *,"rtime",rtime," iret",iret !d
472 call xdrffloat_(ixdrf, rpotE, iret)
473 ! write (iout,*) "rpotE",rpotE," iret",iret !d
475 call xdrffloat_(ixdrf, ruconst, iret)
476 call xdrffloat_(ixdrf, rt_bath, iret)
477 call xdrfint_(ixdrf, nss, iret)
479 call xdrfint_(ixdrf, ihpb(j), iret)
480 call xdrfint_(ixdrf, jhpb(j), iret)
482 call xdrfint_(ixdrf, nprop, iret)
483 if (umbrella(iparm) .or. read_iset(iparm) .or. hamil_rep) &
484 call xdrfint(ixdrf, iset, iret)
486 call xdrffloat_(ixdrf, rprop(i), iret)
489 call xdrffloat(ixdrf, rtime, iret)
490 print *,"rtime",rtime," iret",iret !d
491 call xdrffloat(ixdrf, rpotE, iret)
492 ! write (iout,*) "rpotE",rpotE," iret",iret !d
494 call xdrffloat(ixdrf, ruconst, iret)
495 call xdrffloat(ixdrf, rt_bath, iret)
496 call xdrfint(ixdrf, nss, iret)
499 call xdrfint(ixdrf, ihpb(j), iret)
500 call xdrfint(ixdrf, jhpb(j), iret)
502 call xdrfint(ixdrf, nprop, iret)
503 ! write (iout,*) "nprop",nprop !d
504 if (it.gt.0 .and. nprop.ne.nprop_prev) then
505 write (iout,*) "Warning previous nprop",nprop_prev,&
512 if (umbrella(iparm) .or. read_iset(iparm) .or. hamil_rep) &
513 call xdrfint(ixdrf, iset, iret)
515 call xdrffloat(ixdrf, rprop(i), iret)
520 itraj=mod(it,totraj(iR,iparm))
523 write (iout,*) "ii",ii," itraj",itraj," it",it
525 if (iset.eq.0) iset = 1
528 if (itraj.gt.ntraj) ntraj=itraj
529 nstep(itraj)=nstep(itraj)+1
530 ! rprop(2)=dsqrt(rprop(2))
531 ! rprop(3)=dsqrt(rprop(3))
533 write (iout,*) "umbrella ",umbrella
534 write (iout,*) rtime,rpotE,rt_bath,nss,&
535 (ihpb(j),jhpb(j),j=1,nss),(rprop(j),j=1,nprop)
536 write (iout,*) "nprop",nprop," iset",iset," myparm",myparm
541 #if (defined(AIX) && !defined(JUBL))
542 call xdrf3dfcoord_(ixdrf, xoord, itmp, prec, iret)
544 print *,"before xdrf3dcoord"
545 call xdrf3dfcoord(ixdrf, xoord, itmp, prec, iret)
548 write (iout,'(10f8.3)') ((xoord(j,i),j=1,3),i=1,2*nres+2)
552 if (itmp .ne. nres + nct - nnt + 1) then
553 write (iout,*) "Error: inconsistent sizes",itmp,nres+nct-nnt+1
559 ! write (iout,*) "calling slice" !d
561 islice=slice(nstep(itraj),time,is,ie,ts,te)
562 ! write (iout,*) "islice",islice !d
572 c(j,i+nres+nnt-1)=xoord(j,i+nres)
576 if (islice.gt.0 .and. islice.le.nslice .and. (.not.separate_parset &
577 .or. iset.eq.myparm)) then
579 kk(islice)=kk(islice)+1
580 mm(islice)=mm(islice)+1
581 if (mod(nstep(itraj),isampl(iparm)).eq.0 .and. &
582 conf_check(ll(islice)+1,1)) then
583 if (replica(iparm)) then
584 rt_bath=1.0d0/(rt_bath*1.987D-3)
586 if (abs(real(beta_h(i,iparm))-rt_bath).lt.1.0e-4) then
592 if (i.gt.nT_h(iparm)) then
593 write (iout,*) "Error - temperature of conformation",&
594 ii,1.0d0/(rt_bath*1.987D-3),&
595 " does not match any of the list"
597 1.0d0/(rt_bath*1.987D-3),&
598 (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
601 ! call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
603 kk(islice)=kk(islice)-1
604 mm(islice)=mm(islice)-1
612 jj(islice)=jj(islice)+1
613 if (umbrella(iparm)) then
614 snk(iset,iib,iparm,islice)=snk(iset,iib,iparm,islice)+1
615 else if (hamil_rep) then
616 snk(1,iib,iparm,islice)=snk(1,iib,iparm,islice)+1
618 snk(iR,iib,iparm,islice)=snk(iR,iib,iparm,islice)+1
620 ll(islice)=ll(islice)+1
622 write (iout,*) "Writing conformation, record",ll(islice)
623 write (iout,*) "ib",ib," iib",iib
624 write (iout,*) "ntraj",ntraj," itraj",itraj,&
625 " nstep",nstep(itraj)
626 write (iout,*) "pote",rpotE," time",rtime
627 ! if (replica(iparm)) then
628 ! write (iout,*) "TEMP",1.0d0/(rt_bath*1.987D-3)
629 ! write (iout,*) "TEMP list"
631 ! & (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
633 write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
634 ! write (iout,*) "nres",nres," nnt",nnt," nct",nct," nss",nss
635 ! write (iout,*) "length",nres*4+(nct-nnt+1)*4+4+2*nss*4
638 if (islice.ne.islice1) then
639 ! write (iout,*) "islice",islice," islice1",islice1
641 ! write (iout,*) "Closing file ",
642 ! & bprotfile_temp(:ilen(bprotfile_temp))
643 call opentmp(islice,ientout,bprotfile_temp)
644 ! write (iout,*) "Opening file ",
645 ! & bprotfile_temp(:ilen(bprotfile_temp))
648 if (umbrella(iparm)) then
649 write(ientout,rec=ll(islice)) &
650 ((xoord(l,k),l=1,3),k=1,nres),&
651 ((xoord(l,k),l=1,3),k=nres+1,nres+nct-nnt+1),&
652 nss,(ihpb(k),jhpb(k),k=1,nss),&
653 rpotE+0.0d0,efree,rmsdev,(rprop(i)+0.0d0,i=1,nQ),&
655 else if (hamil_rep) then
656 write(ientout,rec=ll(islice)) &
657 ((xoord(l,k),l=1,3),k=1,nres),&
658 ((xoord(l,k),l=1,3),k=nres+1,nres+nct-nnt+1),&
659 nss,(ihpb(k),jhpb(k),k=1,nss),&
660 rpotE+0.0d0,efree,rmsdev,(rprop(i)+0.0d0,i=1,nQ),&
663 write(ientout,rec=ll(islice)) &
664 ((xoord(l,k),l=1,3),k=1,nres),&
665 ((xoord(l,k),l=1,3),k=nres+1,nres+nct-nnt+1),&
666 nss,(ihpb(k),jhpb(k),k=1,nss),&
667 rpotE+0.0d0,efree,rmsdev,(rprop(i)+0.0d0,i=1,nQ),&
671 call int_from_cart1(.false.)
672 write (iout,*) "Writing conformation, record",ll(islice)
673 write (iout,*) "Cartesian coordinates"
674 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
675 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
676 write (iout,*) "Internal coordinates"
677 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
678 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
679 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
680 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
681 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
682 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
683 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
684 ! write (iout,'(8f10.5)') (rprop(j),j=1,nQ)
685 write (iout,'(16i5)') iscor
695 #if (defined(AIX) && !defined(JUBL))
696 call xdrfclose_(ixdrf, iret)
698 call xdrfclose(ixdrf, iret)
700 write (iout,'(i10," trajectories found in file.")') ntraj+1
701 write (iout,'(a)') "Numbers of steps in trajectories:"
702 write (iout,'(8i10)') (nstep(i),i=0,ntraj)
703 write (iout,*) ii," conformations read from file",&
706 write (iout,*) mm(islice)," conformations read so far, slice",&
708 write (iout,*) ll(islice),&
709 " conformations stored so far, slice",islice
712 print *,"before cxread return"
715 end subroutine cxread
716 !--------------------------------------------------------------------------------
718 !--------------------------------------------------------------------------------
719 subroutine xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
723 ! include "DIMENSIONS"
724 ! include "DIMENSIONS.ZSCOPT"
725 ! include "DIMENSIONS.FREE"
726 use MPI_data, only:nprocs
729 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
730 ! include "COMMON.MPI"
732 integer,parameter :: MaxTraj=2050
733 ! include "COMMON.CHAIN"
734 ! include "COMMON.IOUNITS"
735 ! include "COMMON.PROTFILES"
736 ! include "COMMON.NAMES"
737 ! include "COMMON.VAR"
738 ! include "COMMON.GEO"
739 ! include "COMMON.ENEPS"
740 ! include "COMMON.PROT"
741 ! include "COMMON.INTERACT"
742 ! include "COMMON.FREE"
743 ! include "COMMON.SBRIDGE"
744 ! include "COMMON.OBCINKA"
745 real(kind=4) :: csingle(3,nres*2)
746 character(len=64) :: nazwa,bprotfile_temp
747 integer :: i,j,k,l,ii,jj(nslice),kk(nslice),ll(nslice),&
748 mm(nslice) !(maxslice)
749 integer :: iscor,islice,islice1 !el,slice
750 real(kind=8) :: energ
752 ! external ilen,iroof
753 real(kind=8) :: rmsdev,energia(0:n_ene),efree,eini,temp
754 !el real(kind=8) :: rmsdev,energia(0:max_eneW),efree,eini,temp
755 real(kind=8) :: prop(nQ) !(maxQ)
756 integer :: ntot_all(0:nprocs-1)!(0:maxprocs-1)
757 integer :: iparm,ib,iib,ir,nprop,nthr
758 real(kind=8) :: etot,time,ts(nslice),te(nslice)
759 integer :: is(nslice),ie(nslice),itraj,ntraj,it,iset
760 integer :: nstep(0:MaxTraj-1)
763 call set_slices(is,ie,ts,te,iR,ib,iparm)
773 call opentmp(islice1,ientout,bprotfile_temp)
775 if (replica(iparm)) then
776 if (hamil_rep .or. umbrella(iparm)) then
777 read (ientin,*,end=1112,err=1112) time,eini,&
778 etot,temp,nss,(ihpb(j),jhpb(j),j=1,nss),&
779 nprop,(prop(j),j=1,nprop),iset
781 read (ientin,*,end=1112,err=1112) time,eini,&
782 etot,temp,nss,(ihpb(j),jhpb(j),j=1,nss),&
783 nprop,(prop(j),j=1,nprop)
785 temp=1.0d0/(temp*1.987D-3)
786 ! write (iout,*) time,eini,etot,nss,
787 ! & (ihpb(j),jhpb(j),j=1,nss),(prop(j),j=1,nprop)
790 if (beta_h(i,iparm).eq.temp) then
796 if (i.gt.nT_h(iparm)) then
797 write (iout,*) "Error - temperature of conformation",&
798 ii,1.0d0/(temp*1.987D-3),&
799 " does not match any of the list"
801 1.0d0/(temp*1.987D-3),&
802 (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
805 call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
809 read (ientin,*,end=1112,err=1112) time,eini,&
810 etot,nss,(ihpb(j),jhpb(j),j=1,nss),&
811 nprop,(prop(j),j=1,nprop)
814 itraj=mod(it,totraj(iR,iparm))
815 ! write (*,*) "ii",ii," itraj",itraj
818 if (itraj.gt.ntraj) ntraj=itraj
819 nstep(itraj)=nstep(itraj)+1
820 islice=slice(nstep(itraj),time,is,ie,ts,te)
821 read (ientin,'(8f10.5)',end=1112,err=1112) &
822 ((csingle(l,k),l=1,3),k=1,nres),&
823 ((csingle(l,k+nres),l=1,3),k=nnt,nct)
825 if (islice.gt.0 .and. islice.le.nslice) then
827 kk(islice)=kk(islice)+1
828 mm(islice)=mm(islice)+1
829 if (mod(nstep(itraj),isampl(iparm)).eq.0) then
830 jj(islice)=jj(islice)+1
832 snk(iR,iib,iset,islice)=snk(iR,iib,iset,islice)+1
833 else if (umbrella(iparm)) then
834 snk(iset,iib,iparm,islice)=snk(iset,iib,iparm,islice)+1
836 snk(iR,iib,iparm,islice)=snk(iR,iib,iparm,islice)+1
838 ll(islice)=ll(islice)+1
839 ! write (iout,*) ii,kk,jj,ll,eini,(prop(j),j=1,nprop)
841 ! write (iout,*) "Writing conformation, record",ll(islice)
842 ! write (iout,*) "ib",ib," iib",iib
843 if (replica(iparm)) then
844 write (iout,*) "TEMP",1.0d0/(temp*1.987D-3)
845 write (iout,*) "TEMP list"
847 (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
851 ! write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
852 ! write (iout,*) "nres",nres," nnt",nnt," nct",nct," nss",nss
853 ! write (iout,*) "length",nres*4+(nct-nnt+1)*4+4+2*nss*4
855 if (islice.ne.islice1) then
856 ! write (iout,*) "islice",islice," islice1",islice1
858 ! write (iout,*) "Closing file ",
859 ! & bprotfile_temp(:ilen(bprotfile_temp))
860 call opentmp(islice,ientout,bprotfile_temp)
861 ! write (iout,*) "Opening file ",
862 ! & bprotfile_temp(:ilen(bprotfile_temp))
866 write(ientout,rec=ll(islice)) &
867 ((csingle(l,k),l=1,3),k=1,nres),&
868 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
869 nss,(ihpb(k),jhpb(k),k=1,nss),&
870 eini,efree,rmsdev,(prop(i),i=1,nQ),iR,iib,iparm
877 call int_from_cart1(.false.)
878 write (iout,*) "Writing conformation, record",ll(islice)
879 write (iout,*) "Cartesian coordinates"
880 write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
881 write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
882 write (iout,*) "Internal coordinates"
883 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
884 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
885 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
886 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
887 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
888 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
889 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
890 ! write (iout,'(8f10.5)') (prop(j),j=1,nQ)
891 write (iout,'(16i5)') iscor
899 write (iout,'(i10," trajectories found in file.")') ntraj+1
900 write (iout,'(a)') "Numbers of steps in trajectories:"
901 write (iout,'(8i10)') (nstep(i),i=0,ntraj)
902 write (iout,*) ii," conformations read from file",&
904 write (iout,*) mm(islice)," conformations read so far, slice",&
906 write (iout,*) ll(islice)," conformations stored so far, slice",&
911 !--------------------------------------------------------------------------------
913 !--------------------------------------------------------------------------------
914 subroutine write_dbase(islice,*)
917 use control_data, only:indpdb
919 use conform_compar, only:conf_compar,rmsnat,qwolynes
920 use energy_data, only:icont,ncont,nnt,nct,maxcont!,&
922 ! include "DIMENSIONS"
923 ! include "DIMENSIONS.ZSCOPT"
924 ! include "DIMENSIONS.FREE"
925 ! include "DIMENSIONS.COMPAR"
926 use geometry, only:int_from_cart1
929 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
930 ! include "COMMON.MPI"
932 ! include "COMMON.CONTROL"
933 ! include "COMMON.CHAIN"
934 ! include "COMMON.IOUNITS"
935 ! include "COMMON.PROTFILES"
936 ! include "COMMON.NAMES"
937 ! include "COMMON.VAR"
938 ! include "COMMON.SBRIDGE"
939 ! include "COMMON.GEO"
940 ! include "COMMON.FFIELD"
941 ! include "COMMON.ENEPS"
942 ! include "COMMON.LOCAL"
943 ! include "COMMON.WEIGHTS"
944 ! include "COMMON.INTERACT"
945 ! include "COMMON.FREE"
946 ! include "COMMON.ENERGIES"
947 ! include "COMMON.COMPAR"
948 ! include "COMMON.PROT"
949 ! include "COMMON.CONTACTS1"
950 character(len=64) :: nazwa
951 character(len=80) :: bxname,cxname
952 character(len=64) :: bprotfile_temp
953 character(len=3) :: liczba,licz
954 character(len=2) :: licz2
955 integer :: i,itj,ii,iii,j,k,l
956 integer :: ixdrf,iret
957 integer :: iscor,islice
958 real(kind=8) :: rmsdev,efree,eini,qnat2
959 real(kind=4) :: csingle(3,nres*2)
960 real(kind=8) :: energ
963 ! external ilen,iroof
964 integer :: ir,ib,iparm
965 integer :: isecstr(nres)
967 write (licz2,'(bz,i2.2)') islice
968 call opentmp(islice,ientout,bprotfile_temp)
969 write (iout,*) "bprotfile_temp ",bprotfile_temp
971 if (.not.bxfile .and. .not. cxfile .and. indpdb.eq.0 &
972 .and. ensembles.eq.0) then
973 close(ientout,status="delete")
977 write (liczba,'(bz,i3.3)') me
978 if (bxfile .or. cxfile .or. ensembles.gt.0) then
979 if (.not.separate_parset) then
980 bxname = prefix(:ilen(prefix))//liczba//".bx"
982 write (licz,'(bz,i3.3)') myparm
983 bxname = prefix(:ilen(prefix))//liczba//"_par"//licz//".bx"
986 open (ientin,file=bxname,status="unknown",&
987 form="unformatted",access="direct",recl=lenrec1)
990 if (bxfile .or. cxfile .or. ensembles.gt.0) then
991 if (nslice.eq.1) then
992 bxname = prefix(:ilen(prefix))//".bx"
994 bxname = prefix(:ilen(prefix))// &
995 "_slice_"//licz2//".bx"
997 open (ientin,file=bxname,status="unknown",&
998 form="unformatted",access="direct",recl=lenrec1)
999 write (iout,*) "Calculating energies; writing geometry",&
1000 " and energy components to ",bxname(:ilen(bxname))
1002 #if (defined(AIX) && !defined(JUBL))
1003 call xdrfopen_(ixdrf,cxname, "w", iret)
1005 call xdrfopen(ixdrf,cxname, "w", iret)
1008 write (iout,*) "Error opening cxfile ",cxname(:ilen(cxname))
1014 if (indpdb.gt.0) then
1015 if (nslice.eq.1) then
1017 if (.not.separate_parset) then
1018 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot)) &
1021 write (licz,'(bz,i3.3)') myparm
1022 statname=prefix(:ilen(prefix))//'_par'//licz//'_'// &
1023 pot(:ilen(pot))//liczba//'.stat'
1027 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))//'.stat'
1031 if (.not.separate_parset) then
1032 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))// &
1033 "_slice_"//licz2//liczba//'.stat'
1035 write (licz,'(bz,i3.3)') myparm
1036 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))// &
1037 '_par'//licz//"_slice_"//licz2//liczba//'.stat'
1040 statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot)) &
1041 //"_slice_"//licz2//'.stat'
1044 print *,istat,statname
1045 open(istat,file=statname,status="unknown")
1047 print *,"Tu dochodze"
1054 print *,"before ientout read"
1055 read(ientout,rec=i,err=101) &
1056 ((csingle(l,k),l=1,3),k=1,nres),&
1057 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
1058 nss,(ihpb(k),jhpb(k),k=1,nss),&
1059 eini,efree,rmsdev,(q(k,i),k=1,nQ),iR,ib,iparm
1060 ! write (iout,*) iR,ib,iparm,eini,efree
1066 call int_from_cart1(.false.)
1068 ! write (iout,*) "Calling conf_compar",i
1070 anatemp= 1.0d0/(beta_h(ib,iparm)*1.987D-3)
1071 print *,"before conf_compar"
1072 if (indpdb.gt.0) then
1073 print *,"just before conf_compar",i
1074 print *,icont,ncont,nnt,nct,"maxcont",maxcont
1076 ! call conf_compar(i,.false.,.true.)
1077 ! call conf_compar(i)
1081 print *,"just after conf_compar"
1083 ! call elecont(.false.,ncont,icont,nnt,nct)
1084 ! call secondary2(.false.,.false.,ncont,icont,isecstr)
1086 ! write (iout,*) "Exit conf_compar",i
1088 print *,"before ientin"
1089 if (bxfile .or.cxfile .or. ensembles.gt.0) write (ientin,rec=i) &
1090 ((csingle(l,k),l=1,3),k=1,nres),&
1091 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
1092 nss,(ihpb(k),jhpb(k),k=1,nss),&
1093 ! & potE(i,iparm),-entfac(i),rms_nat,iscore
1094 potE(i,nparmset),-entfac(i),rms_nat,iscore
1095 ! write (iout,'(2i5,3e15.5)') i,me,potE(i,iparm),-entfac(i)
1097 if (cxfile) call cxwrite(ixdrf,csingle,potE(i,nparmset),&
1098 -entfac(i),rms_nat,iscore)
1101 close(ientout,status="delete")
1103 if (bxfile .or. cxfile .or. ensembles.gt.0) close(ientin)
1105 print *,"before MPI_barrier"
1106 call MPI_Barrier(WHAM_COMM,IERROR)
1107 if (me.ne.Master .or. .not.bxfile .and. .not. cxfile &
1108 .and. ensembles.eq.0) return
1110 if (bxfile .or. ensembles.gt.0) then
1111 if (nslice.eq.1) then
1112 if (.not.separate_parset) then
1113 bxname = prefix(:ilen(prefix))//".bx"
1115 write (licz,'(bz,i3.3)') myparm
1116 bxname = prefix(:ilen(prefix))//"_par"//licz//".bx"
1119 if (.not.separate_parset) then
1120 bxname = prefix(:ilen(prefix))//"_slice_"//licz2//".bx"
1122 write (licz,'(bz,i3.3)') myparm
1123 bxname = prefix(:ilen(prefix))//"par_"//licz// &
1124 "_slice_"//licz2//".bx"
1127 open (ientout,file=bxname,status="unknown",&
1128 form="unformatted",access="direct",recl=lenrec1)
1129 write (iout,*) "Master is creating binary database ",&
1130 bxname(:ilen(bxname))
1133 if (nslice.eq.1) then
1134 if (.not.separate_parset) then
1135 cxname = prefix(:ilen(prefix))//".cx"
1137 cxname = prefix(:ilen(prefix))//"_par"//licz//".cx"
1140 if (.not.separate_parset) then
1141 cxname = prefix(:ilen(prefix))// &
1142 "_slice_"//licz2//".cx"
1144 cxname = prefix(:ilen(prefix))//"_par"//licz// &
1145 "_slice_"//licz2//".cx"
1148 #if (defined(AIX) && !defined(JUBL))
1149 call xdrfopen_(ixdrf,cxname, "w", iret)
1151 call xdrfopen(ixdrf,cxname, "w", iret)
1154 write (iout,*) "Error opening cxfile ",cxname(:ilen(cxname))
1159 write (liczba,'(bz,i3.3)') j
1160 if (separate_parset) then
1161 write (licz,'(bz,i3.3)') myparm
1162 bxname = prefix(:ilen(prefix))//liczba//"_par"//licz//".bx"
1164 bxname = prefix(:ilen(prefix))//liczba//".bx"
1166 open (ientin,file=bxname,status="unknown",&
1167 form="unformatted",access="direct",recl=lenrec1)
1168 write (iout,*) "Master is reading conformations from ",&
1169 bxname(:ilen(bxname))
1171 ! write (iout,*) "j",j," indstart",indstart(j)," indend",indend(j)
1173 do i=indstart(j),indend(j)
1175 read(ientin,rec=iii,err=101) &
1176 ((csingle(l,k),l=1,3),k=1,nres),&
1177 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
1178 nss,(ihpb(k),jhpb(k),k=1,nss),&
1179 eini,efree,rmsdev,iscor
1180 if (bxfile .or. ensembles.gt.0) then
1181 write (ientout,rec=i) &
1182 ((csingle(l,k),l=1,3),k=1,nres),&
1183 ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
1184 nss,(ihpb(k),jhpb(k),k=1,nss),&
1185 eini,efree,rmsdev,iscor
1187 print *,"before cxwrite"
1188 if(cxfile)call cxwrite(ixdrf,csingle,eini,efree,rmsdev,iscor)
1189 print *,"after cxwrite"
1196 call int_from_cart1(.false.)
1197 write (iout,'(2i5,3e15.5)') i,iii,eini,efree
1198 write (iout,*) "The Cartesian geometry is:"
1199 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1200 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1201 write (iout,*) "The internal geometry is:"
1202 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
1203 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
1204 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
1205 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
1206 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
1207 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
1208 write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
1209 write (iout,'(f10.5,i5)') rmsdev,iscor
1212 write (iout,*) iii," conformations (from",indstart(j)," to",&
1213 indend(j),") read from ",&
1214 bxname(:ilen(bxname))
1215 close (ientin,status="delete")
1217 if (bxfile .or. cxfile .or. ensembles.gt.0) close (ientout)
1218 #if (defined(AIX) && !defined(JUBL))
1219 if (cxfile) call xdrfclose_(ixdrf,cxname,iret)
1221 if (cxfile) call xdrfclose(ixdrf,cxname,iret)
1225 101 write (iout,*) "Error in scratchfile."
1228 end subroutine write_dbase
1229 !-------------------------------------------------------------------------------
1230 subroutine cxwrite(ixdrf,csingle,eini,efree,rmsdev,iscor)
1232 ! include "DIMENSIONS"
1233 ! include "DIMENSIONS.ZSCOPT"
1234 ! include "DIMENSIONS.FREE"
1235 ! include "DIMENSIONS.COMPAR"
1238 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
1239 ! include "COMMON.MPI"
1241 ! include "COMMON.CONTROL"
1242 ! include "COMMON.CHAIN"
1243 ! include "COMMON.IOUNITS"
1244 ! include "COMMON.PROTFILES"
1245 ! include "COMMON.NAMES"
1246 ! include "COMMON.VAR"
1247 ! include "COMMON.SBRIDGE"
1248 ! include "COMMON.GEO"
1249 ! include "COMMON.FFIELD"
1250 ! include "COMMON.ENEPS"
1251 ! include "COMMON.LOCAL"
1252 ! include "COMMON.WEIGHTS"
1253 ! include "COMMON.INTERACT"
1254 ! include "COMMON.FREE"
1255 ! include "COMMON.ENERGIES"
1256 ! include "COMMON.COMPAR"
1257 ! include "COMMON.PROT"
1258 integer :: i,j,itmp,iscor,iret,ixdrf
1259 real(kind=8) :: rmsdev,efree,eini
1260 real(kind=4) :: csingle(3,nres*2),xoord(3,2*nres+2)
1261 real(kind=4) :: prec
1263 ! write (iout,*) "cxwrite"
1268 xoord(j,i)=csingle(j,i)
1273 xoord(j,nres+i-nnt+1)=csingle(j,i+nres)
1279 ! write (iout,*) "itmp",itmp
1281 #if (defined(AIX) && !defined(JUBL))
1282 call xdrf3dfcoord_(ixdrf, xoord, itmp, prec, iret)
1284 ! write (iout,*) "xdrf3dfcoord"
1286 call xdrfint_(ixdrf, nss, iret)
1287 write (iout,*) "iret",iret
1288 write (iout,*) "nss",nss,i,"TUTU"
1290 call xdrfint_(ixdrf, ihpb(j), iret)
1291 call xdrfint_(ixdrf, jhpb(j), iret)
1292 write(iout,*), ihpb(j),jhpb(j),"TUTU"
1294 call xdrffloat_(ixdrf,real(eini),iret)
1295 call xdrffloat_(ixdrf,real(efree),iret)
1296 write(iout,*) "TUTU", eini
1297 write(iout,*) "TUTU", efree
1298 call xdrffloat_(ixdrf,real(rmsdev),iret)
1299 call xdrfint_(ixdrf,iscor,iret)
1301 call xdrf3dfcoord(ixdrf, xoord, itmp, prec, iret)
1302 write (iout,*) "iret",iret
1303 write (iout,*) "nss",nss,i,"TUTU"
1305 call xdrfint(ixdrf, nss, iret)
1307 call xdrfint(ixdrf, ihpb(j), iret)
1308 call xdrfint(ixdrf, jhpb(j), iret)
1309 write(iout,*), ihpb(j),jhpb(j),"TUTU"
1311 call xdrffloat(ixdrf,real(eini),iret)
1312 call xdrffloat(ixdrf,real(efree),iret)
1313 write(iout,*) "TUTU", eini
1314 write(iout,*) "TUTU", efree
1315 call xdrffloat(ixdrf,real(rmsdev),iret)
1316 call xdrfint(ixdrf,iscor,iret)
1320 end subroutine cxwrite
1321 !-------------------------------------------------------------------------------
1323 !-------------------------------------------------------------------------------
1324 subroutine set_slices(is,ie,ts,te,iR,ib,iparm)
1326 ! include 'DIMENSIONS'
1327 ! include 'DIMENSIONS.ZSCOPT'
1328 ! include 'DIMENSIONS.FREE'
1329 ! include 'COMMON.IOUNITS'
1330 ! include 'COMMON.PROTFILES'
1331 ! include 'COMMON.OBCINKA'
1332 ! include 'COMMON.PROT'
1333 integer :: islice,iR,ib,iparm
1334 integer :: is(MaxSlice),ie(MaxSlice),nrec_slice
1335 real(kind=8) :: ts(MaxSlice),te(MaxSlice),time_slice
1338 if (time_end_collect(iR,ib,iparm).ge.1.0d10) then
1339 ts(islice)=time_start_collect(iR,ib,iparm)
1340 te(islice)=time_end_collect(iR,ib,iparm)
1341 nrec_slice=(rec_end(iR,ib,iparm)- &
1342 rec_start(iR,ib,iparm)+1)/nslice
1343 is(islice)=rec_start(iR,ib,iparm)+(islice-1)*nrec_slice
1344 ie(islice)=rec_start(iR,ib,iparm)+islice*nrec_slice-1
1346 time_slice=(time_end_collect(iR,ib,iparm) &
1347 -time_start_collect(iR,ib,iparm))/nslice
1348 ts(islice)=time_start_collect(iR,ib,iparm)+(islice-1)* &
1350 te(islice)=time_start_collect(iR,ib,iparm)+islice*time_slice
1351 is(islice)=rec_start(iR,ib,iparm)
1352 ie(islice)=rec_end(iR,ib,iparm)
1356 write (iout,*) "nrec_slice",nrec_slice," time_slice",time_slice
1357 write (iout,*) "is",(is(islice),islice=1,nslice)
1358 write (iout,*) "ie",(ie(islice),islice=1,nslice)
1359 write (iout,*) "rec_start",&
1360 rec_start(iR,ib,iparm)," rec_end",rec_end(iR,ib,iparm)
1361 write (iout,*) "ts",(ts(islice),islice=1,nslice)
1362 write (iout,*) "te",(te(islice),islice=1,nslice)
1363 write (iout,*) "time_start",&
1364 time_start_collect(iR,ib,iparm)," time_end",&
1365 time_end_collect(iR,ib,iparm)
1369 end subroutine set_slices
1370 !-----------------------------------------------------------------------------
1371 integer function slice(irecord,time,is,ie,ts,te)
1373 ! include 'DIMENSIONS'
1374 ! include 'DIMENSIONS.ZSCOPT'
1375 ! include 'DIMENSIONS.FREE'
1376 ! include 'COMMON.IOUNITS'
1377 ! include 'COMMON.PROTFILES'
1378 ! include 'COMMON.OBCINKA'
1379 ! include 'COMMON.PROT'
1380 integer :: is(MaxSlice),ie(MaxSlice),nrec_slice
1381 real(kind=8) :: ts(MaxSlice),te(MaxSlice),time_slice
1382 integer :: i,ii,irecord
1383 real(kind=8) :: time
1385 ! write (iout,*) "within slice nslice",nslice
1387 if (irecord.lt.is(1) .or. time.lt.ts(1)) then
1391 do while (ii.le.nslice .and. &
1392 (irecord.lt.is(ii) .or. irecord.gt.ie(ii) .or. &
1393 time.lt.ts(ii) .or. time.gt.te(ii)) )
1394 ! write (iout,*) "ii",ii,time,ts(ii)
1399 ! write (iout,*) "end: ii",ii
1404 !-----------------------------------------------------------------------------
1406 !-----------------------------------------------------------------------------
1407 logical function conf_check(ii,iprint)
1409 use names, only:ntyp1
1411 use energy_data, only:itype,dsc,molnum
1412 use geometry, only:int_from_cart1
1414 ! include "DIMENSIONS"
1415 ! include "DIMENSIONS.ZSCOPT"
1416 ! include "DIMENSIONS.FREE"
1420 ! include "COMMON.MPI"
1422 ! include "COMMON.CHAIN"
1423 ! include "COMMON.IOUNITS"
1424 ! include "COMMON.PROTFILES"
1425 ! include "COMMON.NAMES"
1426 ! include "COMMON.VAR"
1427 ! include "COMMON.SBRIDGE"
1428 ! include "COMMON.GEO"
1429 ! include "COMMON.FFIELD"
1430 ! include "COMMON.ENEPS"
1431 ! include "COMMON.LOCAL"
1432 ! include "COMMON.WEIGHTS"
1433 ! include "COMMON.INTERACT"
1434 ! include "COMMON.FREE"
1435 ! include "COMMON.ENERGIES"
1436 ! include "COMMON.CONTROL"
1437 ! include "COMMON.TORCNSTR"
1441 integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
1443 integer :: j,k,l,ii,itj,iprint,mnum
1444 if (.not. check_conf) then
1448 call int_from_cart1(.false.)
1451 if (mnum.ne.1) cycle
1452 if (mnum.eq.5) cycle
1453 if (itype(j-1,mnum).ne.ntyp1 .and. itype(j,mnum).ne.ntyp1_molec(mnum) .and. &
1454 (vbld(j).lt.2.0d0 .or. vbld(j).gt.5.0d0)) then
1456 write (iout,*) "Bad CA-CA bond length",j," ",vbld(j),&
1457 " for conformation",ii,mnum
1458 if (iprint.gt.1) then
1459 write (iout,*) "The Cartesian geometry is:"
1460 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1461 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1462 write (iout,*) "The internal geometry is:"
1463 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
1464 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
1465 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
1466 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
1467 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
1468 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
1470 if (iprint.gt.0) write (iout,*) &
1471 "This conformation WILL NOT be added to the database."
1478 if (mnum.ne.1) cycle
1480 if (itype(j,mnum).ne.10 .and.itype(j,mnum).ne.ntyp1_molec(mnum) .and. &
1481 (vbld(nres+j)-dsc(iabs(itj))).gt.2.0d0) then
1483 write (iout,*) "Bad CA-SC bond length",j," ",vbld(nres+j),&
1484 " for conformation",ii
1485 if (iprint.gt.1) then
1486 write (iout,*) "The Cartesian geometry is:"
1487 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1488 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1489 write (iout,*) "The internal geometry is:"
1490 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
1491 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
1492 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
1493 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
1494 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
1495 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
1497 if (iprint.gt.0) write (iout,*) &
1498 "This conformation WILL NOT be added to the database."
1504 if (theta(j).le.0.0d0) then
1506 write (iout,*) "Zero theta angle(s) in conformation",ii
1507 if (iprint.gt.1) then
1508 write (iout,*) "The Cartesian geometry is:"
1509 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
1510 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
1511 write (iout,*) "The internal geometry is:"
1512 write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
1513 write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
1514 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
1515 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
1516 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
1517 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
1519 if (iprint.gt.0) write (iout,*) &
1520 "This conformation WILL NOT be added to the database."
1524 if (theta(j).gt.179.97*deg2rad) theta(j)=179.97*deg2rad
1527 ! write (iout,*) "conf_check passed",ii
1529 end function conf_check
1530 !-----------------------------------------------------------------------------
1531 end module io_database