From a4a1fe80279c20e53dedeea567fa8f1b9f33f9ad Mon Sep 17 00:00:00 2001
From: bennibolm <benjamin.bolm@gmx.de>
Date: Thu, 12 Oct 2023 15:30:48 +0200
Subject: [PATCH 1/7] Allow arbitrary nonlinear positivity limiting

---
 src/equations/compressible_euler_2d.jl        |  5 +-
 src/solvers/dgsem_tree/subcell_limiters_2d.jl | 46 +++++++++----------
 2 files changed, 26 insertions(+), 25 deletions(-)

diff --git a/src/equations/compressible_euler_2d.jl b/src/equations/compressible_euler_2d.jl
index d95fbdbcfc9..782c8f253a9 100644
--- a/src/equations/compressible_euler_2d.jl
+++ b/src/equations/compressible_euler_2d.jl
@@ -1572,7 +1572,7 @@ end
 end
 
 # Transformation from conservative variables u to d(p)/d(u)
-@inline function dpdu(u, equations::CompressibleEulerEquations2D)
+@inline function pressure(u, equations::CompressibleEulerEquations2D, derivative::True)
     rho, rho_v1, rho_v2, rho_e = u
 
     v1 = rho_v1 / rho
@@ -1581,6 +1581,9 @@ end
 
     return (equations.gamma - 1.0) * SVector(0.5 * v_square, -v1, -v2, 1.0)
 end
+@inline function pressure(u, equations::CompressibleEulerEquations2D, derivative::False)
+    return pressure(u, equations)
+end
 
 @inline function entropy2cons(w, equations::CompressibleEulerEquations2D)
     # See Hughes, Franca, Mallet (1986) A new finite element formulation for CFD
diff --git a/src/solvers/dgsem_tree/subcell_limiters_2d.jl b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
index a8fcc42c061..c71b3436749 100644
--- a/src/solvers/dgsem_tree/subcell_limiters_2d.jl
+++ b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
@@ -554,27 +554,25 @@ end
 
             # Perform Newton's bisection method to find new alpha
             newton_loops_alpha!(alpha, var_min[i, j, element], u_local, i, j, element,
-                                pressure_goal, pressure_dgoal_dbeta,
-                                pressure_initialCheck, pressure_finalCheck,
-                                dt, mesh, equations, dg, cache, limiter)
+                                variable, dt, mesh, equations, dg, cache, limiter)
         end
     end
 
     return nothing
 end
 
-pressure_goal(bound, u, equations) = bound - pressure(u, equations)
-function pressure_dgoal_dbeta(u, dt, antidiffusive_flux, equations)
-    -dot(dpdu(u, equations), dt * antidiffusive_flux)
+goal_function(variable, bound, u, equations) = bound - variable(u, equations)
+function dgoal_function(variable, u, dt, antidiffusive_flux, equations)
+    -dot(variable(u, equations, true), dt * antidiffusive_flux)
 end
-pressure_initialCheck(bound, goal, newton_abstol) = goal <= 0
-function pressure_finalCheck(bound, goal, newton_abstol)
+
+initialCheck(bound, goal, newton_abstol) = goal <= 0
+function finalCheck(bound, goal, newton_abstol)
     (goal <= eps()) && (goal > -max(newton_abstol, abs(bound) * newton_abstol))
 end
 
 @inline function newton_loops_alpha!(alpha, bound, u, i, j, element,
-                                     goal_fct, dgoal_fct, initialCheck, finalCheck,
-                                     dt, mesh, equations, dg, cache, limiter)
+                                     variable, dt, mesh, equations, dg, cache, limiter)
     @unpack inverse_weights = dg.basis
     @unpack antidiffusive_flux1, antidiffusive_flux2 = cache.antidiffusive_fluxes
     if mesh isa TreeMesh
@@ -589,37 +587,36 @@ end
     antidiffusive_flux = gamma_constant_newton * inverse_jacobian * inverse_weights[i] *
                          get_node_vars(antidiffusive_flux1, equations, dg, i, j,
                                        element)
-    newton_loop!(alpha, bound, u, i, j, element, goal_fct, dgoal_fct, initialCheck,
-                 finalCheck, equations, dt, limiter, antidiffusive_flux)
+    newton_loop!(alpha, bound, u, i, j, element, variable, equations, dt, limiter,
+                 antidiffusive_flux)
 
     # positive xi direction
     antidiffusive_flux = -gamma_constant_newton * inverse_jacobian *
                          inverse_weights[i] *
                          get_node_vars(antidiffusive_flux1, equations, dg, i + 1, j,
                                        element)
-    newton_loop!(alpha, bound, u, i, j, element, goal_fct, dgoal_fct, initialCheck,
-                 finalCheck, equations, dt, limiter, antidiffusive_flux)
+    newton_loop!(alpha, bound, u, i, j, element, variable, equations, dt, limiter,
+                 antidiffusive_flux)
 
     # negative eta direction
     antidiffusive_flux = gamma_constant_newton * inverse_jacobian * inverse_weights[j] *
                          get_node_vars(antidiffusive_flux2, equations, dg, i, j,
                                        element)
-    newton_loop!(alpha, bound, u, i, j, element, goal_fct, dgoal_fct, initialCheck,
-                 finalCheck, equations, dt, limiter, antidiffusive_flux)
+    newton_loop!(alpha, bound, u, i, j, element, variable, equations, dt, limiter,
+                 antidiffusive_flux)
 
     # positive eta direction
     antidiffusive_flux = -gamma_constant_newton * inverse_jacobian *
                          inverse_weights[j] *
                          get_node_vars(antidiffusive_flux2, equations, dg, i, j + 1,
                                        element)
-    newton_loop!(alpha, bound, u, i, j, element, goal_fct, dgoal_fct, initialCheck,
-                 finalCheck, equations, dt, limiter, antidiffusive_flux)
+    newton_loop!(alpha, bound, u, i, j, element, variable, equations, dt, limiter,
+                 antidiffusive_flux)
 
     return nothing
 end
 
-@inline function newton_loop!(alpha, bound, u, i, j, element,
-                              goal_fct, dgoal_fct, initialCheck, finalCheck,
+@inline function newton_loop!(alpha, bound, u, i, j, element, variable,
                               equations, dt, limiter, antidiffusive_flux)
     newton_reltol, newton_abstol = limiter.newton_tolerances
 
@@ -632,7 +629,7 @@ end
 
     # If state is valid, perform initial check and return if correction is not needed
     if isValidState(u_curr, equations)
-        as = goal_fct(bound, u_curr, equations)
+        as = goal_function(variable, bound, u_curr, equations)
 
         initialCheck(bound, as, newton_abstol) && return nothing
     end
@@ -643,7 +640,8 @@ end
 
         # If the state is valid, evaluate d(goal)/d(beta)
         if isValidState(u_curr, equations)
-            dSdbeta = dgoal_fct(u_curr, dt, antidiffusive_flux, equations)
+            dSdbeta = dgoal_function(variable, u_curr, dt, antidiffusive_flux,
+                                     equations)
         else # Otherwise, perform a bisection step
             dSdbeta = 0
         end
@@ -667,7 +665,7 @@ end
             end
 
             # Check new beta for condition and update bounds
-            as = goal_fct(bound, u_curr, equations)
+            as = goal_function(variable, bound, u_curr, equations)
             if initialCheck(bound, as, newton_abstol)
                 # New beta fulfills condition
                 beta_L = beta
@@ -686,7 +684,7 @@ end
             end
 
             # Evaluate goal function
-            as = goal_fct(bound, u_curr, equations)
+            as = goal_function(variable, bound, u_curr, equations)
         end
 
         # Check relative tolerance

From 2eac4cf7444d9b5da4786ac24b6670ef7ed96b03 Mon Sep 17 00:00:00 2001
From: bennibolm <benjamin.bolm@gmx.de>
Date: Thu, 12 Oct 2023 16:13:54 +0200
Subject: [PATCH 2/7] Implement changes for limiting of entropys

---
 src/equations/compressible_euler_2d.jl        |  2 +
 src/solvers/dgsem_tree/subcell_limiters_2d.jl | 63 +++++++++----------
 2 files changed, 31 insertions(+), 34 deletions(-)

diff --git a/src/equations/compressible_euler_2d.jl b/src/equations/compressible_euler_2d.jl
index 782c8f253a9..ddf769049c6 100644
--- a/src/equations/compressible_euler_2d.jl
+++ b/src/equations/compressible_euler_2d.jl
@@ -1544,6 +1544,7 @@ end
 
     return SVector(w1, w2, w3, w4)
 end
+entropy_math(u, equations, derivative::True) = cons2entropy(u, equations)
 
 # Transformation from conservative variables u to entropy vector dSdu, S = -rho*s/(gamma-1), s=ln(p)-gamma*ln(rho)
 @inline function cons2entropy_spec(u, equations::CompressibleEulerEquations2D)
@@ -1570,6 +1571,7 @@ end
 
     return SVector(w1, w2, w3, w4)
 end
+entropy_spec(u, equations, derivative::True) = cons2entropy_spec(u, equations)
 
 # Transformation from conservative variables u to d(p)/d(u)
 @inline function pressure(u, equations::CompressibleEulerEquations2D, derivative::True)
diff --git a/src/solvers/dgsem_tree/subcell_limiters_2d.jl b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
index c71b3436749..9e120db4e33 100644
--- a/src/solvers/dgsem_tree/subcell_limiters_2d.jl
+++ b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
@@ -408,8 +408,8 @@ end
         for j in eachnode(dg), i in eachnode(dg)
             u_local = get_node_vars(u, equations, dg, i, j, element)
             newton_loops_alpha!(alpha, s_min[i, j, element], u_local, i, j, element,
-                                specEntropy_goal, specEntropy_dGoal_dbeta,
-                                specEntropy_initialCheck, standard_finalCheck,
+                                entropy_spec, initial_check_entropy_spec,
+                                final_check_standard,
                                 dt, mesh, equations, dg, cache, limiter)
         end
     end
@@ -417,11 +417,7 @@ end
     return nothing
 end
 
-specEntropy_goal(bound, u, equations) = bound - entropy_spec(u, equations)
-function specEntropy_dGoal_dbeta(u, dt, antidiffusive_flux, equations)
-    -dot(cons2entropy_spec(u, equations), dt * antidiffusive_flux)
-end
-function specEntropy_initialCheck(bound, goal, newton_abstol)
+function initial_check_entropy_spec(bound, goal, newton_abstol)
     goal <= max(newton_abstol, abs(bound) * newton_abstol)
 end
 
@@ -438,8 +434,8 @@ end
         for j in eachnode(dg), i in eachnode(dg)
             u_local = get_node_vars(u, equations, dg, i, j, element)
             newton_loops_alpha!(alpha, s_max[i, j, element], u_local, i, j, element,
-                                mathEntropy_goal, mathEntropy_dGoal_dbeta,
-                                mathEntropy_initialCheck, standard_finalCheck,
+                                entropy_math, initial_check_entropy_math,
+                                final_check_standard,
                                 dt, mesh, equations, dg, cache, limiter)
         end
     end
@@ -447,11 +443,7 @@ end
     return nothing
 end
 
-mathEntropy_goal(bound, u, equations) = bound - entropy_math(u, equations)
-function mathEntropy_dGoal_dbeta(u, dt, antidiffusive_flux, equations)
-    -dot(cons2entropy(u, equations), dt * antidiffusive_flux)
-end
-function mathEntropy_initialCheck(bound, goal, newton_abstol)
+function initial_check_entropy_math(bound, goal, newton_abstol)
     goal >= -max(newton_abstol, abs(bound) * newton_abstol)
 end
 
@@ -554,7 +546,9 @@ end
 
             # Perform Newton's bisection method to find new alpha
             newton_loops_alpha!(alpha, var_min[i, j, element], u_local, i, j, element,
-                                variable, dt, mesh, equations, dg, cache, limiter)
+                                variable, initial_check_nonnegative,
+                                final_check_nonnegative,
+                                dt, mesh, equations, dg, cache, limiter)
         end
     end
 
@@ -563,16 +557,17 @@ end
 
 goal_function(variable, bound, u, equations) = bound - variable(u, equations)
 function dgoal_function(variable, u, dt, antidiffusive_flux, equations)
-    -dot(variable(u, equations, true), dt * antidiffusive_flux)
+    -dot(variable(u, equations, True()), dt * antidiffusive_flux)
 end
 
-initialCheck(bound, goal, newton_abstol) = goal <= 0
-function finalCheck(bound, goal, newton_abstol)
+initial_check_nonnegative(bound, goal, newton_abstol) = goal <= 0
+function final_check_nonnegative(bound, goal, newton_abstol)
     (goal <= eps()) && (goal > -max(newton_abstol, abs(bound) * newton_abstol))
 end
 
-@inline function newton_loops_alpha!(alpha, bound, u, i, j, element,
-                                     variable, dt, mesh, equations, dg, cache, limiter)
+@inline function newton_loops_alpha!(alpha, bound, u, i, j, element, variable,
+                                     initial_check, final_check, dt, mesh, equations,
+                                     dg, cache, limiter)
     @unpack inverse_weights = dg.basis
     @unpack antidiffusive_flux1, antidiffusive_flux2 = cache.antidiffusive_fluxes
     if mesh isa TreeMesh
@@ -587,37 +582,37 @@ end
     antidiffusive_flux = gamma_constant_newton * inverse_jacobian * inverse_weights[i] *
                          get_node_vars(antidiffusive_flux1, equations, dg, i, j,
                                        element)
-    newton_loop!(alpha, bound, u, i, j, element, variable, equations, dt, limiter,
-                 antidiffusive_flux)
+    newton_loop!(alpha, bound, u, i, j, element, variable, initial_check, final_check,
+                 equations, dt, limiter, antidiffusive_flux)
 
     # positive xi direction
     antidiffusive_flux = -gamma_constant_newton * inverse_jacobian *
                          inverse_weights[i] *
                          get_node_vars(antidiffusive_flux1, equations, dg, i + 1, j,
                                        element)
-    newton_loop!(alpha, bound, u, i, j, element, variable, equations, dt, limiter,
-                 antidiffusive_flux)
+    newton_loop!(alpha, bound, u, i, j, element, variable, initial_check, final_check,
+                 equations, dt, limiter, antidiffusive_flux)
 
     # negative eta direction
     antidiffusive_flux = gamma_constant_newton * inverse_jacobian * inverse_weights[j] *
                          get_node_vars(antidiffusive_flux2, equations, dg, i, j,
                                        element)
-    newton_loop!(alpha, bound, u, i, j, element, variable, equations, dt, limiter,
-                 antidiffusive_flux)
+    newton_loop!(alpha, bound, u, i, j, element, variable, initial_check, final_check,
+                 equations, dt, limiter, antidiffusive_flux)
 
     # positive eta direction
     antidiffusive_flux = -gamma_constant_newton * inverse_jacobian *
                          inverse_weights[j] *
                          get_node_vars(antidiffusive_flux2, equations, dg, i, j + 1,
                                        element)
-    newton_loop!(alpha, bound, u, i, j, element, variable, equations, dt, limiter,
-                 antidiffusive_flux)
+    newton_loop!(alpha, bound, u, i, j, element, variable, initial_check, final_check,
+                 equations, dt, limiter, antidiffusive_flux)
 
     return nothing
 end
 
-@inline function newton_loop!(alpha, bound, u, i, j, element, variable,
-                              equations, dt, limiter, antidiffusive_flux)
+@inline function newton_loop!(alpha, bound, u, i, j, element, variable, initial_check,
+                              final_check, equations, dt, limiter, antidiffusive_flux)
     newton_reltol, newton_abstol = limiter.newton_tolerances
 
     beta = 1 - alpha[i, j, element]
@@ -631,7 +626,7 @@ end
     if isValidState(u_curr, equations)
         as = goal_function(variable, bound, u_curr, equations)
 
-        initialCheck(bound, as, newton_abstol) && return nothing
+        initial_check(bound, as, newton_abstol) && return nothing
     end
 
     # Newton iterations
@@ -666,7 +661,7 @@ end
 
             # Check new beta for condition and update bounds
             as = goal_function(variable, bound, u_curr, equations)
-            if initialCheck(bound, as, newton_abstol)
+            if initial_check(bound, as, newton_abstol)
                 # New beta fulfills condition
                 beta_L = beta
             else
@@ -693,7 +688,7 @@ end
         end
 
         # Check absolute tolerance
-        if finalCheck(bound, as, newton_abstol)
+        if final_check(bound, as, newton_abstol)
             break
         end
     end
@@ -708,7 +703,7 @@ end
     return nothing
 end
 
-function standard_finalCheck(bound, goal, newton_abstol)
+function final_check_standard(bound, goal, newton_abstol)
     abs(goal) < max(newton_abstol, abs(bound) * newton_abstol)
 end
 

From 651b7169d6d882f052b45eda7eedfd023d286667 Mon Sep 17 00:00:00 2001
From: bennibolm <benjamin.bolm@gmx.de>
Date: Thu, 12 Oct 2023 16:23:17 +0200
Subject: [PATCH 3/7] Relocate checks and goal functions

---
 src/solvers/dgsem_tree/subcell_limiters_2d.jl | 39 ++++++++++---------
 1 file changed, 21 insertions(+), 18 deletions(-)

diff --git a/src/solvers/dgsem_tree/subcell_limiters_2d.jl b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
index 9e120db4e33..14ffef07fe8 100644
--- a/src/solvers/dgsem_tree/subcell_limiters_2d.jl
+++ b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
@@ -417,10 +417,6 @@ end
     return nothing
 end
 
-function initial_check_entropy_spec(bound, goal, newton_abstol)
-    goal <= max(newton_abstol, abs(bound) * newton_abstol)
-end
-
 @inline function idp_math_entropy!(alpha, limiter, u, t, dt, semi, elements)
     mesh, equations, dg, cache = mesh_equations_solver_cache(semi)
     (; variable_bounds) = limiter.cache.subcell_limiter_coefficients
@@ -443,10 +439,6 @@ end
     return nothing
 end
 
-function initial_check_entropy_math(bound, goal, newton_abstol)
-    goal >= -max(newton_abstol, abs(bound) * newton_abstol)
-end
-
 @inline function idp_positivity!(alpha, limiter, u, dt, semi, elements)
     # Conservative variables
     for variable in limiter.positivity_variables_cons
@@ -555,16 +547,6 @@ end
     return nothing
 end
 
-goal_function(variable, bound, u, equations) = bound - variable(u, equations)
-function dgoal_function(variable, u, dt, antidiffusive_flux, equations)
-    -dot(variable(u, equations, True()), dt * antidiffusive_flux)
-end
-
-initial_check_nonnegative(bound, goal, newton_abstol) = goal <= 0
-function final_check_nonnegative(bound, goal, newton_abstol)
-    (goal <= eps()) && (goal > -max(newton_abstol, abs(bound) * newton_abstol))
-end
-
 @inline function newton_loops_alpha!(alpha, bound, u, i, j, element, variable,
                                      initial_check, final_check, dt, mesh, equations,
                                      dg, cache, limiter)
@@ -703,6 +685,27 @@ end
     return nothing
 end
 
+# Initial checks
+function initial_check_entropy_spec(bound, goal, newton_abstol)
+    goal <= max(newton_abstol, abs(bound) * newton_abstol)
+end
+
+function initial_check_entropy_math(bound, goal, newton_abstol)
+    goal >= -max(newton_abstol, abs(bound) * newton_abstol)
+end
+
+initial_check_nonnegative(bound, goal, newton_abstol) = goal <= 0
+function final_check_nonnegative(bound, goal, newton_abstol)
+    (goal <= eps()) && (goal > -max(newton_abstol, abs(bound) * newton_abstol))
+end
+
+# Goal and d(Goal)d(u) function
+goal_function(variable, bound, u, equations) = bound - variable(u, equations)
+function dgoal_function(variable, u, dt, antidiffusive_flux, equations)
+    -dot(variable(u, equations, True()), dt * antidiffusive_flux)
+end
+
+# Final check
 function final_check_standard(bound, goal, newton_abstol)
     abs(goal) < max(newton_abstol, abs(bound) * newton_abstol)
 end

From 70b4fd846a1bde1ce596acf2d8161a5ed087baef Mon Sep 17 00:00:00 2001
From: bennibolm <benjamin.bolm@gmx.de>
Date: Thu, 12 Oct 2023 16:58:49 +0200
Subject: [PATCH 4/7] Rename functions

---
 src/equations/compressible_euler_2d.jl        | 8 ++++----
 src/solvers/dgsem_tree/subcell_limiters_2d.jl | 8 ++++----
 2 files changed, 8 insertions(+), 8 deletions(-)

diff --git a/src/equations/compressible_euler_2d.jl b/src/equations/compressible_euler_2d.jl
index ddf769049c6..cc7b86ef679 100644
--- a/src/equations/compressible_euler_2d.jl
+++ b/src/equations/compressible_euler_2d.jl
@@ -1574,7 +1574,7 @@ end
 entropy_spec(u, equations, derivative::True) = cons2entropy_spec(u, equations)
 
 # Transformation from conservative variables u to d(p)/d(u)
-@inline function pressure(u, equations::CompressibleEulerEquations2D, derivative::True)
+@inline function dpdu(u, equations::CompressibleEulerEquations2D)
     rho, rho_v1, rho_v2, rho_e = u
 
     v1 = rho_v1 / rho
@@ -1583,8 +1583,8 @@ entropy_spec(u, equations, derivative::True) = cons2entropy_spec(u, equations)
 
     return (equations.gamma - 1.0) * SVector(0.5 * v_square, -v1, -v2, 1.0)
 end
-@inline function pressure(u, equations::CompressibleEulerEquations2D, derivative::False)
-    return pressure(u, equations)
+@inline function pressure(u, equations::CompressibleEulerEquations2D, derivative::True)
+    return dpdu(u, equations)
 end
 
 @inline function entropy2cons(w, equations::CompressibleEulerEquations2D)
@@ -1611,7 +1611,7 @@ end
     return SVector(rho, rho_v1, rho_v2, rho_e)
 end
 
-@inline function isValidState(cons, equations::CompressibleEulerEquations2D)
+@inline function is_valid_state(cons, equations::CompressibleEulerEquations2D)
     p = pressure(cons, equations)
     if cons[1] <= 0.0 || p <= 0.0
         return false
diff --git a/src/solvers/dgsem_tree/subcell_limiters_2d.jl b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
index 14ffef07fe8..00a5996ef6e 100644
--- a/src/solvers/dgsem_tree/subcell_limiters_2d.jl
+++ b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
@@ -605,7 +605,7 @@ end
     u_curr = u + beta * dt * antidiffusive_flux
 
     # If state is valid, perform initial check and return if correction is not needed
-    if isValidState(u_curr, equations)
+    if is_valid_state(u_curr, equations)
         as = goal_function(variable, bound, u_curr, equations)
 
         initial_check(bound, as, newton_abstol) && return nothing
@@ -616,7 +616,7 @@ end
         beta_old = beta
 
         # If the state is valid, evaluate d(goal)/d(beta)
-        if isValidState(u_curr, equations)
+        if is_valid_state(u_curr, equations)
             dSdbeta = dgoal_function(variable, u_curr, dt, antidiffusive_flux,
                                      equations)
         else # Otherwise, perform a bisection step
@@ -636,7 +636,7 @@ end
             u_curr = u + beta * dt * antidiffusive_flux
 
             # If the state is invalid, finish bisection step without checking tolerance and iterate further
-            if !isValidState(u_curr, equations)
+            if !is_valid_state(u_curr, equations)
                 beta_R = beta
                 continue
             end
@@ -655,7 +655,7 @@ end
             u_curr = u + beta * dt * antidiffusive_flux
 
             # If the state is invalid, redefine right bound without checking tolerance and iterate further
-            if !isValidState(u_curr, equations)
+            if !is_valid_state(u_curr, equations)
                 beta_R = beta
                 continue
             end

From 463de3abccbe605612758e8b11254f3a36c72a74 Mon Sep 17 00:00:00 2001
From: bennibolm <benjamin.bolm@gmx.de>
Date: Fri, 13 Oct 2023 13:27:31 +0200
Subject: [PATCH 5/7] Add @inline

---
 src/equations/compressible_euler_2d.jl        | 4 ++--
 src/solvers/dgsem_tree/subcell_limiters_2d.jl | 7 ++++---
 2 files changed, 6 insertions(+), 5 deletions(-)

diff --git a/src/equations/compressible_euler_2d.jl b/src/equations/compressible_euler_2d.jl
index cc7b86ef679..8fd99bbfc6b 100644
--- a/src/equations/compressible_euler_2d.jl
+++ b/src/equations/compressible_euler_2d.jl
@@ -1544,7 +1544,7 @@ end
 
     return SVector(w1, w2, w3, w4)
 end
-entropy_math(u, equations, derivative::True) = cons2entropy(u, equations)
+@inline entropy_math(u, equations, derivative::True) = cons2entropy(u, equations)
 
 # Transformation from conservative variables u to entropy vector dSdu, S = -rho*s/(gamma-1), s=ln(p)-gamma*ln(rho)
 @inline function cons2entropy_spec(u, equations::CompressibleEulerEquations2D)
@@ -1571,7 +1571,7 @@ entropy_math(u, equations, derivative::True) = cons2entropy(u, equations)
 
     return SVector(w1, w2, w3, w4)
 end
-entropy_spec(u, equations, derivative::True) = cons2entropy_spec(u, equations)
+@inline entropy_spec(u, equations, derivative::True) = cons2entropy_spec(u, equations)
 
 # Transformation from conservative variables u to d(p)/d(u)
 @inline function dpdu(u, equations::CompressibleEulerEquations2D)
diff --git a/src/solvers/dgsem_tree/subcell_limiters_2d.jl b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
index 00a5996ef6e..0a778aa8813 100644
--- a/src/solvers/dgsem_tree/subcell_limiters_2d.jl
+++ b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
@@ -695,9 +695,6 @@ function initial_check_entropy_math(bound, goal, newton_abstol)
 end
 
 initial_check_nonnegative(bound, goal, newton_abstol) = goal <= 0
-function final_check_nonnegative(bound, goal, newton_abstol)
-    (goal <= eps()) && (goal > -max(newton_abstol, abs(bound) * newton_abstol))
-end
 
 # Goal and d(Goal)d(u) function
 goal_function(variable, bound, u, equations) = bound - variable(u, equations)
@@ -710,6 +707,10 @@ function final_check_standard(bound, goal, newton_abstol)
     abs(goal) < max(newton_abstol, abs(bound) * newton_abstol)
 end
 
+function final_check_nonnegative(bound, goal, newton_abstol)
+    (goal <= eps()) && (goal > -max(newton_abstol, abs(bound) * newton_abstol))
+end
+
 # this method is used when the limiter is constructed as for shock-capturing volume integrals
 function create_cache(limiter::Type{SubcellLimiterMCL}, equations::AbstractEquations{2},
                       basis::LobattoLegendreBasis, PressurePositivityLimiterKuzmin)

From e4d3a9d6b573908769db0e6c2ed6614a46d3e1d4 Mon Sep 17 00:00:00 2001
From: bennibolm <benjamin.bolm@gmx.de>
Date: Fri, 13 Oct 2023 13:39:17 +0200
Subject: [PATCH 6/7] Add @inline to checks and goal functions; Update test
 required

---
 src/solvers/dgsem_tree/subcell_limiters_2d.jl | 14 +++++++-------
 1 file changed, 7 insertions(+), 7 deletions(-)

diff --git a/src/solvers/dgsem_tree/subcell_limiters_2d.jl b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
index 0a778aa8813..0055b2181b0 100644
--- a/src/solvers/dgsem_tree/subcell_limiters_2d.jl
+++ b/src/solvers/dgsem_tree/subcell_limiters_2d.jl
@@ -686,28 +686,28 @@ end
 end
 
 # Initial checks
-function initial_check_entropy_spec(bound, goal, newton_abstol)
+@inline function initial_check_entropy_spec(bound, goal, newton_abstol)
     goal <= max(newton_abstol, abs(bound) * newton_abstol)
 end
 
-function initial_check_entropy_math(bound, goal, newton_abstol)
+@inline function initial_check_entropy_math(bound, goal, newton_abstol)
     goal >= -max(newton_abstol, abs(bound) * newton_abstol)
 end
 
-initial_check_nonnegative(bound, goal, newton_abstol) = goal <= 0
+@inline initial_check_nonnegative(bound, goal, newton_abstol) = goal <= 0
 
 # Goal and d(Goal)d(u) function
-goal_function(variable, bound, u, equations) = bound - variable(u, equations)
-function dgoal_function(variable, u, dt, antidiffusive_flux, equations)
+@inline goal_function(variable, bound, u, equations) = bound - variable(u, equations)
+@inline function dgoal_function(variable, u, dt, antidiffusive_flux, equations)
     -dot(variable(u, equations, True()), dt * antidiffusive_flux)
 end
 
 # Final check
-function final_check_standard(bound, goal, newton_abstol)
+@inline function final_check_standard(bound, goal, newton_abstol)
     abs(goal) < max(newton_abstol, abs(bound) * newton_abstol)
 end
 
-function final_check_nonnegative(bound, goal, newton_abstol)
+@inline function final_check_nonnegative(bound, goal, newton_abstol)
     (goal <= eps()) && (goal > -max(newton_abstol, abs(bound) * newton_abstol))
 end
 

From d04afadf3d80938dc25c17ba1d7741e81f8ee167 Mon Sep 17 00:00:00 2001
From: bennibolm <benjamin.bolm@gmx.de>
Date: Fri, 13 Oct 2023 15:25:38 +0200
Subject: [PATCH 7/7] Update tests

---
 test/test_structured_2d.jl | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/test/test_structured_2d.jl b/test/test_structured_2d.jl
index be1a1d3138b..f776de954bc 100644
--- a/test/test_structured_2d.jl
+++ b/test/test_structured_2d.jl
@@ -187,8 +187,8 @@ isdir(outdir) && rm(outdir, recursive=true)
 
   @trixi_testset "elixir_euler_source_terms_sc_subcell.jl" begin
     @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_source_terms_sc_subcell.jl"),
-      l2   = [0.008160127272557726, 0.008658253869683077, 0.009351900401871649, 0.02775701488343099],
-      linf = [0.027225608222781528, 0.0407340321806311, 0.0381940733564341, 0.08080650914262844],
+      l2   = [0.00816013114351954, 0.008658251709937477, 0.009351905651482216, 0.027757012781694318],
+      linf = [0.027225615981281148, 0.040734036539016305, 0.0381940733564341, 0.08080650914262844],
       tspan = (0.0, 0.5))
   end
 
@@ -248,7 +248,7 @@ isdir(outdir) && rm(outdir, recursive=true)
   @trixi_testset "elixir_euler_shock_upstream_sc_subcell.jl" begin
     @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_shock_upstream_sc_subcell.jl"),
       l2   = [1.2351468819080416, 1.1269856120551724, 1.7239124305681928, 11.715260007491556],
-      linf = [5.385491808683259, 6.575446013701839, 10.065227889186632, 51.008985921289565],
+      linf = [5.385492532917423, 6.575446146030286, 10.0652310822613, 51.00901293102744],
       cells_per_dimension = (8, 12),
       tspan = (0.0, 0.5))
   end