Skip to content

Commit

Permalink
Fix normalization for Charactersitics and slip wall
Browse files Browse the repository at this point in the history
  • Loading branch information
@bennibolm
bennibolm committed Nov 28, 2023
1 parent 2beb8e5 commit bcfc1e9
Show file tree
Hide file tree
Showing 6 changed files with 64 additions and 49 deletions.
15 changes: 8 additions & 7 deletions examples/structured_2d_dgsem/elixir_euler_double_mach.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -76,22 +76,23 @@ end
@inline function Trixi.get_boundary_outer_state(u_inner, cache, t,
boundary_condition::typeof(boundary_condition_mixed_characteristic_wall),
normal_direction::AbstractVector, direction,
equations,
equations::CompressibleEulerEquations2D,
dg, indices...)
x = Trixi.get_node_coords(cache.elements.node_coordinates, equations, dg, indices...)
if x[1] < 1 / 6 # BoundaryConditionCharacteristic
u_outer = Trixi.characteristic_boundary_value_function(initial_condition_double_mach_reflection,
u_inner,
normal_direction,
direction, x, t,
equations)
normal_direction /
norm(normal_direction),
direction, x, t, equations)

else # if x[1] >= 1 / 6 # boundary_condition_slip_wall
u_rotate = Trixi.rotate_to_x(u_inner, normal_direction, equations)
factor = (normal_direction[1] * u_inner[2] + normal_direction[2] * u_inner[3])
u_normal = (factor / sum(normal_direction .^ 2)) * normal_direction

u_outer = SVector(u_inner[1],
u_inner[2] - 2.0 * u_rotate[2],
u_inner[3] - 2.0 * u_rotate[3],
u_inner[2] - 2.0 * u_normal[1],
u_inner[3] - 2.0 * u_normal[2],
u_inner[4])
end
return u_outer
Expand Down
15 changes: 8 additions & 7 deletions examples/structured_2d_dgsem/elixir_euler_double_mach_MCL.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -76,22 +76,23 @@ end
@inline function Trixi.get_boundary_outer_state(u_inner, cache, t,
boundary_condition::typeof(boundary_condition_mixed_characteristic_wall),
normal_direction::AbstractVector, direction,
equations,
equations::CompressibleEulerEquations2D,
dg, indices...)
x = Trixi.get_node_coords(cache.elements.node_coordinates, equations, dg, indices...)
if x[1] < 1 / 6 # BoundaryConditionCharacteristic
u_outer = Trixi.characteristic_boundary_value_function(initial_condition_double_mach_reflection,
u_inner,
normal_direction,
direction, x, t,
equations)
normal_direction /
norm(normal_direction),
direction, x, t, equations)

else # if x[1] >= 1 / 6 # boundary_condition_slip_wall
u_rotate = Trixi.rotate_to_x(u_inner, normal_direction, equations)
factor = (normal_direction[1] * u_inner[2] + normal_direction[2] * u_inner[3])
u_normal = (factor / sum(normal_direction .^ 2)) * normal_direction

u_outer = SVector(u_inner[1],
u_inner[2] - 2.0 * u_rotate[2],
u_inner[3] - 2.0 * u_rotate[3],
u_inner[2] - 2.0 * u_normal[1],
u_inner[3] - 2.0 * u_normal[2],
u_inner[4])
end
return u_outer
Expand Down
2 changes: 1 addition & 1 deletion src/solvers/dgsem_structured/dg_2d_subcell_limiters.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -109,7 +109,7 @@
return nothing
end

@inline function calc_lambdas_bar_states!(u, t, mesh::StructuredMesh,
@inline function calc_lambdas_bar_states!(u, t, mesh::StructuredMesh{2},
nonconservative_terms, equations, limiter,
dg, cache, boundary_conditions;
calc_bar_states = true)
Expand Down
45 changes: 29 additions & 16 deletions src/solvers/dgsem_tree/dg_2d_subcell_limiters.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1831,29 +1831,25 @@ end
@inline function get_boundary_outer_state(u_inner, cache, t,
boundary_condition::typeof(boundary_condition_slip_wall),
orientation::Integer, direction,
equations, dg, indices...)
equations::CompressibleEulerEquations2D,
dg, indices...)
return SVector(u_inner[1], -u_inner[2], -u_inner[3], u_inner[4])
end

@inline function get_boundary_outer_state(u_inner, cache, t,
boundary_condition::typeof(boundary_condition_slip_wall),
normal_direction::AbstractVector,
direction, equations, dg, indices...)
u_rotate = rotate_to_x(u_inner, normal_direction, equations)
normal_direction::AbstractVector, direction,
equations::CompressibleEulerEquations2D,
dg, indices...)
factor = (normal_direction[1] * u_inner[2] + normal_direction[2] * u_inner[3])
u_normal = (factor / sum(normal_direction .^ 2)) * normal_direction

return SVector(u_inner[1],
u_inner[2] - 2.0 * u_rotate[2],
u_inner[3] - 2.0 * u_rotate[3],
u_inner[2] - 2.0 * u_normal[1],
u_inner[3] - 2.0 * u_normal[2],
u_inner[4])
end

# Default implementation of `get_boundary_outer_state` returns inner value.
@inline function get_boundary_outer_state(u_inner, cache, t, boundary_condition,
orientation_or_normal, direction, equations,
dg, indices...)
return u_inner
end

@inline function get_boundary_outer_state(u_inner, cache, t,
boundary_condition::BoundaryConditionDirichlet,
orientation_or_normal, direction, equations,
Expand All @@ -1868,13 +1864,30 @@ end

@inline function get_boundary_outer_state(u_inner, cache, t,
boundary_condition::BoundaryConditionCharacteristic,
orientation_or_normal, direction, equations,
orientation::Integer, direction, equations,
dg, indices...)
@unpack node_coordinates = cache.elements
(; node_coordinates) = cache.elements

x = get_node_coords(node_coordinates, equations, dg, indices...)
u_outer = boundary_condition.boundary_value_function(boundary_condition.outer_boundary_value_function,
u_inner, orientation_or_normal,
u_inner, orientation,
direction, x, t, equations)

return u_outer
end

@inline function get_boundary_outer_state(u_inner, cache, t,
boundary_condition::BoundaryConditionCharacteristic,
normal_direction::AbstractVector, direction,
equations, dg, indices...)
(; node_coordinates) = cache.elements

x = get_node_coords(node_coordinates, equations, dg, indices...)

u_outer = boundary_condition.boundary_value_function(boundary_condition.outer_boundary_value_function,
u_inner,
normal_direction /
norm(normal_direction),
direction, x, t, equations)

return u_outer
Expand Down
4 changes: 2 additions & 2 deletions src/time_integration/methods_SSP.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -269,8 +269,8 @@ function calc_normal_directions!(container_bar_states, mesh::TreeMesh, equations
nothing
end

function calc_normal_directions!(container_bar_states, mesh::StructuredMesh, equations,
dg, cache)
function calc_normal_directions!(container_bar_states, mesh::StructuredMesh{2},
equations, dg, cache)
(; weights, derivative_matrix) = dg.basis
(; contravariant_vectors) = cache.elements

Expand Down
32 changes: 16 additions & 16 deletions test/test_structured_2d.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -695,16 +695,16 @@ end
@test_trixi_include(joinpath(EXAMPLES_DIR,
"elixir_euler_shock_upstream_sc_subcell.jl"),
l2=[
1.2351650580585378,
1.1268555757641623,
1.723421391322672,
11.715156856828612,
1.2360964282793594,
1.1280398030799386,
1.7267648591561005,
11.7190241252948,
],
linf=[
5.378363125399197,
6.5624165241817245,
10.008377369949162,
51.24533429408666,
5.393863681425582,
6.476830142076269,
10.13110957325047,
49.922573759404,
],
tspan=(0.0, 0.5))
# Ensure that we do not have excessive memory allocations
Expand All @@ -720,16 +720,16 @@ end
@trixi_testset "elixir_euler_shock_upstream_MCL.jl" begin
@test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_shock_upstream_MCL.jl"),
l2=[
1.2607430289878367,
1.1565837325291748,
1.779179030245946,
11.891223800389836,
1.2655219565562812,
1.1981291403309833,
1.7964926709729876,
11.968876601483794,
],
linf=[
5.68876088478026,
8.16555442595031,
10.859100194836678,
50.258224089906975,
5.675285252773072,
10.063582862371597,
12.3837111905078,
51.34700027665198,
],
tspan=(0.0, 0.5))
# Ensure that we do not have excessive memory allocations
Expand Down

0 comments on commit bcfc1e9

Please sign in to comment.