From 0792bd32565204020afad118ac3dbbe357ae0113 Mon Sep 17 00:00:00 2001 From: Andrew Winters Date: Thu, 6 Jul 2023 08:23:57 +0200 Subject: [PATCH] simplify elixir as we can set discontinuous ICs in 1D. Also update beach test values --- .../elixir_shallowwater_beach.jl | 8 ++-- ...ixir_shallowwater_well_balanced_wet_dry.jl | 37 ++----------------- src/solvers/dgsem_tree/indicators_1d.jl | 13 +++---- test/test_tree_1d_shallowwater.jl | 4 +- 4 files changed, 15 insertions(+), 47 deletions(-) diff --git a/examples/tree_1d_dgsem/elixir_shallowwater_beach.jl b/examples/tree_1d_dgsem/elixir_shallowwater_beach.jl index d025430b48..c920333bef 100644 --- a/examples/tree_1d_dgsem/elixir_shallowwater_beach.jl +++ b/examples/tree_1d_dgsem/elixir_shallowwater_beach.jl @@ -11,7 +11,7 @@ equations = ShallowWaterEquations1D(gravity_constant=9.812) initial_condition_beach(x, t, equations:: ShallowWaterEquations1D) Initial condition to simulate a wave running towards a beach and crashing. Difficult test including both wetting and drying in the domain using slip wall boundary conditions. -The water height and speed function used here, are an adaption of the initial condition +The water height and speed function used here, are an adaption of the initial condition found in section 5.2 of the paper: - Andreas Bollermann, Sebastian Noelle, Maria Lukáčová-Medvid’ová (2011) Finite volume evolution Galerkin methods for the shallow water equations with dry beds\n @@ -34,7 +34,7 @@ function initial_condition_beach(x, t, equations:: ShallowWaterEquations1D) v = 0.0 else H = f - v = sqrt(equations.gravity/D) * H + v = sqrt(equations.gravity / D) * H end # It is mandatory to shift the water level at dry areas to make sure the water height h @@ -69,7 +69,7 @@ volume_integral = VolumeIntegralShockCapturingHG(indicator_sc; volume_flux_fv=surface_flux) solver = DGSEM(basis, surface_flux, volume_integral) - + ############################################################################### # Create the TreeMesh for the domain [0, 8] @@ -101,7 +101,7 @@ analysis_callback = AnalysisCallback(semi, interval=analysis_interval, save_anal alive_callback = AliveCallback(analysis_interval=analysis_interval) -save_solution = SaveSolutionCallback(interval=250, +save_solution = SaveSolutionCallback(dt=0.5, save_initial_solution=true, save_final_solution=true) diff --git a/examples/tree_1d_dgsem/elixir_shallowwater_well_balanced_wet_dry.jl b/examples/tree_1d_dgsem/elixir_shallowwater_well_balanced_wet_dry.jl index ac108a8c69..c68c54e5ad 100644 --- a/examples/tree_1d_dgsem/elixir_shallowwater_well_balanced_wet_dry.jl +++ b/examples/tree_1d_dgsem/elixir_shallowwater_well_balanced_wet_dry.jl @@ -65,7 +65,7 @@ volume_integral = VolumeIntegralShockCapturingHG(indicator_sc; volume_flux_fv=surface_flux) solver = DGSEM(basis, surface_flux, volume_integral) - + ############################################################################### # Create the TreeMesh for the domain [0, 1] @@ -85,38 +85,9 @@ semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver) tspan = (0.0, 25.0) ode = semidiscretize(semi, tspan) - -############################################################################### -# Workaround to set a discontinuous water and bottom topography for -# debugging and testing. Essentially, this is a slight augmentation of the -# `compute_coefficients` where the `x` node value passed here is slightly -# perturbed to the left / right in order to set a true discontinuity that avoids -# the doubled value of the LGL nodes at a particular element interface. -# -# Note! The errors from the analysis callback are not important but the error -# for this lake at rest test case `∑|H0-(h+b)|` should be near machine roundoff. - -# point to the data we want to augment -u = Trixi.wrap_array(ode.u0, semi) -# reset the initial condition -for element in eachelement(semi.solver, semi.cache) - for i in eachnode(semi.solver) - x_node = Trixi.get_node_coords(semi.cache.elements.node_coordinates, equations, semi.solver, i, element) - # We know that the discontinuity is on an interface. Slightly augment the x value by a factor - # of unit roundoff to avoid the repeated value from the LGL nodes at the interface. - if i == 1 - x_node = SVector(nextfloat(x_node[1])) - elseif i == nnodes(semi.solver) - x_node = SVector(prevfloat(x_node[1])) - end - u_node = initial_condition_complex_bottom_well_balanced(x_node, first(tspan), equations) - Trixi.set_node_vars!(u, u_node, equations, semi.solver, i, element) - end -end - summary_callback = SummaryCallback() -analysis_interval = 1000 +analysis_interval = 5000 analysis_callback = AnalysisCallback(semi, interval=analysis_interval, save_analysis=false) alive_callback = AliveCallback(analysis_interval=analysis_interval) @@ -125,7 +96,7 @@ save_solution = SaveSolutionCallback(interval=5000, save_initial_solution=true, save_final_solution=true) -stepsize_callback = StepsizeCallback(cfl=1.5) +stepsize_callback = StepsizeCallback(cfl=1.5) callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, save_solution, stepsize_callback) @@ -141,7 +112,7 @@ sol = solve(ode, SSPRK43(stage_limiter!); dt=1.0, summary_callback() # print the timer summary ############################################################################### -# workaround to compute the well-balancedness error for this partiular problem +# Workaround to compute the well-balancedness error for this particular problem # that has two reference water heights. One for a lake to the left of the # discontinuous bottom topography `H0_upper = 2.5` and another for a lake to the # right of the discontinuous bottom topography `H0_lower = 1.5`. diff --git a/src/solvers/dgsem_tree/indicators_1d.jl b/src/solvers/dgsem_tree/indicators_1d.jl index f6c6ea6275..04252b68bd 100644 --- a/src/solvers/dgsem_tree/indicators_1d.jl +++ b/src/solvers/dgsem_tree/indicators_1d.jl @@ -27,15 +27,12 @@ end # Modified indicator for ShallowWaterEquations1D to apply full FV method on cells # containing some "dry" LGL nodes. That is, if an element is partially "wet" then it becomes a # full FV element. -function (indicator_hg::IndicatorHennemannGassnerShallowWater)(u, - mesh::Union{TreeMesh{1}, - StructuredMesh{ - 1 - } - }, +function (indicator_hg::IndicatorHennemannGassnerShallowWater)(u::AbstractArray{<:Any, 3 + }, + mesh, equations::ShallowWaterEquations1D, - dg::DGSEM, - cache; kwargs...) + dg::DGSEM, cache; + kwargs...) @unpack alpha_max, alpha_min, alpha_smooth, variable = indicator_hg @unpack alpha, alpha_tmp, indicator_threaded, modal_threaded = indicator_hg.cache # TODO: Taal refactor, when to `resize!` stuff changed possibly by AMR? diff --git a/test/test_tree_1d_shallowwater.jl b/test/test_tree_1d_shallowwater.jl index e1201d4e2f..573d02f363 100644 --- a/test/test_tree_1d_shallowwater.jl +++ b/test/test_tree_1d_shallowwater.jl @@ -98,8 +98,8 @@ EXAMPLES_DIR = pkgdir(Trixi, "examples", "tree_1d_dgsem") @trixi_testset "elixir_shallowwater_beach.jl" begin @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_shallowwater_beach.jl"), - l2 = [0.17978981087947196, 1.237741563904404, 6.289469352788183e-8], - linf = [0.8460686638571695, 3.374190220014994, 4.452604152049844e-7], + l2 = [0.17979210479598923, 1.2377495706611434, 6.289818963361573e-8], + linf = [0.845938394800688, 3.3740800777086575, 4.4541473087633676e-7], tspan = (0.0, 0.05)) end