HLL MHD Breaking changes for Consistency (#1708)

* Breaking changes HLL MHD * format * format examples * hlle * fix * news, tests, example changes * fmt * remove left-right-biased flux from test * Set version to v0.6.0 * Migrate neural network-based indicators to new repository (#1701) * Remove all neural network indicator stuff from `src/` * Migrate neural network tests * Migrate neural network examples * Migrate test dependencies * Update NEWS.md * Fix typo * Remove Requires.jl-based use of Flux.jl * Fix formatting * Add migration of indicators to section with breaking changes --------- Co-authored-by: Hendrik Ranocha <[email protected]> * fix hlle noncartesian 2d * remove parantheses * correct test vals * Make parabolic terms nonexperimental (#1714) * Make parabolic terms non-experimental * Make NSE a separate item * Add MPI to supported features * Mention that parabolic terms are now officially supported in NEWS.md Co-authored-by: Hendrik Ranocha <[email protected]> * Deprecate some `DGMultiMesh` constructors (#1709) * remove previously deprecated functions * fix typo in NEWS.md about deprecation vs removal * fix literate tutorial * removing other deprecation * format * Revert "fix typo in NEWS.md about deprecation vs removal" This reverts commit 6b03020. * add gradient variable type parameter to `AbstractEquationsParabolic` (#1409) * add gradient variable type parameter * fix parabolic literate test * remove trailing comment * remove unnecessary abstract type * move gradient variable structs * formatting * fix dropped changes * try to fix doc tests * fixing navier stokes 1D * formatting * remove duplicate GradientVariablesPrimitive/Entropy definition * update news * bring downloads back * fix failing test * fmt --------- Co-authored-by: Michael Schlottke-Lakemper <[email protected]> Co-authored-by: Hendrik Ranocha <[email protected]>
trixi-framework · Nov 11, 2023 · 5fd3d73 · 5fd3d73
1 parent b75cab9
commit 5fd3d73
Show file tree

Hide file tree

Showing 17 changed files with 177 additions and 35 deletions.
diff --git a/NEWS.md b/NEWS.md
@@ -10,6 +10,11 @@ for human readability.
 
 #### Changed
 
+- The wave speed estimates for `flux_hll`, `FluxHLL()` are now consistent across equations.
+  In particular, the functions `min_max_speed_naive`, `min_max_speed_einfeldt` are now 
+  conceptually identical across equations.
+  Users, who have been using `flux_hll` for MHD have now to use `flux_hlle` in order to use the
+  Einfeldt wave speed estimate.
 - Parabolic diffusion terms are now officially supported and not marked as experimental
   anymore.
 

diff --git a/examples/p4est_2d_dgsem/elixir_mhd_alfven_wave.jl b/examples/p4est_2d_dgsem/elixir_mhd_alfven_wave.jl
@@ -12,7 +12,9 @@ initial_condition = initial_condition_convergence_test
 
 # Get the DG approximation space
 volume_flux = (flux_central, flux_nonconservative_powell)
-solver = DGSEM(polydeg = 4, surface_flux = (flux_hll, flux_nonconservative_powell),
+solver = DGSEM(polydeg = 4,
+               surface_flux = (flux_hlle,
+                               flux_nonconservative_powell),
                volume_integral = VolumeIntegralFluxDifferencing(volume_flux))
 
 coordinates_min = (0.0, 0.0)

diff --git a/examples/p4est_3d_dgsem/elixir_mhd_alfven_wave_nonconforming.jl b/examples/p4est_3d_dgsem/elixir_mhd_alfven_wave_nonconforming.jl
@@ -10,7 +10,9 @@ equations = IdealGlmMhdEquations3D(5 / 3)
 initial_condition = initial_condition_convergence_test
 
 volume_flux = (flux_hindenlang_gassner, flux_nonconservative_powell)
-solver = DGSEM(polydeg = 3, surface_flux = (flux_hll, flux_nonconservative_powell),
+solver = DGSEM(polydeg = 3,
+               surface_flux = (flux_hlle,
+                               flux_nonconservative_powell),
                volume_integral = VolumeIntegralFluxDifferencing(volume_flux))
 
 coordinates_min = (-1.0, -1.0, -1.0)

diff --git a/examples/structured_3d_dgsem/elixir_mhd_alfven_wave.jl b/examples/structured_3d_dgsem/elixir_mhd_alfven_wave.jl
@@ -10,7 +10,9 @@ equations = IdealGlmMhdEquations3D(5 / 3)
 initial_condition = initial_condition_convergence_test
 
 volume_flux = (flux_central, flux_nonconservative_powell)
-solver = DGSEM(polydeg = 5, surface_flux = (flux_hll, flux_nonconservative_powell),
+solver = DGSEM(polydeg = 5,
+               surface_flux = (flux_hlle,
+                               flux_nonconservative_powell),
                volume_integral = VolumeIntegralFluxDifferencing(volume_flux))
 
 # Create the mesh

diff --git a/examples/t8code_2d_dgsem/elixir_mhd_alfven_wave.jl b/examples/t8code_2d_dgsem/elixir_mhd_alfven_wave.jl
@@ -11,7 +11,7 @@ initial_condition = initial_condition_convergence_test
 
 # Get the DG approximation space
 volume_flux = (flux_central, flux_nonconservative_powell)
-solver = DGSEM(polydeg = 4, surface_flux = (flux_hll, flux_nonconservative_powell),
+solver = DGSEM(polydeg = 4, surface_flux = (flux_hlle, flux_nonconservative_powell),
                volume_integral = VolumeIntegralFluxDifferencing(volume_flux))
 
 coordinates_min = (0.0, 0.0)

diff --git a/examples/tree_1d_dgsem/elixir_mhd_briowu_shock_tube.jl b/examples/tree_1d_dgsem/elixir_mhd_briowu_shock_tube.jl
@@ -32,7 +32,7 @@ initial_condition = initial_condition_briowu_shock_tube
 
 boundary_conditions = BoundaryConditionDirichlet(initial_condition)
 
-surface_flux = flux_hll
+surface_flux = flux_hlle
 volume_flux = flux_derigs_etal
 basis = LobattoLegendreBasis(4)
 

diff --git a/examples/tree_1d_dgsem/elixir_mhd_ryujones_shock_tube.jl b/examples/tree_1d_dgsem/elixir_mhd_ryujones_shock_tube.jl
@@ -39,7 +39,7 @@ initial_condition = initial_condition_ryujones_shock_tube
 
 boundary_conditions = BoundaryConditionDirichlet(initial_condition)
 
-surface_flux = flux_hll
+surface_flux = flux_hlle
 volume_flux = flux_hindenlang_gassner
 basis = LobattoLegendreBasis(3)
 

diff --git a/examples/tree_1d_dgsem/elixir_mhd_shu_osher_shock_tube.jl b/examples/tree_1d_dgsem/elixir_mhd_shu_osher_shock_tube.jl
@@ -61,7 +61,7 @@ initial_condition = initial_condition_shu_osher_shock_tube
 
 boundary_conditions = BoundaryConditionDirichlet(initial_condition)
 
-surface_flux = flux_hll
+surface_flux = flux_hlle
 volume_flux = flux_hindenlang_gassner
 basis = LobattoLegendreBasis(4)
 

diff --git a/examples/tree_2d_dgsem/elixir_mhd_alfven_wave_mortar.jl b/examples/tree_2d_dgsem/elixir_mhd_alfven_wave_mortar.jl
@@ -10,7 +10,9 @@ equations = IdealGlmMhdEquations2D(gamma)
 initial_condition = initial_condition_convergence_test
 
 volume_flux = (flux_hindenlang_gassner, flux_nonconservative_powell)
-solver = DGSEM(polydeg = 3, surface_flux = (flux_hll, flux_nonconservative_powell),
+solver = DGSEM(polydeg = 3,
+               surface_flux = (flux_hlle,
+                               flux_nonconservative_powell),
                volume_integral = VolumeIntegralFluxDifferencing(volume_flux))
 
 coordinates_min = (0.0, 0.0)

diff --git a/examples/tree_3d_dgsem/elixir_mhd_alfven_wave_mortar.jl b/examples/tree_3d_dgsem/elixir_mhd_alfven_wave_mortar.jl
@@ -10,7 +10,9 @@ equations = IdealGlmMhdEquations3D(5 / 3)
 initial_condition = initial_condition_convergence_test
 
 volume_flux = (flux_hindenlang_gassner, flux_nonconservative_powell)
-solver = DGSEM(polydeg = 3, surface_flux = (flux_hll, flux_nonconservative_powell),
+solver = DGSEM(polydeg = 3,
+               surface_flux = (flux_hlle,
+                               flux_nonconservative_powell),
                volume_integral = VolumeIntegralFluxDifferencing(volume_flux))
 
 coordinates_min = (-1.0, -1.0, -1.0)

diff --git a/examples/unstructured_2d_dgsem/elixir_mhd_alfven_wave.jl b/examples/unstructured_2d_dgsem/elixir_mhd_alfven_wave.jl
@@ -11,7 +11,9 @@ equations = IdealGlmMhdEquations2D(gamma)
 initial_condition = initial_condition_convergence_test
 
 volume_flux = (flux_central, flux_nonconservative_powell)
-solver = DGSEM(polydeg = 7, surface_flux = (flux_hll, flux_nonconservative_powell),
+solver = DGSEM(polydeg = 7,
+               surface_flux = (flux_hlle,
+                               flux_nonconservative_powell),
                volume_integral = VolumeIntegralFluxDifferencing(volume_flux))
 
 # Get the unstructured quad mesh from a file (downloads the file if not available locally)

diff --git a/src/equations/ideal_glm_mhd_1d.jl b/src/equations/ideal_glm_mhd_1d.jl
@@ -277,6 +277,22 @@ end
     λ_max = max(abs(v_ll), abs(v_rr)) + max(cf_ll, cf_rr)
 end
 
+# Calculate estimates for minimum and maximum wave speeds for HLL-type fluxes
+@inline function min_max_speed_naive(u_ll, u_rr, orientation::Integer,
+                                     equations::IdealGlmMhdEquations1D)
+    rho_ll, rho_v1_ll, _ = u_ll
+    rho_rr, rho_v1_rr, _ = u_rr
+
+    # Calculate primitive variables
+    v1_ll = rho_v1_ll / rho_ll
+    v1_rr = rho_v1_rr / rho_rr
+
+    λ_min = v1_ll - calc_fast_wavespeed(u_ll, orientation, equations)
+    λ_max = v1_rr + calc_fast_wavespeed(u_rr, orientation, equations)
+
+    return λ_min, λ_max
+end
+
 # More refined estimates for minimum and maximum wave speeds for HLL-type fluxes
 @inline function min_max_speed_davis(u_ll, u_rr, orientation::Integer,
                                      equations::IdealGlmMhdEquations1D)
@@ -298,15 +314,15 @@ end
 end
 
 """
-    min_max_speed_naive(u_ll, u_rr, orientation::Integer, equations::IdealGlmMhdEquations1D)
+    min_max_speed_einfeldt(u_ll, u_rr, orientation::Integer, equations::IdealGlmMhdEquations1D)
 
 Calculate minimum and maximum wave speeds for HLL-type fluxes as in
 - Li (2005)
   An HLLC Riemann solver for magneto-hydrodynamics
   [DOI: 10.1016/j.jcp.2004.08.020](https://doi.org/10.1016/j.jcp.2004.08.020).
 """
-@inline function min_max_speed_naive(u_ll, u_rr, orientation::Integer,
-                                     equations::IdealGlmMhdEquations1D)
+@inline function min_max_speed_einfeldt(u_ll, u_rr, orientation::Integer,
+                                        equations::IdealGlmMhdEquations1D)
     rho_ll, rho_v1_ll, _ = u_ll
     rho_rr, rho_v1_rr, _ = u_rr
 

diff --git a/src/equations/ideal_glm_mhd_2d.jl b/src/equations/ideal_glm_mhd_2d.jl
@@ -775,6 +775,56 @@ end
     return max(abs(v_ll), abs(v_rr)) + max(cf_ll, cf_rr)
 end
 
+# Calculate estimate for minimum and maximum wave speeds for HLL-type fluxes
+@inline function min_max_speed_naive(u_ll, u_rr, orientation::Integer,
+                                     equations::IdealGlmMhdEquations2D)
+    rho_ll, rho_v1_ll, rho_v2_ll, _ = u_ll
+    rho_rr, rho_v1_rr, rho_v2_rr, _ = u_rr
+
+    # Calculate primitive velocity variables
+    v1_ll = rho_v1_ll / rho_ll
+    v2_ll = rho_v2_ll / rho_ll
+
+    v1_rr = rho_v1_rr / rho_rr
+    v2_rr = rho_v2_rr / rho_rr
+
+    # Approximate the left-most and right-most eigenvalues in the Riemann fan
+    if orientation == 1 # x-direction
+        λ_min = v1_ll - calc_fast_wavespeed(u_ll, orientation, equations)
+        λ_max = v1_rr + calc_fast_wavespeed(u_rr, orientation, equations)
+    else # y-direction
+        λ_min = v2_ll - calc_fast_wavespeed(u_ll, orientation, equations)
+        λ_max = v2_rr + calc_fast_wavespeed(u_rr, orientation, equations)
+    end
+
+    return λ_min, λ_max
+end
+
+@inline function min_max_speed_naive(u_ll, u_rr, normal_direction::AbstractVector,
+                                     equations::IdealGlmMhdEquations2D)
+    rho_ll, rho_v1_ll, rho_v2_ll, _ = u_ll
+    rho_rr, rho_v1_rr, rho_v2_rr, _ = u_rr
+
+    # Calculate primitive velocity variables
+    v1_ll = rho_v1_ll / rho_ll
+    v2_ll = rho_v2_ll / rho_ll
+
+    v1_rr = rho_v1_rr / rho_rr
+    v2_rr = rho_v2_rr / rho_rr
+
+    v_normal_ll = (v1_ll * normal_direction[1] + v2_ll * normal_direction[2])
+    v_normal_rr = (v1_rr * normal_direction[1] + v2_rr * normal_direction[2])
+
+    c_f_ll = calc_fast_wavespeed(u_ll, normal_direction, equations)
+    c_f_rr = calc_fast_wavespeed(u_rr, normal_direction, equations)
+
+    # Estimate the min/max eigenvalues in the normal direction
+    λ_min = min(v_normal_ll - c_f_ll, v_normal_rr - c_f_rr)
+    λ_max = max(v_normal_rr + c_f_rr, v_normal_rr + c_f_rr)
+
+    return λ_min, λ_max
+end
+
 # More refined estimates for minimum and maximum wave speeds for HLL-type fluxes
 @inline function min_max_speed_davis(u_ll, u_rr, orientation::Integer,
                                      equations::IdealGlmMhdEquations2D)
@@ -833,15 +883,15 @@ end
 end
 
 """
-    min_max_speed_naive(u_ll, u_rr, orientation::Integer, equations::IdealGlmMhdEquations2D)
+    min_max_speed_einfeldt(u_ll, u_rr, orientation::Integer, equations::IdealGlmMhdEquations2D)
 
 Calculate minimum and maximum wave speeds for HLL-type fluxes as in
 - Li (2005)
   An HLLC Riemann solver for magneto-hydrodynamics
   [DOI: 10.1016/j.jcp.2004.08.020](https://doi.org/10.1016/j.jcp.2004.08.020).
 """
-@inline function min_max_speed_naive(u_ll, u_rr, orientation::Integer,
-                                     equations::IdealGlmMhdEquations2D)
+@inline function min_max_speed_einfeldt(u_ll, u_rr, orientation::Integer,
+                                        equations::IdealGlmMhdEquations2D)
     rho_ll, rho_v1_ll, rho_v2_ll, _ = u_ll
     rho_rr, rho_v1_rr, rho_v2_rr, _ = u_rr
 
@@ -870,8 +920,8 @@ Calculate minimum and maximum wave speeds for HLL-type fluxes as in
     return λ_min, λ_max
 end
 
-@inline function min_max_speed_naive(u_ll, u_rr, normal_direction::AbstractVector,
-                                     equations::IdealGlmMhdEquations2D)
+@inline function min_max_speed_einfeldt(u_ll, u_rr, normal_direction::AbstractVector,
+                                        equations::IdealGlmMhdEquations2D)
     rho_ll, rho_v1_ll, rho_v2_ll, _ = u_ll
     rho_rr, rho_v1_rr, rho_v2_rr, _ = u_rr
 

diff --git a/src/equations/ideal_glm_mhd_3d.jl b/src/equations/ideal_glm_mhd_3d.jl
@@ -670,6 +670,64 @@ end
     return max(abs(v_ll), abs(v_rr)) + max(cf_ll, cf_rr)
 end
 
+# Calculate estimate for minimum and maximum wave speeds for HLL-type fluxes
+@inline function min_max_speed_naive(u_ll, u_rr, orientation::Integer,
+                                     equations::IdealGlmMhdEquations3D)
+    rho_ll, rho_v1_ll, rho_v2_ll, rho_v3_ll, _ = u_ll
+    rho_rr, rho_v1_rr, rho_v2_rr, rho_v3_rr, _ = u_rr
+
+    # Calculate primitive variables and speed of sound
+    v1_ll = rho_v1_ll / rho_ll
+    v2_ll = rho_v2_ll / rho_ll
+    v3_ll = rho_v3_ll / rho_ll
+
+    v1_rr = rho_v1_rr / rho_rr
+    v2_rr = rho_v2_rr / rho_rr
+    v3_rr = rho_v3_rr / rho_rr
+
+    # Approximate the left-most and right-most eigenvalues in the Riemann fan
+    if orientation == 1 # x-direction
+        λ_min = v1_ll - calc_fast_wavespeed(u_ll, orientation, equations)
+        λ_max = v1_rr + calc_fast_wavespeed(u_rr, orientation, equations)
+    elseif orientation == 2 # y-direction
+        λ_min = v2_ll - calc_fast_wavespeed(u_ll, orientation, equations)
+        λ_max = v2_rr + calc_fast_wavespeed(u_rr, orientation, equations)
+    else # z-direction
+        λ_min = v3_ll - calc_fast_wavespeed(u_ll, orientation, equations)
+        λ_max = v3_rr + calc_fast_wavespeed(u_rr, orientation, equations)
+    end
+
+    return λ_min, λ_max
+end
+
+@inline function min_max_speed_naive(u_ll, u_rr, normal_direction::AbstractVector,
+                                     equations::IdealGlmMhdEquations3D)
+    rho_ll, rho_v1_ll, rho_v2_ll, rho_v3_ll, _ = u_ll
+    rho_rr, rho_v1_rr, rho_v2_rr, rho_v3_rr, _ = u_rr
+
+    # Calculate primitive velocity variables
+    v1_ll = rho_v1_ll / rho_ll
+    v2_ll = rho_v2_ll / rho_ll
+    v3_ll = rho_v3_ll / rho_ll
+
+    v1_rr = rho_v1_rr / rho_rr
+    v2_rr = rho_v2_rr / rho_rr
+    v3_rr = rho_v3_rr / rho_rr
+
+    v_normal_ll = (v1_ll * normal_direction[1] +
+                   v2_ll * normal_direction[2] +
+                   v3_ll * normal_direction[3])
+    v_normal_rr = (v1_rr * normal_direction[1] +
+                   v2_rr * normal_direction[2] +
+                   v3_rr * normal_direction[3])
+
+    # Estimate the min/max eigenvalues in the normal direction
+    λ_min = v_normal_ll - calc_fast_wavespeed(u_ll, normal_direction, equations)
+    λ_max = v_normal_rr + calc_fast_wavespeed(u_rr, normal_direction, equations)
+
+    return λ_min, λ_max
+end
+
 # More refined estimates for minimum and maximum wave speeds for HLL-type fluxes
 @inline function min_max_speed_davis(u_ll, u_rr, orientation::Integer,
                                      equations::IdealGlmMhdEquations3D)
@@ -742,15 +800,15 @@ end
 end
 
 """
-    min_max_speed_naive(u_ll, u_rr, orientation_or_normal_direction, equations::IdealGlmMhdEquations3D)
+    min_max_speed_einfeldt(u_ll, u_rr, orientation_or_normal_direction, equations::IdealGlmMhdEquations3D)
 
 Calculate minimum and maximum wave speeds for HLL-type fluxes as in
 - Li (2005)
   An HLLC Riemann solver for magneto-hydrodynamics
   [DOI: 10.1016/j.jcp.2004.08.020](https://doi.org/10.1016/j.jcp.2004.08.020)
 """
-@inline function min_max_speed_naive(u_ll, u_rr, orientation::Integer,
-                                     equations::IdealGlmMhdEquations3D)
+@inline function min_max_speed_einfeldt(u_ll, u_rr, orientation::Integer,
+                                        equations::IdealGlmMhdEquations3D)
     rho_ll, rho_v1_ll, rho_v2_ll, rho_v3_ll, _ = u_ll
     rho_rr, rho_v1_rr, rho_v2_rr, rho_v3_rr, _ = u_rr
 
@@ -787,8 +845,8 @@ Calculate minimum and maximum wave speeds for HLL-type fluxes as in
     return λ_min, λ_max
 end
 
-@inline function min_max_speed_naive(u_ll, u_rr, normal_direction::AbstractVector,
-                                     equations::IdealGlmMhdEquations3D)
+@inline function min_max_speed_einfeldt(u_ll, u_rr, normal_direction::AbstractVector,
+                                        equations::IdealGlmMhdEquations3D)
     rho_ll, rho_v1_ll, rho_v2_ll, rho_v3_ll, _ = u_ll
     rho_rr, rho_v1_rr, rho_v2_rr, rho_v3_rr, _ = u_rr
 

diff --git a/test/test_tree_2d_mhd.jl b/test/test_tree_2d_mhd.jl
@@ -208,7 +208,8 @@ end
                             0.04879208429337193,
                         ],
                         tspan=(0.0, 0.06),
-                        surface_flux=(flux_hll, flux_nonconservative_powell))
+                        surface_flux=(flux_hlle,
+                                      flux_nonconservative_powell))
     # Ensure that we do not have excessive memory allocations
     # (e.g., from type instabilities)
     let

diff --git a/test/test_tree_3d_mhd.jl b/test/test_tree_3d_mhd.jl
@@ -240,7 +240,8 @@ end
                                       1000
                             end
                         end,
-                        surface_flux=(flux_hll, flux_nonconservative_powell),
+                        surface_flux=(flux_hlle,
+                                      flux_nonconservative_powell),
                         volume_flux=(flux_central, flux_nonconservative_powell),
                         coordinates_min=(0.0, 0.0, 0.0),
                         coordinates_max=(1.0, 1.0, 1.0),