anchor.diff:            
             if (hse_bc_const == "density") then
                 
  -             ! constant density
  +             ! constant density in the ghost cells
                 
                ! zero gradient to rho, u
                ic = U%grid%lo
  
positive.code:     endif
      
  
      !-------------------------------------------------------------------------    
      ! lower boundary (-x)
      !-------------------------------------------------------------------------
      select case (U%grid%xlboundary)
  
negative.code:              
            else if (hse_bc_const == "temperature") then
               
               ! constant temperature (or internal energy) in the ghost
               ! cells
  
               ic = U%grid%lo
               
               econst = (U%data(ic,iuener) - &
                    0.5_dp_t*U%data(ic,iumomx)**2/U%data(ic,iudens))/U%data(ic,iudens)
               
               ! zero gradient for u
               vel = min(0.0_dp_t, U%data(ic,iumomx)/U%data(ic,iudens))
               
               ! get the pressure above (we already know that it has the
               ! internal energy econst)
               call eos(eos_input_e, p_above, econst, U%data(i+1, iudens))
  
               ! p via HSE.  Here we explicitly make use of the
               ! gamma-law nature of the EOS.  Otherwise, we would need
               ! to iterate to simultaneously satisfy HSE + the EOS.
               p = (p_above - 0.5_dp_t*U%data(i+1,iudens)*grav*U%grid%dx) / &
                    (1.0_dp_t + 0.5_dp_t*grav*U%grid%dx/(econst*(gamma - 1.0_dp_t)))
  
               ! now find the density that corresponds to this p, e
               call eos(eos_input_pe, p, econst, dens)
               
               U%data(i,iudens) = dens
               U%data(i,iumomx) = dens*vel
  
               ! compute E
               U%data(i,iuener) = dens*econst + &
                    0.5_dp_t*U%data(i,iumomx)**2/dens
  