--- /dev/null
+ module io_database
+!-----------------------------------------------------------------------------
+ use names
+ use wham_data
+ use io_units
+ use io_base, only:ilen
+ use energy_data, only:nnt,nct,nss,ihpb,jhpb,iset
+ use geometry_data, only:nres,c
+#ifdef MPI
+ use MPI_data
+! include "COMMON.MPI"
+#endif
+
+ implicit none
+!-----------------------------------------------------------------------------
+!
+!
+!-----------------------------------------------------------------------------
+ contains
+!-----------------------------------------------------------------------------
+! readrtns.F
+!-------------------------------------------------------------------------------
+ subroutine opentmp(islice,iunit,bprotfile_temp)
+! implicit none
+! include "DIMENSIONS"
+! include "DIMENSIONS.ZSCOPT"
+! include "DIMENSIONS.FREE"
+! use MPI_data, only:me
+#ifdef MPI
+ include "mpif.h"
+ integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+! include "COMMON.MPI"
+#endif
+! include "COMMON.IOUNITS"
+! include "COMMON.PROTFILES"
+! include "COMMON.PROT"
+! include "COMMON.FREE"
+ character(len=64) :: bprotfile_temp
+ character(len=3) :: liczba,liczba2
+ character(len=2) :: liczba1
+ integer :: iunit,islice
+! integer ilen,iroof
+! external ilen,iroof
+! logical :: lerr
+! integer :: lenrec,lenrec2
+
+!el
+! lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
+! lenrec=lenrec2+8
+ write (liczba1,'(bz,i2.2)') islice
+#ifdef MPI
+ write (liczba,'(bz,i3.3)') me
+!#ifdef MPI
+! write (iout,*) "separate_parset ",separate_parset,
+! & " myparm",myparm
+ if (separate_parset) then
+ write (liczba2,'(bz,i3.3)') myparm
+ bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
+ prefix(:ilen(prefix))//liczba//"_"//liczba2//".xbin.tmp"//liczba1
+ open (iunit,file=bprotfile_temp,status="unknown",&
+ form="unformatted",access="direct",recl=lenrec)
+ else
+ bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
+ prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
+ open (iunit,file=bprotfile_temp,status="unknown",&
+ form="unformatted",access="direct",recl=lenrec)
+ endif
+#else
+ bprotfile_temp = scratchdir(:ilen(scratchdir))// &
+ "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
+ open (iunit,file=bprotfile_temp,status="unknown",&
+ form="unformatted",access="direct",recl=lenrec)
+#endif
+! write (iout,*) "OpenTmp iunit",iunit," bprotfile_temp",
+! & bprotfile_temp
+! call flush(iout)
+ return
+ end subroutine opentmp
+!-------------------------------------------------------------------------------
+ subroutine read_database(*)
+
+! use energy_data, only:nct,nnt,nss
+! implicit none
+! include "DIMENSIONS"
+! include "DIMENSIONS.ZSCOPT"
+! include "DIMENSIONS.FREE"
+ use MPI_data, only:me,nprocs
+#ifdef MPI
+ include "mpif.h"
+ integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+! include "COMMON.MPI"
+#endif
+! include "COMMON.CHAIN"
+! include "COMMON.IOUNITS"
+! include "COMMON.PROTFILES"
+! include "COMMON.NAMES"
+! include "COMMON.VAR"
+! include "COMMON.GEO"
+! include "COMMON.ENEPS"
+! include "COMMON.PROT"
+! include "COMMON.INTERACT"
+! include "COMMON.FREE"
+! include "COMMON.SBRIDGE"
+! include "COMMON.OBCINKA"
+ real(kind=4) :: csingle(3,nres*2) !(3,maxres2)
+ character(len=64) :: nazwa,bprotfile_temp
+ character(len=3) :: liczba
+ character(len=2) :: liczba1
+ integer :: i,j,ii,jj(nslice),k,kk(nslice),l,&
+ ll(nslice),mm(nslice),if
+ integer :: nrec,nlines,iscor,iunit,islice
+ real(kind=8) :: energ
+! integer ilen,iroof
+! external ilen,iroof
+ real(kind=8) :: rmsdev,energia(0:n_ene),efree,eini,temp
+!el real(kind=8) :: rmsdev,energia(0:max_eneW),efree,eini,temp
+ real(kind=8) :: prop(nQ) !(maxQ)
+ integer :: ntot_all(nslice,0:nprocs-1)!(maxslice,0:maxprocs-1)
+ integer :: iparm,ib,iib,ir,nprop,nthr,npars
+ real(kind=8) :: etot,time
+ integer :: ixdrf,iret
+ logical :: lerr,linit
+
+ lenrec1=12*(nres+nct-nnt+1)+4*(2*nss+2)+24
+ lenrec2=12*(nres+nct-nnt+1)+4*(2*nss+2)+24+8*nQ
+ lenrec=lenrec2+8
+ write (iout,*) "lenrec",lenrec," lenrec1",lenrec1,&
+ " lenrec2",lenrec2
+
+ do i=1,nQ
+ prop(i)=0.0d0
+ enddo
+ do islice=1,nslice
+ ll(islice)=0
+ mm(islice)=0
+ enddo
+ write (iout,*) "nparmset",nparmset
+ if (hamil_rep) then
+ npars=1
+ else
+ npars=nparmset
+ endif
+ do iparm=1,npars
+
+ if (replica(iparm)) then
+ nthr = 1
+ else
+ nthr = nT_h(iparm)
+ endif
+
+ do ib=1,nthr
+ do iR=1,nRR(ib,iparm)
+
+ write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
+ do islice=1,nslice
+ jj(islice)=0
+ kk(islice)=0
+ enddo
+
+ IF (NFILE_BIN(iR,ib,iparm).GT.0) THEN
+! Read conformations from binary DA files (one per batch) and write them to
+! a binary DA scratchfile.
+ write (liczba,'(bz,i3.3)') me
+ do if=1,nfile_bin(iR,ib,iparm)
+ nazwa=protfiles(if,1,iR,ib,iparm) &
+ (:ilen(protfiles(if,1,iR,ib,iparm)))//".bx"
+ open (ientin,file=nazwa,status="old",form="unformatted",&
+ access="direct",recl=lenrec2,err=1111)
+ ii=0
+ do islice=1,nslice
+ call opentmp(islice,ientout,bprotfile_temp)
+ call bxread(nazwa,islice,ii,jj(islice),kk(islice),ll(islice),&
+ mm(islice),iR,ib,iparm)
+ close(ientout)
+ enddo
+ close(ientin)
+ enddo
+ ENDIF ! NFILE_BIN>0
+!
+ IF (NFILE_ASC(iR,ib,iparm).GT.0) THEN
+! Read conformations from multiple ASCII int files and write them to a binary
+! DA scratchfile.
+ do if=1,nfile_asc(iR,ib,iparm)
+ nazwa=protfiles(if,2,iR,ib,iparm) &
+ (:ilen(protfiles(if,2,iR,ib,iparm)))//".x"
+ open(unit=ientin,file=nazwa,status='old',err=1111)
+ write(iout,*) "reading ",nazwa(:ilen(nazwa))
+ ii=0
+ call xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
+ enddo ! if
+ ENDIF
+ IF (NFILE_CX(iR,ib,iparm).gt.0) THEN
+! Read conformations from cx files and write them to a binary
+! DA scratchfile.
+ do if=1,nfile_cx(iR,ib,iparm)
+ nazwa=protfiles(if,2,iR,ib,iparm) &
+ (:ilen(protfiles(if,2,iR,ib,iparm)))//".cx"
+ write(iout,*) "reading ",nazwa(:ilen(nazwa))
+ ii=0
+ print *,"Calling cxread"
+ call cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,&
+ *1111)
+write(iout,*)"after call cxread"
+ close(ientout)
+ write (iout,*) "exit cxread"
+ call flush(iout)
+ enddo
+ ENDIF
+write(iout,*)"*********************in read database"
+
+ do islice=1,nslice
+! stot(islice)=0
+ stot(islice)=stot(islice)+jj(islice)
+ enddo
+
+ enddo
+ enddo
+ write (iout,*) "IPARM",iparm
+ enddo
+
+ if (nslice.eq.1) then
+#ifdef MPI
+ write (liczba,'(bz,i3.3)') me
+ bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
+ prefix(:ilen(prefix))//liczba//".xbin.tmp"
+#else
+ bprotfile_temp = scratchdir(:ilen(scratchdir))// &
+ "/"//prefix(:ilen(prefix))//".xbin.tmp"
+#endif
+ write(iout,*) mm(1)," conformations read",ll(1),&
+ " conformations written to ",&
+ bprotfile_temp(:ilen(bprotfile_temp))
+ else
+ do islice=1,nslice
+ write (liczba1,'(bz,i2.2)') islice
+#ifdef MPI
+ write (liczba,'(bz,i3.3)') me
+ bprotfile_temp = scratchdir(:ilen(scratchdir))//"/"// &
+ prefix(:ilen(prefix))//liczba//".xbin.tmp"//liczba1
+#else
+ bprotfile_temp = scratchdir(:ilen(scratchdir))// &
+ "/"//prefix(:ilen(prefix))//".xbin.tmp"//liczba1
+#endif
+ write(iout,*) mm(islice)," conformations read",ll(islice),&
+ " conformations written to ",&
+ bprotfile_temp(:ilen(bprotfile_temp))
+ enddo
+ endif
+
+#ifdef MPI
+! Check if everyone has the same number of conformations
+ call MPI_Allgather(stot(1),nslice,MPI_INTEGER,&
+ ntot_all(1,0),nslice,MPI_INTEGER,MPI_Comm_World,IERROR)
+ lerr=.false.
+ do i=0,nprocs-1
+ if (i.ne.me) then
+ do islice=1,nslice
+ if (stot(islice).ne.ntot_all(islice,i)) then
+ write (iout,*) "Number of conformations at processor",i,&
+ " differs from that at processor",me,&
+ stot(islice),ntot_all(islice,i)," slice",islice
+ lerr = .true.
+ endif
+ enddo
+ endif
+ enddo
+ if (lerr) then
+ write (iout,*)
+ write (iout,*) "Numbers of conformations read by processors"
+ write (iout,*)
+ do i=0,nprocs-1
+ write (iout,'(8i10)') i,(ntot_all(islice,i),islice=1,nslice)
+ enddo
+ write (iout,*) "Calculation terminated."
+ call flush(iout)
+ return 1
+ endif
+ do islice=1,nslice
+ ntot(islice)=stot(islice)
+ enddo
+write(iout,*) "end of read database"
+ return
+#endif
+ 1111 write(iout,*) "Error opening coordinate file ",nazwa(:ilen(nazwa))
+ call flush(iout)
+ return 1
+ end subroutine read_database
+!--------------------------------------------------------------------------------
+ integer function iroof(n,m)
+ integer :: n,m,ii
+ ii = n/m
+ if (ii*m .lt. n) ii=ii+1
+ iroof = ii
+ return
+ end function iroof
+!--------------------------------------------------------------------------------
+! bxread.F
+!--------------------------------------------------------------------------------
+ subroutine bxread(nazwa,islice,ii,jj,kk,ll,mm,iR,ib,iparm)
+! implicit none
+! include "DIMENSIONS"
+! include "DIMENSIONS.ZSCOPT"
+! include "DIMENSIONS.FREE"
+! use energy_data, only:nnt,nct,nss,ihpb,jhpbi
+ use MPI_data, only:nprocs
+#ifdef MPI
+ include "mpif.h"
+ integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+! include "COMMON.MPI"
+#endif
+! include "COMMON.CHAIN"
+! include "COMMON.IOUNITS"
+! include "COMMON.PROTFILES"
+! include "COMMON.NAMES"
+! include "COMMON.VAR"
+! include "COMMON.GEO"
+! include "COMMON.ENEPS"
+! include "COMMON.PROT"
+! include "COMMON.INTERACT"
+! include "COMMON.FREE"
+! include "COMMON.SBRIDGE"
+ real(kind=4) :: csingle(3,nres*2) !(3,maxres2)
+ character(len=64) :: nazwa,bprotfile_temp
+ character(len=3) :: liczba
+ integer :: i,is,ie,j,ii,jj,k,kk,l,ll,mm,if
+ integer :: nrec,nlines,iscor,islice
+ real(kind=8) :: energ
+! integer ilen,iroof
+! external ilen,iroof
+ real(kind=8) :: rmsdev,energia(0:n_ene),efree,eini,temp
+!el real(kind=8) :: rmsdev,energia(0:max_eneW),efree,eini,temp
+ real(kind=8) :: prop(nQ) !(maxQ)
+ integer :: ntot_all(0:nprocs-1)!(0:maxprocs-1)
+ integer :: iparm,ib,iib,ir,nprop,nthr,nrec_slice
+ real(kind=8) :: etot,time
+ logical :: lerr
+ nrec_slice=(rec_end(iR,ib,iparm)-rec_start(iR,ib,iparm)+1)/nslice
+ is=rec_start(iR,ib,iparm)+(islice-1)*nrec_slice
+ ie=rec_start(iR,ib,iparm)+islice*nrec_slice-1
+ write (iout,*) "bxread: islice",islice," nslice",nslice,&
+ " nrec_slice",nrec_slice
+ write (iout,*) "is",is," ie",ie,"rec_start",&
+ rec_start(iR,ib,iparm)," rec_end",rec_end(iR,ib,iparm)
+ do i=is,ie
+ read(ientin,rec=i+1,err=101) &
+ ((csingle(l,k),l=1,3),k=1,nres),&
+ ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+ eini,efree,rmsdev,(prop(j),j=1,nQ),iscor
+ ii=ii+1
+ kk=kk+1
+ if (mod(kk,isampl(iparm)).eq.0) then
+ jj=jj+1
+ write(ientout,rec=jj) &
+ ((csingle(l,k),l=1,3),k=1,nres),&
+ ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+ eini,efree,rmsdev,(prop(j),j=1,nQ),iR,ib,iparm
+#ifdef DEBUG
+ do i=1,2*nres
+ do j=1,3
+ c(j,i)=csingle(j,i)
+ enddo
+ enddo
+ call int_from_cart1(.false.)
+ write (iout,*) "Writing conformation, record",jj
+ write (iout,*) "Cartesian coordinates"
+ write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
+ write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
+ write (iout,*) "Internal coordinates"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
+ write (iout,'(f10.5,i5)') rmsdev,iscor
+#endif
+ endif
+ enddo
+ 101 continue
+ close(ientin)
+ write (iout,*) ii," conformations read from DA file ",&
+ nazwa(:ilen(nazwa))
+ write (iout,*) kk," conformations read so far, slice",islice
+ write (iout,*) jj," conformations stored so far, slice",islice
+
+ return
+ end subroutine bxread
+!--------------------------------------------------------------------------------
+! cxread.F
+!--------------------------------------------------------------------------------
+ subroutine cxread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm,*)
+
+#define DEBUG
+#ifdef DEBUG
+ use geometry, only:int_from_cart1
+ use geometry_data, only:vbld,rad2deg,theta,phi,alph,omeg
+ integer :: iscor
+#endif
+! implicit real*8 (a-h,o-z)
+! include 'DIMENSIONS'
+! include 'DIMENSIONS.ZSCOPT'
+! include 'DIMENSIONS.FREE'
+ integer,parameter :: MaxTraj=2050
+! include 'COMMON.CHAIN'
+! include 'COMMON.INTERACT'
+! include 'COMMON.NAMES'
+! include 'COMMON.IOUNITS'
+! include 'COMMON.HEADER'
+! include 'COMMON.SBRIDGE'
+! include 'COMMON.PROTFILES'
+! include 'COMMON.OBCINKA'
+! include 'COMMON.FREE'
+! include 'COMMON.VAR'
+! include 'COMMON.GEO'
+! include 'COMMON.PROT'
+ character(len=64) :: nazwa,bprotfile_temp
+ real(kind=4) :: rtime,rpotE,ruconst,rt_bath,rprop(nQ) !(2000) !(maxQ)
+ real(kind=8) :: time
+ integer :: iret,itmp,itraj,ntraj
+ real(kind=4) :: xoord(3,2*nres+2),prec
+ integer :: nstep(0:MaxTraj-1)
+! integer ilen
+! external ilen
+ integer :: ii,jj(nslice),kk(nslice),ll(nslice),mm(nslice) !(maxslice)
+ integer :: is(nSlice),ie(nSlice),nrec_slice
+ real(kind=8) :: ts(nSlice),te(nSlice),time_slice
+ integer :: iR,ib,iparm,i,j,it,islice,nprop_prev
+ integer :: k,l,iib,islice1,nprop
+ real(kind=8) :: efree,rmsdev
+ integer :: ixdrf
+!el integer :: slice
+! logical :: conf_check
+! ixdrf=0
+! nprop=0
+
+! ruconst=0.0d0
+! rtime=0.0d0
+! rpotE=0.0d0
+! rt_bath=0.0d0
+
+ call set_slices(is,ie,ts,te,iR,ib,iparm)
+ nprop_prev=0
+ do i=1,nQ
+ rprop(i)=0.0d0
+ enddo
+ do i=0,MaxTraj-1
+ nstep(i)=0
+ enddo
+ ntraj=0
+ it=0
+ iret=1
+#if (defined(AIX) && !defined(JUBL))
+ call xdrfopen_(ixdrf,nazwa, "r", iret)
+#else
+ call xdrfopen(ixdrf,nazwa, "r", iret)
+#endif
+ if (iret.eq.0) return 1
+
+ islice1=1
+ call opentmp(islice1,ientout,bprotfile_temp)
+ print *,"bumbum" !d
+ do while (iret.gt.0)
+
+#if (defined(AIX) && !defined(JUBL))
+ call xdrffloat_(ixdrf, rtime, iret)
+ print *,"rtime",rtime," iret",iret !d
+ call xdrffloat_(ixdrf, rpotE, iret)
+ write (iout,*) "rpotE",rpotE," iret",iret !d
+ call flush(iout)
+ call xdrffloat_(ixdrf, ruconst, iret)
+ call xdrffloat_(ixdrf, rt_bath, iret)
+ call xdrfint_(ixdrf, nss, iret)
+ do j=1,nss
+ call xdrfint_(ixdrf, ihpb(j), iret)
+ call xdrfint_(ixdrf, jhpb(j), iret)
+ enddo
+ call xdrfint_(ixdrf, nprop, iret)
+ if (umbrella(iparm) .or. read_iset(iparm) .or. hamil_rep) &
+ call xdrfint(ixdrf, iset, iret)
+ do i=1,nprop
+ call xdrffloat_(ixdrf, rprop(i), iret)
+ enddo
+#else
+ call xdrffloat(ixdrf, rtime, iret)
+ call xdrffloat(ixdrf, rpotE, iret)
+ write (iout,*) "rpotE",rpotE," iret",iret !d
+ call flush(iout)
+ call xdrffloat(ixdrf, ruconst, iret)
+ call xdrffloat(ixdrf, rt_bath, iret)
+ call xdrfint(ixdrf, nss, iret)
+ do j=1,nss
+ call xdrfint(ixdrf, ihpb(j), iret)
+ call xdrfint(ixdrf, jhpb(j), iret)
+ enddo
+ call xdrfint(ixdrf, nprop, iret)
+ write (iout,*) "nprop",nprop !d
+ if (it.gt.0 .and. nprop.ne.nprop_prev) then
+ write (iout,*) "Warning previous nprop",nprop_prev,&
+ " current",nprop
+ nprop=nprop_prev
+ else
+ nprop_prev=nprop
+ endif
+ call flush(iout)
+ if (umbrella(iparm) .or. read_iset(iparm) .or. hamil_rep) &
+ call xdrfint(ixdrf, iset, iret)
+ do i=1,nprop
+ call xdrffloat(ixdrf, rprop(i), iret)
+ enddo
+#endif
+ if (iret.eq.0) exit
+ itraj=mod(it,totraj(iR,iparm))
+#define DEBUG
+#ifdef DEBUG
+ write (iout,*) "ii",ii," itraj",itraj," it",it
+#endif
+ if (iset.eq.0) iset = 1
+ call flush(iout)
+ it=it+1
+ if (itraj.gt.ntraj) ntraj=itraj
+ nstep(itraj)=nstep(itraj)+1
+! rprop(2)=dsqrt(rprop(2))
+! rprop(3)=dsqrt(rprop(3))
+#ifdef DEBUG
+ write (iout,*) "umbrella ",umbrella
+ write (iout,*) rtime,rpotE,rt_bath,nss,&
+ (ihpb(j),jhpb(j),j=1,nss),(rprop(j),j=1,nprop)
+ write (iout,*) "nprop",nprop," iset",iset," myparm",myparm
+ call flush(iout)
+#endif
+ prec=10000.0
+ itmp=0
+#if (defined(AIX) && !defined(JUBL))
+ call xdrf3dfcoord_(ixdrf, xoord, itmp, prec, iret)
+#else
+ call xdrf3dfcoord(ixdrf, xoord, itmp, prec, iret)
+#endif
+#ifdef DEBUG
+ write (iout,'(10f8.3)') ((xoord(j,i),j=1,3),i=1,2*nres+2)
+#endif
+#undef DEBUG
+ if (iret.eq.0) exit
+ if (itmp .ne. nres + nct - nnt + 1) then
+ write (iout,*) "Error: inconsistent sizes",itmp,nres+nct-nnt+1
+ call flush(iout)
+ exit
+ endif
+
+ time=rtime
+ write (iout,*) "calling slice" !d
+ call flush(iout) !d
+ islice=slice(nstep(itraj),time,is,ie,ts,te)
+ write (iout,*) "islice",islice !d
+ call flush(iout) !d
+
+ do i=1,nres
+ do j=1,3
+ c(j,i)=xoord(j,i)
+ enddo
+ enddo
+ do i=1,nct-nnt+1
+ do j=1,3
+ c(j,i+nres+nnt-1)=xoord(j,i+nres)
+ enddo
+ enddo
+
+ if (islice.gt.0 .and. islice.le.nslice .and. (.not.separate_parset &
+ .or. iset.eq.myparm)) then
+ ii=ii+1
+ kk(islice)=kk(islice)+1
+ mm(islice)=mm(islice)+1
+ if (mod(nstep(itraj),isampl(iparm)).eq.0 .and. &
+ conf_check(ll(islice)+1,1)) then
+ if (replica(iparm)) then
+ rt_bath=1.0d0/(rt_bath*1.987D-3)
+ do i=1,nT_h(iparm)
+ if (abs(real(beta_h(i,iparm))-rt_bath).lt.1.0e-4) then
+ iib = i
+ goto 22
+ endif
+ enddo
+ 22 continue
+ if (i.gt.nT_h(iparm)) then
+ write (iout,*) "Error - temperature of conformation",&
+ ii,1.0d0/(rt_bath*1.987D-3),&
+ " does not match any of the list"
+ write (iout,*) &
+ 1.0d0/(rt_bath*1.987D-3),&
+ (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
+ call flush(iout)
+! exit
+! call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
+ ii=ii-1
+ kk(islice)=kk(islice)-1
+ mm(islice)=mm(islice)-1
+ goto 112
+ endif
+ else
+ iib = ib
+ endif
+
+ efree=0.0d0
+ jj(islice)=jj(islice)+1
+ if (umbrella(iparm)) then
+ snk(iset,iib,iparm,islice)=snk(iset,iib,iparm,islice)+1
+ else if (hamil_rep) then
+ snk(1,iib,iparm,islice)=snk(1,iib,iparm,islice)+1
+ else
+ snk(iR,iib,iparm,islice)=snk(iR,iib,iparm,islice)+1
+ endif
+ ll(islice)=ll(islice)+1
+#ifdef DEBUG
+ write (iout,*) "Writing conformation, record",ll(islice)
+ write (iout,*) "ib",ib," iib",iib
+ write (iout,*) "ntraj",ntraj," itraj",itraj,&
+ " nstep",nstep(itraj)
+ write (iout,*) "pote",rpotE," time",rtime
+! if (replica(iparm)) then
+! write (iout,*) "TEMP",1.0d0/(rt_bath*1.987D-3)
+! write (iout,*) "TEMP list"
+! write (iout,*)
+! & (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
+! endif
+ write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
+! write (iout,*) "nres",nres," nnt",nnt," nct",nct," nss",nss
+! write (iout,*) "length",nres*4+(nct-nnt+1)*4+4+2*nss*4
+ call flush(iout)
+#endif
+ if (islice.ne.islice1) then
+! write (iout,*) "islice",islice," islice1",islice1
+ close(ientout)
+! write (iout,*) "Closing file ",
+! & bprotfile_temp(:ilen(bprotfile_temp))
+ call opentmp(islice,ientout,bprotfile_temp)
+! write (iout,*) "Opening file ",
+! & bprotfile_temp(:ilen(bprotfile_temp))
+ islice1=islice
+ endif
+ if (umbrella(iparm)) then
+ write(ientout,rec=ll(islice)) &
+ ((xoord(l,k),l=1,3),k=1,nres),&
+ ((xoord(l,k),l=1,3),k=nres+1,nres+nct-nnt+1),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+ rpotE+0.0d0,efree,rmsdev,(rprop(i)+0.0d0,i=1,nQ),&
+ iset,iib,iparm
+ else if (hamil_rep) then
+ write(ientout,rec=ll(islice)) &
+ ((xoord(l,k),l=1,3),k=1,nres),&
+ ((xoord(l,k),l=1,3),k=nres+1,nres+nct-nnt+1),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+ rpotE+0.0d0,efree,rmsdev,(rprop(i)+0.0d0,i=1,nQ),&
+ iR,iib,iset
+ else
+ write(ientout,rec=ll(islice)) &
+ ((xoord(l,k),l=1,3),k=1,nres),&
+ ((xoord(l,k),l=1,3),k=nres+1,nres+nct-nnt+1),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+ rpotE+0.0d0,efree,rmsdev,(rprop(i)+0.0d0,i=1,nQ),&
+ iR,iib,iparm
+ endif
+#ifdef DEBUG
+ call int_from_cart1(.false.)
+ write (iout,*) "Writing conformation, record",ll(islice)
+ write (iout,*) "Cartesian coordinates"
+ write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
+ write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
+ write (iout,*) "Internal coordinates"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
+! write (iout,'(8f10.5)') (rprop(j),j=1,nQ)
+ write (iout,'(16i5)') iscor
+ call flush(iout)
+#endif
+ endif
+ endif
+
+ 112 continue
+
+ enddo
+ close(ientout)
+#if (defined(AIX) && !defined(JUBL))
+ call xdrfclose_(ixdrf, iret)
+#else
+ call xdrfclose(ixdrf, iret)
+#endif
+ write (iout,'(i10," trajectories found in file.")') ntraj+1
+ write (iout,'(a)') "Numbers of steps in trajectories:"
+ write (iout,'(8i10)') (nstep(i),i=0,ntraj)
+ write (iout,*) ii," conformations read from file",&
+ nazwa(:ilen(nazwa))
+ do islice=1,nslice
+ write (iout,*) mm(islice)," conformations read so far, slice",&
+ islice
+ write (iout,*) ll(islice),&
+ " conformations stored so far, slice",islice
+ enddo
+ call flush(iout)
+#undef DEBUG
+ return
+ end subroutine cxread
+!--------------------------------------------------------------------------------
+! xread.F
+!--------------------------------------------------------------------------------
+ subroutine xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
+
+ use geometry_data
+! implicit none
+! include "DIMENSIONS"
+! include "DIMENSIONS.ZSCOPT"
+! include "DIMENSIONS.FREE"
+ use MPI_data, only:nprocs
+#ifdef MPI
+ include "mpif.h"
+ integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+! include "COMMON.MPI"
+#endif
+ integer,parameter :: MaxTraj=2050
+! include "COMMON.CHAIN"
+! include "COMMON.IOUNITS"
+! include "COMMON.PROTFILES"
+! include "COMMON.NAMES"
+! include "COMMON.VAR"
+! include "COMMON.GEO"
+! include "COMMON.ENEPS"
+! include "COMMON.PROT"
+! include "COMMON.INTERACT"
+! include "COMMON.FREE"
+! include "COMMON.SBRIDGE"
+! include "COMMON.OBCINKA"
+ real(kind=4) :: csingle(3,nres*2)
+ character(len=64) :: nazwa,bprotfile_temp
+ integer :: i,j,k,l,ii,jj(nslice),kk(nslice),ll(nslice),&
+ mm(nslice) !(maxslice)
+ integer :: iscor,islice,islice1 !el,slice
+ real(kind=8) :: energ
+! integer ilen,iroof
+! external ilen,iroof
+ real(kind=8) :: rmsdev,energia(0:n_ene),efree,eini,temp
+!el real(kind=8) :: rmsdev,energia(0:max_eneW),efree,eini,temp
+ real(kind=8) :: prop(nQ) !(maxQ)
+ integer :: ntot_all(0:nprocs-1)!(0:maxprocs-1)
+ integer :: iparm,ib,iib,ir,nprop,nthr
+ real(kind=8) :: etot,time,ts(nslice),te(nslice)
+ integer :: is(nslice),ie(nslice),itraj,ntraj,it,iset
+ integer :: nstep(0:MaxTraj-1)
+ logical :: lerr
+
+ call set_slices(is,ie,ts,te,iR,ib,iparm)
+ do i=1,nQ
+ prop(i)=0.0d0
+ enddo
+ do i=0,MaxTraj-1
+ nstep(i)=0
+ enddo
+ ntraj=0
+ it=0
+ islice1=1
+ call opentmp(islice1,ientout,bprotfile_temp)
+ do while (.true.)
+ if (replica(iparm)) then
+ if (hamil_rep .or. umbrella(iparm)) then
+ read (ientin,*,end=1112,err=1112) time,eini,&
+ etot,temp,nss,(ihpb(j),jhpb(j),j=1,nss),&
+ nprop,(prop(j),j=1,nprop),iset
+ else
+ read (ientin,*,end=1112,err=1112) time,eini,&
+ etot,temp,nss,(ihpb(j),jhpb(j),j=1,nss),&
+ nprop,(prop(j),j=1,nprop)
+ endif
+ temp=1.0d0/(temp*1.987D-3)
+! write (iout,*) time,eini,etot,nss,
+! & (ihpb(j),jhpb(j),j=1,nss),(prop(j),j=1,nprop)
+! call flush(iout)
+ do i=1,nT_h(iparm)
+ if (beta_h(i,iparm).eq.temp) then
+ iib = i
+ goto 22
+ endif
+ enddo
+ 22 continue
+ if (i.gt.nT_h(iparm)) then
+ write (iout,*) "Error - temperature of conformation",&
+ ii,1.0d0/(temp*1.987D-3),&
+ " does not match any of the list"
+ write (iout,*) &
+ 1.0d0/(temp*1.987D-3),&
+ (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
+ call flush(iout)
+#ifdef MPI
+ call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
+#endif
+ endif
+ else
+ read (ientin,*,end=1112,err=1112) time,eini,&
+ etot,nss,(ihpb(j),jhpb(j),j=1,nss),&
+ nprop,(prop(j),j=1,nprop)
+ iib = ib
+ endif
+ itraj=mod(it,totraj(iR,iparm))
+! write (*,*) "ii",ii," itraj",itraj
+! call flush(iout)
+ it=it+1
+ if (itraj.gt.ntraj) ntraj=itraj
+ nstep(itraj)=nstep(itraj)+1
+ islice=slice(nstep(itraj),time,is,ie,ts,te)
+ read (ientin,'(8f10.5)',end=1112,err=1112) &
+ ((csingle(l,k),l=1,3),k=1,nres),&
+ ((csingle(l,k+nres),l=1,3),k=nnt,nct)
+ efree=0.0d0
+ if (islice.gt.0 .and. islice.le.nslice) then
+ ii=ii+1
+ kk(islice)=kk(islice)+1
+ mm(islice)=mm(islice)+1
+ if (mod(nstep(itraj),isampl(iparm)).eq.0) then
+ jj(islice)=jj(islice)+1
+ if (hamil_rep) then
+ snk(iR,iib,iset,islice)=snk(iR,iib,iset,islice)+1
+ else if (umbrella(iparm)) then
+ snk(iset,iib,iparm,islice)=snk(iset,iib,iparm,islice)+1
+ else
+ snk(iR,iib,iparm,islice)=snk(iR,iib,iparm,islice)+1
+ endif
+ ll(islice)=ll(islice)+1
+! write (iout,*) ii,kk,jj,ll,eini,(prop(j),j=1,nprop)
+#ifdef DEBUG
+! write (iout,*) "Writing conformation, record",ll(islice)
+! write (iout,*) "ib",ib," iib",iib
+ if (replica(iparm)) then
+ write (iout,*) "TEMP",1.0d0/(temp*1.987D-3)
+ write (iout,*) "TEMP list"
+ write (iout,*) &
+ (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
+ endif
+ call flush(iout)
+#endif
+! write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
+! write (iout,*) "nres",nres," nnt",nnt," nct",nct," nss",nss
+! write (iout,*) "length",nres*4+(nct-nnt+1)*4+4+2*nss*4
+! call flush(iout)
+ if (islice.ne.islice1) then
+! write (iout,*) "islice",islice," islice1",islice1
+ close(ientout)
+! write (iout,*) "Closing file ",
+! & bprotfile_temp(:ilen(bprotfile_temp))
+ call opentmp(islice,ientout,bprotfile_temp)
+! write (iout,*) "Opening file ",
+! & bprotfile_temp(:ilen(bprotfile_temp))
+! call flush(iout)
+ islice1=islice
+ endif
+ write(ientout,rec=ll(islice)) &
+ ((csingle(l,k),l=1,3),k=1,nres),&
+ ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+ eini,efree,rmsdev,(prop(i),i=1,nQ),iR,iib,iparm
+#ifdef DEBUG
+ do i=1,2*nres
+ do j=1,3
+ c(j,i)=csingle(j,i)
+ enddo
+ enddo
+ call int_from_cart1(.false.)
+ write (iout,*) "Writing conformation, record",ll(islice)
+ write (iout,*) "Cartesian coordinates"
+ write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
+ write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
+ write (iout,*) "Internal coordinates"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
+! write (iout,'(8f10.5)') (prop(j),j=1,nQ)
+ write (iout,'(16i5)') iscor
+ call flush(iout)
+#endif
+ endif
+ endif
+ enddo
+ 1112 continue
+ close(ientout)
+ write (iout,'(i10," trajectories found in file.")') ntraj+1
+ write (iout,'(a)') "Numbers of steps in trajectories:"
+ write (iout,'(8i10)') (nstep(i),i=0,ntraj)
+ write (iout,*) ii," conformations read from file",&
+ nazwa(:ilen(nazwa))
+ write (iout,*) mm(islice)," conformations read so far, slice",&
+ islice
+ write (iout,*) ll(islice)," conformations stored so far, slice",&
+ islice
+ call flush(iout)
+ return
+ end subroutine xread
+!--------------------------------------------------------------------------------
+! enecalc1.F
+!--------------------------------------------------------------------------------
+ subroutine write_dbase(islice,*)
+
+ use geometry_data
+ use control_data, only:indpdb
+ use w_compar_data
+ use conform_compar, only:conf_compar
+! implicit none
+! include "DIMENSIONS"
+! include "DIMENSIONS.ZSCOPT"
+! include "DIMENSIONS.FREE"
+! include "DIMENSIONS.COMPAR"
+ use geometry, only:int_from_cart1
+#ifdef MPI
+ include "mpif.h"
+ integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+! include "COMMON.MPI"
+#endif
+! include "COMMON.CONTROL"
+! include "COMMON.CHAIN"
+! include "COMMON.IOUNITS"
+! include "COMMON.PROTFILES"
+! include "COMMON.NAMES"
+! include "COMMON.VAR"
+! include "COMMON.SBRIDGE"
+! include "COMMON.GEO"
+! include "COMMON.FFIELD"
+! include "COMMON.ENEPS"
+! include "COMMON.LOCAL"
+! include "COMMON.WEIGHTS"
+! include "COMMON.INTERACT"
+! include "COMMON.FREE"
+! include "COMMON.ENERGIES"
+! include "COMMON.COMPAR"
+! include "COMMON.PROT"
+! include "COMMON.CONTACTS1"
+ character(len=64) :: nazwa
+ character(len=80) :: bxname,cxname
+ character(len=64) :: bprotfile_temp
+ character(len=3) :: liczba,licz
+ character(len=2) :: licz2
+ integer :: i,itj,ii,iii,j,k,l
+ integer :: ixdrf,iret
+ integer :: iscor,islice
+ real(kind=8) :: rmsdev,efree,eini
+ real(kind=4) :: csingle(3,nres*2)
+ real(kind=8) :: energ
+! integer ilen,iroof
+! external ilen,iroof
+ integer :: ir,ib,iparm
+ integer :: isecstr(nres)
+ write (licz2,'(bz,i2.2)') islice
+ call opentmp(islice,ientout,bprotfile_temp)
+ write (iout,*) "bprotfile_temp ",bprotfile_temp
+ call flush(iout)
+ if (.not.bxfile .and. .not. cxfile .and. indpdb.eq.0 &
+ .and. ensembles.eq.0) then
+ close(ientout,status="delete")
+ return
+ endif
+#ifdef MPI
+ write (liczba,'(bz,i3.3)') me
+ if (bxfile .or. cxfile .or. ensembles.gt.0) then
+ if (.not.separate_parset) then
+ bxname = prefix(:ilen(prefix))//liczba//".bx"
+ else
+ write (licz,'(bz,i3.3)') myparm
+ bxname = prefix(:ilen(prefix))//liczba//"_par"//licz//".bx"
+ endif
+ open (ientin,file=bxname,status="unknown",&
+ form="unformatted",access="direct",recl=lenrec1)
+ endif
+#else
+ if (bxfile .or. cxfile .or. ensembles.gt.0) then
+ if (nslice.eq.1) then
+ bxname = prefix(:ilen(prefix))//".bx"
+ else
+ bxname = prefix(:ilen(prefix))// &
+ "_slice_"//licz2//".bx"
+ endif
+ open (ientin,file=bxname,status="unknown",&
+ form="unformatted",access="direct",recl=lenrec1)
+ write (iout,*) "Calculating energies; writing geometry",&
+ " and energy components to ",bxname(:ilen(bxname))
+ endif
+#if (defined(AIX) && !defined(JUBL))
+ call xdrfopen_(ixdrf,cxname, "w", iret)
+#else
+ call xdrfopen(ixdrf,cxname, "w", iret)
+#endif
+ if (iret.eq.0) then
+ write (iout,*) "Error opening cxfile ",cxname(:ilen(cxname))
+ cxfile=.false.
+ endif
+!el endif
+#endif
+ if (indpdb.gt.0) then
+ if (nslice.eq.1) then
+#ifdef MPI
+ if (.not.separate_parset) then
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot)) &
+ //liczba//'.stat'
+ else
+ write (licz,'(bz,i3.3)') myparm
+ statname=prefix(:ilen(prefix))//'_par'//licz//'_'// &
+ pot(:ilen(pot))//liczba//'.stat'
+ endif
+
+#else
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))//'.stat'
+#endif
+ else
+#ifdef MPI
+ if (.not.separate_parset) then
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))// &
+ "_slice_"//licz2//liczba//'.stat'
+ else
+ write (licz,'(bz,i3.3)') myparm
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))// &
+ '_par'//licz//"_slice_"//licz2//liczba//'.stat'
+ endif
+#else
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot)) &
+ //"_slice_"//licz2//'.stat'
+#endif
+ endif
+ open(istat,file=statname,status="unknown")
+ endif
+
+#ifdef MPI
+ do i=1,scount(me)
+#else
+ do i=1,ntot(islice)
+#endif
+ read(ientout,rec=i,err=101) &
+ ((csingle(l,k),l=1,3),k=1,nres),&
+ ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+ eini,efree,rmsdev,(q(k,i),k=1,nQ),iR,ib,iparm
+! write (iout,*) iR,ib,iparm,eini,efree
+ do j=1,2*nres
+ do k=1,3
+ c(k,j)=csingle(k,j)
+ enddo
+ enddo
+ call int_from_cart1(.false.)
+ iscore=0
+! write (iout,*) "Calling conf_compar",i
+! call flush(iout)
+ anatemp= 1.0d0/(beta_h(ib,iparm)*1.987D-3)
+ if (indpdb.gt.0) then
+ call conf_compar(i,.false.,.true.)
+! else
+! call elecont(.false.,ncont,icont,nnt,nct)
+! call secondary2(.false.,.false.,ncont,icont,isecstr)
+ endif
+! write (iout,*) "Exit conf_compar",i
+! call flush(iout)
+ if (bxfile .or.cxfile .or. ensembles.gt.0) write (ientin,rec=i) &
+ ((csingle(l,k),l=1,3),k=1,nres),&
+ ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+! & potE(i,iparm),-entfac(i),rms_nat,iscore
+ potE(i,nparmset),-entfac(i),rms_nat,iscore
+! write (iout,'(2i5,3e15.5)') i,me,potE(i,iparm),-entfac(i)
+#ifndef MPI
+ if (cxfile) call cxwrite(ixdrf,csingle,potE(i,nparmset),&
+ -entfac(i),rms_nat,iscore)
+#endif
+ enddo
+ close(ientout,status="delete")
+ close(istat)
+ if (bxfile .or. cxfile .or. ensembles.gt.0) close(ientin)
+#ifdef MPI
+ call MPI_Barrier(WHAM_COMM,IERROR)
+ if (me.ne.Master .or. .not.bxfile .and. .not. cxfile &
+ .and. ensembles.eq.0) return
+ write (iout,*)
+ if (bxfile .or. ensembles.gt.0) then
+ if (nslice.eq.1) then
+ if (.not.separate_parset) then
+ bxname = prefix(:ilen(prefix))//".bx"
+ else
+ write (licz,'(bz,i3.3)') myparm
+ bxname = prefix(:ilen(prefix))//"_par"//licz//".bx"
+ endif
+ else
+ if (.not.separate_parset) then
+ bxname = prefix(:ilen(prefix))//"_slice_"//licz2//".bx"
+ else
+ write (licz,'(bz,i3.3)') myparm
+ bxname = prefix(:ilen(prefix))//"par_"//licz// &
+ "_slice_"//licz2//".bx"
+ endif
+ endif
+ open (ientout,file=bxname,status="unknown",&
+ form="unformatted",access="direct",recl=lenrec1)
+ write (iout,*) "Master is creating binary database ",&
+ bxname(:ilen(bxname))
+ endif
+ if (cxfile) then
+ if (nslice.eq.1) then
+ if (.not.separate_parset) then
+ cxname = prefix(:ilen(prefix))//".cx"
+ else
+ cxname = prefix(:ilen(prefix))//"_par"//licz//".cx"
+ endif
+ else
+ if (.not.separate_parset) then
+ cxname = prefix(:ilen(prefix))// &
+ "_slice_"//licz2//".cx"
+ else
+ cxname = prefix(:ilen(prefix))//"_par"//licz// &
+ "_slice_"//licz2//".cx"
+ endif
+ endif
+#if (defined(AIX) && !defined(JUBL))
+ call xdrfopen_(ixdrf,cxname, "w", iret)
+#else
+ call xdrfopen(ixdrf,cxname, "w", iret)
+#endif
+ if (iret.eq.0) then
+ write (iout,*) "Error opening cxfile ",cxname(:ilen(cxname))
+ cxfile=.false.
+ endif
+ endif
+ do j=0,nprocs-1
+ write (liczba,'(bz,i3.3)') j
+ if (separate_parset) then
+ write (licz,'(bz,i3.3)') myparm
+ bxname = prefix(:ilen(prefix))//liczba//"_par"//licz//".bx"
+ else
+ bxname = prefix(:ilen(prefix))//liczba//".bx"
+ endif
+ open (ientin,file=bxname,status="unknown",&
+ form="unformatted",access="direct",recl=lenrec1)
+ write (iout,*) "Master is reading conformations from ",&
+ bxname(:ilen(bxname))
+ iii = 0
+! write (iout,*) "j",j," indstart",indstart(j)," indend",indend(j)
+! call flush(iout)
+ do i=indstart(j),indend(j)
+ iii = iii+1
+ read(ientin,rec=iii,err=101) &
+ ((csingle(l,k),l=1,3),k=1,nres),&
+ ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+ eini,efree,rmsdev,iscor
+ if (bxfile .or. ensembles.gt.0) then
+ write (ientout,rec=i) &
+ ((csingle(l,k),l=1,3),k=1,nres),&
+ ((csingle(l,k+nres),l=1,3),k=nnt,nct),&
+ nss,(ihpb(k),jhpb(k),k=1,nss),&
+ eini,efree,rmsdev,iscor
+ endif
+ if(cxfile)call cxwrite(ixdrf,csingle,eini,efree,rmsdev,iscor)
+#ifdef DEBUG
+ do k=1,2*nres
+ do l=1,3
+ c(l,k)=csingle(l,k)
+ enddo
+ enddo
+ call int_from_cart1(.false.)
+ write (iout,'(2i5,3e15.5)') i,iii,eini,efree
+ write (iout,*) "The Cartesian geometry is:"
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,*) "The internal geometry is:"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
+ write (iout,'(f10.5,i5)') rmsdev,iscor
+#endif
+ enddo ! i
+ write (iout,*) iii," conformations (from",indstart(j)," to",&
+ indend(j),") read from ",&
+ bxname(:ilen(bxname))
+ close (ientin,status="delete")
+ enddo ! j
+ if (bxfile .or. cxfile .or. ensembles.gt.0) close (ientout)
+#if (defined(AIX) && !defined(JUBL))
+ if (cxfile) call xdrfclose_(ixdrf,cxname,iret)
+#else
+ if (cxfile) call xdrfclose(ixdrf,cxname,iret)
+#endif
+#endif
+ return
+ 101 write (iout,*) "Error in scratchfile."
+ call flush(iout)
+ return 1
+ end subroutine write_dbase
+!-------------------------------------------------------------------------------
+ subroutine cxwrite(ixdrf,csingle,eini,efree,rmsdev,iscor)
+! implicit none
+! include "DIMENSIONS"
+! include "DIMENSIONS.ZSCOPT"
+! include "DIMENSIONS.FREE"
+! include "DIMENSIONS.COMPAR"
+#ifdef MPI
+ include "mpif.h"
+ integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+! include "COMMON.MPI"
+#endif
+! include "COMMON.CONTROL"
+! include "COMMON.CHAIN"
+! include "COMMON.IOUNITS"
+! include "COMMON.PROTFILES"
+! include "COMMON.NAMES"
+! include "COMMON.VAR"
+! include "COMMON.SBRIDGE"
+! include "COMMON.GEO"
+! include "COMMON.FFIELD"
+! include "COMMON.ENEPS"
+! include "COMMON.LOCAL"
+! include "COMMON.WEIGHTS"
+! include "COMMON.INTERACT"
+! include "COMMON.FREE"
+! include "COMMON.ENERGIES"
+! include "COMMON.COMPAR"
+! include "COMMON.PROT"
+ integer :: i,j,itmp,iscor,iret,ixdrf
+ real(kind=8) :: rmsdev,efree,eini
+ real(kind=4) :: csingle(3,nres*2),xoord(3,2*nres+2)
+ real(kind=4) :: prec
+
+! write (iout,*) "cxwrite"
+! call flush(iout)
+ prec=10000.0
+ do i=1,nres
+ do j=1,3
+ xoord(j,i)=csingle(j,i)
+ enddo
+ enddo
+ do i=nnt,nct
+ do j=1,3
+ xoord(j,nres+i-nnt+1)=csingle(j,i+nres)
+ enddo
+ enddo
+
+ itmp=nres+nct-nnt+1
+
+! write (iout,*) "itmp",itmp
+! call flush(iout)
+#if (defined(AIX) && !defined(JUBL))
+ call xdrf3dfcoord_(ixdrf, xoord, itmp, prec, iret)
+
+! write (iout,*) "xdrf3dfcoord"
+! call flush(iout)
+ call xdrfint_(ixdrf, nss, iret)
+ do j=1,nss
+ call xdrfint_(ixdrf, ihpb(j), iret)
+ call xdrfint_(ixdrf, jhpb(j), iret)
+ enddo
+ call xdrffloat_(ixdrf,real(eini),iret)
+ call xdrffloat_(ixdrf,real(efree),iret)
+ call xdrffloat_(ixdrf,real(rmsdev),iret)
+ call xdrfint_(ixdrf,iscor,iret)
+#else
+ call xdrf3dfcoord(ixdrf, xoord, itmp, prec, iret)
+
+ call xdrfint(ixdrf, nss, iret)
+ do j=1,nss
+ call xdrfint(ixdrf, ihpb(j), iret)
+ call xdrfint(ixdrf, jhpb(j), iret)
+ enddo
+ call xdrffloat(ixdrf,real(eini),iret)
+ call xdrffloat(ixdrf,real(efree),iret)
+ call xdrffloat(ixdrf,real(rmsdev),iret)
+ call xdrfint(ixdrf,iscor,iret)
+#endif
+
+ return
+ end subroutine cxwrite
+!-------------------------------------------------------------------------------
+! slices.F
+!-------------------------------------------------------------------------------
+ subroutine set_slices(is,ie,ts,te,iR,ib,iparm)
+! implicit none
+! include 'DIMENSIONS'
+! include 'DIMENSIONS.ZSCOPT'
+! include 'DIMENSIONS.FREE'
+! include 'COMMON.IOUNITS'
+! include 'COMMON.PROTFILES'
+! include 'COMMON.OBCINKA'
+! include 'COMMON.PROT'
+ integer :: islice,iR,ib,iparm
+ integer :: is(MaxSlice),ie(MaxSlice),nrec_slice
+ real(kind=8) :: ts(MaxSlice),te(MaxSlice),time_slice
+
+ do islice=1,nslice
+ if (time_end_collect(iR,ib,iparm).ge.1.0d10) then
+ ts(islice)=time_start_collect(iR,ib,iparm)
+ te(islice)=time_end_collect(iR,ib,iparm)
+ nrec_slice=(rec_end(iR,ib,iparm)- &
+ rec_start(iR,ib,iparm)+1)/nslice
+ is(islice)=rec_start(iR,ib,iparm)+(islice-1)*nrec_slice
+ ie(islice)=rec_start(iR,ib,iparm)+islice*nrec_slice-1
+ else
+ time_slice=(time_end_collect(iR,ib,iparm) &
+ -time_start_collect(iR,ib,iparm))/nslice
+ ts(islice)=time_start_collect(iR,ib,iparm)+(islice-1)* &
+ time_slice
+ te(islice)=time_start_collect(iR,ib,iparm)+islice*time_slice
+ is(islice)=rec_start(iR,ib,iparm)
+ ie(islice)=rec_end(iR,ib,iparm)
+ endif
+ enddo
+
+ write (iout,*) "nrec_slice",nrec_slice," time_slice",time_slice
+ write (iout,*) "is",(is(islice),islice=1,nslice)
+ write (iout,*) "ie",(ie(islice),islice=1,nslice)
+ write (iout,*) "rec_start",&
+ rec_start(iR,ib,iparm)," rec_end",rec_end(iR,ib,iparm)
+ write (iout,*) "ts",(ts(islice),islice=1,nslice)
+ write (iout,*) "te",(te(islice),islice=1,nslice)
+ write (iout,*) "time_start",&
+ time_start_collect(iR,ib,iparm)," time_end",&
+ time_end_collect(iR,ib,iparm)
+ call flush(iout)
+
+ return
+ end subroutine set_slices
+!-----------------------------------------------------------------------------
+ integer function slice(irecord,time,is,ie,ts,te)
+! implicit none
+! include 'DIMENSIONS'
+! include 'DIMENSIONS.ZSCOPT'
+! include 'DIMENSIONS.FREE'
+! include 'COMMON.IOUNITS'
+! include 'COMMON.PROTFILES'
+! include 'COMMON.OBCINKA'
+! include 'COMMON.PROT'
+ integer :: is(MaxSlice),ie(MaxSlice),nrec_slice
+ real(kind=8) :: ts(MaxSlice),te(MaxSlice),time_slice
+ integer :: i,ii,irecord
+ real(kind=8) :: time
+
+! write (iout,*) "within slice nslice",nslice
+! call flush(iout)
+ if (irecord.lt.is(1) .or. time.lt.ts(1)) then
+ ii=0
+ else
+ ii=1
+ do while (ii.le.nslice .and. &
+ (irecord.lt.is(ii) .or. irecord.gt.ie(ii) .or. &
+ time.lt.ts(ii) .or. time.gt.te(ii)) )
+! write (iout,*) "ii",ii,time,ts(ii)
+! call flush(iout)
+ ii=ii+1
+ enddo
+ endif
+! write (iout,*) "end: ii",ii
+! call flush(iout)
+ slice=ii
+ return
+ end function slice
+!-----------------------------------------------------------------------------
+! enecalc1.F
+!-----------------------------------------------------------------------------
+ logical function conf_check(ii,iprint)
+
+ use names, only:ntyp1
+ use geometry_data
+ use energy_data, only:itype,dsc
+ use geometry, only:int_from_cart1
+! use
+! include "DIMENSIONS"
+! include "DIMENSIONS.ZSCOPT"
+! include "DIMENSIONS.FREE"
+!#ifdef MPI
+! use MPI_data
+! include "mpif.h"
+! include "COMMON.MPI"
+!#endif
+! include "COMMON.CHAIN"
+! include "COMMON.IOUNITS"
+! include "COMMON.PROTFILES"
+! include "COMMON.NAMES"
+! include "COMMON.VAR"
+! include "COMMON.SBRIDGE"
+! include "COMMON.GEO"
+! include "COMMON.FFIELD"
+! include "COMMON.ENEPS"
+! include "COMMON.LOCAL"
+! include "COMMON.WEIGHTS"
+! include "COMMON.INTERACT"
+! include "COMMON.FREE"
+! include "COMMON.ENERGIES"
+! include "COMMON.CONTROL"
+! include "COMMON.TORCNSTR"
+! implicit none
+#ifdef MPI
+ include "mpif.h"
+ integer :: IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+#endif
+ integer :: j,k,l,ii,itj,iprint
+ if (.not. check_conf) then
+ conf_check=.true.
+ return
+ endif
+ call int_from_cart1(.false.)
+ do j=nnt+1,nct
+ if (itype(j-1).ne.ntyp1 .and. itype(j).ne.ntyp1 .and. &
+ (vbld(j).lt.2.0d0 .or. vbld(j).gt.5.0d0)) then
+ if (iprint.gt.0) &
+ write (iout,*) "Bad CA-CA bond length",j," ",vbld(j),&
+ " for conformation",ii
+ if (iprint.gt.1) then
+ write (iout,*) "The Cartesian geometry is:"
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,*) "The internal geometry is:"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ endif
+ if (iprint.gt.0) write (iout,*) &
+ "This conformation WILL NOT be added to the database."
+ conf_check=.false.
+ return
+ endif
+ enddo
+ do j=nnt,nct
+ itj=itype(j)
+ if (itype(j).ne.10 .and.itype(j).ne.ntyp1 .and. &
+ (vbld(nres+j)-dsc(iabs(itj))).gt.2.0d0) then
+ if (iprint.gt.0) &
+ write (iout,*) "Bad CA-SC bond length",j," ",vbld(nres+j),&
+ " for conformation",ii
+ if (iprint.gt.1) then
+ write (iout,*) "The Cartesian geometry is:"
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,*) "The internal geometry is:"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ endif
+ if (iprint.gt.0) write (iout,*) &
+ "This conformation WILL NOT be added to the database."
+ conf_check=.false.
+ return
+ endif
+ enddo
+ do j=3,nres
+ if (theta(j).le.0.0d0) then
+ if (iprint.gt.0) &
+ write (iout,*) "Zero theta angle(s) in conformation",ii
+ if (iprint.gt.1) then
+ write (iout,*) "The Cartesian geometry is:"
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,*) "The internal geometry is:"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ endif
+ if (iprint.gt.0) write (iout,*) &
+ "This conformation WILL NOT be added to the database."
+ conf_check=.false.
+ return
+ endif
+ if (theta(j).gt.179.97*deg2rad) theta(j)=179.97*deg2rad
+ enddo
+ conf_check=.true.
+! write (iout,*) "conf_check passed",ii
+ return
+ end function conf_check
+!-----------------------------------------------------------------------------
+ end module io_database
projects / unres4.git / blobdiff