Fix flux calculation for MHD in 2D (thanks to @andrewwinters5000!)

trixi-framework · Oct 28, 2020 · c6f4cf4 · c6f4cf4
1 parent f25c1a1
commit c6f4cf4
Showing 1 changed file with 4 additions and 3 deletions.
diff --git a/src/equations/2d/ideal_glm_mhd.jl b/src/equations/2d/ideal_glm_mhd.jl
@@ -158,14 +158,15 @@ end
   v1 = rho_v1/rho
   v2 = rho_v2/rho
   v3 = rho_v3/rho
+  kin_en = 0.5*rho*(v1^2 + v2^2 + v3^2)
   mag_en = 0.5*(B1^2 + B2^2 + B3^2)
-  p = (equation.gamma - 1) * (rho_e - 0.5*rho*(v1^2 + v2^2 + v3^2) - mag_en - 0.5*psi^2)
+  p = (equation.gamma - 1) * (rho_e - kin_en - mag_en - 0.5*psi^2)
   if orientation == 1
     f1 = rho_v1
     f2 = rho_v1*v1 + p + mag_en - B1^2
     f3 = rho_v1*v2 - B1*B2
     f4 = rho_v1*v3 - B1*B3
-    f5 = (rho_e + p + mag_en)*v1 - B1*(v1*B1 + v2*B2 + v3*B3) + equation.c_h*psi*B1
+    f5 = (kin_en + equation.gamma*p/(equation.gamma - 1) + 2.0*mag_en)*v1 - B1*(v1*B1 + v2*B2 + v3*B3) + equation.c_h*psi*B1
     f6 = equation.c_h*psi
     f7 = v1*B2 - v2*B1
     f8 = v1*B3 - v3*B1
@@ -175,7 +176,7 @@ end
     f2 = rho_v2*v1 - B1*B2
     f3 = rho_v2*v2 + p + mag_en - B2^2
     f4 = rho_v2*v3 - B2*B3
-    f5 = (rho_e + p + mag_en)*v2 - B2*(v1*B1 + v2*B2 + v3*B3) + equation.c_h*psi*B2
+    f5 = (kin_en + equation.gamma*p/(equation.gamma - 1) + 2.0*mag_en)*v2 - B2*(v1*B1 + v2*B2 + v3*B3) + equation.c_h*psi*B2
     f6 = v2*B1 - v1*B2
     f7 = equation.c_h*psi
     f8 = v2*B3 - v3*B2