compute_roe_state Subroutine

   subroutine compute_roe_state(eos_left, state_left, eos_right, state_right, r_d, r_u, r_e, r_a)
   !< Evaluate the intermediate state from the known states U1,U4 using the Roe linearization.
   class(eos_object),                      intent(in)    :: eos_left     !< Equation of state for left state.
   class(conservative_object),             intent(in)    :: state_left   !< Left Riemann state.
   class(eos_object),                      intent(in)    :: eos_right    !< Equation of state for right state.
   class(conservative_object),             intent(in)    :: state_right  !< Right Riemann state.
   real(R8P),                              intent(out)   :: r_d          !< Roe intermediate state density.
   type(vector),                           intent(out)   :: r_u          !< Roe intermediate state velocity vector..
   real(R8P),                              intent(out)   :: r_e          !< Roe intermediate state enthalpy.
   real(R8P),                              intent(out)   :: r_a          !< Roe intermediate state sound speed.
   real(R8P)                                             :: x, omx       !< x = sqrt(r1)/(sqrt(r1)+sqrt(r4)),  omx = 1-x
   real(R8P)                                             :: cp, cv       !< Roe intermediate state Cp and Cv.
   type(conservative_compressible), pointer              :: state_left_  !< Left Riemann state, local variable.
   type(conservative_compressible), pointer              :: state_right_ !< Right Riemann state, local variable.
   type(eos_compressible),          pointer              :: eos_left_    !< Left Riemann state, local variable.
   type(eos_compressible),          pointer              :: eos_right_   !< Right Riemann state, local variable.

   state_left_ => conservative_compressible_pointer(to=state_left)
   state_right_ => conservative_compressible_pointer(to=state_right)
   eos_left_ => eos_compressible_pointer(to=eos_left)
   eos_right_ => eos_compressible_pointer(to=eos_right)
   x = sqrt(state_left_%density)/(sqrt(state_left_%density)+sqrt(state_right_%density)) ; omx = 1._R8P - x

   r_d   = sqrt(state_left_%density*state_right_%density)
   r_u%x = state_left_%momentum%x / state_left_%density * x + state_right_%momentum%x / state_right_%density * omx
   r_u%y = state_left_%momentum%y / state_left_%density * x + state_right_%momentum%y / state_right_%density * omx
   r_u%z = state_left_%momentum%z / state_left_%density * x + state_right_%momentum%z / state_right_%density * omx
   r_e   = (state_left_%energy  + state_left_%pressure(eos_left))   / state_left_%density  * x + &
           (state_right_%energy + state_right_%pressure(eos_right)) / state_right_%density * omx
   cp    = eos_left_%cp_ * x + eos_right_%cp_ * omx
   cv    = eos_left_%cv_ * x + eos_right_%cv_ * omx
   r_a   = sqrt((cp/cv - 1._R8P)*(r_e - 0.5_R8P * r_u%sq_norm()))

   endsubroutine compute_roe_state