gkowal/amun-code
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

SCHEMES: Make lmean() numerically symmetric.

Browse Source 
Signed-off-by: Grzegorz Kowal <grzegorz@amuncode.org>
This commit is contained in:
Grzegorz Kowal 2022-01-03 11:38:08 -03:00
parent e6cca5468e
commit 894a3362a2
1 changed files with 13 additions and 6 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

19
sources/schemes.F90
View File

@ -4110,19 +4110,26 @@ module schemes
real(kind=8), intent(in) :: l, r
real(kind=8) :: u, rr
real(kind=8) :: u, d, s
real(kind=8), parameter :: c1 = 2.0d+00, c2 = 2.0d+00 / 3.0d+00, &
c3 = 4.0d-01, c4 = 2.0d+00 / 7.0d+00
!-------------------------------------------------------------------------------
!
rr = r * r
u = (l * (l - 2.0d+00 * r) + rr) / (l * (l + 2.0d+00 * r) + rr)
if (u < 1.0d-04) then
lmean = (l + r) / (c1 + u * (c2 + u * (c3 + u * c4)))
d = r - l
s = r + l
u = d / s
if (u < 1.0d-02) then
u = u**2
lmean = s / (c1 + u * (c2 + u * (c3 + u * c4)))
else
lmean = (r - l) / log(r / l)
if (d >= 0.0d+00) then
s = log(r / l)
else
s = log(l / r)
end if
lmean = abs(d) / s
end if
!-------------------------------------------------------------------------------