diff --git a/examples/p4est_2d_dgsem/elixir_euler_sedov_blast_wave_sc_subcell.jl b/examples/p4est_2d_dgsem/elixir_euler_sedov_blast_wave_sc_subcell.jl index 6fcb3b0fa92..f346f6ad59f 100644 --- a/examples/p4est_2d_dgsem/elixir_euler_sedov_blast_wave_sc_subcell.jl +++ b/examples/p4est_2d_dgsem/elixir_euler_sedov_blast_wave_sc_subcell.jl @@ -86,10 +86,14 @@ save_solution = SaveSolutionCallback(interval = 300, stepsize_callback = StepsizeCallback(cfl = 0.5) +limiting_analysis_callback = LimitingAnalysisCallback(output_directory = "out", + interval = 1) + callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, save_solution, + limiting_analysis_callback, stepsize_callback) ############################################################################### diff --git a/examples/structured_2d_dgsem/elixir_euler_double_mach.jl b/examples/structured_2d_dgsem/elixir_euler_double_mach.jl index 33be05c615c..12bc936b1cc 100644 --- a/examples/structured_2d_dgsem/elixir_euler_double_mach.jl +++ b/examples/structured_2d_dgsem/elixir_euler_double_mach.jl @@ -153,9 +153,13 @@ save_solution = SaveSolutionCallback(interval = 1000, stepsize_callback = StepsizeCallback(cfl = 0.9) +limiting_analysis_callback = LimitingAnalysisCallback(output_directory = "out", + interval = 1) + callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, stepsize_callback, + limiting_analysis_callback, save_solution) ############################################################################### diff --git a/examples/structured_2d_dgsem/elixir_euler_double_mach_MCL.jl b/examples/structured_2d_dgsem/elixir_euler_double_mach_MCL.jl index 68290ac31ff..9c2a03a26f0 100644 --- a/examples/structured_2d_dgsem/elixir_euler_double_mach_MCL.jl +++ b/examples/structured_2d_dgsem/elixir_euler_double_mach_MCL.jl @@ -155,9 +155,13 @@ save_solution = SaveSolutionCallback(interval = 1000, stepsize_callback = StepsizeCallback(cfl = 0.9) +limiting_analysis_callback = LimitingAnalysisCallback(output_directory = "out", + interval = 1) + callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, stepsize_callback, + limiting_analysis_callback, save_solution) ############################################################################### diff --git a/examples/tree_2d_dgsem/elixir_euler_blast_wave_MCL.jl b/examples/tree_2d_dgsem/elixir_euler_blast_wave_MCL.jl index ba5f27d9f9f..d05cd8d1172 100644 --- a/examples/tree_2d_dgsem/elixir_euler_blast_wave_MCL.jl +++ b/examples/tree_2d_dgsem/elixir_euler_blast_wave_MCL.jl @@ -82,9 +82,13 @@ save_solution = SaveSolutionCallback(interval = 500, stepsize_callback = StepsizeCallback(cfl = 0.9) +limiting_analysis_callback = LimitingAnalysisCallback(output_directory = "out", + interval = 1) + callbacks = CallbackSet(summary_callback, analysis_callback, alive_callback, save_solution, + limiting_analysis_callback, stepsize_callback) ############################################################################### diff --git a/src/callbacks_step/limiting_analysis.jl b/src/callbacks_step/limiting_analysis.jl index 5bbf8c8b801..5dce041b549 100644 --- a/src/callbacks_step/limiting_analysis.jl +++ b/src/callbacks_step/limiting_analysis.jl @@ -145,7 +145,7 @@ end @unpack output_directory = limiting_analysis_callback @unpack alpha = limiter.cache.subcell_limiter_coefficients - alpha_avg = analyze_coefficient_IDP(mesh, equations, dg, cache, limiter) + alpha_avg = analyze_coefficient(mesh, equations, dg, cache, limiter) open("$output_directory/alphas.txt", "a") do f println(f, iter, ", ", time, ", ", maximum(alpha), ", ", alpha_avg) @@ -156,6 +156,8 @@ end dg, cache, limiter::SubcellLimiterMCL, time, iter) + @assert limiter.Plotting "Parameter `Plotting` needs to be activated for analysis of limiting factor with `LimitingAnalysisCallback`" + @unpack output_directory = limiting_analysis_callback @unpack weights = dg.basis @unpack alpha, alpha_pressure, alpha_entropy, @@ -163,8 +165,8 @@ end n_vars = nvariables(equations) - alpha_min_avg, alpha_mean_avg = analyze_coefficient_MCL(mesh, equations, dg, cache, - limiter) + alpha_min_avg, alpha_mean_avg = analyze_coefficient(mesh, equations, dg, cache, + limiter) open("$output_directory/alphas_min.txt", "a") do f print(f, iter, ", ", time) diff --git a/src/callbacks_step/limiting_analysis_2d.jl b/src/callbacks_step/limiting_analysis_2d.jl index ef547c01967..2b77257228d 100644 --- a/src/callbacks_step/limiting_analysis_2d.jl +++ b/src/callbacks_step/limiting_analysis_2d.jl @@ -5,7 +5,8 @@ @muladd begin #! format: noindent -function analyze_coefficient_IDP(mesh::TreeMesh2D, equations, dg, cache, limiter) +function analyze_coefficient(mesh::TreeMesh2D, equations, dg, cache, + limiter::SubcellLimiterIDP) @unpack weights = dg.basis @unpack alpha = limiter.cache.subcell_limiter_coefficients @@ -22,7 +23,9 @@ function analyze_coefficient_IDP(mesh::TreeMesh2D, equations, dg, cache, limiter return alpha_avg / total_volume end -function analyze_coefficient_IDP(mesh::StructuredMesh{2}, equations, dg, cache, limiter) +function analyze_coefficient(mesh::Union{StructuredMesh{2}, P4estMesh{2}}, + equations, dg, cache, + limiter::SubcellLimiterIDP) @unpack weights = dg.basis @unpack alpha = limiter.cache.subcell_limiter_coefficients @@ -39,7 +42,8 @@ function analyze_coefficient_IDP(mesh::StructuredMesh{2}, equations, dg, cache, return alpha_avg / total_volume end -function analyze_coefficient_MCL(mesh::TreeMesh2D, equations, dg, cache, limiter) +function analyze_coefficient(mesh::TreeMesh2D, equations, dg, cache, + limiter::SubcellLimiterMCL) @unpack weights = dg.basis @unpack alpha, alpha_mean, alpha_pressure, alpha_mean_pressure, alpha_entropy, alpha_mean_entropy = limiter.cache.subcell_limiter_coefficients @@ -83,8 +87,9 @@ function analyze_coefficient_MCL(mesh::TreeMesh2D, equations, dg, cache, limiter return alpha_avg ./ total_volume, alpha_mean_avg ./ total_volume end -function analyze_coefficient_MCL(mesh::StructuredMesh{2}, equations, dg, cache, - limiter) +function analyze_coefficient(mesh::Union{StructuredMesh{2}, P4estMesh{2}}, + equations, dg, cache, + limiter::SubcellLimiterMCL) @unpack weights = dg.basis @unpack alpha, alpha_mean, alpha_pressure, alpha_mean_pressure, alpha_entropy, alpha_mean_entropy = limiter.cache.subcell_limiter_coefficients diff --git a/test/test_p4est_2d.jl b/test/test_p4est_2d.jl index a8c20ee3bd1..36d1a3d6672 100644 --- a/test/test_p4est_2d.jl +++ b/test/test_p4est_2d.jl @@ -335,6 +335,7 @@ end end @trixi_testset "elixir_euler_sedov_blast_wave_sc_subcell.jl" begin + rm(joinpath("out", "alphas.txt"), force = true) @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_sedov_blast_wave_sc_subcell.jl"), l2=[ @@ -350,6 +351,20 @@ end 6.268843623142863 ], tspan=(0.0, 0.3)) + lines = readlines(joinpath("out", "alphas.txt")) + @test lines[1] == + "# iter, simu_time, alpha_max, alpha_avg" + cmd = string(Base.julia_cmd()) + coverage = occursin("--code-coverage", cmd) && + !occursin("--code-coverage=none", cmd) + if coverage + # Run with coverage takes 6 time steps. + @test startswith(lines[end], "6, 0.005") + @test occursin(r"1.0, 0.968", lines[end]) + else + # Run without coverage takes 128 time steps. + @test startswith(lines[end], "128, 0.1, 1.0, 0.902") + end # Ensure that we do not have excessive memory allocations # (e.g., from type instabilities) let diff --git a/test/test_structured_2d.jl b/test/test_structured_2d.jl index 8254712cf3f..7168ebe20d1 100644 --- a/test/test_structured_2d.jl +++ b/test/test_structured_2d.jl @@ -719,6 +719,7 @@ end end @trixi_testset "elixir_euler_double_mach.jl" begin + rm(joinpath("out", "alphas.txt"), force = true) @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_double_mach.jl"), l2=[ 0.87417841433288, @@ -734,6 +735,20 @@ end ], initial_refinement_level=2, tspan=(0.0, 0.05)) + lines = readlines(joinpath("out", "alphas.txt")) + @test lines[1] == + "# iter, simu_time, alpha_max, alpha_avg" + cmd = string(Base.julia_cmd()) + coverage = occursin("--code-coverage", cmd) && + !occursin("--code-coverage=none", cmd) + if coverage + # Run with coverage takes 1 time steps. + @test startswith(lines[end], "1, 0.0002") + @test occursin(r"1.0, 0.9809", lines[end]) + else + # Run without coverage takes 193 time steps. + @test startswith(lines[end], "193, 0.05, 1.0, 0.3160") + end # Ensure that we do not have excessive memory allocations # (e.g., from type instabilities) let @@ -749,6 +764,8 @@ end end @trixi_testset "elixir_euler_double_mach_MCL.jl" begin + rm(joinpath("out", "alphas_mean.txt"), force = true) + rm(joinpath("out", "alphas_min.txt"), force = true) @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_double_mach_MCL.jl"), l2=[ 0.8887316108902574, @@ -764,6 +781,32 @@ end ], initial_refinement_level=2, tspan=(0.0, 0.05)) + lines = readlines(joinpath("out", "alphas_mean.txt")) + @test lines[1] == + "# iter, simu_time, alpha_min_rho, alpha_avg_rho, alpha_min_rho_v1, alpha_avg_rho_v1, alpha_min_rho_v2, alpha_avg_rho_v2, alpha_min_rho_e, alpha_avg_rho_e" + cmd = string(Base.julia_cmd()) + coverage = occursin("--code-coverage", cmd) && + !occursin("--code-coverage=none", cmd) + if coverage + # Run with coverage takes 1 time steps. + @test startswith(lines[end], "1, 0.0002") + else + # Run without coverage takes 191 time steps. + @test startswith(lines[end], "191, 0.05, 3.7017") + end + lines = readlines(joinpath("out", "alphas_min.txt")) + @test lines[1] == + "# iter, simu_time, alpha_min_rho, alpha_avg_rho, alpha_min_rho_v1, alpha_avg_rho_v1, alpha_min_rho_v2, alpha_avg_rho_v2, alpha_min_rho_e, alpha_avg_rho_e" + cmd = string(Base.julia_cmd()) + coverage = occursin("--code-coverage", cmd) && + !occursin("--code-coverage=none", cmd) + if coverage + # Run with coverage takes 1 time steps. + @test startswith(lines[end], "1, 0.0002") # TODO + else + # Run without coverage takes 191 time steps. + @test startswith(lines[end], "191, 0.05, -0.0, 0.7216") + end # Ensure that we do not have excessive memory allocations # (e.g., from type instabilities) let diff --git a/test/test_tree_2d_euler.jl b/test/test_tree_2d_euler.jl index 65448da0e3e..63328d6fdb9 100644 --- a/test/test_tree_2d_euler.jl +++ b/test/test_tree_2d_euler.jl @@ -532,6 +532,7 @@ end @trixi_testset "elixir_euler_sedov_blast_wave_sc_subcell.jl" begin rm(joinpath("out", "deviations.txt"), force = true) + rm(joinpath("out", "alphas.txt"), force = true) @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_sedov_blast_wave_sc_subcell.jl"), l2=[ @@ -563,6 +564,20 @@ end # Run without coverage takes 381 time steps. @test startswith(lines[end], "381") end + lines = readlines(joinpath("out", "alphas.txt")) + @test lines[1] == + "# iter, simu_time, alpha_max, alpha_avg" + cmd = string(Base.julia_cmd()) + coverage = occursin("--code-coverage", cmd) && + !occursin("--code-coverage=none", cmd) + if coverage + # Run with coverage takes 6 time steps. + @test startswith(lines[end], "6, 0.014") + @test occursin(r"1.0, 0.953", lines[end]) + else + # Run without coverage takes 381 time steps. + @test startswith(lines[end], "381, 1.0, 1.0, 0.544") + end # Ensure that we do not have excessive memory allocations # (e.g., from type instabilities) let @@ -579,6 +594,8 @@ end @trixi_testset "elixir_euler_sedov_blast_wave_MCL.jl" begin rm(joinpath("out", "deviations.txt"), force = true) + rm(joinpath("out", "alphas_mean.txt"), force = true) + rm(joinpath("out", "alphas_min.txt"), force = true) @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_sedov_blast_wave_MCL.jl"), l2=[ 0.4740321851943766, @@ -609,6 +626,32 @@ end # Run without coverage takes 349 time steps. @test startswith(lines[end], "349") end + lines = readlines(joinpath("out", "alphas_mean.txt")) + @test lines[1] == + "# iter, simu_time, alpha_min_rho, alpha_avg_rho, alpha_min_rho_v1, alpha_avg_rho_v1, alpha_min_rho_v2, alpha_avg_rho_v2, alpha_min_rho_e, alpha_avg_rho_e, alpha_min_pressure, alpha_avg_pressure, alpha_min_entropy, alpha_avg_entropy" + cmd = string(Base.julia_cmd()) + coverage = occursin("--code-coverage", cmd) && + !occursin("--code-coverage=none", cmd) + if coverage + # Run with coverage takes 6 time steps. + @test startswith(lines[end], "6, 0.014") + else + # Run without coverage takes 349 time steps. + @test startswith(lines[end], "349, 1.0, 0.0002") + end + lines = readlines(joinpath("out", "alphas_min.txt")) + @test lines[1] == + "# iter, simu_time, alpha_min_rho, alpha_avg_rho, alpha_min_rho_v1, alpha_avg_rho_v1, alpha_min_rho_v2, alpha_avg_rho_v2, alpha_min_rho_e, alpha_avg_rho_e, alpha_min_pressure, alpha_avg_pressure, alpha_min_entropy, alpha_avg_entropy" + cmd = string(Base.julia_cmd()) + coverage = occursin("--code-coverage", cmd) && + !occursin("--code-coverage=none", cmd) + if coverage + # Run with coverage takes 6 time steps. + @test startswith(lines[end], "6, 0.014") + else + # Run without coverage takes 349 time steps. + @test startswith(lines[end], "349, 1.0, -0.0, 0.773") + end # Ensure that we do not have excessive memory allocations # (e.g., from type instabilities) let