Assure conservation for SSP scheme (#1640)

* Add denominator variable for SSP scheme * Fix format * Implement suggestions
trixi-framework · Sep 16, 2023 · 73384ac · 73384ac
1 parent 6069149
commit 73384ac
Showing 1 changed file with 11 additions and 6 deletions.
diff --git a/src/time_integration/methods_SSP.jl b/src/time_integration/methods_SSP.jl
@@ -24,14 +24,16 @@ The third-order SSP Runge-Kutta method of Shu and Osher.
     This is an experimental feature and may change in future releases.
 """
 struct SimpleSSPRK33{StageCallbacks} <: SimpleAlgorithmSSP
-    a::SVector{3, Float64}
-    b::SVector{3, Float64}
+    numerator_a::SVector{3, Float64}
+    numerator_b::SVector{3, Float64}
+    denominator::SVector{3, Float64}
     c::SVector{3, Float64}
     stage_callbacks::StageCallbacks
 
     function SimpleSSPRK33(; stage_callbacks = ())
-        a = SVector(0.0, 3 / 4, 1 / 3)
-        b = SVector(1.0, 1 / 4, 2 / 3)
+        numerator_a = SVector(0.0, 3.0, 1.0) # a = numerator_a / denominator
+        numerator_b = SVector(1.0, 1.0, 2.0) # b = numerator_b / denominator
+        denominator = SVector(1.0, 4.0, 3.0)
         c = SVector(0.0, 1.0, 1 / 2)
 
         # Butcher tableau
@@ -42,7 +44,8 @@ struct SimpleSSPRK33{StageCallbacks} <: SimpleAlgorithmSSP
         # --------------------
         #   b | 1/6  1/6  2/3
 
-        new{typeof(stage_callbacks)}(a, b, c, stage_callbacks)
+        new{typeof(stage_callbacks)}(numerator_a, numerator_b, denominator, c,
+                                     stage_callbacks)
     end
 end
 
@@ -166,7 +169,9 @@ function solve!(integrator::SimpleIntegratorSSP)
             end
 
             # perform convex combination
-            @. integrator.u = alg.a[stage] * integrator.r0 + alg.b[stage] * integrator.u
+            @. integrator.u = (alg.numerator_a[stage] * integrator.r0 +
+                               alg.numerator_b[stage] * integrator.u) /
+                              alg.denominator[stage]
         end
 
         integrator.iter += 1