Check BCs for periodicity for periodic Tree & Structured meshes (#1860)

* Check BCs for periodicity for periodic meshes

* default case for periodic bcs

* fmt

* specialize

* error ant fmt

* isperiodic TreeMesh

* avoid if

* shorten

* shorten

* Make meshes non-periodic

* fix rti

* shorten dispatch

---------

Co-authored-by: Hendrik Ranocha <[email protected]>

examples/structured_2d_dgsem/elixir_advection_coupled.jl

@@ -53,7 +53,8 @@ cells_per_dimension = (8, 8)
 coordinates_min1 = (-1.0, 0.0) # minimum coordinates (min(x), min(y))
 coordinates_max1 = (0.0, 1.0) # maximum coordinates (max(x), max(y))
 
-mesh1 = StructuredMesh(cells_per_dimension, coordinates_min1, coordinates_max1)
+mesh1 = StructuredMesh(cells_per_dimension, coordinates_min1, coordinates_max1,
+                       periodicity = false)
 
 # Define the coupling functions
 coupling_function12 = (x, u, equations_other, equations_own) -> u
@@ -84,7 +85,8 @@ semi1 = SemidiscretizationHyperbolic(mesh1, equations, initial_condition_converg
 coordinates_min2 = (0.0, 0.0) # minimum coordinates (min(x), min(y))
 coordinates_max2 = (1.0, 1.0) # maximum coordinates (max(x), max(y))
 
-mesh2 = StructuredMesh(cells_per_dimension, coordinates_min2, coordinates_max2)
+mesh2 = StructuredMesh(cells_per_dimension, coordinates_min2, coordinates_max2,
+                       periodicity = false)
 
 # Define the coupling functions
 coupling_function21 = (x, u, equations_other, equations_own) -> u
@@ -115,7 +117,8 @@ semi2 = SemidiscretizationHyperbolic(mesh2, equations, initial_condition_converg
 coordinates_min3 = (-1.0, -1.0) # minimum coordinates (min(x), min(y))
 coordinates_max3 = (0.0, 0.0) # maximum coordinates (max(x), max(y))
 
-mesh3 = StructuredMesh(cells_per_dimension, coordinates_min3, coordinates_max3)
+mesh3 = StructuredMesh(cells_per_dimension, coordinates_min3, coordinates_max3,
+                       periodicity = false)
 
 # Define the coupling functions
 coupling_function34 = (x, u, equations_other, equations_own) -> u
@@ -146,7 +149,8 @@ semi3 = SemidiscretizationHyperbolic(mesh3, equations, initial_condition_converg
 coordinates_min4 = (0.0, -1.0) # minimum coordinates (min(x), min(y))
 coordinates_max4 = (1.0, 0.0) # maximum coordinates (max(x), max(y))
 
-mesh4 = StructuredMesh(cells_per_dimension, coordinates_min4, coordinates_max4)
+mesh4 = StructuredMesh(cells_per_dimension, coordinates_min4, coordinates_max4,
+                       periodicity = false)
 
 # Define the coupling functions
 coupling_function43 = (x, u, equations_other, equations_own) -> u

examples/structured_2d_dgsem/elixir_euler_rayleigh_taylor_instability.jl

@@ -69,7 +69,7 @@ mapping(xi, eta) = SVector(0.25 * 0.5 * (1.0 + xi), 0.5 * (1.0 + eta))
 
 num_elements_per_dimension = 32
 cells_per_dimension = (num_elements_per_dimension, num_elements_per_dimension * 4)
-mesh = StructuredMesh(cells_per_dimension, mapping)
+mesh = StructuredMesh(cells_per_dimension, mapping, periodicity = false)
 
 initial_condition = initial_condition_rayleigh_taylor_instability
 boundary_conditions = (x_neg = boundary_condition_slip_wall,

examples/structured_2d_dgsem/elixir_euler_warm_bubble.jl

@@ -100,7 +100,8 @@ coordinates_min = (0.0, 0.0)
 coordinates_max = (20_000.0, 10_000.0)
 
 cells_per_dimension = (64, 32)
-mesh = StructuredMesh(cells_per_dimension, coordinates_min, coordinates_max)
+mesh = StructuredMesh(cells_per_dimension, coordinates_min, coordinates_max,
+                      periodicity = (true, false))
 
 semi = SemidiscretizationHyperbolic(mesh, equations, warm_bubble_setup, solver,
                                     source_terms = warm_bubble_setup,

src/meshes/tree_mesh.jl

@@ -228,5 +228,8 @@ function total_volume(mesh::TreeMesh)
     return mesh.tree.length_level_0^ndims(mesh)
 end
 
+isperiodic(mesh::TreeMesh) = isperiodic(mesh.tree)
+isperiodic(mesh::TreeMesh, dimension) = isperiodic(mesh.tree, dimension)
+
 include("parallel_tree_mesh.jl")
 end # @muladd

src/semidiscretization/semidiscretization_hyperbolic.jl

@@ -72,6 +72,8 @@ function SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver
     _boundary_conditions = digest_boundary_conditions(boundary_conditions, mesh, solver,
                                                       cache)
 
+    check_periodicity_mesh_boundary_conditions(mesh, _boundary_conditions)
+
     SemidiscretizationHyperbolic{typeof(mesh), typeof(equations),
                                  typeof(initial_condition),
                                  typeof(_boundary_conditions), typeof(source_terms),
@@ -210,6 +212,74 @@ function digest_boundary_conditions(boundary_conditions::AbstractArray, mesh, so
     throw(ArgumentError("Please use a (named) tuple instead of an (abstract) array to supply multiple boundary conditions (to improve performance)."))
 end
 
+# No checks for these meshes yet available
+function check_periodicity_mesh_boundary_conditions(mesh::Union{P4estMesh,
+                                                                UnstructuredMesh2D,
+                                                                T8codeMesh,
+                                                                DGMultiMesh},
+                                                    boundary_conditions)
+end
+
+# No actions needed for periodic boundary conditions
+function check_periodicity_mesh_boundary_conditions(mesh::Union{TreeMesh,
+                                                                StructuredMesh},
+                                                    boundary_conditions::BoundaryConditionPeriodic)
+end
+
+function check_periodicity_mesh_boundary_conditions_x(mesh, x_neg, x_pos)
+    if isperiodic(mesh, 1) &&
+       (x_neg != BoundaryConditionPeriodic() ||
+        x_pos != BoundaryConditionPeriodic())
+        @error "For periodic mesh non-periodic boundary conditions in x-direction are supplied."
+    end
+end
+
+function check_periodicity_mesh_boundary_conditions_y(mesh, y_neg, y_pos)
+    if isperiodic(mesh, 2) &&
+       (y_neg != BoundaryConditionPeriodic() ||
+        y_pos != BoundaryConditionPeriodic())
+        @error "For periodic mesh non-periodic boundary conditions in y-direction are supplied."
+    end
+end
+
+function check_periodicity_mesh_boundary_conditions_z(mesh, z_neg, z_pos)
+    if isperiodic(mesh, 3) &&
+       (z_neg != BoundaryConditionPeriodic() ||
+        z_pos != BoundaryConditionPeriodic())
+        @error "For periodic mesh non-periodic boundary conditions in z-direction are supplied."
+    end
+end
+
+function check_periodicity_mesh_boundary_conditions(mesh::Union{TreeMesh{1},
+                                                                StructuredMesh{1}},
+                                                    boundary_conditions::Union{NamedTuple,
+                                                                               Tuple})
+    check_periodicity_mesh_boundary_conditions_x(mesh, boundary_conditions[1],
+                                                 boundary_conditions[2])
+end
+
+function check_periodicity_mesh_boundary_conditions(mesh::Union{TreeMesh{2},
+                                                                StructuredMesh{2}},
+                                                    boundary_conditions::Union{NamedTuple,
+                                                                               Tuple})
+    check_periodicity_mesh_boundary_conditions_x(mesh, boundary_conditions[1],
+                                                 boundary_conditions[2])
+    check_periodicity_mesh_boundary_conditions_y(mesh, boundary_conditions[3],
+                                                 boundary_conditions[4])
+end
+
+function check_periodicity_mesh_boundary_conditions(mesh::Union{TreeMesh{3},
+                                                                StructuredMesh{3}},
+                                                    boundary_conditions::Union{NamedTuple,
+                                                                               Tuple})
+    check_periodicity_mesh_boundary_conditions_x(mesh, boundary_conditions[1],
+                                                 boundary_conditions[2])
+    check_periodicity_mesh_boundary_conditions_y(mesh, boundary_conditions[3],
+                                                 boundary_conditions[4])
+    check_periodicity_mesh_boundary_conditions_z(mesh, boundary_conditions[5],
+                                                 boundary_conditions[6])
+end
+
 function Base.show(io::IO, semi::SemidiscretizationHyperbolic)
     @nospecialize semi # reduce precompilation time
 

src/semidiscretization/semidiscretization_hyperbolic_parabolic.jl

@@ -136,6 +136,8 @@ function SemidiscretizationHyperbolicParabolic(mesh, equations, equations_parabo
     _boundary_conditions_parabolic = digest_boundary_conditions(boundary_conditions_parabolic,
                                                                 mesh, solver, cache)
 
+    check_periodicity_mesh_boundary_conditions(mesh, _boundary_conditions)
+
     cache_parabolic = (;
                        create_cache_parabolic(mesh, equations, equations_parabolic,
                                               solver, solver_parabolic, RealT,