source/wham/src/xread.F

   1       subroutine xread(nazwa,ii,jj,kk,ll,mm,iR,ib,iparm)
   2       implicit none
   3       include "DIMENSIONS"
   4       include "DIMENSIONS.ZSCOPT"
   5       include "DIMENSIONS.FREE"
   6       integer MaxTraj
   7       parameter (MaxTraj=2050)
   8 #ifdef MPI
   9       include "mpif.h"
  10       integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
  11       include "COMMON.MPI"
  12 #endif
  13       include "COMMON.CHAIN"
  14       include "COMMON.IOUNITS"
  15       include "COMMON.PROTFILES"
  16       include "COMMON.NAMES"
  17       include "COMMON.VAR"
  18       include "COMMON.GEO"
  19       include "COMMON.ENEPS"
  20       include "COMMON.PROT"
  21       include "COMMON.INTERACT"
  22       include "COMMON.FREE"
  23       include "COMMON.SBRIDGE"
  24       include "COMMON.OBCINKA"
  25       real*4 csingle(3,maxres2)
  26       character*64 nazwa,bprotfile_temp
  27       integer i,j,k,l,ii,jj(maxslice),kk(maxslice),ll(maxslice),
  28      &  mm(maxslice)
  29       integer iscor,islice,islice1,slice
  30       double precision energ
  31       integer ilen,iroof
  32       external ilen,iroof
  33       double precision rmsdev,energia(0:max_ene),efree,eini,temp
  34       double precision prop(maxQ)
  35       integer ntot_all(0:maxprocs-1)
  36       integer iparm,ib,iib,ir,nprop,nthr
  37       double precision etot,time,ts(maxslice),te(maxslice)
  38       integer is(maxslice),ie(maxslice),itraj,ntraj,it,iset
  39       integer nstep(0:MaxTraj-1)
  40       logical lerr
  41
  42       call set_slices(is,ie,ts,te,iR,ib,iparm)
  43       do i=1,nQ
  44         prop(i)=0.0d0
  45       enddo
  46       do i=0,MaxTraj-1
  47         nstep(i)=0
  48       enddo
  49       ntraj=0
  50       it=0
  51       islice1=1
  52       call opentmp(islice1,ientout,bprotfile_temp)
  53       do while (.true.)
  54         if (replica(iparm)) then
  55           if (hamil_rep .or. umbrella(iparm)) then
  56           read (ientin,*,end=1112,err=1112) time,eini,
  57      &      etot,temp,nss,(ihpb(j),jhpb(j),j=1,nss),
  58      &      nprop,(prop(j),j=1,nprop),iset
  59           else
  60           read (ientin,*,end=1112,err=1112) time,eini,
  61      &      etot,temp,nss,(ihpb(j),jhpb(j),j=1,nss),
  62      &      nprop,(prop(j),j=1,nprop)
  63           endif
  64           temp=1.0d0/(temp*1.987D-3)
  65 c           write (iout,*) time,eini,etot,nss,
  66 c     &     (ihpb(j),jhpb(j),j=1,nss),(prop(j),j=1,nprop)
  67 c           call flush(iout)
  68            do i=1,nT_h(iparm)
  69              if (beta_h(i,iparm).eq.temp) then
  70                iib = i
  71                goto 22
  72              endif
  73            enddo
  74   22       continue
  75            if (i.gt.nT_h(iparm)) then
  76              write (iout,*) "Error - temperature of conformation",
  77      &       ii,1.0d0/(temp*1.987D-3),
  78      &       " does not match any of the list"
  79              write (iout,*)
  80      &        1.0d0/(temp*1.987D-3),
  81      &        (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
  82              call flush(iout)
  83              call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE)
  84            endif
  85         else
  86            read (ientin,*,end=1112,err=1112) time,eini,
  87      &       etot,nss,(ihpb(j),jhpb(j),j=1,nss),
  88      &       nprop,(prop(j),j=1,nprop)
  89              iib = ib
  90         endif
  91         itraj=mod(it,totraj(iR,iparm))
  92 c        write (*,*) "ii",ii," itraj",itraj
  93 c        call flush(iout)
  94         it=it+1
  95         if (itraj.gt.ntraj) ntraj=itraj
  96         nstep(itraj)=nstep(itraj)+1
  97         islice=slice(nstep(itraj),time,is,ie,ts,te)
  98         read (ientin,'(8f10.5)',end=1112,err=1112)
  99      &    ((csingle(l,k),l=1,3),k=1,nres),
 100      &    ((csingle(l,k+nres),l=1,3),k=nnt,nct)
 101         efree=0.0d0
 102         if (islice.gt.0 .and. islice.le.nslice) then
 103         ii=ii+1
 104         kk(islice)=kk(islice)+1
 105         mm(islice)=mm(islice)+1
 106         if (mod(nstep(itraj),isampl(iparm)).eq.0) then
 107         jj(islice)=jj(islice)+1
 108         if (hamil_rep) then
 109           snk(iR,iib,iset,islice)=snk(iR,iib,iset,islice)+1
 110         else if (umbrella(iparm)) then
 111           snk(iset,iib,iparm,islice)=snk(iset,iib,iparm,islice)+1
 112         else
 113           snk(iR,iib,iparm,islice)=snk(iR,iib,iparm,islice)+1
 114         endif
 115         ll(islice)=ll(islice)+1
 116 c         write (iout,*) ii,kk,jj,ll,eini,(prop(j),j=1,nprop)
 117 #ifdef DEBUG
 118 c        write (iout,*) "Writing conformation, record",ll(islice)
 119 c        write (iout,*) "ib",ib," iib",iib
 120          if (replica(iparm)) then
 121            write (iout,*) "TEMP",1.0d0/(temp*1.987D-3)
 122            write (iout,*) "TEMP list"
 123            write (iout,*)
 124      &      (1.0d0/(beta_h(i,iparm)*1.987D-3),i=1,nT_h(iparm))
 125          endif
 126          call flush(iout)
 127 #endif
 128 c         write (iout,*) "iparm",iparm," ib",ib," iR",iR," nQ",nQ
 129 c         write (iout,*) "nres",nres," nnt",nnt," nct",nct," nss",nss
 130 c         write (iout,*) "length",nres*4+(nct-nnt+1)*4+4+2*nss*4
 131 c         call flush(iout)
 132          if (islice.ne.islice1) then
 133 c             write (iout,*) "islice",islice," islice1",islice1
 134              close(ientout)
 135 c             write (iout,*) "Closing file ",
 136 c     &             bprotfile_temp(:ilen(bprotfile_temp))
 137              call opentmp(islice,ientout,bprotfile_temp)
 138 c             write (iout,*) "Opening file ",
 139 c     &             bprotfile_temp(:ilen(bprotfile_temp))
 140 c             call flush(iout)
 141              islice1=islice
 142          endif
 143          write(ientout,rec=ll(islice))
 144      &     ((csingle(l,k),l=1,3),k=1,nres),
 145      &     ((csingle(l,k+nres),l=1,3),k=nnt,nct),
 146      &     nss,(ihpb(k),jhpb(k),k=1,nss),
 147      &     eini,efree,rmsdev,(prop(i),i=1,nQ),iR,iib,iparm
 148 #ifdef DEBUG
 149          do i=1,2*nres
 150            do j=1,3
 151              c(j,i)=csingle(j,i)
 152            enddo
 153          enddo
 154          call int_from_cart1(.false.)
 155          write (iout,*) "Writing conformation, record",ll(islice)
 156          write (iout,*) "Cartesian coordinates"
 157          write (iout,'(8f10.5)') ((c(j,i),j=1,3),i=1,nres)
 158          write (iout,'(8f10.5)') ((c(j,i+nres),j=1,3),i=nnt,nct)
 159          write (iout,*) "Internal coordinates"
 160          write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
 161          write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
 162          write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
 163          write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
 164          write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
 165          write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
 166          write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
 167 c         write (iout,'(8f10.5)') (prop(j),j=1,nQ)
 168          write (iout,'(16i5)') iscor
 169          call flush(iout)
 170 #endif
 171          endif
 172          endif
 173        enddo
 174  1112  continue
 175        close(ientout)
 176        write (iout,'(i10," trajectories found in file.")') ntraj+1
 177        write (iout,'(a)') "Numbers of steps in trajectories:"
 178        write (iout,'(8i10)') (nstep(i),i=0,ntraj)
 179        write (iout,*) ii," conformations read from file",
 180      &   nazwa(:ilen(nazwa))
 181        write (iout,*) mm(islice)," conformations read so far, slice",
 182      &    islice
 183        write (iout,*) ll(islice)," conformations stored so far, slice",
 184      &   islice
 185        call flush(iout)
 186        return
 187        end