trixi-framework · DanielDoehring · Jun 20, 2023 · Jun 20, 2023 · Jun 21, 2023 · Jun 21, 2023
diff --git a/examples/dgmulti_1d/elixir_euler_fdsbp_periodic.jl b/examples/dgmulti_1d/elixir_euler_fdsbp_periodic.jl
@@ -4,7 +4,7 @@ using Trixi, OrdinaryDiffEq
 dg = DGMulti(element_type = Line(),
              approximation_type = periodic_derivative_operator(
                derivative_order=1, accuracy_order=4, xmin=0.0, xmax=1.0, N=50),
-             surface_flux = flux_hll,
+             surface_flux = FluxHLL(min_max_speed_naive),
              volume_integral = VolumeIntegralWeakForm())
 
 equations = CompressibleEulerEquations1D(1.4)

diff --git a/examples/dgmulti_2d/elixir_euler_bilinear.jl b/examples/dgmulti_2d/elixir_euler_bilinear.jl
@@ -2,7 +2,7 @@
 using Trixi, OrdinaryDiffEq
 
 dg = DGMulti(polydeg = 3, element_type = Quad(), approximation_type = SBP(),
-             surface_integral = SurfaceIntegralWeakForm(FluxHLL()),
+             surface_integral = SurfaceIntegralWeakForm(FluxHLL(min_max_speed_naive)),
              volume_integral = VolumeIntegralFluxDifferencing(flux_ranocha))
 
 equations = CompressibleEulerEquations2D(1.4)

diff --git a/examples/dgmulti_2d/elixir_euler_curved.jl b/examples/dgmulti_2d/elixir_euler_curved.jl
@@ -2,7 +2,7 @@
 using Trixi, OrdinaryDiffEq
 
 dg = DGMulti(polydeg = 3, element_type = Quad(), approximation_type = SBP(),
-             surface_integral = SurfaceIntegralWeakForm(FluxHLL()),
+             surface_integral = SurfaceIntegralWeakForm(FluxHLL(min_max_speed_naive)),
              volume_integral = VolumeIntegralFluxDifferencing(flux_ranocha))
 
 equations = CompressibleEulerEquations2D(1.4)

diff --git a/examples/dgmulti_2d/elixir_euler_fdsbp_periodic.jl b/examples/dgmulti_2d/elixir_euler_fdsbp_periodic.jl
@@ -4,7 +4,7 @@ using Trixi, OrdinaryDiffEq
 dg = DGMulti(element_type = Quad(),
              approximation_type = periodic_derivative_operator(
                derivative_order=1, accuracy_order=4, xmin=0.0, xmax=1.0, N=50),
-             surface_flux = flux_hll,
+             surface_flux = FluxHLL(min_max_speed_naive),
              volume_integral = VolumeIntegralWeakForm())
 
 equations = CompressibleEulerEquations2D(1.4)

diff --git a/examples/dgmulti_2d/elixir_euler_rayleigh_taylor_instability.jl b/examples/dgmulti_2d/elixir_euler_rayleigh_taylor_instability.jl
@@ -61,7 +61,7 @@ end
 
 # numerical parameters
 dg = DGMulti(polydeg = 3, element_type = Quad(), approximation_type = Polynomial(),
-             surface_integral = SurfaceIntegralWeakForm(flux_hll),
+             surface_integral = SurfaceIntegralWeakForm(FluxHLL(min_max_speed_naive)),
              volume_integral = VolumeIntegralFluxDifferencing(flux_ranocha))
 
 num_elements = 16

diff --git a/examples/dgmulti_2d/elixir_euler_triangulate_pkg_mesh.jl b/examples/dgmulti_2d/elixir_euler_triangulate_pkg_mesh.jl
@@ -1,7 +1,7 @@
 using Trixi, OrdinaryDiffEq
 
 dg = DGMulti(polydeg = 3, element_type = Tri(),
-             surface_integral = SurfaceIntegralWeakForm(FluxHLL()),
+             surface_integral = SurfaceIntegralWeakForm(FluxHLL(min_max_speed_naive)),
              volume_integral = VolumeIntegralWeakForm())
 
 equations = CompressibleEulerEquations2D(1.4)

diff --git a/examples/dgmulti_2d/elixir_euler_weakform.jl b/examples/dgmulti_2d/elixir_euler_weakform.jl
@@ -2,7 +2,7 @@
 using Trixi, OrdinaryDiffEq
 
 dg = DGMulti(polydeg = 3, element_type = Tri(), approximation_type = Polynomial(),
-             surface_integral = SurfaceIntegralWeakForm(FluxHLL()),
+             surface_integral = SurfaceIntegralWeakForm(FluxHLL(min_max_speed_naive)),
              volume_integral = VolumeIntegralWeakForm())
 
 equations = CompressibleEulerEquations2D(1.4)

diff --git a/examples/dgmulti_2d/elixir_euler_weakform_periodic.jl b/examples/dgmulti_2d/elixir_euler_weakform_periodic.jl
@@ -2,7 +2,7 @@
 using Trixi, OrdinaryDiffEq
 
 dg = DGMulti(polydeg = 3, element_type = Tri(), approximation_type = Polynomial(),
-             surface_integral = SurfaceIntegralWeakForm(FluxHLL()),
+             surface_integral = SurfaceIntegralWeakForm(FluxHLL(min_max_speed_naive)),
              volume_integral = VolumeIntegralWeakForm())
 
 equations = CompressibleEulerEquations2D(1.4)

diff --git a/examples/dgmulti_3d/elixir_euler_curved.jl b/examples/dgmulti_3d/elixir_euler_curved.jl
@@ -2,7 +2,7 @@
 using Trixi, OrdinaryDiffEq
 
 dg = DGMulti(polydeg = 3, element_type = Hex(), approximation_type=SBP(),
-             surface_integral = SurfaceIntegralWeakForm(FluxHLL()),
+             surface_integral = SurfaceIntegralWeakForm(FluxHLL(min_max_speed_naive)),
              volume_integral = VolumeIntegralFluxDifferencing(flux_ranocha))
 
 equations = CompressibleEulerEquations3D(1.4)

diff --git a/examples/dgmulti_3d/elixir_euler_weakform.jl b/examples/dgmulti_3d/elixir_euler_weakform.jl
@@ -2,7 +2,7 @@
 using Trixi, OrdinaryDiffEq
 
 dg = DGMulti(polydeg = 3, element_type = Tet(),
-             surface_integral = SurfaceIntegralWeakForm(FluxHLL()),
+             surface_integral = SurfaceIntegralWeakForm(FluxHLL(min_max_speed_naive)),
              volume_integral = VolumeIntegralWeakForm())
 
 equations = CompressibleEulerEquations3D(1.4)

diff --git a/examples/dgmulti_3d/elixir_euler_weakform_periodic.jl b/examples/dgmulti_3d/elixir_euler_weakform_periodic.jl
@@ -2,7 +2,7 @@
 using Trixi, OrdinaryDiffEq
 
 dg = DGMulti(polydeg = 3, element_type = Tet(), approximation_type = Polynomial(),
-             surface_integral = SurfaceIntegralWeakForm(FluxHLL()),
+             surface_integral = SurfaceIntegralWeakForm(FluxHLL(min_max_speed_naive)),
              volume_integral = VolumeIntegralWeakForm())
 
 equations = CompressibleEulerEquations3D(1.4)

diff --git a/examples/p4est_2d_dgsem/elixir_eulergravity_convergence.jl b/examples/p4est_2d_dgsem/elixir_eulergravity_convergence.jl
@@ -12,7 +12,7 @@ gamma = 2.0
 equations_euler = CompressibleEulerEquations2D(gamma)
 
 polydeg = 3
-solver_euler = DGSEM(polydeg, flux_hll)
+solver_euler = DGSEM(polydeg, FluxHLL(min_max_speed_naive))
 
 coordinates_min = (0.0, 0.0)
 coordinates_max = (2.0, 2.0)

diff --git a/examples/p4est_2d_dgsem/elixir_mhd_alfven_wave.jl b/examples/p4est_2d_dgsem/elixir_mhd_alfven_wave.jl
@@ -13,7 +13,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=(FluxHLL(min_max_speed_einfeldt), flux_nonconservative_powell),
                volume_integral=VolumeIntegralFluxDifferencing(volume_flux))
 
 coordinates_min = (0.0      , 0.0      )

diff --git a/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonconforming_earth.jl b/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonconforming_earth.jl
@@ -68,7 +68,7 @@ boundary_conditions = Dict(
   :outside => boundary_condition
 )
 
-surface_flux = flux_hll
+surface_flux = FluxHLL(min_max_speed_naive)
 # Note that a free stream is not preserved if N < 2 * N_geo, where N is the
 # polydeg of the solver and N_geo is the polydeg of the mesh.
 # However, the FSP error is negligible in this example.

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,7 @@ 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=(FluxHLL(min_max_speed_einfeldt), flux_nonconservative_powell),
                volume_integral=VolumeIntegralFluxDifferencing(volume_flux))
 
 coordinates_min = (-1.0, -1.0, -1.0)

diff --git a/examples/paper_self_gravitating_gas_dynamics/elixir_euler_convergence.jl b/examples/paper_self_gravitating_gas_dynamics/elixir_euler_convergence.jl
@@ -8,7 +8,7 @@ equations = CompressibleEulerEquations2D(2.0)
 
 initial_condition = initial_condition_eoc_test_coupled_euler_gravity
 
-solver = DGSEM(polydeg=3, surface_flux=flux_hll)
+solver = DGSEM(polydeg=3, surface_flux=FluxHLL(min_max_speed_naive))
 
 coordinates_min = (0.0, 0.0)
 coordinates_max = (2.0, 2.0)

diff --git a/examples/paper_self_gravitating_gas_dynamics/elixir_eulergravity_convergence.jl b/examples/paper_self_gravitating_gas_dynamics/elixir_eulergravity_convergence.jl
@@ -12,7 +12,7 @@ gamma = 2.0
 equations_euler = CompressibleEulerEquations2D(gamma)
 
 polydeg = 3
-solver_euler = DGSEM(polydeg, flux_hll)
+solver_euler = DGSEM(polydeg, FluxHLL(min_max_speed_naive))
 
 coordinates_min = (0.0, 0.0)
 coordinates_max = (2.0, 2.0)

diff --git a/examples/paper_self_gravitating_gas_dynamics/elixir_eulergravity_jeans_instability.jl b/examples/paper_self_gravitating_gas_dynamics/elixir_eulergravity_jeans_instability.jl
@@ -68,7 +68,7 @@ gamma = 5/3
 equations_euler = CompressibleEulerEquations2D(gamma)
 
 polydeg = 3
-solver_euler = DGSEM(polydeg, flux_hll)
+solver_euler = DGSEM(polydeg, FluxHLL(min_max_speed_naive))
 
 coordinates_min = (0.0, 0.0)
 coordinates_max = (1.0, 1.0)

diff --git a/examples/paper_self_gravitating_gas_dynamics/elixir_eulergravity_sedov_blast_wave.jl b/examples/paper_self_gravitating_gas_dynamics/elixir_eulergravity_sedov_blast_wave.jl
@@ -85,7 +85,7 @@ function boundary_condition_sedov_self_gravity(u_inner, orientation, direction,
 end
 boundary_conditions = boundary_condition_sedov_self_gravity
 
-surface_flux = flux_hll
+surface_flux = FluxHLL(min_max_speed_naive)
 volume_flux  = flux_chandrashekar
 polydeg = 3
 basis = LobattoLegendreBasis(polydeg)

diff --git a/examples/structured_1d_dgsem/elixir_euler_sedov_hll.jl b/examples/structured_1d_dgsem/elixir_euler_sedov_hll.jl
@@ -0,0 +1,93 @@
+
+using OrdinaryDiffEq
+using Trixi
+
+###############################################################################
+# semidiscretization of the compressible Euler equations
+
+equations = CompressibleEulerEquations1D(1.4)
+
+"""
+    initial_condition_sedov_blast_wave(x, t, equations::CompressibleEulerEquations1D)
+
+The Sedov blast wave setup based on Flash
+- http://flash.uchicago.edu/site/flashcode/user_support/flash_ug_devel/node184.html#SECTION010114000000000000000
+"""
+function initial_condition_sedov_blast_wave(x, t, equations::CompressibleEulerEquations1D)
+  # Set up polar coordinates
+  inicenter = SVector(0.0)
+  x_norm = x[1] - inicenter[1]
+  r = abs(x_norm)
+
+  # Setup based on http://flash.uchicago.edu/site/flashcode/user_support/flash_ug_devel/node184.html#SECTION010114000000000000000
+  r0 = 0.21875 # = 3.5 * smallest dx (for domain length=4 and max-ref=6)
+  # r0 = 0.5 # = more reasonable setup
+  E = 1.0
+  p0_inner = 6 * (equations.gamma - 1) * E / (3 * pi * r0)
+  p0_outer = 1.0e-5 # = true Sedov setup
+
+  # Calculate primitive variables
+  rho = 1.0
+  v1  = 0.0
+  p   = r > r0 ? p0_outer : p0_inner
+
+  return prim2cons(SVector(rho, v1, p), equations)
+end
+
+initial_condition = initial_condition_sedov_blast_wave
+surface_flux = flux_hll
+volume_flux  = flux_ranocha
+basis = LobattoLegendreBasis(3)
+shock_indicator_variable = density_pressure
+indicator_sc = IndicatorHennemannGassner(equations, basis,
+                                         alpha_max=1.0,
+                                         alpha_min=0.001,
+                                         alpha_smooth=true,
+                                         variable=shock_indicator_variable)
+volume_integral = VolumeIntegralShockCapturingHG(indicator_sc;
+                                                 volume_flux_dg=volume_flux,
+                                                 volume_flux_fv=surface_flux)
+solver = DGSEM(basis, surface_flux, volume_integral)
+
+cells_per_dimension = (64,)
+coordinates_min = (-2.0,)
+coordinates_max = ( 2.0,)
+mesh = StructuredMesh(cells_per_dimension, coordinates_min, coordinates_max, periodicity=true)
+
+semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver)
+
+
+###############################################################################
+# ODE solvers, callbacks etc.
+
+tspan = (0.0, 12.5)
+ode = semidiscretize(semi, tspan)
+
+summary_callback = SummaryCallback()
+
+analysis_interval = 1000
+analysis_callback = AnalysisCallback(semi, interval=analysis_interval)
+
+alive_callback = AliveCallback(analysis_interval=analysis_interval)
+
+save_solution = SaveSolutionCallback(interval=100,
+                                     save_initial_solution=true,
+                                     save_final_solution=true,
+                                     solution_variables=cons2prim)
+
+stepsize_callback = StepsizeCallback(cfl=0.5)
+
+callbacks = CallbackSet(summary_callback, 
+                        analysis_callback, 
+                        alive_callback,
+                        save_solution,
+                        stepsize_callback)
+
+
+###############################################################################
+# run the simulation
+
+sol = solve(ode, CarpenterKennedy2N54(williamson_condition=false),
+            dt=stepsize_callback(ode), # solve needs some value here but it will be overwritten by the stepsize_callback
+            save_everystep=false, callback=callbacks);
+summary_callback() # print the timer summary
diff --git a/examples/structured_1d_dgsem/elixir_euler_source_terms.jl b/examples/structured_1d_dgsem/elixir_euler_source_terms.jl
@@ -12,7 +12,7 @@ equations = CompressibleEulerEquations1D(1.4)
 initial_condition = initial_condition_convergence_test
 
 # Note that the expected EOC of 5 is not reached with this flux.
-# Using flux_hll instead yields the expected EOC.
+# Using FluxHLL(min_max_speed_naive) instead yields the expected EOC.
 solver = DGSEM(polydeg=4, surface_flux=flux_lax_friedrichs)
 
 coordinates_min = (0.0,)

diff --git a/examples/structured_2d_dgsem/elixir_euler_rayleigh_taylor_instability.jl b/examples/structured_2d_dgsem/elixir_euler_rayleigh_taylor_instability.jl
@@ -61,7 +61,7 @@ end
 
 # numerical parameters
 volume_flux = flux_ranocha
-solver = DGSEM(polydeg=3, surface_flux=flux_hll,
+solver = DGSEM(polydeg=3, surface_flux=FluxHLL(min_max_speed_naive),
                volume_integral=VolumeIntegralFluxDifferencing(volume_flux))
 
 # The domain is [0, 0.25] x [0, 1]

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,7 @@ 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=(FluxHLL(min_max_speed_einfeldt), flux_nonconservative_powell),
                volume_integral=VolumeIntegralFluxDifferencing(volume_flux))
 
 # Create the mesh

diff --git a/examples/tree_1d_dgsem/elixir_euler_source_terms.jl b/examples/tree_1d_dgsem/elixir_euler_source_terms.jl
@@ -10,7 +10,7 @@ equations = CompressibleEulerEquations1D(1.4)
 initial_condition = initial_condition_convergence_test
 
 # Note that the expected EOC of 5 is not reached with this flux.
-# Using flux_hll instead yields the expected EOC.
+# Using FluxHLL(min_max_speed_naive) instead yields the expected EOC.
 solver = DGSEM(polydeg=4, surface_flux=flux_lax_friedrichs)
 
 coordinates_min = 0.0

diff --git a/examples/tree_1d_dgsem/elixir_eulergravity_convergence.jl b/examples/tree_1d_dgsem/elixir_eulergravity_convergence.jl
@@ -10,7 +10,7 @@ equations_euler = CompressibleEulerEquations1D(gamma)
 initial_condition = initial_condition_eoc_test_coupled_euler_gravity
 
 polydeg = 3
-solver_euler = DGSEM(polydeg, flux_hll)
+solver_euler = DGSEM(polydeg, FluxHLL(min_max_speed_naive))
 
 coordinates_min = 0.0
 coordinates_max = 2.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 = FluxHLL(min_max_speed_einfeldt)
 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 = FluxHLL(min_max_speed_einfeldt)
 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
@@ -60,7 +60,7 @@ initial_condition = initial_condition_shu_osher_shock_tube
 
 boundary_conditions = BoundaryConditionDirichlet(initial_condition)
 
-surface_flux = flux_hll
+surface_flux = FluxHLL(min_max_speed_einfeldt)
 volume_flux  = flux_hindenlang_gassner
 basis = LobattoLegendreBasis(4)
 

diff --git a/examples/tree_1d_dgsem/elixir_shallowwater_well_balanced_nonperiodic.jl b/examples/tree_1d_dgsem/elixir_shallowwater_well_balanced_nonperiodic.jl
@@ -26,7 +26,7 @@ boundary_condition = BoundaryConditionDirichlet(initial_condition)
 # Get the DG approximation space
 
 volume_flux = (flux_wintermeyer_etal, flux_nonconservative_wintermeyer_etal)
-solver = DGSEM(polydeg=4, surface_flux=(flux_hll, flux_nonconservative_fjordholm_etal),
+solver = DGSEM(polydeg=4, surface_flux=(FluxHLL(min_max_speed_naive), flux_nonconservative_fjordholm_etal),
                volume_integral=VolumeIntegralFluxDifferencing(volume_flux))
 
 ###############################################################################

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,7 @@ 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=(FluxHLL(min_max_speed_einfeldt), flux_nonconservative_powell),
                volume_integral=VolumeIntegralFluxDifferencing(volume_flux))
 
 coordinates_min = (0.0, 0.0)

diff --git a/examples/tree_3d_dgsem/elixir_euler_convergence.jl b/examples/tree_3d_dgsem/elixir_euler_convergence.jl
@@ -9,7 +9,7 @@ equations = CompressibleEulerEquations3D(2.0)
 
 initial_condition = initial_condition_eoc_test_coupled_euler_gravity
 
-solver = DGSEM(polydeg=3, surface_flux=flux_hll,
+solver = DGSEM(polydeg=3, surface_flux=FluxHLL(min_max_speed_naive),
                volume_integral=VolumeIntegralWeakForm())
 
 coordinates_min = (0.0, 0.0, 0.0)

diff --git a/examples/tree_3d_dgsem/elixir_euler_sedov_blast_wave.jl b/examples/tree_3d_dgsem/elixir_euler_sedov_blast_wave.jl
@@ -87,7 +87,7 @@ function boundary_condition_sedov_self_gravity(u_inner, orientation, direction,
 end
 boundary_conditions = boundary_condition_sedov_self_gravity
 
-surface_flux = flux_hll
+surface_flux = FluxHLL(min_max_speed_naive)
 volume_flux = flux_ranocha
 polydeg = 3
 basis = LobattoLegendreBasis(polydeg)

diff --git a/examples/tree_3d_dgsem/elixir_eulergravity_convergence.jl b/examples/tree_3d_dgsem/elixir_eulergravity_convergence.jl
@@ -10,7 +10,7 @@ equations_euler = CompressibleEulerEquations3D(gamma)
 initial_condition = initial_condition_eoc_test_coupled_euler_gravity
 
 polydeg = 3
-solver_euler = DGSEM(polydeg, flux_hll)
+solver_euler = DGSEM(polydeg, FluxHLL(min_max_speed_naive))
 
 coordinates_min = (0.0, 0.0, 0.0)
 coordinates_max = (2.0, 2.0, 2.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,7 @@ 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=(FluxHLL(min_max_speed_einfeldt), flux_nonconservative_powell),
                volume_integral=VolumeIntegralFluxDifferencing(volume_flux))
 
 coordinates_min = (-1.0, -1.0, -1.0)