diff --git a/docs/src/callbacks.md b/docs/src/callbacks.md index 7f44dfd5925..f018bcf7c39 100644 --- a/docs/src/callbacks.md +++ b/docs/src/callbacks.md @@ -30,6 +30,11 @@ An example elixir using AMR can be found at [`examples/tree_2d_dgsem/elixir_adve The [`AnalysisCallback`](@ref) can be used to analyze the numerical solution, e.g. calculate errors or user-specified integrals, and print the results to the screen. The results can also be saved in a file. An example can be found at [`examples/tree_2d_dgsem/elixir_euler_vortex.jl`](https://github.com/trixi-framework/Trixi.jl/blob/main/examples/tree_2d_dgsem/elixir_euler_vortex.jl). +Note that the errors (e.g. `L2 error` or `Linf error`) are computed with respect to the initial condition. +The percentage of the simulation time refers to the ratio of the current time and the final time, i.e. it does +not consider the maximal number of iterations. So the simulation could finish before 100% are reached. +Note that, e.g., due to AMR or smaller time step sizes, the simulation can actually take longer than +the percentage indicates. In [Performance metrics of the `AnalysisCallback`](@ref performance-metrics) you can find a detailed description of the different performance metrics the `AnalysisCallback` computes. diff --git a/src/callbacks_step/alive.jl b/src/callbacks_step/alive.jl index eeacd9681d8..9df7181521e 100644 --- a/src/callbacks_step/alive.jl +++ b/src/callbacks_step/alive.jl @@ -86,9 +86,15 @@ function (alive_callback::AliveCallback)(integrator) println("─"^100) println() elseif mpi_isroot() + t = integrator.t + t_initial = first(integrator.sol.prob.tspan) + t_final = last(integrator.sol.prob.tspan) + sim_time_percentage = (t - t_initial) / (t_final - t_initial) * 100 runtime_absolute = 1.0e-9 * (time_ns() - alive_callback.start_time) - @printf("#timesteps: %6d │ Δt: %.4e │ sim. time: %.4e │ run time: %.4e s\n", - integrator.stats.naccept, integrator.dt, integrator.t, runtime_absolute) + println(rpad(@sprintf("#timesteps: %6d │ Δt: %.4e │ sim. time: %.4e (%5.3f%%)", + integrator.stats.naccept, integrator.dt, t, + sim_time_percentage), 71) * + @sprintf("│ run time: %.4e s", runtime_absolute)) end # avoid re-evaluating possible FSAL stages diff --git a/src/callbacks_step/analysis.jl b/src/callbacks_step/analysis.jl index a3e49007fb6..e5b4a01a885 100644 --- a/src/callbacks_step/analysis.jl +++ b/src/callbacks_step/analysis.jl @@ -232,6 +232,12 @@ function (analysis_callback::AnalysisCallback)(u_ode, du_ode, integrator, semi) @unpack dt, t = integrator iter = integrator.stats.naccept + # Compute the percentage of the simulation that is done + t = integrator.t + t_initial = first(integrator.sol.prob.tspan) + t_final = last(integrator.sol.prob.tspan) + sim_time_percentage = (t - t_initial) / (t_final - t_initial) * 100 + # Record performance measurements and compute performance index (PID) runtime_since_last_analysis = 1.0e-9 * (time_ns() - analysis_callback.start_time_last_analysis) @@ -291,8 +297,8 @@ function (analysis_callback::AnalysisCallback)(u_ode, du_ode, integrator, semi) " " * " └── GC time: " * @sprintf("%10.8e s (%5.3f%%)", gc_time_absolute, gc_time_percentage)) - mpi_println(" sim. time: " * @sprintf("%10.8e", t) * - " " * + mpi_println(rpad(" sim. time: " * + @sprintf("%10.8e (%5.3f%%)", t, sim_time_percentage), 46) * " time/DOF/rhs!: " * @sprintf("%10.8e s", runtime_relative)) mpi_println(" " * " " * " " * diff --git a/src/solvers/dgsem_tree/containers_2d.jl b/src/solvers/dgsem_tree/containers_2d.jl index 5f5807a181b..9e9fe88c15b 100644 --- a/src/solvers/dgsem_tree/containers_2d.jl +++ b/src/solvers/dgsem_tree/containers_2d.jl @@ -1325,16 +1325,16 @@ mutable struct ContainerSubcellLimiterIDP2D{uEltype <: Real} alpha::Array{uEltype, 3} # [i, j, element] alpha1::Array{uEltype, 3} alpha2::Array{uEltype, 3} - variable_bounds::Vector{Array{uEltype, 3}} + variable_bounds::Dict{Symbol, Array{uEltype, 3}} # internal `resize!`able storage _alpha::Vector{uEltype} _alpha1::Vector{uEltype} _alpha2::Vector{uEltype} - _variable_bounds::Vector{Vector{uEltype}} + _variable_bounds::Dict{Symbol, Vector{uEltype}} end function ContainerSubcellLimiterIDP2D{uEltype}(capacity::Integer, n_nodes, - length) where {uEltype <: Real} + bound_keys) where {uEltype <: Real} nan_uEltype = convert(uEltype, NaN) # Initialize fields with defaults @@ -1345,12 +1345,12 @@ function ContainerSubcellLimiterIDP2D{uEltype}(capacity::Integer, n_nodes, _alpha2 = fill(nan_uEltype, n_nodes * (n_nodes + 1) * capacity) alpha2 = unsafe_wrap(Array, pointer(_alpha2), (n_nodes, n_nodes + 1, capacity)) - _variable_bounds = Vector{Vector{uEltype}}(undef, length) - variable_bounds = Vector{Array{uEltype, 3}}(undef, length) - for i in 1:length - _variable_bounds[i] = fill(nan_uEltype, n_nodes * n_nodes * capacity) - variable_bounds[i] = unsafe_wrap(Array, pointer(_variable_bounds[i]), - (n_nodes, n_nodes, capacity)) + _variable_bounds = Dict{Symbol, Vector{uEltype}}() + variable_bounds = Dict{Symbol, Array{uEltype, 3}}() + for key in bound_keys + _variable_bounds[key] = fill(nan_uEltype, n_nodes * n_nodes * capacity) + variable_bounds[key] = unsafe_wrap(Array, pointer(_variable_bounds[key]), + (n_nodes, n_nodes, capacity)) end return ContainerSubcellLimiterIDP2D{uEltype}(alpha, alpha1, alpha2, @@ -1369,7 +1369,7 @@ nnodes(container::ContainerSubcellLimiterIDP2D) = size(container.alpha, 1) function Base.resize!(container::ContainerSubcellLimiterIDP2D, capacity) n_nodes = nnodes(container) - @unpack _alpha, _alpha1, _alpha2 = container + (; _alpha, _alpha1, _alpha2) = container resize!(_alpha, n_nodes * n_nodes * capacity) container.alpha = unsafe_wrap(Array, pointer(_alpha), (n_nodes, n_nodes, capacity)) container.alpha .= convert(eltype(container.alpha), NaN) @@ -1380,11 +1380,12 @@ function Base.resize!(container::ContainerSubcellLimiterIDP2D, capacity) container.alpha2 = unsafe_wrap(Array, pointer(_alpha2), (n_nodes, n_nodes + 1, capacity)) - @unpack _variable_bounds = container - for i in 1:length(_variable_bounds) - resize!(_variable_bounds[i], n_nodes * n_nodes * capacity) - container.variable_bounds[i] = unsafe_wrap(Array, pointer(_variable_bounds[i]), - (n_nodes, n_nodes, capacity)) + (; _variable_bounds) = container + for (key, _) in _variable_bounds + resize!(_variable_bounds[key], n_nodes * n_nodes * capacity) + container.variable_bounds[key] = unsafe_wrap(Array, + pointer(_variable_bounds[key]), + (n_nodes, n_nodes, capacity)) end return nothing diff --git a/src/solvers/dgsem_tree/subcell_limiters.jl b/src/solvers/dgsem_tree/subcell_limiters.jl index e1f3c71f4ea..6197ca6b85d 100644 --- a/src/solvers/dgsem_tree/subcell_limiters.jl +++ b/src/solvers/dgsem_tree/subcell_limiters.jl @@ -52,9 +52,10 @@ function SubcellLimiterIDP(equations::AbstractEquations, basis; positivity_variables_cons = [], positivity_correction_factor = 0.1) positivity = (length(positivity_variables_cons) > 0) - number_bounds = length(positivity_variables_cons) - cache = create_cache(SubcellLimiterIDP, equations, basis, number_bounds) + bound_keys = Tuple(Symbol("$(i)_min") for i in positivity_variables_cons) + + cache = create_cache(SubcellLimiterIDP, equations, basis, bound_keys) SubcellLimiterIDP{typeof(positivity_correction_factor), typeof(cache)}(positivity, positivity_variables_cons, diff --git a/src/solvers/dgsem_tree/subcell_limiters_2d.jl b/src/solvers/dgsem_tree/subcell_limiters_2d.jl index 748e42966a5..4497217fb56 100644 --- a/src/solvers/dgsem_tree/subcell_limiters_2d.jl +++ b/src/solvers/dgsem_tree/subcell_limiters_2d.jl @@ -7,15 +7,13 @@ # this method is used when the limiter is constructed as for shock-capturing volume integrals function create_cache(limiter::Type{SubcellLimiterIDP}, equations::AbstractEquations{2}, - basis::LobattoLegendreBasis, number_bounds) + basis::LobattoLegendreBasis, bound_keys) subcell_limiter_coefficients = Trixi.ContainerSubcellLimiterIDP2D{real(basis) }(0, nnodes(basis), - number_bounds) + bound_keys) - cache = (; subcell_limiter_coefficients) - - return cache + return (; subcell_limiter_coefficients) end function (limiter::SubcellLimiterIDP)(u::AbstractArray{<:Any, 4}, semi, dg::DGSEM, t, @@ -25,8 +23,7 @@ function (limiter::SubcellLimiterIDP)(u::AbstractArray{<:Any, 4}, semi, dg::DGSE alpha .= zero(eltype(alpha)) if limiter.positivity - @trixi_timeit timer() "positivity" idp_positivity!(alpha, limiter, u, dt, - semi) + @trixi_timeit timer() "positivity" idp_positivity!(alpha, limiter, u, dt, semi) end # Calculate alpha1 and alpha2 @@ -49,22 +46,21 @@ end @inline function idp_positivity!(alpha, limiter, u, dt, semi) # Conservative variables - for (index, variable) in enumerate(limiter.positivity_variables_cons) - idp_positivity!(alpha, limiter, u, dt, semi, variable, index) + for variable in limiter.positivity_variables_cons + idp_positivity!(alpha, limiter, u, dt, semi, variable) end return nothing end -@inline function idp_positivity!(alpha, limiter, u, dt, semi, variable, index) +@inline function idp_positivity!(alpha, limiter, u, dt, semi, variable) mesh, equations, dg, cache = mesh_equations_solver_cache(semi) - @unpack antidiffusive_flux1, antidiffusive_flux2 = cache.antidiffusive_fluxes - @unpack inverse_weights = dg.basis - @unpack positivity_correction_factor = limiter - - @unpack variable_bounds = limiter.cache.subcell_limiter_coefficients + (; antidiffusive_flux1, antidiffusive_flux2) = cache.antidiffusive_fluxes + (; inverse_weights) = dg.basis + (; positivity_correction_factor) = limiter - var_min = variable_bounds[index] + (; variable_bounds) = limiter.cache.subcell_limiter_coefficients + var_min = variable_bounds[Symbol("$(variable)_min")] @threaded for element in eachelement(dg, semi.cache) inverse_jacobian = cache.elements.inverse_jacobian[element]