Skip to content

Commit

Permalink
update 1D convergence test
Browse files Browse the repository at this point in the history
  • Loading branch information
patrickersing committed Apr 26, 2024
1 parent 7f00ff7 commit 737ebcc
Show file tree
Hide file tree
Showing 5 changed files with 63 additions and 56 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ using TrixiShallowWater
###############################################################################
# Semidiscretization of the multilayer shallow water equations with three layers

equations = ShallowWaterMultiLayerEquations1D(gravity_constant = 10.0,
equations = ShallowWaterMultiLayerEquations1D(gravity_constant = 1.1,
rhos = (0.9, 1.0, 1.1))

initial_condition = initial_condition_convergence_test
Expand All @@ -25,7 +25,7 @@ solver = DGSEM(polydeg = 3,
coordinates_min = 0.0
coordinates_max = sqrt(2.0)
mesh = TreeMesh(coordinates_min, coordinates_max,
initial_refinement_level = 4,
initial_refinement_level = 2,
n_cells_max = 10_000,
periodicity = true)

Expand All @@ -50,12 +50,16 @@ save_solution = SaveSolutionCallback(interval = 500,
save_initial_solution = true,
save_final_solution = true)

callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, save_solution)
stepsize_callback = StepsizeCallback(cfl = 1.0)

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

###############################################################################
# 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,
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
Expand Up @@ -58,7 +58,6 @@ callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, sav
###############################################################################
# run the simulation

# use a Runge-Kutta method with automatic (error based) time step size control
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);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -52,12 +52,12 @@ save_solution = SaveSolutionCallback(interval = 500,

stepsize_callback = StepsizeCallback(cfl = 0.9)

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

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

# use a Runge-Kutta method with automatic (error based) time step size control
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);
Expand Down
44 changes: 24 additions & 20 deletions src/equations/shallow_water_multilayer_1d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -138,15 +138,15 @@ function Trixi.initial_condition_convergence_test(x, t,
end

# Some constants are chosen such that the function is periodic on the domain [0,sqrt(2)]
ω = 2.0 * pi * sqrt(2.0)
ω = pi * sqrt(2.0)

H1 = 4.0 + 0.1 * cos* x[1] + t)
H2 = 2.0 + 0.1 * sin* x[1] + t)
H3 = 1.5 + 0.1 * cos* x[1] + t)
v1 = 0.9
v2 = 1.0
v3 = 1.1
b = 1.0 + 0.1 * cos(2.0 * ω * x[1])
b = 1.0 + 0.1 * cos* x[1])

return prim2cons(SVector(H1, H2, H3, v1, v2, v3, b), equations)
end
Expand All @@ -163,39 +163,43 @@ in non-periodic domains).
equations::ShallowWaterMultiLayerEquations1D)
# Same settings as in `initial_condition_convergence_test`. Some derivative simplify because
# this manufactured solution velocity is taken to be constant
ω = 2 * pi * sqrt(2.0)
ω = pi * sqrt(2.0)
g = equations.gravity

du1 = -0.1 * sin(t + x[1] * ω) - 0.1 * cos(t + x[1] * ω) +
0.9 * (-0.1 * sin(t + x[1] * ω) * ω - 0.1 * cos(t + x[1] * ω) * ω)

du2 = 0.1 * sin(t + x[1] * ω) + 0.1 * cos(t + x[1] * ω) +
0.1 * sin(t + x[1] * ω) * ω +
0.1 * cos(t + x[1] * ω) * ω

du3 = -0.1 * sin(t + x[1] * ω) +
1.1 * (-0.1 * sin(t + x[1] * ω) * ω + 0.2 * sin(2.0 * x[1] * ω) * ω)
1.1 * (0.1 * sin(x[1] * ω) * ω - 0.1 * sin(t + x[1] * ω) * ω)

du4 = 0.9 * (-0.1 * sin(t + x[1] * ω) - 0.1 * cos(t + x[1] * ω)) +
0.81 * (-0.1 * sin(t + x[1] * ω) * ω - 0.1 * cos(t + x[1] * ω) * ω) +
10.0 * (2.0 - 0.1 * sin(t + x[1] * ω) + 0.1 * cos(t + x[1] * ω)) *
g * (2.0 - 0.1 * sin(t + x[1] * ω) + 0.1 * cos(t + x[1] * ω)) *
(-0.1 * sin(t + x[1] * ω) * ω - 0.1 * cos(t + x[1] * ω) * ω) +
(2.0 - 0.1 * sin(t + x[1] * ω) + 0.1 * cos(t + x[1] * ω)) *
0.1 * g * (2.0 - 0.1 * sin(t + x[1] * ω) + 0.1 * cos(t + x[1] * ω)) *
cos(t + x[1] * ω) * ω

du5 = 0.1 * sin(t + x[1] * ω) + 0.1 * cos(t + x[1] * ω) +
0.1 * sin(t + x[1] * ω) * ω +
0.1 * cos(t + x[1] * ω) * ω +
10.0 * (0.5 + 0.1 * sin(t + x[1] * ω) - 0.1 * cos(t + x[1] * ω)) *
g * (0.5 + 0.1 * sin(t + x[1] * ω) - 0.1 * cos(t + x[1] * ω)) *
(0.1 * sin(t + x[1] * ω) * ω + 0.1 * cos(t + x[1] * ω) * ω) +
g * (0.5 + 0.1 * sin(t + x[1] * ω) - 0.1 * cos(t + x[1] * ω)) *
(-0.1 * sin(t + x[1] * ω) * ω +
0.9 * (-0.1 * sin(t + x[1] * ω) * ω - 0.1 * cos(t + x[1] * ω) * ω)) +
10.0 * (0.5 + 0.1 * sin(t + x[1] * ω) - 0.1 * cos(t + x[1] * ω)) *
(0.1 * sin(t + x[1] * ω) * ω + 0.1 * cos(t + x[1] * ω) * ω)
du6 = -0.11 * sin(t + x[1] * ω) +
1.21 * (-0.1 * sin(t + x[1] * ω) * ω + 0.2 * sin(2.0x[1] * ω) * ω) +
10.0 * (0.5 - 0.1 * cos(2.0 * x[1] * ω) + 0.1 * cos(t + x[1] * ω)) *
(-0.1 * sin(t + x[1] * ω) * ω + 0.2 * sin(2.0 * x[1] * ω) * ω) +
10.0 * (0.5 - 0.1 * cos(2.0 * x[1] * ω) + 0.1 * cos(t + x[1] * ω)) *
(0.8181818181818181 *
(-0.1 * sin(t + x[1] * ω) * ω - 0.1 * cos(t + x[1] * ω) * ω) +
0.9090909090909091 *
(0.1 * sin(t + x[1] * ω) * ω + 0.1 * cos(t + x[1] * ω) * ω) -
0.2 * sin(2.0 * x[1] * ω) * ω)
0.9 * (-0.1 * sin(t + x[1] * ω) * ω - 0.1 * cos(t + x[1] * ω) * ω))

du6 = -0.11000000000000001 * sin(t + x[1] * ω) +
1.2100000000000002 * (0.1 * sin(x[1] * ω) * ω - 0.1 * sin(t + x[1] * ω) * ω) +
g * (0.5 - 0.1 * cos(x[1] * ω) + 0.1 * cos(t + x[1] * ω)) *
(-0.1 * sin(x[1] * ω) * ω +
1.0 / 1.1 * (0.1 * sin(t + x[1] * ω) * ω + 0.1 * cos(t + x[1] * ω) * ω) +
0.9 / 1.1 * (-0.1 * sin(t + x[1] * ω) * ω - 0.1 * cos(t + x[1] * ω) * ω)) +
g * (0.5 - 0.1 * cos(x[1] * ω) + 0.1 * cos(t + x[1] * ω)) *
(0.1 * sin(x[1] * ω) * ω - 0.1 * sin(t + x[1] * ω) * ω)

return SVector(du1, du2, du3, du4, du5, du6, zero(eltype(u)))
end
Expand Down
56 changes: 28 additions & 28 deletions test/test_tree_1d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -480,22 +480,22 @@ end # 2LSWE
@test_trixi_include(joinpath(EXAMPLES_DIR,
"elixir_shallowwater_multilayer_convergence.jl"),
l2=[
0.0024190350732389444,
0.0011553525849735248,
0.001427790649735925,
0.0006301962973580928,
0.00047816381506058254,
0.0005518957180338281,
0.00013912637728693124,
0.0012190787287581911,
0.0007900811065315603,
0.00024265888138436312,
0.00044695521604641624,
0.0009421189161161149,
0.0001681269882698457,
0.000139126377287091,
],
linf=[
0.011265033642474886,
0.005221927831678352,
0.007355554947688914,
0.0021140505291696865,
0.001562655691777548,
0.002409077213895827,
0.0002187445586192549,
0.0035625489595796367,
0.0025099269565310167,
0.000677952137713711,
0.0019817824651813254,
0.002749046992753801,
0.0006393958093637853,
0.00021874455861881081,
],
tspan=(0.0, 0.25))
# Ensure that we do not have excessive memory allocations
Expand All @@ -512,22 +512,22 @@ end # 2LSWE
@test_trixi_include(joinpath(EXAMPLES_DIR,
"elixir_shallowwater_multilayer_convergence.jl"),
l2=[
0.0019005751369889307,
0.0006416934635148739,
0.001323748067661235,
0.0002514848406032196,
0.00019551914601560465,
0.00029396199135642574,
0.00013912637728693124,
0.0007413378395106811,
0.0005405760373610748,
8.772497321219398e-5,
0.0006969915241438401,
0.0004924779641660143,
0.00011816455915431224,
0.000139126377287091,
],
linf=[
0.005947927726240643,
0.0028032845734302647,
0.003944642659854947,
0.0009866355929912807,
0.0009301080749861135,
0.0012826136652659414,
0.0002187445586192549,
0.002657085141947846,
0.0028102251075707296,
0.00022855706065982861,
0.003310637567035979,
0.002490976099376596,
0.0005211086991627756,
0.00021874455861881081,
],
surface_flux=(flux_lax_friedrichs,
flux_nonconservative_ersing_etal),
Expand Down

0 comments on commit 737ebcc

Please sign in to comment.