Skip to content

Commit

Permalink
add tests for the new solver
Browse files Browse the repository at this point in the history
  • Loading branch information
andrewwinters5000 committed Mar 4, 2024
1 parent f1cabf8 commit a8cf691
Show file tree
Hide file tree
Showing 3 changed files with 235 additions and 0 deletions.
82 changes: 82 additions & 0 deletions examples/unstructured_2d_fdsbp/elixir_euler_free_stream_upwind.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,82 @@

using OrdinaryDiffEq
using Trixi

###############################################################################
# semidiscretization of the compressible Euler equations

equations = CompressibleEulerEquations2D(1.4)

initial_condition = initial_condition_constant

# Boundary conditions for free-stream preservation test
boundary_condition_free_stream = BoundaryConditionDirichlet(initial_condition)

boundary_conditions = Dict(:outerCircle => boundary_condition_free_stream,
:cone1 => boundary_condition_free_stream,
:cone2 => boundary_condition_free_stream,
:iceCream => boundary_condition_free_stream)

###############################################################################
# Get the Upwind FDSBP approximation space

# Note, one must set `xmin=-1` and `xmax=1` due to the reuse
# of interpolation routines from `calc_node_coordinates!` to create
# the physical coordinates in the mappings.
D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017,
derivative_order=1,
accuracy_order=8,
xmin=-1.0, xmax=1.0,
N=17)

flux_splitting = splitting_vanleer_haenel
solver = FDSBP(D_upw,
surface_integral = SurfaceIntegralStrongForm(FluxUpwind(flux_splitting)),
volume_integral = VolumeIntegralUpwind(flux_splitting))

###############################################################################
# Get the curved quad mesh from a file (downloads the file if not available locally)

# Mesh with second order boundary polynomials requires an upwind SBP operator
# with (at least) 4th order boundary closure to guarantee the approximation is free-stream
mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/ec9a345f09199ebe471d35d5c1e4e08f/raw/15975943d8642e42f8292235314b6f1b30aa860d/mesh_inner_outer_boundaries.mesh",
joinpath(@__DIR__, "mesh_inner_outer_boundaries.mesh"))

mesh = UnstructuredMesh2D(mesh_file)

###############################################################################
# create the semi discretization object

semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver,
boundary_conditions = boundary_conditions)

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

tspan = (0.0, 5.0)
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 = 1000,
save_initial_solution = true,
save_final_solution = true)

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

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

# set small tolerances for the free-stream preservation test
sol = solve(ode, SSPRK43(), abstol = 1.0e-12, reltol = 1.0e-12,
save_everystep = false, callback = callbacks)

summary_callback() # print the timer summary
83 changes: 83 additions & 0 deletions examples/unstructured_2d_fdsbp/elixir_euler_source_terms_upwind.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,83 @@

using OrdinaryDiffEq
using Trixi

###############################################################################
# semidiscretization of the compressible Euler equations

equations = CompressibleEulerEquations2D(1.4)

initial_condition = initial_condition_convergence_test

source_term = source_terms_convergence_test

boundary_condition_eoc = BoundaryConditionDirichlet(initial_condition)

boundary_conditions = Dict( :Top => boundary_condition_eoc,
:Bottom => boundary_condition_eoc,
:Right => boundary_condition_eoc,
:Left => boundary_condition_eoc )

###############################################################################
# Get the Upwind FDSBP approximation space

# Note, one must set `xmin=-1` and `xmax=1` due to the reuse
# of interpolation routines from `calc_node_coordinates!` to create
# the physical coordinates in the mappings.
D_upw = upwind_operators(SummationByPartsOperators.Mattsson2017,
derivative_order=1,
accuracy_order=4,
xmin=-1.0, xmax=1.0,
N=9)

flux_splitting = splitting_steger_warming
solver = FDSBP(D_upw,
surface_integral = SurfaceIntegralStrongForm(FluxUpwind(flux_splitting)),
volume_integral = VolumeIntegralUpwind(flux_splitting))

###############################################################################
# Get the curved quad mesh from a file (downloads the file if not available locally)

# Mesh with first order boundary polynomials requires an upwind SBP operator
# with (at least) 2nd order boundary closure to guarantee the approximation is free-stream
mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/a4f4743008bf3233957a9ea6ac7a62e0/raw/8b36cc6649153fe0a5723b200368a210a1d74eaf/mesh_refined_box.mesh",
joinpath(@__DIR__, "mesh_refined_box.mesh"))

mesh = UnstructuredMesh2D(mesh_file)

###############################################################################
# create the semi discretization object

semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver,
source_terms = source_term,
boundary_conditions = boundary_conditions)

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

tspan = (0.0, 1.0)
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 = 1000,
save_initial_solution = true,
save_final_solution = true)

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

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

sol = solve(ode, SSPRK43(), abstol=1.0e-6, reltol=1.0e-6,
save_everystep = false, callback = callbacks)

summary_callback() # print the timer summary
70 changes: 70 additions & 0 deletions test/test_unstructured_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -610,6 +610,76 @@ end
@test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000
end
end

@trixi_testset "FDSBP (upwind): elixir_euler_source_terms_upwind.jl" begin
@test_trixi_include(joinpath(pkgdir(Trixi, "examples", "unstructured_2d_fdsbp"),
"elixir_euler_source_terms_upwind.jl"),
l2=[4.085391175504837e-5,
7.19179253772227e-5,
7.191792537723135e-5,
0.00021775241532855398],
linf=[0.0004054489124620808,
0.0006164432358217731,
0.0006164432358186644,
0.001363103391379461],
tspan=(0.0, 0.05),
atol=1.0e-10)
# Ensure that we do not have excessive memory allocations
# (e.g., from type instabilities)
let
t = sol.t[end]
u_ode = sol.u[end]
du_ode = similar(u_ode)
@test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000
end
end

@trixi_testset "FDSBP (upwind): elixir_euler_source_terms_upwind.jl with LF splitting" begin
@test_trixi_include(joinpath(pkgdir(Trixi, "examples", "unstructured_2d_fdsbp"),
"elixir_euler_source_terms_upwind.jl"),
l2=[3.8300267071890586e-5,
5.295846741663533e-5,
5.295846741663526e-5,
0.00017564759295593478],
linf=[0.00018810716496542312,
0.0003794187430412599,
0.0003794187430412599,
0.0009632958510650269],
tspan=(0.0, 0.025),
flux_splitting=splitting_lax_friedrichs,
atol=1.0e-10)
# Ensure that we do not have excessive memory allocations
# (e.g., from type instabilities)
let
t = sol.t[end]
u_ode = sol.u[end]
du_ode = similar(u_ode)
@test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000
end
end

@trixi_testset "FDSBP (upwind): elixir_euler_free_stream_upwind.jl" begin
@test_trixi_include(joinpath(pkgdir(Trixi, "examples", "unstructured_2d_fdsbp"),
"elixir_euler_free_stream_upwind.jl"),
l2=[3.2114065566681054e-14,
2.132488788134846e-14,
2.106144937311659e-14,
8.609642264224197e-13],
linf=[3.354871935812298e-11,
7.006478730531285e-12,
1.148153794261475e-11,
9.041265514042607e-10],
tspan=(0.0, 0.05),
atol=1.0e-10)
# Ensure that we do not have excessive memory allocations
# (e.g., from type instabilities)
let
t = sol.t[end]
u_ode = sol.u[end]
du_ode = similar(u_ode)
@test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000
end
end
end

# Clean up afterwards: delete Trixi.jl output directory
Expand Down

0 comments on commit a8cf691

Please sign in to comment.