        program pimc
c C. Pierleoni & S. Ciuchi, Phys. Dep. University of Aquila (ITALY)
c Path Integral MC code for Frohlich model hamiltonian in
c electron-phonon interaction problem
c put global comments about major version updates here
c 10/12/1997 possibility of many slices moves by driver omove 
c 04/02/1998 Path reconstruction by Levy flight algorithm: driver levy
c 04/02/1998 Ideal path sampling by Levy flight algorithm: driver levyid
c            Path reconstruction and CBMC : driver cbmc
c 09/06/1998 many particles version (boltzmannions) for bipolaron
c 23/09/1998 improved fourier regularized action (phonon2)
c 20/10/1998 improved fourier regularized action with gaussian int. (phonon3)
c 04/01/2000 new trial action for cumulant approx. (not for nparts>1)
c 07/04/2000 npol2a.cp 2dimensional extension
c 07/08/2000 npol2b.cp tabulation of action and potential functions
        implicit none
        include 'pimc.par'
        integer il,nb,l,js,jsp,i,is,k,ik
        integer n,m,js1,is1,it
        real*8 yself
        real*8 phonon,gaussn,coulomb,phonon2,phonon3
        real*4 etime,dtime,tim(2),te
        include 'pimc.cm'
c setting up the system and the run (read the input file)
        call setsys 
c
c    ** Zero average value **
        call zeroav(0)
c
c compute initial action
        call inaction
        ttanal=0.
        te=etime(tim)
c    ** Carry on the run **
        do nb=nin+1,nout
         ntry = 0
         nacc = 0
         ntrydsp = 0
         naccdsp = 0
         call zeroav(1)
         if (nstart.eq.-1) call readpc
         if (nstart.ge.0) then
c         if (nalg.eq.0) call omove
          if (nalg.eq.1) then
           if (ncumtr.eq.0) call levy
          endif
c         if (nalg.eq.2) call levyid
          if (nalg.eq.3) call cbmc
         endif
         if (np.eq.1.and.nacc.eq.0.and.nstart.ge.0) then
          write (*,1166) nb
 1166     format ('block #',i3,': no accepted moves !')
          stop
         endif
c update averages
         call analbl
         anorm(1)=anorm(1)+1.d0
c spill out data
         if(nalg.ne.2) call spill(nb)
c print out data
         call sprint(nb)
        enddo        ! over blocks
        tt=etime(tim)-te
        write(2,*)
        write(2,*)'Run completed:'
        write(2,*)'startup time (sec.) = ',te
        write(2,*)'run time (sec.) = ',tt
        write(2,*)'analysis time (sec.) = ',ttanal
        write(2,*)'total time (sec.) = ',tt+te
        call rmaget(0)
        stop
        end
