From 53dfc00487c4527dbad287db70ea222fb6459e43 Mon Sep 17 00:00:00 2001 From: patrickersing Date: Tue, 26 Mar 2024 16:38:17 +0100 Subject: [PATCH] add changes from code review --- ...r_shallowwater_multilayer_well_balanced.jl | 2 +- src/equations/shallow_water_multilayer_1d.jl | 31 ++++++++++++------- src/equations/shallow_water_two_layer_1d.jl | 4 +-- src/equations/shallow_water_two_layer_2d.jl | 4 +-- src/equations/shallow_water_wet_dry_1d.jl | 2 +- src/equations/shallow_water_wet_dry_2d.jl | 2 +- 6 files changed, 27 insertions(+), 18 deletions(-) diff --git a/examples/tree_1d_dgsem/elixir_shallowwater_multilayer_well_balanced.jl b/examples/tree_1d_dgsem/elixir_shallowwater_multilayer_well_balanced.jl index d539461..4bc9c63 100644 --- a/examples/tree_1d_dgsem/elixir_shallowwater_multilayer_well_balanced.jl +++ b/examples/tree_1d_dgsem/elixir_shallowwater_multilayer_well_balanced.jl @@ -12,7 +12,7 @@ equations = ShallowWaterMultiLayerEquations1D(gravity_constant = 1.0, H0 = 0.7, """ initial_condition_fjordholm_well_balanced(x, t, equations::ShallowWaterMultiLayerEquations1D) -Initial condition to test well balanced with a bottom topography from Fjordholm +Initial condition to test well balanced with a bottom topography adapted from Fjordholm - Ulrik Skre Fjordholm (2012) Energy conservative and stable schemes for the two-layer shallow water equations. [DOI: 10.1142/9789814417099_0039](https://doi.org/10.1142/9789814417099_0039) diff --git a/src/equations/shallow_water_multilayer_1d.jl b/src/equations/shallow_water_multilayer_1d.jl index 4861297..75e041c 100644 --- a/src/equations/shallow_water_multilayer_1d.jl +++ b/src/equations/shallow_water_multilayer_1d.jl @@ -20,10 +20,13 @@ Multi-Layer Shallow Water equations (MLSWE) in one space dimension. The equation where ``m = 1, 2, ..., M`` is the layer index and the unknown variables are the water height ``h`` and the velocity ``v``. Furthermore, ``g`` denotes the gravitational constant, ``b(x)`` the bottom -topography and ``\rho`` the layer density, that must be chosen such that -``\rho_1 < \rho_2 < ... < \rho_M``, to make sure that layers are ordered from top to bottom, with +topography and ``\rho_m`` the m-th layer density, that must be chosen such that +``\rho_1 < \rho_2 < ... < \rho_M``, to ensure that different layers are ordered from top to bottom, with increasing density. +We use a specific formulation of the system, where the pressure term is reformulated as a +nonconservative term, which has some benefits for the design of well-balanced schemes. + The additional quantity ``H_0`` is also available to store a reference value for the total water height that is useful to set initial conditions or test the "lake-at-rest" well-balancedness. @@ -129,6 +132,11 @@ A smooth initial condition for a three-layer configuration used for convergence """ function Trixi.initial_condition_convergence_test(x, t, equations::ShallowWaterMultiLayerEquations1D) + # Check if the equation has been set up with three layers + if nlayers(equations) != 3 + throw(ArgumentError("This initial condition is only valid for three layers")) + end + # Some constants are chosen such that the function is periodic on the domain [0,sqrt(2)] ω = 2.0 * pi * sqrt(2.0) @@ -226,7 +234,7 @@ For details see Section 9.2.5 of the book: return f end -# Calculate 1D flux for a single point +# Calculate 1D advective portion of the flux for a single point @inline function Trixi.flux(u, orientation::Integer, equations::ShallowWaterMultiLayerEquations1D) # Extract waterheight and momentum and compute velocities @@ -260,7 +268,7 @@ In the two-layer setting this combination is equivalent to the fluxes in: - Patrick Ersing, Andrew R. Winters (2023) An entropy stable discontinuous Galerkin method for the two-layer shallow water equations on curvilinear meshes - [DOI: 10.48550/arXiv.2306.12699](https://doi.org/10.48550/arXiv.2306.12699) + [DOI: 10.1007/s10915-024-02451-2](https://doi.org/10.1007/s10915-024-02451-2) """ @inline function Trixi.flux_nonconservative_ersing_etal(u_ll, u_rr, orientation::Integer, @@ -279,8 +287,9 @@ In the two-layer setting this combination is equivalent to the fluxes in: # Initialize flux vector f = zero(MVector{2 * nlayers(equations) + 1, real(equations)}) - # Compute the nonconservative flux in each layer (0, ..., 0, f_hv[1], ..., f_hv[NLAYERS], 0) - # with f_hv[i] = gh[i] * (b + ∑h[k] + ∑σ[k]h[k])_x. + # Compute the nonconservative flux in each layer (0, ..., 0, f_hv[1], ..., f_hv[NLAYERS], 0), + # where f_hv[i] = gh[i] * (b + ∑h[k] + ∑σ[k]h[k])_x and σ[k] = ρ[k] / ρ[i] denotes the density + # ratio of different layers for i in eachlayer(equations) f_hv = g * h_ll[i] * b_jump for j in eachlayer(equations) @@ -313,7 +322,7 @@ In the two-layer setting this combination is equivalent to the fluxes in: - Patrick Ersing, Andrew R. Winters (2023) An entropy stable discontinuous Galerkin method for the two-layer shallow water equations on curvilinear meshes - [DOI: 10.48550/arXiv.2306.12699](https://doi.org/10.48550/arXiv.2306.12699) + [DOI: 10.1007/s10915-024-02451-2](https://doi.org/10.1007/s10915-024-02451-2) """ @inline function flux_ersing_etal(u_ll, u_rr, orientation::Integer, @@ -371,7 +380,6 @@ end # Calculate approximation for maximum wave speed for local Lax-Friedrichs-type dissipation as the # maximum velocity magnitude plus the maximum speed of sound. This function uses approximate # eigenvalues as there is no simple way to calculate them analytically. - @inline function Trixi.max_abs_speed_naive(u_ll, u_rr, orientation::Integer, equations::ShallowWaterMultiLayerEquations1D) @@ -437,8 +445,8 @@ end end # Convert conservative variables to entropy variables -# Note, only the first four are the entropy variables, the fifth entry still just carries the -# bottom topography values for convenience +# Note, only the first `2 * nlayers(equations)` are the entropy variables +# The last entry still just carries the bottom topography values for convenience @inline function Trixi.cons2entropy(u, equations::ShallowWaterMultiLayerEquations1D) # Extract conservative variables and compute velocity h = waterheight(u, equations) @@ -451,7 +459,8 @@ end # Calculate entropy variables in each layer for i in eachlayer(equations) - # Compute w1[i] = ρ[i]g * (b + ∑h[k] + ∑σ[k]h[k]) + # Compute w1[i] = ρ[i]g * (b + ∑h[k] + ∑σ[k]h[k]), where σ[k] = ρ[k] / ρ[i] denotes the + # density ratio of different layers w1 = equations.rhos[i] * (g * b - 0.5 * v[i]^2) for j in eachlayer(equations) if j < i diff --git a/src/equations/shallow_water_two_layer_1d.jl b/src/equations/shallow_water_two_layer_1d.jl index cd114b4..3f127ce 100644 --- a/src/equations/shallow_water_two_layer_1d.jl +++ b/src/equations/shallow_water_two_layer_1d.jl @@ -227,7 +227,7 @@ For further details see: - Patrick Ersing, Andrew R. Winters (2023) An entropy stable discontinuous Galerkin method for the two-layer shallow water equations on curvilinear meshes - [DOI: 10.48550/arXiv.2306.12699](https://doi.org/10.48550/arXiv.2306.12699) + [DOI: 10.1007/s10915-024-02451-2](https://doi.org/10.1007/s10915-024-02451-2) """ @inline function Trixi.flux_nonconservative_ersing_etal(u_ll, u_rr, orientation::Integer, @@ -308,7 +308,7 @@ For further details see: - Patrick Ersing, Andrew R. Winters (2023) An entropy stable discontinuous Galerkin method for the two-layer shallow water equations on curvilinear meshes - [DOI: 10.48550/arXiv.2306.12699](https://doi.org/10.48550/arXiv.2306.12699) + [DOI: 10.1007/s10915-024-02451-2](https://doi.org/10.1007/s10915-024-02451-2) """ @inline function flux_es_ersing_etal(u_ll, u_rr, orientation::Integer, diff --git a/src/equations/shallow_water_two_layer_2d.jl b/src/equations/shallow_water_two_layer_2d.jl index ac0df7d..87c9dbe 100644 --- a/src/equations/shallow_water_two_layer_2d.jl +++ b/src/equations/shallow_water_two_layer_2d.jl @@ -340,7 +340,7 @@ For further details see: - Patrick Ersing, Andrew R. Winters (2023) An entropy stable discontinuous Galerkin method for the two-layer shallow water equations on curvilinear meshes - [DOI: 10.48550/arXiv.2306.12699](https://doi.org/10.48550/arXiv.2306.12699) + [DOI: 10.1007/s10915-024-02451-2](https://doi.org/10.1007/s10915-024-02451-2) """ @inline function Trixi.flux_nonconservative_ersing_etal(u_ll, u_rr, orientation::Integer, @@ -514,7 +514,7 @@ For further details see: - Patrick Ersing, Andrew R. Winters (2023) An entropy stable discontinuous Galerkin method for the two-layer shallow water equations on curvilinear meshes - [DOI: 10.48550/arXiv.2306.12699](https://doi.org/10.48550/arXiv.2306.12699) + [DOI: 10.1007/s10915-024-02451-2](https://doi.org/10.1007/s10915-024-02451-2) """ @inline function flux_es_ersing_etal(u_ll, u_rr, orientation_or_normal_direction, diff --git a/src/equations/shallow_water_wet_dry_1d.jl b/src/equations/shallow_water_wet_dry_1d.jl index f3eaa69..6111e4a 100644 --- a/src/equations/shallow_water_wet_dry_1d.jl +++ b/src/equations/shallow_water_wet_dry_1d.jl @@ -333,7 +333,7 @@ For further details see: - Patrick Ersing, Andrew R. Winters (2023) An entropy stable discontinuous Galerkin method for the two-layer shallow water equations on curvilinear meshes - [DOI: 10.48550/arXiv.2306.12699](https://doi.org/10.48550/arXiv.2306.12699) + [DOI: 10.1007/s10915-024-02451-2](https://doi.org/10.1007/s10915-024-02451-2) """ @inline function Trixi.flux_nonconservative_ersing_etal(u_ll, u_rr, orientation::Integer, diff --git a/src/equations/shallow_water_wet_dry_2d.jl b/src/equations/shallow_water_wet_dry_2d.jl index 7762168..5c4aa23 100644 --- a/src/equations/shallow_water_wet_dry_2d.jl +++ b/src/equations/shallow_water_wet_dry_2d.jl @@ -559,7 +559,7 @@ For further details see: - Patrick Ersing, Andrew R. Winters (2023) An entropy stable discontinuous Galerkin method for the two-layer shallow water equations on curvilinear meshes - [DOI: 10.48550/arXiv.2306.12699](https://doi.org/10.48550/arXiv.2306.12699) + [DOI: 10.1007/s10915-024-02451-2](https://doi.org/10.1007/s10915-024-02451-2) """ @inline function Trixi.flux_nonconservative_ersing_etal(u_ll, u_rr, orientation::Integer,