From 5fd3d73630fd973f9216fab3c40c1c52509cca68 Mon Sep 17 00:00:00 2001 From: Daniel Doehring Date: Fri, 10 Nov 2023 16:02:11 +0100 Subject: [PATCH] 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 * 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 * 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 6b03020a8c881a86550484891a0f53bca018e447. * 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 Co-authored-by: Hendrik Ranocha --- NEWS.md | 5 ++ .../p4est_2d_dgsem/elixir_mhd_alfven_wave.jl | 4 +- .../elixir_mhd_alfven_wave_nonconforming.jl | 4 +- .../elixir_mhd_alfven_wave.jl | 4 +- .../t8code_2d_dgsem/elixir_mhd_alfven_wave.jl | 2 +- .../elixir_mhd_briowu_shock_tube.jl | 2 +- .../elixir_mhd_ryujones_shock_tube.jl | 2 +- .../elixir_mhd_shu_osher_shock_tube.jl | 2 +- .../elixir_mhd_alfven_wave_mortar.jl | 4 +- .../elixir_mhd_alfven_wave_mortar.jl | 4 +- .../elixir_mhd_alfven_wave.jl | 4 +- src/equations/ideal_glm_mhd_1d.jl | 22 +++++- src/equations/ideal_glm_mhd_2d.jl | 60 ++++++++++++++-- src/equations/ideal_glm_mhd_3d.jl | 68 +++++++++++++++++-- test/test_tree_2d_mhd.jl | 3 +- test/test_tree_3d_mhd.jl | 3 +- test/test_unit.jl | 19 +++--- 17 files changed, 177 insertions(+), 35 deletions(-) diff --git a/NEWS.md b/NEWS.md index 5d258fa65bb..54fbb90b8fc 100644 --- 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 index 93db7bfdbaf..b2b402a25f6 100644 --- 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 index 6a62368ef99..12ddf9e4a5f 100644 --- 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 index 6eb35078ef4..03f50ff3e55 100644 --- 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 index 8040f79cafd..1e2362a123c 100644 --- 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 index bb294af68cb..4d537508b47 100644 --- 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 index b6f856fbc64..a7d7689a806 100644 --- 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 index 8c0317277b0..689537ebdd4 100644 --- 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 index 229360f266e..0200b591844 100644 --- 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 index 3ce166a7fa7..2fa22d535d7 100644 --- 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 index fdafed98c8d..3ed3e828ca8 100644 --- 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 index 7e5c94c7bc3..85030e8a5ad 100644 --- 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 index e8de0cedde1..43d1991e34b 100644 --- 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 index 09990837706..321e501b051 100644 --- 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 index bd6a95bba50..953c077c0a3 100644 --- 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 index 7ce5ef1d18f..e75685f0b43 100644 --- 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), diff --git a/test/test_unit.jl b/test/test_unit.jl index 609100793ba..1d768c0df69 100644 --- a/test/test_unit.jl +++ b/test/test_unit.jl @@ -952,14 +952,12 @@ end end @timed_testset "Consistency check for HLLE flux: MHD" begin - # Note: min_max_speed_naive for MHD is essentially min_max_speed_einfeldt - equations = IdealGlmMhdEquations1D(1.4) u_values = [SVector(1.0, 0.4, -0.5, 0.1, 1.0, 0.1, -0.2, 0.1), SVector(1.5, -0.2, 0.1, 0.2, 5.0, -0.1, 0.1, 0.2)] for u in u_values - @test flux_hll(u, u, 1, equations) ≈ flux(u, 1, equations) + @test flux_hlle(u, u, 1, equations) ≈ flux(u, 1, equations) end equations = IdealGlmMhdEquations2D(1.4, 5.0) #= c_h =# @@ -973,11 +971,11 @@ end SVector(1.5, -0.2, 0.1, 0.2, 5.0, -0.1, 0.1, 0.2, 0.2)] for u in u_values, orientation in orientations - @test flux_hll(u, u, orientation, equations) ≈ flux(u, orientation, equations) + @test flux_hlle(u, u, orientation, equations) ≈ flux(u, orientation, equations) end for u in u_values, normal_direction in normal_directions - @test flux_hll(u, u, normal_direction, equations) ≈ + @test flux_hlle(u, u, normal_direction, equations) ≈ flux(u, normal_direction, equations) end @@ -993,11 +991,11 @@ end SVector(1.5, -0.2, 0.1, 0.2, 5.0, -0.1, 0.1, 0.2, 0.2)] for u in u_values, orientation in orientations - @test flux_hll(u, u, orientation, equations) ≈ flux(u, orientation, equations) + @test flux_hlle(u, u, orientation, equations) ≈ flux(u, orientation, equations) end for u in u_values, normal_direction in normal_directions - @test flux_hll(u, u, normal_direction, equations) ≈ + @test flux_hlle(u, u, normal_direction, equations) ≈ flux(u, normal_direction, equations) end end @@ -1244,8 +1242,8 @@ end fluxes = [ flux_central, flux_hindenlang_gassner, - flux_hll, FluxHLL(min_max_speed_davis), + flux_hlle, ] for f_std in fluxes @@ -1271,8 +1269,8 @@ end fluxes = [ flux_central, flux_hindenlang_gassner, - flux_hll, FluxHLL(min_max_speed_davis), + flux_hlle, ] for f_std in fluxes @@ -1322,7 +1320,8 @@ end dg = DGMulti(polydeg = 1, element_type = Line(), approximation_type = Polynomial(), surface_integral = SurfaceIntegralWeakForm(flux_central), volume_integral = VolumeIntegralFluxDifferencing(flux_central)) - mesh = DGMultiMesh(dg, cells_per_dimension = (1,), periodicity = false) + cells_per_dimension = (1,) + mesh = DGMultiMesh(dg, cells_per_dimension, periodicity = false) @test mesh.boundary_faces[:entire_boundary] == [1, 2] end