Skip to content
Original file line number Diff line number Diff line change
Expand Up @@ -59,11 +59,11 @@ coordinates_max = (1.0, 1.0, 1.0)
mesh = TreeMesh(coordinates_min, coordinates_max,
initial_refinement_level = 3,
n_cells_max = 100_000,
periodicity = true)
periodicity = false)

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

###############################################################################
# ODE solvers, callbacks etc.
Expand Down
5 changes: 2 additions & 3 deletions src/solvers/dgsem_p4est/subcell_limiters_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,6 @@
function calc_bounds_twosided_interface!(var_min, var_max, variable,
u, t, semi, mesh::P4estMesh{2}, equations)
_, _, dg, cache = mesh_equations_solver_cache(semi)
(; boundary_conditions) = semi

(; neighbor_ids, node_indices) = cache.interfaces
index_range = eachnode(dg)
Expand Down Expand Up @@ -69,6 +68,7 @@ function calc_bounds_twosided_interface!(var_min, var_max, variable,
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_twosided_boundary!(var_min, var_max, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)
Expand Down Expand Up @@ -136,11 +136,9 @@ end
function calc_bounds_onesided_interface!(var_minmax, minmax, variable, u, t, semi,
mesh::P4estMesh{2})
_, equations, dg, cache = mesh_equations_solver_cache(semi)
(; boundary_conditions) = semi

(; neighbor_ids, node_indices) = cache.interfaces
index_range = eachnode(dg)
index_end = last(index_range)

# Calc bounds at interfaces and periodic boundaries
for interface in eachinterface(dg, cache)
Expand Down Expand Up @@ -192,6 +190,7 @@ function calc_bounds_onesided_interface!(var_minmax, minmax, variable, u, t, sem
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_onesided_boundary!(var_minmax, minmax, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)
Expand Down
4 changes: 2 additions & 2 deletions src/solvers/dgsem_p4est/subcell_limiters_3d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,6 @@
function calc_bounds_twosided_interface!(var_min, var_max, variable,
u, t, semi, mesh::P4estMesh{3}, equations)
_, _, dg, cache = mesh_equations_solver_cache(semi)
(; boundary_conditions) = semi

(; neighbor_ids, node_indices) = cache.interfaces
index_range = eachnode(dg)
Expand Down Expand Up @@ -96,6 +95,7 @@ function calc_bounds_twosided_interface!(var_min, var_max, variable,
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_twosided_boundary!(var_min, var_max, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)
Expand Down Expand Up @@ -180,7 +180,6 @@ end
function calc_bounds_onesided_interface!(var_minmax, minmax, variable, u, t, semi,
mesh::P4estMesh{3})
_, equations, dg, cache = mesh_equations_solver_cache(semi)
(; boundary_conditions) = semi

(; neighbor_ids, node_indices) = cache.interfaces
index_range = eachnode(dg)
Expand Down Expand Up @@ -260,6 +259,7 @@ function calc_bounds_onesided_interface!(var_minmax, minmax, variable, u, t, sem
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_onesided_boundary!(var_minmax, minmax, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)
Expand Down
32 changes: 28 additions & 4 deletions src/solvers/dgsem_structured/subcell_limiters_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,8 +8,6 @@
function calc_bounds_twosided_interface!(var_min, var_max, variable,
u, t, semi, mesh::StructuredMesh{2}, equations)
_, _, dg, cache = mesh_equations_solver_cache(semi)
(; boundary_conditions) = semi
(; contravariant_vectors) = cache.elements

# Calc bounds at interfaces and periodic boundaries
for element in eachelement(dg, cache)
Expand Down Expand Up @@ -48,9 +46,23 @@ function calc_bounds_twosided_interface!(var_min, var_max, variable,
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_twosided_boundary!(var_min, var_max, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)

return nothing
end

@inline function calc_bounds_twosided_boundary!(var_min, var_max, variable, u, t,
boundary_conditions,
mesh::StructuredMesh{2},
equations, dg, cache)
if isperiodic(mesh)
return nothing
end
(; contravariant_vectors) = cache.elements

linear_indices = LinearIndices(size(mesh))
if !isperiodic(mesh, 1)
# - xi direction
Expand Down Expand Up @@ -135,8 +147,6 @@ end
function calc_bounds_onesided_interface!(var_minmax, minmax, variable, u, t, semi,
mesh::StructuredMesh{2})
_, equations, dg, cache = mesh_equations_solver_cache(semi)
(; boundary_conditions) = semi
(; contravariant_vectors) = cache.elements

# Calc bounds at interfaces and periodic boundaries
for element in eachelement(dg, cache)
Expand Down Expand Up @@ -172,9 +182,23 @@ function calc_bounds_onesided_interface!(var_minmax, minmax, variable, u, t, sem
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_onesided_boundary!(var_minmax, minmax, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)

return nothing
end

@inline function calc_bounds_onesided_boundary!(var_minmax, minmax, variable, u, t,
boundary_conditions,
mesh::StructuredMesh{2},
equations, dg, cache)
if isperiodic(mesh)
return nothing
end
(; contravariant_vectors) = cache.elements

linear_indices = LinearIndices(size(mesh))
if !isperiodic(mesh, 1)
# - xi direction
Expand Down
28 changes: 26 additions & 2 deletions src/solvers/dgsem_tree/subcell_limiters_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,7 @@ end
u, t, semi, mesh::TreeMesh2D,
equations)
_, _, dg, cache = mesh_equations_solver_cache(semi)
(; boundary_conditions) = semi

# Calc bounds at interfaces and periodic boundaries
for interface in eachinterface(dg, cache)
# Get neighboring element ids
Expand Down Expand Up @@ -93,6 +93,18 @@ end
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_twosided_boundary!(var_min, var_max, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)

return nothing
end

@inline function calc_bounds_twosided_boundary!(var_min, var_max, variable, u, t,
boundary_conditions,
mesh::TreeMesh{2}, equations,
dg, cache)
for boundary in eachboundary(dg, cache)
element = cache.boundaries.neighbor_ids[boundary]
orientation = cache.boundaries.orientations[boundary]
Expand Down Expand Up @@ -179,7 +191,7 @@ end
@inline function calc_bounds_onesided_interface!(var_minmax, min_or_max, variable, u, t,
semi, mesh::TreeMesh2D)
_, equations, dg, cache = mesh_equations_solver_cache(semi)
(; boundary_conditions) = semi

# Calc bounds at interfaces and periodic boundaries
for interface in eachinterface(dg, cache)
# Get neighboring element ids
Expand Down Expand Up @@ -214,6 +226,18 @@ end
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_onesided_boundary!(var_minmax, min_or_max, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)

return nothing
end

@inline function calc_bounds_onesided_boundary!(var_minmax, min_or_max, variable, u, t,
boundary_conditions,
mesh::TreeMesh{2}, equations,
dg, cache)
for boundary in eachboundary(dg, cache)
element = cache.boundaries.neighbor_ids[boundary]
orientation = cache.boundaries.orientations[boundary]
Expand Down
115 changes: 115 additions & 0 deletions src/solvers/dgsem_tree/subcell_limiters_3d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -109,6 +109,64 @@ end
end
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_twosided_boundary!(var_min, var_max, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)

return nothing
end

@inline function calc_bounds_twosided_boundary!(var_min, var_max, variable,
u, t, boundary_conditions,
mesh::TreeMesh{3}, equations,
dg, cache)
for boundary in eachboundary(dg, cache)
element = cache.boundaries.neighbor_ids[boundary]
orientation = cache.boundaries.orientations[boundary]
neighbor_side = cache.boundaries.neighbor_sides[boundary]

for j in eachnode(dg), i in eachnode(dg)
# Define node indices and boundary index based on the orientation and neighbor_side
if neighbor_side == 2 # Element is on the right, boundary on the left
if orientation == 1 # boundary in x-direction
node_index = (1, i, j)
boundary_index = 1
elseif orientation == 2 # boundary in y-direction
node_index = (i, 1, j)
boundary_index = 3
else # orientation == 3 # boundary in z-direction
node_index = (i, j, 1)
boundary_index = 5
end
else # Element is on the left, boundary on the right
if orientation == 1 # boundary in x-direction
node_index = (nnodes(dg), i, j)
boundary_index = 2
elseif orientation == 2 # boundary in y-direction
node_index = (i, nnodes(dg), j)
boundary_index = 4
else # orientation == 3 # boundary in z-direction
node_index = (i, j, nnodes(dg))
boundary_index = 6
end
end
u_inner = get_node_vars(u, equations, dg, node_index..., element)
u_outer = get_boundary_outer_state(u_inner, t,
boundary_conditions[boundary_index],
orientation, boundary_index,
mesh, equations, dg, cache,
node_index..., element)
var_outer = u_outer[variable]

var_min[node_index..., element] = min(var_min[node_index..., element],
var_outer)
var_max[node_index..., element] = max(var_max[node_index..., element],
var_outer)
end
end

return nothing
end

Expand Down Expand Up @@ -211,6 +269,63 @@ end
end
end

# Calc bounds at physical boundaries
(; boundary_conditions) = semi
calc_bounds_onesided_boundary!(var_minmax, min_or_max, variable, u, t,
boundary_conditions,
mesh, equations, dg, cache)

return nothing
end

@inline function calc_bounds_onesided_boundary!(var_minmax, min_or_max, variable,
u, t, boundary_conditions,
mesh::TreeMesh{3}, equations,
dg, cache)
for boundary in eachboundary(dg, cache)
element = cache.boundaries.neighbor_ids[boundary]
orientation = cache.boundaries.orientations[boundary]
neighbor_side = cache.boundaries.neighbor_sides[boundary]

for j in eachnode(dg), i in eachnode(dg)
# Define node indices and boundary index based on the orientation and neighbor_side
if neighbor_side == 2 # Element is on the right, boundary on the left
if orientation == 1 # boundary in x-direction
node_index = (1, i, j)
boundary_index = 1
elseif orientation == 2 # boundary in y-direction
node_index = (i, 1, j)
boundary_index = 3
else # orientation == 3 # boundary in z-direction
node_index = (i, j, 1)
boundary_index = 5
end
else # Element is on the left, boundary on the right
if orientation == 1 # boundary in x-direction
node_index = (nnodes(dg), i, j)
boundary_index = 2
elseif orientation == 2 # boundary in y-direction
node_index = (i, nnodes(dg), j)
boundary_index = 4
else # orientation == 3 # boundary in z-direction
node_index = (i, j, nnodes(dg))
boundary_index = 6
end
end
u_inner = get_node_vars(u, equations, dg, node_index..., element)
u_outer = get_boundary_outer_state(u_inner, t,
boundary_conditions[boundary_index],
orientation, boundary_index,
mesh, equations, dg, cache,
node_index..., element)
var_outer = variable(u_outer, equations)

var_minmax[node_index..., element] = min_or_max(var_minmax[node_index...,
element],
var_outer)
end
end

return nothing
end

Expand Down
Loading