exchanged old IC of KHI to new KHI setup (#730)

* exchanged old IC of KHI to new KHI setup * Update src/equations/compressible_euler_2d.jl Co-authored-by: Hendrik Ranocha <ranocha@users.noreply.github.com> * Update src/equations/compressible_euler_2d.jl Co-authored-by: Hendrik Ranocha <ranocha@users.noreply.github.com> * some small changes Co-authored-by: Hendrik Ranocha <ranocha@users.noreply.github.com>
trixi-framework · Jul 22, 2021 · 39c01b8 · 39c01b8
1 parent 50ba5ff
commit 39c01b8
Show file tree

Hide file tree

Showing 5 changed files with 27 additions and 64 deletions.
diff --git a/examples/tree_2d_dgsem/elixir_euler_kelvin_helmholtz_instability.jl b/examples/tree_2d_dgsem/elixir_euler_kelvin_helmholtz_instability.jl
@@ -1,5 +1,4 @@
 
-using Random: seed!
 using OrdinaryDiffEq
 using Trixi
 
@@ -8,7 +7,6 @@ using Trixi
 gamma = 1.4
 equations = CompressibleEulerEquations2D(gamma)
 
-seed!(0)
 initial_condition = initial_condition_khi
 
 surface_flux = flux_lax_friedrichs
@@ -25,8 +23,8 @@ volume_integral = VolumeIntegralShockCapturingHG(indicator_sc;
                                                  volume_flux_fv=surface_flux)
 solver = DGSEM(basis, surface_flux, volume_integral)
 
-coordinates_min = (-0.5, -0.5)
-coordinates_max = ( 0.5,  0.5)
+coordinates_min = (-1.0, -1.0)
+coordinates_max = ( 1.0,  1.0)
 mesh = TreeMesh(coordinates_min, coordinates_max,
                 initial_refinement_level=5,
                 n_cells_max=100_000)
@@ -35,7 +33,7 @@ semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver)
 ###############################################################################
 # ODE solvers, callbacks etc.
 
-tspan = (0.0, 5.0)
+tspan = (0.0, 3.0)
 ode = semidiscretize(semi, tspan)
 
 summary_callback = SummaryCallback()

diff --git a/examples/tree_2d_dgsem/elixir_euler_kelvin_helmholtz_instability_amr.jl b/examples/tree_2d_dgsem/elixir_euler_kelvin_helmholtz_instability_amr.jl
@@ -1,5 +1,4 @@
 
-using Random: seed!
 using OrdinaryDiffEq
 using Trixi
 
@@ -8,7 +7,6 @@ using Trixi
 gamma = 1.4
 equations = CompressibleEulerEquations2D(gamma)
 
-seed!(0)
 initial_condition = initial_condition_khi
 
 surface_flux = flux_lax_friedrichs
@@ -25,8 +23,8 @@ volume_integral = VolumeIntegralShockCapturingHG(indicator_sc;
                                                  volume_flux_fv=surface_flux)
 solver = DGSEM(basis, surface_flux, volume_integral)
 
-coordinates_min = (-0.5, -0.5)
-coordinates_max = ( 0.5,  0.5)
+coordinates_min = (-1.0, -1.0)
+coordinates_max = ( 1.0,  1.0)
 mesh = TreeMesh(coordinates_min, coordinates_max,
                 initial_refinement_level=5,
                 n_cells_max=100_000)
@@ -38,7 +36,7 @@ semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver)
 ###############################################################################
 # ODE solvers, callbacks etc.
 
-tspan = (0.0, 5.0)
+tspan = (0.0, 3.0)
 ode = semidiscretize(semi, tspan)
 
 summary_callback = SummaryCallback()

diff --git a/src/equations/compressible_euler_2d.jl b/src/equations/compressible_euler_2d.jl
@@ -348,36 +348,26 @@ end
 """
     initial_condition_khi(x, t, equations::CompressibleEulerEquations2D)
 
-The classical Kelvin-Helmholtz instability based on
-- https://rsaa.anu.edu.au/research/established-projects/fyris/2-d-kelvin-helmholtz-test
+A version of the classical Kelvin-Helmholtz instability based on
+- Andrés M. Rueda-Ramírez, Gregor J. Gassner (2021)
+  A Subcell Finite Volume Positivity-Preserving Limiter for DGSEM Discretizations 
+  of the Euler Equations
+  [arXiv: 2102.06017](https://arxiv.org/abs/2102.06017)
 """
 function initial_condition_khi(x, t, equations::CompressibleEulerEquations2D)
-  # https://rsaa.anu.edu.au/research/established-projects/fyris/2-d-kelvin-helmholtz-test
   # change discontinuity to tanh
   # typical resolution 128^2, 256^2
-  # domain size is [-0.5,0.5]^2
-  dens0 = 1.0 # outside density
-  dens1 = 2.0 # inside density
-  velx0 = -0.5 # outside velocity
-  velx1 = 0.5 # inside velocity
-  slope = 50 # used for tanh instead of discontinuous initial condition
-  # pressure equilibrium
-  p     = 2.5
-  # density
-  rho = dens0 + (dens1-dens0) * 0.5*(1+(tanh(slope*(x[2]+0.25)) - (tanh(slope*(x[2]-0.25)) + 1)))
-  if iszero(t) # initial condition
-    #  y velocity v2 is only white noise
-    v2  = 0.01*(rand(Float64,1)[1]-0.5)
-    #  x velocity is also augmented with noise
-    v1 = velx0 + (velx1-velx0) * 0.5*(1+(tanh(slope*(x[2]+0.25)) - (tanh(slope*(x[2]-0.25)) + 1)))+0.01*(rand(Float64,1)[1]-0.5)
-  else # background values to compute reference values for CI
-    v2 = 0.0
-    v1 = velx0 + (velx1-velx0) * 0.5*(1+(tanh(slope*(x[2]+0.25)) - (tanh(slope*(x[2]-0.25)) + 1)))
-  end
+  # domain size is [-1,+1]^2
+  slope = 15
+  amplitude = 0.02
+  B = tanh(slope * x[2] + 7.5) - tanh(slope * x[2] - 7.5)
+  rho = 0.5 + 0.75 * B
+  v1 = 0.5 * (B - 1)
+  v2 = 0.1 * sin(2 * pi * x[1])
+  p = 1.0
   return prim2cons(SVector(rho, v1, v2, p), equations)
 end
 
-
 """
     initial_condition_blob(x, t, equations::CompressibleEulerEquations2D)
 

diff --git a/test/Project.toml b/test/Project.toml
@@ -6,7 +6,6 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 MPI = "da04e1cc-30fd-572f-bb4f-1f8673147195"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
-Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [compat]

diff --git a/test/test_examples_2d_euler.jl b/test/test_examples_2d_euler.jl
@@ -96,39 +96,17 @@ EXAMPLES_DIR = joinpath(pathof(Trixi) |> dirname |> dirname, "examples", "tree_2
   end
 
   @trixi_testset "elixir_euler_kelvin_helmholtz_instability.jl" begin
-    if Threads.nthreads() == 1
-      # This example uses random numbers to generate the initial condition.
-      # Hence, we can only check "errors" if everything is made reproducible.
-      # However, that's not enough to ensure reproducibility since the stream
-      # of random numbers is not guaranteed to be the same across different
-      # minor versions of Julia.
-      # See https://github.com/trixi-framework/Trixi.jl/issues/232#issuecomment-709738400
-      @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_kelvin_helmholtz_instability.jl"),
-        l2   = [0.002046113073936985, 0.002862623943300569, 0.001971116879236713, 0.004816623657677065],
-        linf = [0.024375050653856478, 0.01803061241763637, 0.009938942915093363, 0.02097211774984231],
-        tspan = (0.0, 0.2))
-    else
-      @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_kelvin_helmholtz_instability.jl"),
-        tspan = (0.0, 0.2))
-    end
+    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_kelvin_helmholtz_instability.jl"),
+      l2   = [5.56898597e-02,   3.29845866e-02,   5.22436730e-02,   8.00923511e-02],
+      linf = [2.40499700e-01,   1.66109782e-01,   1.23559478e-01,   2.69558145e-01],
+      tspan = (0.0, 0.2))
   end
 
   @trixi_testset "elixir_euler_kelvin_helmholtz_instability_amr.jl" begin
-    if Threads.nthreads() == 1
-      # This example uses random numbers to generate the initial condition.
-      # Hence, we can only check "errors" if everything is made reproducible.
-      # However, that's not enough to ensure reproducibility since the stream
-      # of random numbers is not guaranteed to be the same across different
-      # minor versions of Julia.
-      # See https://github.com/trixi-framework/Trixi.jl/issues/232#issuecomment-709738400
-      @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_kelvin_helmholtz_instability_amr.jl"),
-        l2   = [0.0015773943210539198, 0.0023012411276933705, 0.0013793453564426015, 0.003239622733845669],
-        linf = [0.025791632274568155, 0.015170620391733003, 0.00958159135384009, 0.013265593557157729],
-        tspan = (0.0, 0.2))
-    else
-      @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_kelvin_helmholtz_instability_amr.jl"),
-        tspan = (0.0, 0.2))
-    end
+    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_kelvin_helmholtz_instability_amr.jl"),
+      l2   = [5.56928413e-02,   3.31135409e-02,   5.22350998e-02,   8.00669862e-02],
+      linf = [2.53988861e-01,   1.74418201e-01,   1.23234549e-01,   2.69116662e-01],
+      tspan = (0.0, 0.2))
   end
 
   @trixi_testset "elixir_euler_vortex.jl" begin