Numerically stable logarithmic mean (Ismail-Roe 2009 / Ranocha 2018).
_log = (a - b) / (ln(a) - ln(b)) for a /= b _log = (a + b) / 2 for a == b
Falls back to a Taylor series in u = f^2 with f = (a-b)/(a+b) when a and b are close, to avoid 0/0.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=prec), | intent(in) | :: | a | |||
| real(kind=prec), | intent(in) | :: | b |
pure function log_mean(a,b) result(am)
!! Numerically stable logarithmic mean (Ismail-Roe 2009 / Ranocha 2018).
!!
!! <a>_log = (a - b) / (ln(a) - ln(b)) for a /= b
!! <a>_log = (a + b) / 2 for a == b
!!
!! Falls back to a Taylor series in u = f^2 with f = (a-b)/(a+b)
!! when a and b are close, to avoid 0/0.
real(prec),intent(in) :: a,b
real(prec) :: am
! Local
real(prec) :: zeta,f,u,F_F
real(prec),parameter :: eps_lm = 1.0e-2_prec
zeta = a/b
f = (zeta-1.0_prec)/(zeta+1.0_prec)
u = f*f
if(u < eps_lm) then
F_F = 1.0_prec+u/3.0_prec+u*u/5.0_prec+u*u*u/7.0_prec
else
F_F = log(zeta)/(2.0_prec*f)
endif
am = (a+b)/(2.0_prec*F_F)
endfunction log_mean