Skip to content

Commit

Permalink
add tests and elixirs
Browse files Browse the repository at this point in the history
  • Loading branch information
patrickersing committed Dec 18, 2023
1 parent c9cb2ea commit 6f23516
Show file tree
Hide file tree
Showing 9 changed files with 1,043 additions and 3 deletions.
93 changes: 93 additions & 0 deletions examples/tree_1d_dgsem/elixir_shallowwater_ec.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,93 @@

using OrdinaryDiffEq
using Trixi
using TrixiShallowWater

###############################################################################
# Semidiscretization of the shallow water equations with a discontinuous
# bottom topography function

equations = ShallowWaterEquationsWetDry1D(gravity_constant = 9.81)

# Initial condition with a truly discontinuous water height, velocity, and bottom
# topography function as an academic testcase for entropy conservation.
# The errors from the analysis callback are not important but `∑∂S/∂U ⋅ Uₜ` should
# be around machine roundoff.
# Works as intended for TreeMesh1D with `initial_refinement_level=4`. If the mesh
# refinement level is changed the initial condition below may need changed as well to
# ensure that the discontinuities lie on an element interface.
function initial_condition_ec_discontinuous_bottom(x, t,
equations::ShallowWaterEquationsWetDry1D)
# Set the background values
H = 4.25
v = 0.0
b = sin(x[1]) # arbitrary continuous function

# Setup the discontinuous water height and velocity
if x[1] >= 0.125 && x[1] <= 0.25
H = 5.0
v = 0.1882
end

# Setup a discontinuous bottom topography
if x[1] >= -0.25 && x[1] <= -0.125
b = 2.0 + 0.5 * sin(2.0 * pi * x[1])
end

return prim2cons(SVector(H, v, b), equations)
end

initial_condition = initial_condition_ec_discontinuous_bottom

###############################################################################
# Get the DG approximation space

volume_flux = (flux_wintermeyer_etal, flux_nonconservative_wintermeyer_etal)
solver = DGSEM(polydeg = 4,
surface_flux = (flux_fjordholm_etal, flux_nonconservative_fjordholm_etal),
volume_integral = VolumeIntegralFluxDifferencing(volume_flux))

###############################################################################
# Get the TreeMesh and setup a periodic mesh

coordinates_min = -1.0
coordinates_max = 1.0
mesh = TreeMesh(coordinates_min, coordinates_max,
initial_refinement_level = 4,
n_cells_max = 10_000)

# Create the semi discretization object
semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver)

###############################################################################
# ODE solver

tspan = (0.0, 2.0)
ode = semidiscretize(semi, tspan)

###############################################################################
# Callbacks

summary_callback = SummaryCallback()

analysis_interval = 100
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)

stepsize_callback = StepsizeCallback(cfl = 3.0)

callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, save_solution,
stepsize_callback)

###############################################################################
# run the simulation

sol = solve(ode, CarpenterKennedy2N54(williamson_condition = false),
dt = 1.0, # 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
60 changes: 60 additions & 0 deletions examples/tree_1d_dgsem/elixir_shallowwater_source_terms.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,60 @@

using OrdinaryDiffEq
using Trixi
using TrixiShallowWater

###############################################################################
# Semidiscretization of the shallow water equations

equations = ShallowWaterEquationsWetDry1D(gravity_constant = 9.81)

initial_condition = initial_condition_convergence_test

###############################################################################
# Get the DG approximation space

volume_flux = (flux_wintermeyer_etal, flux_nonconservative_wintermeyer_etal)
solver = DGSEM(polydeg = 3,
surface_flux = (flux_lax_friedrichs, flux_nonconservative_fjordholm_etal),
volume_integral = VolumeIntegralFluxDifferencing(volume_flux))

###############################################################################
# Get the TreeMesh and setup a periodic mesh

coordinates_min = 0.0
coordinates_max = sqrt(2.0)
mesh = TreeMesh(coordinates_min, coordinates_max,
initial_refinement_level = 3,
n_cells_max = 10_000,
periodicity = true)

# create the semi discretization object
semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver,
source_terms = source_terms_convergence_test)

###############################################################################
# ODE solvers, callbacks etc.

tspan = (0.0, 1.0)
ode = semidiscretize(semi, tspan)

summary_callback = SummaryCallback()

analysis_interval = 500
analysis_callback = AnalysisCallback(semi, interval = analysis_interval)

alive_callback = AliveCallback(analysis_interval = analysis_interval)

save_solution = SaveSolutionCallback(interval = 200,
save_initial_solution = true,
save_final_solution = true)

callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, save_solution)

###############################################################################
# run the simulation

# use a Runge-Kutta method with automatic (error based) time step size control
sol = solve(ode, RDPK3SpFSAL49(); abstol = 1.0e-8, reltol = 1.0e-8,
ode_default_options()..., callback = callbacks);
summary_callback() # print the timer summary
Original file line number Diff line number Diff line change
@@ -0,0 +1,63 @@

using OrdinaryDiffEq
using Trixi
using TrixiShallowWater

###############################################################################
# Semidiscretization of the shallow water equations

equations = ShallowWaterEquationsWetDry1D(gravity_constant = 9.81)

initial_condition = initial_condition_convergence_test

boundary_condition = BoundaryConditionDirichlet(initial_condition)

###############################################################################
# Get the DG approximation space

volume_flux = (flux_wintermeyer_etal, flux_nonconservative_wintermeyer_etal)
surface_flux = (flux_lax_friedrichs, flux_nonconservative_fjordholm_etal)
solver = DGSEM(polydeg = 3, surface_flux = surface_flux,
volume_integral = VolumeIntegralFluxDifferencing(volume_flux))

###############################################################################
# Get the TreeMesh and setup a periodic mesh

coordinates_min = 0.0
coordinates_max = sqrt(2.0)
mesh = TreeMesh(coordinates_min, coordinates_max,
initial_refinement_level = 3,
n_cells_max = 10_000,
periodicity = false)

# create the semi discretization object
semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver,
boundary_conditions = boundary_condition,
source_terms = source_terms_convergence_test)

###############################################################################
# ODE solvers, callbacks etc.

tspan = (0.0, 1.0)
ode = semidiscretize(semi, tspan)

summary_callback = SummaryCallback()

analysis_interval = 500
analysis_callback = AnalysisCallback(semi, interval = analysis_interval)

alive_callback = AliveCallback(analysis_interval = analysis_interval)

save_solution = SaveSolutionCallback(interval = 200,
save_initial_solution = true,
save_final_solution = true)

callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, save_solution)

###############################################################################
# run the simulation

# use a Runge-Kutta method with automatic (error based) time step size control
sol = solve(ode, RDPK3SpFSAL49(); abstol = 1.0e-8, reltol = 1.0e-8,
ode_default_options()..., callback = callbacks);
summary_callback() # print the timer summary
88 changes: 88 additions & 0 deletions examples/tree_1d_dgsem/elixir_shallowwater_well_balanced.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,88 @@

using OrdinaryDiffEq
using Trixi
using TrixiShallowWater

###############################################################################
# semidiscretization of the shallow water equations with a discontinuous
# bottom topography function

equations = ShallowWaterEquationsWetDry1D(gravity_constant = 9.81, H0 = 3.25)

# Setup a truly discontinuous bottom topography function for this academic
# testcase of well-balancedness. The errors from the analysis callback are
# not important but the error for this lake-at-rest test case
# `∑|H0-(h+b)|` should be around machine roundoff.
# Works as intended for TreeMesh1D with `initial_refinement_level=3`. If the mesh
# refinement level is changed the initial condition below may need changed as well to
# ensure that the discontinuities lie on an element interface.
function initial_condition_discontinuous_well_balancedness(x, t,
equations::ShallowWaterEquationsWetDry1D)
# Set the background values
H = equations.H0
v = 0.0
b = 0.0

# Setup a discontinuous bottom topography
if x[1] >= 0.5 && x[1] <= 0.75
b = 2.0 + 0.5 * sin(2.0 * pi * x[1])
end

return prim2cons(SVector(H, v, b), equations)
end

initial_condition = initial_condition_discontinuous_well_balancedness

###############################################################################
# Get the DG approximation space

volume_flux = (flux_wintermeyer_etal, flux_nonconservative_wintermeyer_etal)
surface_flux = (flux_fjordholm_etal, flux_nonconservative_fjordholm_etal)
solver = DGSEM(polydeg = 4, surface_flux = surface_flux,
volume_integral = VolumeIntegralFluxDifferencing(volume_flux))

###############################################################################
# Get the TreeMesh and setup a periodic mesh

coordinates_min = -1.0
coordinates_max = 1.0
mesh = TreeMesh(coordinates_min, coordinates_max,
initial_refinement_level = 3,
n_cells_max = 10_000)

# Create the semi discretization object
semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver)

###############################################################################
# ODE solver

tspan = (0.0, 100.0)
ode = semidiscretize(semi, tspan)

###############################################################################
# Callbacks

summary_callback = SummaryCallback()

analysis_interval = 1000
analysis_callback = AnalysisCallback(semi, interval = analysis_interval,
extra_analysis_integrals = (lake_at_rest_error,))

alive_callback = AliveCallback(analysis_interval = analysis_interval)

save_solution = SaveSolutionCallback(interval = 1000,
save_initial_solution = true,
save_final_solution = true)

stepsize_callback = StepsizeCallback(cfl = 3.0)

callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, save_solution,
stepsize_callback)

###############################################################################
# run the simulation

sol = solve(ode, CarpenterKennedy2N54(williamson_condition = false),
dt = 1.0, # 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
Original file line number Diff line number Diff line change
@@ -0,0 +1,77 @@

using OrdinaryDiffEq
using Trixi
using TrixiShallowWater

###############################################################################
# Semidiscretization of the shallow water equations

equations = ShallowWaterEquationsWetDry1D(gravity_constant = 1.0, H0 = 3.0)

# An initial condition with constant total water height and zero velocities to test well-balancedness.
function initial_condition_well_balancedness(x, t, equations::ShallowWaterEquationsWetDry1D)
# Set the background values
H = equations.H0
v = 0.0

b = (1.5 / exp(0.5 * ((x[1] - 1.0)^2)) + 0.75 / exp(0.5 * ((x[1] + 1.0)^2)))

return prim2cons(SVector(H, v, b), equations)
end

initial_condition = initial_condition_well_balancedness

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),
volume_integral = VolumeIntegralFluxDifferencing(volume_flux))

###############################################################################
# Get the TreeMesh and setup a periodic mesh

coordinates_min = 0.0
coordinates_max = sqrt(2.0)
mesh = TreeMesh(coordinates_min, coordinates_max,
initial_refinement_level = 3,
n_cells_max = 10_000,
periodicity = false)

# create the semi discretization object
semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver,
boundary_conditions = boundary_condition)

###############################################################################
# ODE solvers, callbacks etc.

tspan = (0.0, 100.0)
ode = semidiscretize(semi, tspan)

summary_callback = SummaryCallback()

analysis_interval = 1000
analysis_callback = AnalysisCallback(semi, interval = analysis_interval,
save_analysis = true,
extra_analysis_integrals = (lake_at_rest_error,))

alive_callback = AliveCallback(analysis_interval = analysis_interval)

save_solution = SaveSolutionCallback(interval = 1000,
save_initial_solution = true,
save_final_solution = true)

stepsize_callback = StepsizeCallback(cfl = 1.0)

callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, save_solution,
stepsize_callback)

###############################################################################
# run the simulation

sol = solve(ode, CarpenterKennedy2N54(williamson_condition = false),
dt = 1.0, # 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
4 changes: 1 addition & 3 deletions test/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -11,9 +11,7 @@ const TRIXI_NTHREADS = clamp(Sys.CPU_THREADS, 2, 3)

@time @testset "TrixiShallowWater.jl tests" begin
@time if TRIXI_TEST == "all"
@test TrixiShallowWater.foo() == true
@test TrixiShallowWater.bar() == false
@test TrixiShallowWater.baz() isa String
include("test_tree_1d_shallowwater_wet_dry.jl")
end

@time if TRIXI_TEST == "all" || TRIXI_TEST == "upstream"
Expand Down
Loading

0 comments on commit 6f23516

Please sign in to comment.