Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

HLLE CEE 2D3D NonCartesian Meshes #1692

Merged
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
25 commits
Select commit Hold shift + click to select a range
487d08d
HLLE CEE 2D3D NonCartesian Meshes
DanielDoehring Oct 27, 2023
653e2cc
format
DanielDoehring Oct 27, 2023
dfc538b
Merge branch 'main' into HLLE_CEE_2D3D_NonCartesian
DanielDoehring Oct 27, 2023
fd8a922
hlle via hll
DanielDoehring Oct 30, 2023
8dbe9b3
format test
DanielDoehring Oct 30, 2023
55c4326
Merge branch 'HLLE_CEE_2D3D_NonCartesian' of github.com:DanielDoehrin…
DanielDoehring Oct 30, 2023
0f8dd92
format test
DanielDoehring Oct 30, 2023
ccf55a5
Merge branch 'main' into HLLE_CEE_2D3D_NonCartesian
DanielDoehring Oct 30, 2023
aca8cb7
format
DanielDoehring Oct 30, 2023
6e10fa5
do not export hlle
DanielDoehring Oct 30, 2023
7d089e5
Correct test vals
DanielDoehring Oct 31, 2023
c34dc6e
Merge branch 'main' into HLLE_CEE_2D3D_NonCartesian
DanielDoehring Oct 31, 2023
c60f434
test values CI
DanielDoehring Oct 31, 2023
05ddef4
Merge branch 'HLLE_CEE_2D3D_NonCartesian' of github.com:DanielDoehrin…
DanielDoehring Oct 31, 2023
c11d008
Update src/equations/compressible_euler_2d.jl
DanielDoehring Oct 31, 2023
ad17e0e
Update src/equations/compressible_euler_1d.jl
DanielDoehring Oct 31, 2023
2bb19d0
Update src/equations/compressible_euler_2d.jl
DanielDoehring Oct 31, 2023
3825ce1
Update src/equations/compressible_euler_3d.jl
DanielDoehring Oct 31, 2023
2a08efa
Update src/equations/compressible_euler_3d.jl
DanielDoehring Oct 31, 2023
a0c0221
apply suggestions
DanielDoehring Oct 31, 2023
fb3dfc9
additional sentence
DanielDoehring Oct 31, 2023
81bd303
Fix typo
DanielDoehring Oct 31, 2023
d35503f
typos
DanielDoehring Oct 31, 2023
1fc62b9
correct test vals
DanielDoehring Oct 31, 2023
dcb5185
Merge branch 'main' into HLLE_CEE_2D3D_NonCartesian
ranocha Oct 31, 2023
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
36 changes: 4 additions & 32 deletions src/equations/compressible_euler_1d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -808,7 +808,7 @@ function flux_hllc(u_ll, u_rr, orientation::Integer,
end

"""
flux_hlle(u_ll, u_rr, orientation, equations::CompressibleEulerEquations1D)
min_max_speed_einfeldt(u_ll, u_rr, orientation, equations::CompressibleEulerEquations1D)

Computes the HLLE (Harten-Lax-van Leer-Einfeldt) flux for the compressible Euler equations.
Special estimates of the signal velocites and linearization of the Riemann problem developed
Expand All @@ -825,8 +825,8 @@ Compactly summarized:
Numerical methods for conservation laws and related equations.
[Link](https://metaphor.ethz.ch/x/2019/hs/401-4671-00L/literature/mishra_hyperbolic_pdes.pdf)
"""
function flux_hlle(u_ll, u_rr, orientation::Integer,
equations::CompressibleEulerEquations1D)
@inline function min_max_speed_einfeldt(u_ll, u_rr, orientation::Integer,
equations::CompressibleEulerEquations1D)
# Calculate primitive variables, enthalpy and speed of sound
rho_ll, v_ll, p_ll = cons2prim(u_ll, equations)
rho_rr, v_rr, p_rr = cons2prim(u_rr, equations)
Expand Down Expand Up @@ -858,35 +858,7 @@ function flux_hlle(u_ll, u_rr, orientation::Integer,
SsL = min(v_roe - c_roe, v_ll - beta * c_ll, zero(v_roe))
SsR = max(v_roe + c_roe, v_rr + beta * c_rr, zero(v_roe))

if SsL >= 0.0 && SsR > 0.0
# Positive supersonic speed
f_ll = flux(u_ll, orientation, equations)

f1 = f_ll[1]
f2 = f_ll[2]
f3 = f_ll[3]
elseif SsR <= 0.0 && SsL < 0.0
# Negative supersonic speed
f_rr = flux(u_rr, orientation, equations)

f1 = f_rr[1]
f2 = f_rr[2]
f3 = f_rr[3]
else
# Subsonic case
# Compute left and right fluxes
f_ll = flux(u_ll, orientation, equations)
f_rr = flux(u_rr, orientation, equations)

f1 = (SsR * f_ll[1] - SsL * f_rr[1] + SsL * SsR * (u_rr[1] - u_ll[1])) /
(SsR - SsL)
f2 = (SsR * f_ll[2] - SsL * f_rr[2] + SsL * SsR * (u_rr[2] - u_ll[2])) /
(SsR - SsL)
f3 = (SsR * f_ll[3] - SsL * f_rr[3] + SsL * SsR * (u_rr[3] - u_ll[3])) /
(SsR - SsL)
end

return SVector(f1, f2, f3)
return SsL, SsR
end

@inline function max_abs_speeds(u, equations::CompressibleEulerEquations1D)
Expand Down
92 changes: 59 additions & 33 deletions src/equations/compressible_euler_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1253,7 +1253,7 @@ function flux_hllc(u_ll, u_rr, orientation::Integer,
end

"""
flux_hlle(u_ll, u_rr, orientation, equations::CompressibleEulerEquations2D)
min_max_speed_einfeldt(u_ll, u_rr, orientation, equations::CompressibleEulerEquations2D)

Computes the HLLE (Harten-Lax-van Leer-Einfeldt) flux for the compressible Euler equations.
Special estimates of the signal velocites and linearization of the Riemann problem developed
Expand All @@ -1267,8 +1267,8 @@ of the numerical flux.
On Godunov-type methods near low densities.
[DOI: 10.1016/0021-9991(91)90211-3](https://doi.org/10.1016/0021-9991(91)90211-3)
"""
function flux_hlle(u_ll, u_rr, orientation::Integer,
equations::CompressibleEulerEquations2D)
@inline function min_max_speed_einfeldt(u_ll, u_rr, orientation::Integer,
equations::CompressibleEulerEquations2D)
# Calculate primitive variables, enthalpy and speed of sound
rho_ll, v1_ll, v2_ll, p_ll = cons2prim(u_ll, equations)
rho_rr, v1_rr, v2_rr, p_rr = cons2prim(u_rr, equations)
Expand Down Expand Up @@ -1306,39 +1306,65 @@ function flux_hlle(u_ll, u_rr, orientation::Integer,
SsR = max(v2_roe + c_roe, v2_rr + beta * c_rr, zero(v2_roe))
end

if SsL >= 0.0 && SsR > 0.0
# Positive supersonic speed
f_ll = flux(u_ll, orientation, equations)
return SsL, SsR
end

f1 = f_ll[1]
f2 = f_ll[2]
f3 = f_ll[3]
f4 = f_ll[4]
elseif SsR <= 0.0 && SsL < 0.0
# Negative supersonic speed
f_rr = flux(u_rr, orientation, equations)
"""
min_max_speed_einfeldt(u_ll, u_rr, normal_direction, equations::CompressibleEulerEquations2D)

f1 = f_rr[1]
f2 = f_rr[2]
f3 = f_rr[3]
f4 = f_rr[4]
else
# Subsonic case
# Compute left and right fluxes
f_ll = flux(u_ll, orientation, equations)
f_rr = flux(u_rr, orientation, equations)

f1 = (SsR * f_ll[1] - SsL * f_rr[1] + SsL * SsR * (u_rr[1] - u_ll[1])) /
(SsR - SsL)
f2 = (SsR * f_ll[2] - SsL * f_rr[2] + SsL * SsR * (u_rr[2] - u_ll[2])) /
(SsR - SsL)
f3 = (SsR * f_ll[3] - SsL * f_rr[3] + SsL * SsR * (u_rr[3] - u_ll[3])) /
(SsR - SsL)
f4 = (SsR * f_ll[4] - SsL * f_rr[4] + SsL * SsR * (u_rr[4] - u_ll[4])) /
(SsR - SsL)
end
Computes the HLLE (Harten-Lax-van Leer-Einfeldt) flux for the compressible Euler equations.
Special estimates of the signal velocites and linearization of the Riemann problem developed
by Einfeldt to ensure that the internal energy and density remain positive during the computation
of the numerical flux.

return SVector(f1, f2, f3, f4)
- Bernd Einfeldt (1988)
On Godunov-type methods for gas dynamics.
[DOI: 10.1137/0725021](https://doi.org/10.1137/0725021)
- Bernd Einfeldt, Claus-Dieter Munz, Philip L. Roe and Björn Sjögreen (1991)
On Godunov-type methods near low densities.
[DOI: 10.1016/0021-9991(91)90211-3](https://doi.org/10.1016/0021-9991(91)90211-3)
"""
@inline function min_max_speed_einfeldt(u_ll, u_rr, normal_direction::AbstractVector,
equations::CompressibleEulerEquations2D)
# Calculate primitive variables, enthalpy and speed of sound
rho_ll, v1_ll, v2_ll, p_ll = cons2prim(u_ll, equations)
rho_rr, v1_rr, v2_rr, p_rr = cons2prim(u_rr, equations)

v_dot_n_ll = v1_ll * normal_direction[1] + v2_ll * normal_direction[2]
v_dot_n_rr = v1_rr * normal_direction[1] + v2_rr * normal_direction[2]

norm_ = norm(normal_direction)

# `u_ll[4]` is total energy `rho_e_ll` on the left
H_ll = (u_ll[4] + p_ll) / rho_ll
c_ll = sqrt(equations.gamma * p_ll / rho_ll) * norm_

# `u_rr[4]` is total energy `rho_e_rr` on the right
H_rr = (u_rr[4] + p_rr) / rho_rr
c_rr = sqrt(equations.gamma * p_rr / rho_rr) * norm_

# Compute Roe averages
sqrt_rho_ll = sqrt(rho_ll)
sqrt_rho_rr = sqrt(rho_rr)
inv_sum_sqrt_rho = inv(sqrt_rho_ll + sqrt_rho_rr)

v1_roe = (sqrt_rho_ll * v1_ll + sqrt_rho_rr * v1_rr) * inv_sum_sqrt_rho
v2_roe = (sqrt_rho_ll * v2_ll + sqrt_rho_rr * v2_rr) * inv_sum_sqrt_rho
v_roe = v1_roe * normal_direction[1] + v2_roe * normal_direction[2]
v_roe_mag = (v1_roe * normal_direction[1])^2 + (v2_roe * normal_direction[2])^2

H_roe = (sqrt_rho_ll * H_ll + sqrt_rho_rr * H_rr) * inv_sum_sqrt_rho
c_roe = sqrt((equations.gamma - 1) * (H_roe - 0.5 * v_roe_mag)) * norm_

# Compute convenience constant for positivity preservation, see
# https://doi.org/10.1016/0021-9991(91)90211-3
beta = sqrt(0.5 * (equations.gamma - 1) / equations.gamma)

# Estimate the edges of the Riemann fan (with positivity conservation)
SsL = min(v_roe - c_roe, v_dot_n_ll - beta * c_ll, zero(v_roe))
SsR = max(v_roe + c_roe, v_dot_n_rr + beta * c_rr, zero(v_roe))

return SsL, SsR
end

@inline function max_abs_speeds(u, equations::CompressibleEulerEquations2D)
Expand Down
100 changes: 63 additions & 37 deletions src/equations/compressible_euler_3d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1322,7 +1322,7 @@ function flux_hllc(u_ll, u_rr, orientation::Integer,
end

"""
flux_hlle(u_ll, u_rr, orientation, equations::CompressibleEulerEquations3D)
min_max_speed_einfeldt(u_ll, u_rr, orientation, equations::CompressibleEulerEquations3D)

Computes the HLLE (Harten-Lax-van Leer-Einfeldt) flux for the compressible Euler equations.
Special estimates of the signal velocites and linearization of the Riemann problem developed
Expand All @@ -1336,8 +1336,8 @@ of the numerical flux.
On Godunov-type methods near low densities.
[DOI: 10.1016/0021-9991(91)90211-3](https://doi.org/10.1016/0021-9991(91)90211-3)
"""
function flux_hlle(u_ll, u_rr, orientation::Integer,
equations::CompressibleEulerEquations3D)
@inline function min_max_speed_einfeldt(u_ll, u_rr, orientation::Integer,
equations::CompressibleEulerEquations3D)
# Calculate primitive variables, enthalpy and speed of sound
rho_ll, v1_ll, v2_ll, v3_ll, p_ll = cons2prim(u_ll, equations)
rho_rr, v1_rr, v2_rr, v3_rr, p_rr = cons2prim(u_rr, equations)
Expand Down Expand Up @@ -1379,43 +1379,69 @@ function flux_hlle(u_ll, u_rr, orientation::Integer,
SsR = max(v3_roe + c_roe, v3_rr + beta * c_rr, zero(v3_roe))
end

if SsL >= 0.0 && SsR > 0.0
# Positive supersonic speed
f_ll = flux(u_ll, orientation, equations)
return SsL, SsR
end

f1 = f_ll[1]
f2 = f_ll[2]
f3 = f_ll[3]
f4 = f_ll[4]
f5 = f_ll[5]
elseif SsR <= 0.0 && SsL < 0.0
# Negative supersonic speed
f_rr = flux(u_rr, orientation, equations)
"""
min_max_speed_einfeldt(u_ll, u_rr, normal_direction, equations::CompressibleEulerEquations3D)

f1 = f_rr[1]
f2 = f_rr[2]
f3 = f_rr[3]
f4 = f_rr[4]
f5 = f_rr[5]
else
# Subsonic case
# Compute left and right fluxes
f_ll = flux(u_ll, orientation, equations)
f_rr = flux(u_rr, orientation, equations)

f1 = (SsR * f_ll[1] - SsL * f_rr[1] + SsL * SsR * (u_rr[1] - u_ll[1])) /
(SsR - SsL)
f2 = (SsR * f_ll[2] - SsL * f_rr[2] + SsL * SsR * (u_rr[2] - u_ll[2])) /
(SsR - SsL)
f3 = (SsR * f_ll[3] - SsL * f_rr[3] + SsL * SsR * (u_rr[3] - u_ll[3])) /
(SsR - SsL)
f4 = (SsR * f_ll[4] - SsL * f_rr[4] + SsL * SsR * (u_rr[4] - u_ll[4])) /
(SsR - SsL)
f5 = (SsR * f_ll[5] - SsL * f_rr[5] + SsL * SsR * (u_rr[5] - u_ll[5])) /
(SsR - SsL)
end
Computes the HLLE (Harten-Lax-van Leer-Einfeldt) flux for the compressible Euler equations.
Special estimates of the signal velocites and linearization of the Riemann problem developed
by Einfeldt to ensure that the internal energy and density remain positive during the computation
of the numerical flux.

return SVector(f1, f2, f3, f4, f5)
- Bernd Einfeldt (1988)
On Godunov-type methods for gas dynamics.
[DOI: 10.1137/0725021](https://doi.org/10.1137/0725021)
- Bernd Einfeldt, Claus-Dieter Munz, Philip L. Roe and Björn Sjögreen (1991)
On Godunov-type methods near low densities.
[DOI: 10.1016/0021-9991(91)90211-3](https://doi.org/10.1016/0021-9991(91)90211-3)
"""
@inline function min_max_speed_einfeldt(u_ll, u_rr, normal_direction::AbstractVector,
equations::CompressibleEulerEquations3D)
# Calculate primitive variables, enthalpy and speed of sound
rho_ll, v1_ll, v2_ll, v3_ll, p_ll = cons2prim(u_ll, equations)
rho_rr, v1_rr, v2_rr, v3_rr, p_rr = cons2prim(u_rr, equations)

v_dot_n_ll = v1_ll * normal_direction[1] + v2_ll * normal_direction[2] +
v3_ll * normal_direction[3]
v_dot_n_rr = v1_rr * normal_direction[1] + v2_rr * normal_direction[2] +
v3_rr * normal_direction[3]

norm_ = norm(normal_direction)

# `u_ll[5]` is total energy `rho_e_ll` on the left
H_ll = (u_ll[5] + p_ll) / rho_ll
c_ll = sqrt(equations.gamma * p_ll / rho_ll) * norm_

# `u_rr[5]` is total energy `rho_e_rr` on the right
H_rr = (u_rr[5] + p_rr) / rho_rr
c_rr = sqrt(equations.gamma * p_rr / rho_rr) * norm_

# Compute Roe averages
sqrt_rho_ll = sqrt(rho_ll)
sqrt_rho_rr = sqrt(rho_rr)
inv_sum_sqrt_rho = inv(sqrt_rho_ll + sqrt_rho_rr)

v1_roe = (sqrt_rho_ll * v1_ll + sqrt_rho_rr * v1_rr) * inv_sum_sqrt_rho
v2_roe = (sqrt_rho_ll * v2_ll + sqrt_rho_rr * v2_rr) * inv_sum_sqrt_rho
v3_roe = (sqrt_rho_ll * v3_ll + sqrt_rho_rr * v3_rr) * inv_sum_sqrt_rho
v_roe = v1_roe * normal_direction[1] + v2_roe * normal_direction[2] +
v3_roe * normal_direction[3]
v_roe_mag = v1_roe^2 + v2_roe^2 + v3_roe^2

H_roe = (sqrt_rho_ll * H_ll + sqrt_rho_rr * H_rr) * inv_sum_sqrt_rho
c_roe = sqrt((equations.gamma - 1) * (H_roe - 0.5 * v_roe_mag)) * norm_

# Compute convenience constant for positivity preservation, see
# https://doi.org/10.1016/0021-9991(91)90211-3
beta = sqrt(0.5 * (equations.gamma - 1) / equations.gamma)

# Estimate the edges of the Riemann fan (with positivity conservation)
SsL = min(v_roe - c_roe, v_dot_n_ll - beta * c_ll, zero(v_roe))
SsR = max(v_roe + c_roe, v_dot_n_rr + beta * c_rr, zero(v_roe))

return SsL, SsR
end

@inline function max_abs_speeds(u, equations::CompressibleEulerEquations3D)
Expand Down
8 changes: 8 additions & 0 deletions src/equations/numerical_fluxes.jl
Original file line number Diff line number Diff line change
Expand Up @@ -304,6 +304,14 @@ See [`FluxHLL`](@ref).
"""
const flux_hll = FluxHLL()

"""
flux_hlle

See [`min_max_speed_einfeldt`](@ref).
This is a [`FluxHLL`](@ref)-type two-wave solver with special estimates of the wave speeds.
"""
const flux_hlle = FluxHLL(min_max_speed_einfeldt)

# TODO: TrixiShallowWater: move the chen_noelle flux structure to the new package

# An empty version of the `min_max_speed_chen_noelle` function is declared here
Expand Down
26 changes: 26 additions & 0 deletions test/test_p4est_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -200,6 +200,32 @@ end
end
end

@trixi_testset "elixir_euler_sedov.jl (HLLE)" begin
@test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_sedov.jl"),
l2=[
0.411541263324004,
0.2558929632770186,
0.2558929632770193,
1.298715766843915,
],
linf=[
1.3457201726152221,
1.3138961427140758,
1.313896142714079,
6.293305112638921,
],
surface_flux=flux_hlle,
tspan=(0.0, 0.3))
# Ensure that we do not have excessive memory allocations
# (e.g., from type instabilities)
let
t = sol.t[end]
u_ode = sol.u[end]
du_ode = similar(u_ode)
@test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000
end
end

@trixi_testset "elixir_euler_blast_wave_amr.jl" begin
@test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_blast_wave_amr.jl"),
l2=[
Expand Down
28 changes: 28 additions & 0 deletions test/test_p4est_3d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -287,6 +287,34 @@ end
end
end

@trixi_testset "elixir_euler_sedov.jl (HLLE)" begin
@test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_sedov.jl"),
l2=[
0.09946224487902565,
0.04863386374672001,
0.048633863746720116,
0.04863386374672032,
0.3751015774232693,
],
linf=[
0.789241521871487,
0.42046970270100276,
0.42046970270100276,
0.4204697027010028,
4.730877375538398,
],
tspan=(0.0, 0.3),
surface_flux=flux_hlle)
# Ensure that we do not have excessive memory allocations
# (e.g., from type instabilities)
let
t = sol.t[end]
u_ode = sol.u[end]
du_ode = similar(u_ode)
@test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000
end
end

@trixi_testset "elixir_euler_source_terms_nonconforming_earth.jl" begin
@test_trixi_include(joinpath(EXAMPLES_DIR,
"elixir_euler_source_terms_nonconforming_earth.jl"),
Expand Down
Loading
Loading