ouch #2

Author: benegee

examples/moist_euler/moist_bubble/elixir_rainy_euler_moist_bubble.jl

@@ -3,8 +3,8 @@ using Trixi
 using TrixiAtmo
 using TrixiAtmo: source_terms_no_phase_change, saturation_residual,
                  saturation_residual_jacobian, NonlinearSolveDG,
-                 cons2eq_pot_temp, flux_LMARS, flux_chandrashekar,
-                 flux_ec_rain, boundary_condition_simple_slip_wall
+                 cons2eq_pot_temp, flux_chandrashekar,
+                 boundary_condition_simple_slip_wall
 using NLsolve: nlsolve
 
 # Moist bubble test case from paper:

examples/moist_euler/rainy_bubble/elixir_rainy_euler_rainy_bubble.jl

@@ -4,7 +4,7 @@ using TrixiAtmo
 using TrixiAtmo: source_terms_rainy, initial_condition_bubble_rainy,
                  saturation_residual, saturation_residual_jacobian,
                  cons2eq_pot_temp, saturation_vapour_pressure,
-                 flux_ec_rain, boundary_condition_simple_slip_wall,
+                 boundary_condition_simple_slip_wall,
                  generate_hydrostatic_residual, generate_perturbation_residual
 using NLsolve: nlsolve
 

examples/moist_euler/rainy_bubble/elixir_rainy_euler_rainy_bubble_diffusion.jl

@@ -4,11 +4,9 @@ using TrixiAtmo
 using TrixiAtmo: source_terms_rainy, initial_condition_bubble_rainy,
                  saturation_residual, saturation_residual_jacobian,
                  cons2eq_pot_temp, saturation_vapour_pressure,
-                 flux_chandrashekar, flux_LMARS,
                  source_terms_no_phase_change,
                  boundary_condition_laplace,
-                 boundary_condition_simple_slip_wall,
-                 flux_ec_rain
+                 boundary_condition_simple_slip_wall
 using NLsolve: nlsolve
 
 # domain
@@ -55,7 +53,7 @@ mesh = TreeMesh(coordinates_min, coordinates_max,
 
 semi = SemidiscretizationHyperbolicParabolic(mesh, (equations, equations_parabolic),
                                              initial_condition_rainy, solver;
-                                             source_terms_rainy,
+                                             source_terms = source_terms_rainy,
                                              boundary_conditions = (boundary_conditions,
                                                                     boundary_conditions_parabolic))
 

src/TrixiAtmo.jl

@@ -64,7 +64,8 @@ export flux_nonconservative_zeros, flux_nonconservative_ec,
        source_terms_coriolis, source_terms_coriolis_lagrange_multiplier,
        flux_tec, flux_etec, flux_nonconservative_souza_etal,
        flux_nonconservative_artiano_etal,
-       flux_nonconservative_waruzewski_etal, flux_zero
+       flux_nonconservative_waruzewski_etal, flux_zero,
+       flux_ec_rain
 
 export source_terms_lagrange_multiplier, clean_solution_lagrange_multiplier!
 

src/equations/reference_data.jl

@@ -2,7 +2,7 @@
 #! format: noindent
 
 # Physical constants in SI units (reference values from the Williamson et al. test suite)
-const EARTH_RADIUS = 6.371229 # 6.37122e6  # m
+const EARTH_RADIUS = 6.371229e6 # 6.37122e6  # m
 const EARTH_GRAVITATIONAL_ACCELERATION = 9.81 # 9.80616  # m/s²
 const EARTH_ROTATION_RATE = 7.29212e-5 # 7.292e-5  # rad/s
 const SECONDS_PER_DAY = 8.64e4

test/test_unit.jl

@@ -523,4 +523,23 @@ end
           flux(u, aux_vars, orientation, equations)
 end
 
+@testset "Consistency check for EC flux with Rainy Euler" begin
+    # Set up equations and dummy conservative variables state
+    equations = CompressibleRainyEulerEquations2D()
+    # Example state vector (ρ_d, ρ_m, ρ_r, ρu, ρv, ρe, ρq_v, ρq_c, T)
+    u = SVector(1.0, 0.2, 0.1, 0.5, -0.4, 2.2, 0.1, 0.1, 300)
+
+    normal_directions = [SVector(1.0, 0.0),
+        SVector(0.0, 1.0),
+        SVector(0.5, -0.5),
+        SVector(-1.2, 0.3)]
+
+    for normal_direction in normal_directions
+        equal = flux_ec_rain(u, u, normal_direction, equations) .≈
+                flux(u, normal_direction, equations)
+        # TODO
+        expected = [true, true, true, true, true, false, true, true, true]
+        @test equal == expected
+    end
+end
 end