@testsnippet TreeMesh3DEuler begin
    EXAMPLES_DIR = joinpath(examples_dir(), "tree_3d_dgsem")
end

@testitem "TreeMesh3D Euler: elixir_euler_source_terms.jl" setup=[Setup, TreeMesh3DEuler] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_source_terms.jl"),
                        l2=[
                            0.010385936842224346,
                            0.009776048833895767,
                            0.00977604883389591,
                            0.009776048833895733,
                            0.01506687097416608
                        ],
                        linf=[
                            0.03285848350791731,
                            0.0321792316408982,
                            0.032179231640894645,
                            0.032179231640895534,
                            0.0655408023333299
                        ])
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
    # Extra test to make sure the "TimeSeriesCallback" made correct data.
    # Extracts data at all points from the first step of the time series and compares it to the
    # exact solution and an interpolated reference solution
    point_data = [getindex(time_series.affect!.point_data[i], 1:5) for i in 1:3]
    exact_data = [initial_condition_convergence_test(time_series.affect!.point_coordinates[:,
                                                                                           i],
                                                     time_series.affect!.time[1],
                                                     equations) for i in 1:3]
    ref_data = [
        [
            1.951156832316166,
            1.952073047561595,
            1.9520730475615966,
            1.9520730475615953,
            3.814390510967551
        ],
        [
            2.0506452262144363,
            2.050727319703708,
            2.0507273197037073,
            2.0507273197037077,
            4.203653999433724
        ],
        [
            2.046982357537558,
            2.0463728824399654,
            2.0463728824399654,
            2.0463728824399645,
            4.190033459318115
        ]]
    @test point_data≈exact_data atol=1e-1
    @test point_data ≈ ref_data
end

@testitem "TreeMesh3D Euler: elixir_euler_convergence_pure_fv.jl" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_convergence_pure_fv.jl"),
                        l2=[
                            0.037182410351406,
                            0.032062252638283974,
                            0.032062252638283974,
                            0.03206225263828395,
                            0.12228177813586687
                        ],
                        linf=[
                            0.0693648413632646,
                            0.0622101894740843,
                            0.06221018947408474,
                            0.062210189474084965,
                            0.24196451799555962
                        ])
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_source_terms.jl with split_form" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_source_terms.jl"),
                        l2=[
                            0.010385936842223388,
                            0.009776048833894784,
                            0.009776048833894784,
                            0.009776048833894765,
                            0.015066870974164096
                        ],
                        linf=[
                            0.03285848350791687,
                            0.032179231640897754,
                            0.0321792316408942,
                            0.0321792316408982,
                            0.06554080233333615
                        ],
                        volume_integral=VolumeIntegralFluxDifferencing(flux_central))
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_convergence.jl" setup=[Setup, TreeMesh3DEuler] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_convergence.jl"),
                        l2=[
                            0.0003637241020254673, 0.00039555708663848046,
                            0.00039555708663832644, 0.0003955570866385083,
                            0.0007811613481643962
                        ],
                        linf=[
                            0.0024000660244567484, 0.002963541002521053,
                            0.0029635410025201647, 0.002963541002522385,
                            0.007191437359379549
                        ])
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_mortar.jl" setup=[Setup, TreeMesh3DEuler] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_mortar.jl"),
                        l2=[
                            0.0019428114665068841,
                            0.0018659907926698422,
                            0.0018659907926698589,
                            0.0018659907926698747,
                            0.0034549095578444056
                        ],
                        linf=[
                            0.011355360771142298,
                            0.011526889155693887,
                            0.011526889155689002,
                            0.011526889155701436,
                            0.02299726519821288
                        ])
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_amr.jl" setup=[Setup, TreeMesh3DEuler] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_amr.jl"),
                        l2=[
                            0.0038281920613404716,
                            0.003828192061340465,
                            0.0038281920613404694,
                            0.0038281920613404672,
                            0.005742288092010652
                        ],
                        linf=[
                            0.07390396464027349,
                            0.07390396464027305,
                            0.07390396464027305,
                            0.07390396464027305,
                            0.11085594696041134
                        ],
                        tspan=(0.0, 0.1))
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_taylor_green_vortex.jl" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_taylor_green_vortex.jl"),
                        l2=[
                            0.00034949871748737876,
                            0.03133384111621587,
                            0.03133384111621582,
                            0.04378599329988925,
                            0.015796137903453026
                        ],
                        linf=[
                            0.0013935237751798724,
                            0.0724080091006194,
                            0.07240800910061806,
                            0.12795921224174792,
                            0.07677156293692633
                        ],
                        tspan=(0.0, 0.5))
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)

    _, energy_spectrum = @inferred compute_kinetic_energy_spectrum(sol)
    @test energy_spectrum[1:6]≈[
        7.770967315541477e-35,
        8.808347042794627e-33,
        0.12307995980341144,
        0.001913513967052981,
        1.424996168880913e-5,
        3.8935136212048387e-7
    ] rtol=1.0e-12
end

@testitem "TreeMesh3D Euler: elixir_euler_shockcapturing.jl" setup=[Setup, TreeMesh3DEuler] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_shockcapturing.jl"),
                        l2=[
                            0.02570137197844877,
                            0.016179934130642552,
                            0.01617993413064253,
                            0.016172648598753545,
                            0.09261669328795467
                        ],
                        linf=[
                            0.3954458125573179,
                            0.26876916180359345,
                            0.26876916180359345,
                            0.26933123042178553,
                            1.3724137121660251
                        ])
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_shockcapturing_amr.jl" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    # Note: This setup does not make much practical sense. It is only added to exercise the
    # `sedov_self_gravity` AMR indicator, which in its original configuration is too expensive for
    # CI testing
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_shockcapturing_amr.jl"),
                        l2=[
                            0.02217299067704248,
                            0.012771561294571411,
                            0.01277156129457143,
                            0.012770635779336643,
                            0.08091898488262424
                        ],
                        linf=[
                            0.4047819603427084,
                            0.27493532130155474,
                            0.2749353213015551,
                            0.2749304638368023,
                            1.4053942765487641
                        ],
                        maxiters=10)
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_density_pulse.jl" setup=[Setup, TreeMesh3DEuler] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_density_pulse.jl"),
                        l2=[
                            0.057196526814004715,
                            0.057196526814004715,
                            0.05719652681400473,
                            0.057196526814004736,
                            0.08579479022100575
                        ],
                        linf=[
                            0.27415246703018203,
                            0.2741524670301829,
                            0.2741524670301827,
                            0.27415246703018226,
                            0.41122870054527816
                        ])
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_ec.jl" setup=[Setup, TreeMesh3DEuler] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_ec.jl"),
                        l2=[
                            0.02526341317987378,
                            0.016632068583699623,
                            0.016632068583699623,
                            0.01662548715216875,
                            0.0913477018048886
                        ],
                        linf=[
                            0.4372549540810414,
                            0.28613118232798984,
                            0.28613118232799006,
                            0.28796686065271876,
                            1.5072828647309124
                        ])
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_ec.jl with initial_condition=initial_condition_constant" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_ec.jl"),
                        l2=[
                            4.183721551616214e-16,
                            6.059779958716338e-16,
                            4.916596221090319e-16,
                            9.739943366304456e-16,
                            3.7485908743251566e-15
                        ],
                        linf=[
                            2.4424906541753444e-15,
                            3.733124920302089e-15,
                            4.440892098500626e-15,
                            5.329070518200751e-15,
                            2.4868995751603507e-14
                        ],
                        initial_condition=initial_condition_constant)
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_ec.jl with flux_chandrashekar" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_ec.jl"),
                        l2=[
                            0.025265721172813106,
                            0.016649800693500427,
                            0.01664980069350042,
                            0.01664379306708522,
                            0.09137248646784184
                        ],
                        linf=[
                            0.4373399329742198,
                            0.28434487167605427,
                            0.28434487167605427,
                            0.28522678968890774,
                            1.532471676033761
                        ],
                        surface_flux=flux_chandrashekar, volume_flux=flux_chandrashekar)
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_ec.jl with flux_kennedy_gruber" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_ec.jl"),
                        l2=[
                            0.025280033869871984,
                            0.016675487948639846,
                            0.016675487948639853,
                            0.016668992714991282,
                            0.091455613470441
                        ],
                        linf=[
                            0.43348628145015766,
                            0.28853549062014217,
                            0.28853549062014217,
                            0.2903943042772536,
                            1.5236557526482426
                        ],
                        surface_flux=flux_kennedy_gruber,
                        volume_flux=flux_kennedy_gruber)
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_ec.jl with flux_shima_etal" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_ec.jl"),
                        l2=[
                            0.025261716925811403,
                            0.016637655557848952,
                            0.01663765555784895,
                            0.01663105921013437,
                            0.09136239054024566
                        ],
                        linf=[
                            0.43692416928732536,
                            0.28622033209064734,
                            0.28622033209064746,
                            0.2881197143457632,
                            1.506534270303663
                        ],
                        surface_flux=flux_shima_etal, volume_flux=flux_shima_etal)
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_laplace_diffusion.jl" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_laplace_diffusion.jl"),
                        l2=[
                            0.013299230512542162,
                            0.0073025819009651,
                            0.007302581900965106,
                            0.007300042097573285,
                            0.04888085245959731
                        ],
                        linf=[
                            0.31714843611640464,
                            0.23586839231625517,
                            0.23586839231625506,
                            0.23698123351744782,
                            1.1174271158464726
                        ])
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_blob_amr.jl" setup=[Setup, TreeMesh3DEuler] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_blob_amr.jl"),
                        l2=[
                            0.04867856452253151,
                            0.2640486962336911,
                            0.0354927658652858,
                            0.03549276586528571,
                            1.0777274757408568
                        ],
                        linf=[
                            9.558543313792217,
                            49.4518309553356,
                            10.319859082570309,
                            10.319859082570487,
                            195.1066220797401
                        ],
                        tspan=(0.0, 0.2))
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_sedov_blast_wave.jl" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_sedov_blast_wave.jl"),
                        l2=[
                            0.0007127163978031706,
                            0.0023166296394624025,
                            0.002316629639462401,
                            0.0023166296394624038,
                            0.010200581509653256
                        ],
                        linf=[
                            0.06344190883105805,
                            0.6292607955969378,
                            0.6292607955969377,
                            0.6292607955969377,
                            2.397746252817731
                        ],
                        maxiters=5, max_level=6,
                        surface_flux=FluxHLL(min_max_speed_naive))
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_sedov_blast_wave.jl (HLLE)" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_sedov_blast_wave.jl"),
                        l2=[
                            0.0007871241159752619,
                            0.0037168004033428146,
                            0.0037168004033428094,
                            0.0037168004033428514,
                            0.011119869089205635
                        ],
                        linf=[
                            0.13982864363612468,
                            0.786004687738243,
                            0.786004687738243,
                            0.7860046877382431,
                            1.7082524045150382
                        ],
                        tspan=(0.0, 0.01),
                        surface_flux=flux_hlle)
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_sedov_scO2.jl" setup=[Setup, TreeMesh3DEuler] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_sedov_scO2.jl"),
                        l2=[
                            0.00700044366656764,
                            0.018962720270366887,
                            0.01896272027036692,
                            0.018962720270366863,
                            0.21537167091542966
                        ],
                        linf=[
                            0.6040987203417612,
                            2.26892456902628,
                            2.2689245690262805,
                            2.2689245690262805,
                            14.089789344124894
                        ],
                        tspan=(0.0, 0.01))
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_sedov_blast_weak_form_sc.jl" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR,
                                 "elixir_euler_sedov_blast_weak_form_sc.jl"),
                        l2=[
                            0.006050851473525136,
                            0.017742478887918647,
                            0.017742478887918682,
                            0.01774247888791881,
                            0.21366241418560217
                        ],
                        linf=[
                            0.5104160120781451,
                            1.9379474897743383,
                            1.9379474897743383,
                            1.9379474897743383,
                            12.67096713137104
                        ],
                        tspan=(0.0, 0.01))
    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 1000)
end

@testitem "TreeMesh3D Euler: elixir_euler_sedov_blast_wave_sc_subcell.jl" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR,
                                 "elixir_euler_sedov_blast_wave_sc_subcell.jl"),
                        l2=[
                            0.24806841083939926,
                            0.07001337223874464,
                            0.07001337223806398,
                            0.0700133722383429,
                            0.3620366037665587
                        ],
                        linf=[
                            0.9384071822566761,
                            0.573009568617271,
                            0.5730095685845291,
                            0.5730095686063774,
                            4.861205850307592
                        ],
                        tspan=(0.0, 0.5))
    limiter = semi.solver.volume_integral.limiter
    deviations = collect(values(limiter.cache.idp_bounds_delta_global))
    @test all(isfinite, deviations)
    @test maximum(deviations) <= 3.0e-13

    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    # Larger values for allowed allocations due to usage of custom
    # integrator which are not *recorded* for the methods from
    # OrdinaryDiffEq.jl
    # Corresponding issue: https://github.com/trixi-framework/Trixi.jl/issues/1877
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 15_000)
end

@testitem "TreeMesh3D Euler: elixir_euler_sedov_blast_wave_sc_subcell.jl (Adaptive Vol Int.)" setup=[
    Setup,
    TreeMesh3DEuler
] tags=[:tree_part4] begin
    @test_trixi_include(joinpath(EXAMPLES_DIR,
                                 "elixir_euler_sedov_blast_wave_sc_subcell.jl"),
                        # Adaptive volume integral selects based on the heuristic (!) a priori `indicator`
                        # if the stabilized volume integral should be employed or if the default one is deemed sufficient.
                        volume_integral=VolumeIntegralAdaptive(indicator = IndicatorHennemannGassner(equations,
                                                                                                     basis,
                                                                                                     alpha_max = 0.1, # irrelevant, only `alpha_min` is used for limiting activation
                                                                                                     alpha_min = 0.001, # governs when subcell limiting is considered
                                                                                                     alpha_smooth = true,
                                                                                                     variable = density_pressure),
                                                               volume_integral_default = VolumeIntegralWeakForm(),
                                                               volume_integral_stabilized = VolumeIntegralSubcellLimiting(limiter_idp;
                                                                                                                          volume_flux_dg = volume_flux,
                                                                                                                          volume_flux_fv = surface_flux)),
                        l2=[
                            0.24806841083830014,
                            0.07001337223848285,
                            0.0700133722385661,
                            0.07001337223834266,
                            0.362036603766589
                        ],
                        linf=[
                            0.9384071822788941,
                            0.5730095686470306,
                            0.5730095679841436,
                            0.5730095679943239,
                            4.861205850307726
                        ],
                        tspan=(0.0, 0.5))
    limiter = semi.solver.volume_integral.volume_integral_stabilized.limiter
    deviations = collect(values(limiter.cache.idp_bounds_delta_global))
    @test all(isfinite, deviations)
    @test maximum(deviations) <= 3.0e-13

    # Ensure that we do not have excessive memory allocations
    # (e.g., from type instabilities)
    # Larger values for allowed allocations due to usage of custom
    # integrator which are not *recorded* for the methods from
    # OrdinaryDiffEq.jl
    # Corresponding issue: https://github.com/trixi-framework/Trixi.jl/issues/1877
    @test_allocations(Trixi.rhs_hyperbolic!, semi, sol, 15_000)
end
