estimate_speeds_toro_pvrs Subroutine

private subroutine estimate_speeds_toro_pvrs(uL, uR, xC, flux_dim, sL, sR)

Wave speed estimates for HLL/HLLC using Toro (2010) PVRS pressure estimate Implements Eq. (10.67)-(10.69)

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(in) :: uL(nw_phys)
real(kind=dp), intent(in) :: uR(nw_phys)
real(kind=dp), intent(in) :: xC(ndim)
integer, intent(in) :: flux_dim
real(kind=dp), intent(out) :: sL
real(kind=dp), intent(out) :: sR

Called by

proc~~estimate_speeds_toro_pvrs~2~~CalledByGraph proc~estimate_speeds_toro_pvrs~2 estimate_speeds_toro_pvrs proc~riemann_hllc_prim~2 riemann_hllc_prim proc~riemann_hllc_prim~2->proc~estimate_speeds_toro_pvrs~2 proc~reconflux_muscl_hllc_prim~2 reconflux_muscl_hllc_prim proc~reconflux_muscl_hllc_prim~2->proc~riemann_hllc_prim~2 proc~finite_volume_local_muscl_hllc~2 finite_volume_local_muscl_hllc proc~finite_volume_local_muscl_hllc~2->proc~reconflux_muscl_hllc_prim~2 proc~finite_volume_local~2 finite_volume_local proc~finite_volume_local~2->proc~finite_volume_local_muscl_hllc~2

Source Code

subroutine estimate_speeds_toro_pvrs(uL, uR, xC, flux_dim, sL, sR)
  
  real(dp), intent(in)  :: uL(nw_phys), uR(nw_phys)
  real(dp), intent(in)  :: xC(ndim)
  integer,  intent(in)  :: flux_dim
  real(dp), intent(out) :: sL, sR

  real(dp) :: rhoL, rhoR, pL, pR, unL, unR
  real(dp) :: aL, aR, abar, rhobar
  real(dp) :: ppvrs, pstar
  real(dp) :: qL, qR, pratio
  real(dp), parameter :: tiny = 1e-30_dp

  rhoL = uL(iw_rho);  rhoR = uR(iw_rho)
  pL   = uL(iw_e);    pR   = uR(iw_e)
  unL  = uL(iw_mom(flux_dim))
  unR  = uR(iw_mom(flux_dim))

  ! sound speeds
  aL = sqrt(hd_gamma * pL / max(rhoL, tiny))
  aR = sqrt(hd_gamma * pR / max(rhoR, tiny))

  ! PVRS pressure estimate (Eq. 10.67)
  rhobar = 0.5_dp*(rhoL + rhoR)
  abar   = 0.5_dp*(aL + aR)
  ppvrs  = 0.5_dp*(pL + pR) - 0.5_dp*(unR - unL)*rhobar*abar
  pstar  = max(0._dp, ppvrs)

  ! (Eq. 10.69)
  if (pstar <= pL) then
    qL = 1._dp
  else
    pratio = pstar / max(pL, tiny)
    qL = sqrt(1._dp + 0.5_dp*(hd_gamma + 1._dp)/hd_gamma * (pratio - 1._dp))
  end if

  if (pstar <= pR) then
    qR = 1._dp
  else
    pratio = pstar / max(pR, tiny)
    qR = sqrt(1._dp + 0.5_dp*(hd_gamma + 1._dp)/hd_gamma * (pratio - 1._dp))
  end if

  sL = unL - aL*qL
  sR = unR + aR*qR
end subroutine estimate_speeds_toro_pvrs