subroutine entmcm
C Does modified entropic sampling in the space of minima.
      implicit real*8 (a-h,o-z)
      include 'DIMENSIONS'
      include 'COMMON.IOUNITS'
#ifdef MPL
      include 'COMMON.INFO'
#endif
      include 'COMMON.GEO'
      include 'COMMON.CHAIN'
      include 'COMMON.MCM'
      include 'COMMON.MCE'
      include 'COMMON.CONTACTS'
      include 'COMMON.CONTROL'
      include 'COMMON.VAR'
      include 'COMMON.INTERACT'
      include 'COMMON.THREAD'
      include 'COMMON.NAMES'
      logical accepted,not_done,over,ovrtim,error,lprint
      integer MoveType,nbond
      integer conf_comp
      double precision RandOrPert
      double precision varia(maxvar),elowest,ehighest,eold
      double precision przes(3),obr(3,3)
      double precision varold(maxvar)
      logical non_conv
      double precision energia(0:n_ene),energia_ave(0:n_ene)
C
cd    write (iout,*) 'print_mc=',print_mc
      WhatsUp=0
      maxtrial_iter=50
c---------------------------------------------------------------------------
C Initialize counters.
c---------------------------------------------------------------------------
C Total number of generated confs.
      ngen=0
C Total number of moves. In general this won't be equal to the number of
C attempted moves, because we may want to reject some "bad" confs just by
C overlap check.
      nmove=0
C Total number of shift (nbond_move(1)), spike, crankshaft, three-bond,...
C motions.
      do i=1,nres
        nbond_move(i)=0
      enddo
C Initialize total and accepted number of moves of various kind.
      do i=0,MaxMoveType
        moves(i)=0
        moves_acc(i)=0
      enddo
C Total number of energy evaluations.
      neneval=0
      nfun=0
      indminn=-max_ene
      indmaxx=max_ene
      delte=0.5D0
      facee=1.0D0/(maxacc*delte)
      conste=dlog(facee)
C Read entropy from previous simulations. 
      if (ent_read) then
        read (ientin,*) indminn,indmaxx,emin,emax 
        print *,'indminn=',indminn,' indmaxx=',indmaxx,' emin=',emin,
     &          ' emax=',emax
        do i=-max_ene,max_ene
          entropy(i)=(emin+i*delte)*betbol
        enddo
        read (ientin,*) (ijunk,ejunk,entropy(i),i=indminn,indmaxx)
        indmin=indminn
        indmax=indmaxx
        write (iout,*) 'indminn=',indminn,' indmaxx=',indmaxx,
     &                 ' emin=',emin,' emax=',emax
        write (iout,'(/a)') 'Initial entropy'
        do i=indminn,indmaxx
          write (iout,'(i5,2f10.5)') i,emin+i*delte,entropy(i)
        enddo
      endif ! ent_read
C Read the pool of conformations
      call read_pool
C----------------------------------------------------------------------------
C Entropy-sampling simulations with continually updated entropy
C Loop thru simulations
C----------------------------------------------------------------------------
      DO ISWEEP=1,NSWEEP
C----------------------------------------------------------------------------
C Take a conformation from the pool
C----------------------------------------------------------------------------
      if (npool.gt.0) then
        ii=iran_num(1,npool)
        do i=1,nvar
          varia(i)=xpool(i,ii)
        enddo
        write (iout,*) 'Took conformation',ii,' from the pool energy=',
     &               epool(ii)
        call var_to_geom(nvar,varia)
C Print internal coordinates of the initial conformation
        call intout
      else
        call gen_rand_conf(1,*20)
      endif
C----------------------------------------------------------------------------
C Compute and print initial energies.
C----------------------------------------------------------------------------
      nsave=0
#ifdef MPL
      if (MyID.eq.MasterID) then
        do i=1,nctasks
          nsave_part(i)=0
        enddo
      endif
#endif
      Kwita=0
      WhatsUp=0
      write (iout,'(/80(1h*)/a,i2/80(1h*)/)') 'MCE iteration #',isweep
      write (iout,'(/80(1h*)/a)') 'Initial energies:'
      call chainbuild
      call etotal(energia(0))
      etot = energia(0)
      call enerprint(energia(0))
C Minimize the energy of the first conformation.
      if (minim) then
        call geom_to_var(nvar,varia)
        call minimize(etot,varia,iretcode,nfun)
        call etotal(energia(0))
        etot = energia(0)
        write (iout,'(/80(1h*)/a/80(1h*))') 
     &    'Results of the first energy minimization:'
        call enerprint(energia(0))
      endif
      if (refstr) then
       kkk=1
        call fitsq(rms,c(1,nstart_seq),cref(1,nstart_sup,kkk),
     &nsup,przes,
     &             obr,non_conv)
        rms=dsqrt(rms)
        call contact(.false.,ncont,icont,co)
        frac=contact_fract(ncont,ncont_ref,icont,icont_ref)
        write (iout,'(a,f8.3,a,f8.3,a,f8.3)') 
     &    'RMS deviation from the reference structure:',rms,
     &    ' % of native contacts:',frac*100,' contact order:',co
        write (istat,'(i5,11(1pe14.5))') 0,
     &   (energia(print_order(i)),i=1,nprint_ene),etot,rms,frac,co
      else
        write (istat,'(i5,9(1pe14.5))') 0,
     &   (energia(print_order(i)),i=1,nprint_ene),etot
      endif
      close(istat)
      neneval=neneval+nfun+1
      if (.not. ent_read) then
C Initialize the entropy array
        do i=-max_ene,max_ene
         emin=etot
C Uncomment the line below for actual entropic sampling (start with uniform
C energy distribution).
c        entropy(i)=0.0D0
C Uncomment the line below for multicanonical sampling (start with Boltzmann
C distribution).
         entropy(i)=(emin+i*delte)*betbol 
        enddo
        emax=10000000.0D0
        emin=etot
        write (iout,'(/a)') 'Initial entropy'
        do i=indminn,indmaxx
          write (iout,'(i5,2f10.5)') i,emin+i*delte,entropy(i)
        enddo
      endif ! ent_read
#ifdef MPL
      call recv_stop_sig(Kwita)
      if (whatsup.eq.1) then
        call send_stop_sig(-2)
        not_done=.false.
      else if (whatsup.le.-2) then
        not_done=.false.
      else if (whatsup.eq.2) then
        not_done=.false.
      else 
        not_done=.true.
      endif
#else
      not_done = (iretcode.ne.11)
#endif 
      write (iout,'(/80(1h*)/20x,a/80(1h*))')
     &    'Enter Monte Carlo procedure.'
      close(igeom)
      call briefout(0,etot)
      do i=1,nvar
        varold(i)=varia(i)
      enddo
      eold=etot
      indeold=(eold-emin)/delte
      deix=eold-(emin+indeold*delte)
      dent=entropy(indeold+1)-entropy(indeold)
cd    write (iout,*) 'indeold=',indeold,' deix=',deix,' dent=',dent
cd    write (*,*) 'Processor',MyID,' indeold=',indeold,' deix=',deix,
cd   & ' dent=',dent
      sold=entropy(indeold)+(dent/delte)*deix
      elowest=etot
      write (iout,*) 'eold=',eold,' sold=',sold,' elowest=',etot
      write (*,*) 'Processor',MyID,' eold=',eold,' sold=',sold,
     & ' elowest=',etot
      if (minim) call zapis(varia,etot)
      nminima(1)=1.0D0
C NACC is the counter for the accepted conformations of a given processor
      nacc=0
C NACC_TOT counts the total number of accepted conformations
      nacc_tot=0
#ifdef MPL
      if (MyID.eq.MasterID) then
        call receive_MCM_info
      else
        call send_MCM_info(2)
      endif
#endif
      do iene=indminn,indmaxx
        nhist(iene)=0.0D0
      enddo
      do i=2,maxsave
        nminima(i)=0.0D0
      enddo
C Main loop.
c----------------------------------------------------------------------------
      elowest=1.0D10
      ehighest=-1.0D10
      it=0
      do while (not_done)
        it=it+1
        if (print_mc.gt.0) write (iout,'(80(1h*)/20x,a,i7)')
     &                             'Beginning iteration #',it
C Initialize local counter.
        ntrial=0 ! # of generated non-overlapping confs.
        noverlap=0 ! # of overlapping confs.
        accepted=.false.
        do while (.not. accepted .and. WhatsUp.eq.0 .and. Kwita.eq.0)
          ntrial=ntrial+1
C Retrieve the angles of previously accepted conformation
          do j=1,nvar
            varia(j)=varold(j)
          enddo
cd        write (iout,'(a)') 'Old variables:'
cd        write (iout,'(10f8.1)') (rad2deg*varia(i),i=1,nvar)
          call var_to_geom(nvar,varia)
C Rebuild the chain.
          call chainbuild
          MoveType=0
          nbond=0
          lprint=.true.
C Decide whether to generate a random conformation or perturb the old one
          RandOrPert=ran_number(0.0D0,1.0D0)
          if (RandOrPert.gt.RanFract) then
            if (print_mc.gt.0) 
     &        write (iout,'(a)') 'Perturbation-generated conformation.'
            call perturb(error,lprint,MoveType,nbond,1.0D0)
            if (error) goto 20
            if (MoveType.lt.1 .or. MoveType.gt.MaxMoveType) then
              write (iout,'(/a,i7,a/)') 'Error - unknown MoveType=',
     &           MoveType,' returned from PERTURB.'
              goto 20
            endif
            call chainbuild
          else
            MoveType=0
            moves(0)=moves(0)+1
            nstart_grow=iran_num(3,nres)
            if (print_mc.gt.0) 
     &        write (iout,'(2a,i3)') 'Random-generated conformation',
     &        ' - chain regrown from residue',nstart_grow
            call gen_rand_conf(nstart_grow,*30)
          endif
          call geom_to_var(nvar,varia)
cd        write (iout,'(a)') 'New variables:'
cd        write (iout,'(10f8.1)') (rad2deg*varia(i),i=1,nvar)
          ngen=ngen+1
          if (print_mc.gt.0) write (iout,'(a,i5,a,i10,a,i10)') 
     &   'Processor',MyId,' trial move',ntrial,' total generated:',ngen
          if (print_mc.gt.0) write (*,'(a,i5,a,i10,a,i10)') 
     &   'Processor',MyId,' trial move',ntrial,' total generated:',ngen
          call etotal(energia(0))
          etot = energia(0)
c         call enerprint(energia(0))
c         write (iout,'(2(a,1pe14.5))') 'Etot=',Etot,' Elowest=',Elowest
          if (etot-elowest.gt.overlap_cut) then
            write (iout,'(a,i5,a,1pe14.5)')  'Iteration',it,
     &      ' Overlap detected in the current conf.; energy is',etot
            neneval=neneval+1 
            accepted=.false.
            noverlap=noverlap+1
            if (noverlap.gt.maxoverlap) then
              write (iout,'(a)') 'Too many overlapping confs.'
              goto 20
            endif
          else
            if (minim) then
              call minimize(etot,varia,iretcode,nfun)
cd            write (iout,'(a)') 'Variables after minimization:'
cd            write (iout,'(10f8.1)') (rad2deg*varia(i),i=1,nvar)
              call etotal(energia(0))
              etot = energia(0)
              neneval=neneval+nfun+1
            endif
            if (print_mc.gt.2) then
              write (iout,'(a)') 'Total energies of trial conf:'
              call enerprint(energia(0))
            else if (print_mc.eq.1) then
               write (iout,'(a,i6,a,1pe16.6)') 
     &         'Trial conformation:',ngen,' energy:',etot
            endif 
C--------------------------------------------------------------------------
C... Acceptance test
C--------------------------------------------------------------------------
            accepted=.false.
            if (WhatsUp.eq.0) 
     &        call accepting(etot,eold,scur,sold,varia,varold,
     &                       accepted)
            if (accepted) then
              nacc=nacc+1
              nacc_tot=nacc_tot+1
              if (elowest.gt.etot) elowest=etot
              if (ehighest.lt.etot) ehighest=etot
              moves_acc(MoveType)=moves_acc(MoveType)+1
              if (MoveType.eq.1) then
                nbond_acc(nbond)=nbond_acc(nbond)+1
              endif
C Check against conformation repetitions.
              irep=conf_comp(varia,etot)
#if defined(AIX) || defined(PGI)
              open (istat,file=statname,position='append')
#else
              open (istat,file=statname,access='append')
#endif
              if (refstr) then
                kkk=1
                call fitsq(rms,c(1,nstart_seq),cref(1,nstart_sup,kkk),
     &                     nsup,
     &                      przes,obr,non_conv)
                rms=dsqrt(rms)
                call contact(.false.,ncont,icont,co)
                frac=contact_fract(ncont,ncont_ref,icont,icont_ref)
                if (print_mc.gt.0) 
     &          write (iout,'(a,f8.3,a,f8.3,a,f8.3)') 
     &          'RMS deviation from the reference structure:',rms,
     &          ' % of native contacts:',frac*100,' contact order:',co
                if (print_stat)
     &          write (istat,'(i5,11(1pe14.5))') it,
     &           (energia(print_order(i)),i=1,nprint_ene),etot,
     &           rms,frac,co
              elseif (print_stat) then
                write (istat,'(i5,10(1pe14.5))') it,
     &             (energia(print_order(i)),i=1,nprint_ene),etot
              endif  
              close(istat)
              if (print_mc.gt.1) 
     &          call statprint(nacc,nfun,iretcode,etot,elowest)
C Print internal coordinates.
              if (print_int) call briefout(nacc,etot)
#ifdef MPL
              if (MyID.ne.MasterID) then
                call recv_stop_sig(Kwita)
cd              print *,'Processor:',MyID,' STOP=',Kwita
                if (irep.eq.0) then
                  call send_MCM_info(2)
                else
                  call send_MCM_info(1)
                endif
              endif
#endif
C Store the accepted conf. and its energy.
              eold=etot
              sold=scur
              do i=1,nvar
                varold(i)=varia(i)
              enddo
              if (irep.eq.0) then
                irep=nsave+1
cd              write (iout,*) 'Accepted conformation:'
cd              write (iout,*) (rad2deg*varia(i),i=1,nphi)
                if (minim) call zapis(varia,etot)
                do i=1,n_ene
                  ener(i,nsave)=energia(i)
                enddo
                ener(n_ene+1,nsave)=etot
                ener(n_ene+2,nsave)=frac
              endif
              nminima(irep)=nminima(irep)+1.0D0
c             print *,'irep=',irep,' nminima=',nminima(irep)
#ifdef MPL
              if (Kwita.eq.0) call recv_stop_sig(kwita)
#endif
            endif ! accepted
          endif ! overlap
#ifdef MPL
          if (MyID.eq.MasterID) then
            call receive_MCM_info
            if (nacc_tot.ge.maxacc) accepted=.true.
          endif
#endif
          if (ntrial.gt.maxtrial_iter .and. npool.gt.0) then
C Take a conformation from the pool
            ii=iran_num(1,npool)
            do i=1,nvar
              varia(i)=xpool(i,ii)
            enddo
            write (iout,*) 'Iteration',it,' max. # of trials exceeded.'
            write (iout,*) 
     &     'Take conformation',ii,' from the pool energy=',epool(ii)
            if (print_mc.gt.2)
     &      write (iout,'(10f8.3)') (rad2deg*varia(i),i=1,nvar)
            ntrial=0
         endif ! (ntrial.gt.maxtrial_iter .and. npool.gt.0)
   30    continue
        enddo ! accepted
#ifdef MPL
        if (MyID.eq.MasterID) then
          call receive_MCM_info
        endif
        if (Kwita.eq.0) call recv_stop_sig(kwita)
#endif
        if (ovrtim()) WhatsUp=-1
cd      write (iout,*) 'WhatsUp=',WhatsUp,' Kwita=',Kwita
        not_done = (nacc_tot.lt.maxacc) .and. (WhatsUp.eq.0) 
     &         .and. (Kwita.eq.0)
cd      write (iout,*) 'not_done=',not_done
#ifdef MPL
        if (Kwita.lt.0) then
          print *,'Processor',MyID,
     &    ' has received STOP signal =',Kwita,' in EntSamp.'
cd        print *,'not_done=',not_done
          if (Kwita.lt.-1) WhatsUp=Kwita
        else if (nacc_tot.ge.maxacc) then
          print *,'Processor',MyID,' calls send_stop_sig,',
     &     ' because a sufficient # of confs. have been collected.'
cd        print *,'not_done=',not_done
          call send_stop_sig(-1)
        else if (WhatsUp.eq.-1) then
          print *,'Processor',MyID,
     &               ' calls send_stop_sig because of timeout.'
cd        print *,'not_done=',not_done
          call send_stop_sig(-2)
        endif
#endif
      enddo ! not_done

C-----------------------------------------------------------------
C... Construct energy histogram & update entropy
C-----------------------------------------------------------------
      go to 21
   20 WhatsUp=-3
#ifdef MPL
      write (iout,*) 'Processor',MyID,
     &       ' is broadcasting ERROR-STOP signal.'
      write (*,*) 'Processor',MyID,
     &       ' is broadcasting ERROR-STOP signal.'
      call send_stop_sig(-3)
#endif
   21 continue
#ifdef MPL
      if (MyID.eq.MasterID) then
c       call receive_MCM_results
        call receive_energies
#endif
      do i=1,nsave
        if (esave(i).lt.elowest) elowest=esave(i)
        if (esave(i).gt.ehighest) ehighest=esave(i)
      enddo
      write (iout,'(a,i10)') '# of accepted confs:',nacc_tot
      write (iout,'(a,f10.5,a,f10.5)') 'Lowest energy:',elowest,
     & ' Highest energy',ehighest
      if (isweep.eq.1 .and. .not.ent_read) then
        emin=elowest
        emax=ehighest
        write (iout,*) 'EMAX=',emax
        indminn=0
        indmaxx=(ehighest-emin)/delte
        indmin=indminn
        indmax=indmaxx
        do i=-max_ene,max_ene
          entropy(i)=(emin+i*delte)*betbol
        enddo
        ent_read=.true.
      else
        indmin=(elowest-emin)/delte
        indmax=(ehighest-emin)/delte
        if (indmin.lt.indminn) indminn=indmin
        if (indmax.gt.indmaxx) indmaxx=indmax
      endif
      write(iout,*)'indminn=',indminn,' indmaxx=',indmaxx
C Construct energy histogram
      do i=1,nsave
        inde=(esave(i)-emin)/delte
        nhist(inde)=nhist(inde)+nminima(i)
      enddo
C Update entropy (density of states)
      do i=indmin,indmax
        if (nhist(i).gt.0) then
          entropy(i)=entropy(i)+dlog(nhist(i)+0.0D0)
        endif
      enddo
Cd    do i=indmaxx+1
Cd      entropy(i)=1.0D+10
Cd    enddo
      write (iout,'(/80(1h*)/a,i2/80(1h*)/)') 
     &      'End of macroiteration',isweep
      write (iout,'(a,f10.5,a,f10.5)') 'Elowest=',elowest,
     &      ' Ehighest=',ehighest
      write (iout,'(a)') 'Frequecies of minima'
      do i=1,nsave
        write (iout,'(i5,f5.0,f10.5)') i,nminima(i),esave(i)
      enddo
      write (iout,'(/a)') 'Energy histogram'
      do i=indminn,indmaxx
        write (iout,'(i5,2f10.5)') i,emin+i*delte,nhist(i)
      enddo
      write (iout,'(/a)') 'Entropy'
      do i=indminn,indmaxx
        write (iout,'(i5,2f10.5)') i,emin+i*delte,entropy(i)
      enddo
C-----------------------------------------------------------------
C... End of energy histogram construction
C-----------------------------------------------------------------
#ifdef MPL
        entropy(-max_ene-4)=dfloat(indminn)
        entropy(-max_ene-3)=dfloat(indmaxx)
        entropy(-max_ene-2)=emin
        entropy(-max_ene-1)=emax
        call send_MCM_update
cd      print *,entname,ientout
        open (ientout,file=entname,status='unknown')
        write (ientout,'(2i5,2e25.17)') indminn,indmaxx,emin,emax
        do i=indminn,indmaxx
          write (ientout,'(i5,f10.5,f20.15)') i,emin+i*delte,entropy(i)
        enddo
        close(ientout)
      else
        write (iout,'(a)') 'Frequecies of minima'
        do i=1,nsave
          write (iout,'(i5,f5.0,f10.5)') i,nminima(i),esave(i)
        enddo
c       call send_MCM_results
        call send_energies
        call receive_MCM_update
        indminn=entropy(-max_ene-4)
        indmaxx=entropy(-max_ene-3)
        emin=entropy(-max_ene-2)
        emax=entropy(-max_ene-1)
        write (iout,*) 'Received from master:'
        write (iout,*) 'indminn=',indminn,' indmaxx=',indmaxx,
     &                 ' emin=',emin,' emax=',emax
        write (iout,'(/a)') 'Entropy'
        do i=indminn,indmaxx
          write (iout,'(i5,2f10.5)') i,emin+i*delte,entropy(i)
        enddo
      endif
      if (WhatsUp.lt.-1) return
#else
      if (ovrtim() .or. WhatsUp.lt.0) return
#endif

      write (iout,'(/80(1h*)/20x,a)') 'Summary run statistics:'
      call statprint(nacc,nfun,iretcode,etot,elowest)
      write (iout,'(a)') 
     & 'Statistics of multiple-bond motions. Total motions:' 
      write (iout,'(16i5)') (nbond_move(i),i=1,Nbm)
      write (iout,'(a)') 'Accepted motions:'
      write (iout,'(16i5)') (nbond_acc(i),i=1,Nbm)
      write (iout,'(a,i10)') 'Number of chain regrowths:',nregrow
      write (iout,'(a,i10)') 'Accepted chain regrowths:',nregrow_acc

C---------------------------------------------------------------------------
      ENDDO ! ISWEEP
C---------------------------------------------------------------------------

      runtime=tcpu()

      if (isweep.eq.nsweep .and. it.ge.maxacc)
     &write (iout,'(/80(1h*)/20x,a/80(1h*)/)') 'All iterations done.'
      return
      end
c------------------------------------------------------------------------------
      subroutine accepting(ecur,eold,scur,sold,x,xold,accepted)
      implicit real*8 (a-h,o-z)
      include 'DIMENSIONS'
      include 'COMMON.MCM'
      include 'COMMON.MCE'
      include 'COMMON.IOUNITS'
      include 'COMMON.VAR'
#ifdef MPL
      include 'COMMON.INFO'
#endif
      include 'COMMON.GEO'
      double precision ecur,eold,xx,ran_number,bol
      double precision x(maxvar),xold(maxvar)
      double precision tole /1.0D-1/, tola /5.0D0/
      logical accepted
C Check if the conformation is similar.
cd    write (iout,*) 'Enter ACCEPTING'
cd    write (iout,*) 'Old PHI angles:'
cd    write (iout,*) (rad2deg*xold(i),i=1,nphi)
cd    write (iout,*) 'Current angles'
cd    write (iout,*) (rad2deg*x(i),i=1,nphi)
cd    ddif=dif_ang(nphi,x,xold)
cd    write (iout,*) 'Angle norm:',ddif
cd    write (iout,*) 'ecur=',ecur,' emax=',emax
      if (ecur.gt.emax) then
        accepted=.false.
        if (print_mc.gt.0)
     & write (iout,'(a)') 'Conformation rejected as too high in energy'
        return
      else if (dabs(ecur-eold).lt.tole .and. 
     &       dif_ang(nphi,x,xold).lt.tola) then
        accepted=.false.
        if (print_mc.gt.0)
     & write (iout,'(a)') 'Conformation rejected as too similar'
        return
      endif
C Else evaluate the entropy of the conf and compare it with that of the previous
C one.
      indecur=(ecur-emin)/delte
      if (iabs(indecur).gt.max_ene) then
        write (iout,'(a,2i5)') 
     &   'Accepting: Index out of range:',indecur
        scur=1000.0D0 
      else if (indecur.eq.indmaxx) then
        scur=entropy(indecur)
        if (print_mc.gt.0) write (iout,*)'Energy boundary reached',
     &            indmaxx,indecur,entropy(indecur)
      else
        deix=ecur-(emin+indecur*delte)
        dent=entropy(indecur+1)-entropy(indecur)
        scur=entropy(indecur)+(dent/delte)*deix
      endif
cd    print *,'Processor',MyID,' ecur=',ecur,' indecur=',indecur,
cd   & ' scur=',scur,' eold=',eold,' sold=',sold
cd    print *,'deix=',deix,' dent=',dent,' delte=',delte
      if (print_mc.gt.1) then
        write(iout,*)'ecur=',ecur,' indecur=',indecur,' scur=',scur
        write(iout,*)'eold=',eold,' sold=',sold
      endif
      if (scur.le.sold) then
        accepted=.true.
      else
C Else carry out acceptance test
        xx=ran_number(0.0D0,1.0D0) 
        xxh=scur-sold
        if (xxh.gt.50.0D0) then
          bol=0.0D0
        else
          bol=exp(-xxh)
        endif
        if (bol.gt.xx) then
          accepted=.true. 
          if (print_mc.gt.0) write (iout,'(a)') 
     &    'Conformation accepted.'
        else
          accepted=.false.
          if (print_mc.gt.0) write (iout,'(a)') 
     & 'Conformation rejected.'
        endif
      endif
      return
      end 
c-----------------------------------------------------------------------------
      subroutine read_pool
      implicit real*8 (a-h,o-z)
      include 'DIMENSIONS'
      include 'COMMON.IOUNITS'
      include 'COMMON.GEO'
      include 'COMMON.MCM'
      include 'COMMON.MCE'
      include 'COMMON.VAR'
      double precision varia(maxvar)
      print '(a)','Call READ_POOL'
      do npool=1,max_pool
        print *,'i=',i
        read (intin,'(i5,f10.5)',end=10,err=10) iconf,epool(npool)
        if (epool(npool).eq.0.0D0) goto 10
        call read_angles(intin,*10)
        call geom_to_var(nvar,xpool(1,npool))
      enddo
      goto 11
   10 npool=npool-1
   11 write (iout,'(a,i5)') 'Number of pool conformations:',npool
      if (print_mc.gt.2) then
      do i=1,npool
        write (iout,'(a,i5,a,1pe14.5)') 'Pool conformation',i,' energy',
     &    epool(i)
        write (iout,'(10f8.3)') (rad2deg*xpool(j,i),j=1,nvar)
      enddo
      endif ! (print_mc.gt.2)
      return
      end