diff --git a/Project.toml b/Project.toml
index 6266fe97bb5..e0ce86c23ba 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "Trixi"
 uuid = "a7f1ee26-1774-49b1-8366-f1abc58fbfcb"
 authors = ["Michael Schlottke-Lakemper <m.schlottke-lakemper@hlrs.de>", "Gregor Gassner <ggassner@uni-koeln.de>", "Hendrik Ranocha <mail@ranocha.de>", "Andrew R. Winters <andrew.ross.winters@liu.se>", "Jesse Chan <jesse.chan@rice.edu>"]
-version = "0.9.6-DEV"
+version = "0.9.7-DEV"
 
 [deps]
 Accessors = "7d9f7c33-5ae7-4f3b-8dc6-eff91059b697"
diff --git a/src/solvers/dgmulti/sbp.jl b/src/solvers/dgmulti/sbp.jl
index 232555e18b5..76444f64899 100644
--- a/src/solvers/dgmulti/sbp.jl
+++ b/src/solvers/dgmulti/sbp.jl
@@ -221,7 +221,10 @@ function calc_volume_integral!(du, u, mesh::DGMultiMesh,
                 for id in nzrange(A_base, i)
                     j = rows[id]
                     u_j = u[j]
-                    A_ij = vals[id]
+
+                    # we use the negative of A_ij since A is skew-symmetric, 
+                    # and we are accessing the transpose of A. 
+                    A_ij = -vals[id]
                     AF_ij = 2 * A_ij *
                             volume_flux(u_i, u_j, normal_direction, equations)
                     du_i = du_i + AF_ij
diff --git a/test/test_threaded.jl b/test/test_threaded.jl
index a3d52c1923f..83544aa0072 100644
--- a/test/test_threaded.jl
+++ b/test/test_threaded.jl
@@ -442,7 +442,7 @@ end
         end
     end
 
-    @trixi_testset "elixir_euler_fdsbp_periodic.jl" begin
+    @trixi_testset "elixir_euler_fdsbp_periodic.jl (2D)" begin
         @test_trixi_include(joinpath(examples_dir(), "dgmulti_2d",
                                      "elixir_euler_fdsbp_periodic.jl"),
                             l2=[
@@ -467,6 +467,33 @@ end
             @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 5000
         end
     end
+
+    @trixi_testset "elixir_euler_fdsbp_periodic.jl (3D)" begin
+        @test_trixi_include(joinpath(examples_dir(),
+                                     "dgmulti_3d/elixir_euler_fdsbp_periodic.jl"),
+                            l2=[
+                                7.561896970325353e-5,
+                                6.884047859361093e-5,
+                                6.884047859363204e-5,
+                                6.884047859361148e-5,
+                                0.000201107274617457
+                            ],
+                            linf=[
+                                0.0001337520020225913,
+                                0.00011571467799287305,
+                                0.0001157146779990903,
+                                0.00011571467799376123,
+                                0.0003446082308800058
+                            ])
+        # Ensure that we do not have excessive memory allocations
+        # (e.g., from type instabilities)
+        let
+            t = sol.t[end]
+            u_ode = sol.u[end]
+            du_ode = similar(u_ode)
+            @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000
+        end
+    end
 end
 end