anchor.diff:        double precision qvec(meqn),flux(meqn)
         double precision auxvec_center(maux), auxvec_edge(maux)
   
  -      integer i,j,k,m
  +      integer i,j,k,m, idir, iface
   
         do j = 1,my
   c        # left side k=0 = in; k=1 = out
           idir = 0
  +        iface = 0
            do k = 0,1
               do m = 1,maux
   c              # Cell centered values;  Each cell-centered value
  
positive.code:             do m = 1,meqn
                 qvec(m) = q(1-k,j,m)
              enddo
              call rpn2_cons(meqn,maux,idir,qvec,
       &               auxvec_center,auxvec_edge,flux)         
              do m = 1,meqn
                 flux0(j,m,k) = flux0(j,m,k) + dt*el0(j)*flux(m)
  
negative.code:       do i = 1,mx
  c        # bottom side 0 = in; 1 = out
           idir = 1
           iface = 0
           do k = 0,1
              do m = 1,maux
  c              # Cell centered values                
  