1d/equations/euler/rp/rec1eu.f

c
c ===================================================================
      subroutine rec1(maxm,meqn,mwaves,mbc,mx,q,method,mthlim,ql,qr)
c ===================================================================
c
c     # Reconstruction in primitive variables for the one-dimensional
c     # Euler equations. This reconstruction is not conservative!
c 
      implicit double precision (a-h,o-z)
c
      dimension    q(1-mbc:maxm+mbc, meqn)
      dimension   ql(1-mbc:maxm+mbc, meqn)
      dimension   qr(1-mbc:maxm+mbc, meqn)
      dimension method(7),mthlim(mwaves)
      common /param/  gamma,gamma1
c
      mlim = 0
      do 90 mw=1,mwaves
	 if (mthlim(mw) .gt. 0) then
            mlim = mthlim(mw)
            goto 95
         endif
 90   continue
 95   continue
c
c     Linear interpolation: om=0.d0, 2nd order: om!=0.d0,
c     3rd order: om=1.d0/3.d0
c
      om = 1.d0/3.d0     
      do 110 i=2-mbc,mx+mbc-1
c
c     # Reconstruction
         rho  = q(i  ,1)
         rho1 = q(i-1,1)
         rho2 = q(i+1,1)   
         u  = q(i  ,2)/rho
         u1 = q(i-1,2)/rho1
         u2 = q(i+1,2)/rho2
         p  = gamma1*(q(i  ,3) - 0.5d0*rho *u**2 )
         p1 = gamma1*(q(i-1,3) - 0.5d0*rho1*u1**2)
         p2 = gamma1*(q(i+1,3) - 0.5d0*rho2*u2**2)
c     
         call reclim(rho,rho1,rho2,mlim,om,rhol,rhor)  
         call reclim(u,u1,u2,mlim,om,ul,ur)  
         call reclim(p,p1,p2,mlim,om,pl,pr)  
c     
         ql(i,1) = rhol
         qr(i,1) = rhor
         ql(i,2) = ul*rhol
         qr(i,2) = ur*rhor
         ql(i,3) = pl/gamma1+0.5d0*rhol*ul**2
         qr(i,3) = pr/gamma1+0.5d0*rhor*ur**2
c     
 110  continue
c     
      return
      end
c
c