Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Subcell limiting: Revise order of bounds using a Dict #1649

Merged
31 changes: 16 additions & 15 deletions src/solvers/dgsem_tree/containers_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand All @@ -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,
Expand All @@ -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))
resize!(_alpha1, (n_nodes + 1) * n_nodes * capacity)
Expand All @@ -1379,11 +1379,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
Expand Down
5 changes: 3 additions & 2 deletions src/solvers/dgsem_tree/subcell_limiters.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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,
Expand Down
31 changes: 13 additions & 18 deletions src/solvers/dgsem_tree/subcell_limiters_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,17 +6,14 @@
#! format: noindent

# this method is used when the limiter is constructed as for shock-capturing volume integrals
function create_cache(indicator::Type{SubcellLimiterIDP},
equations::AbstractEquations{2},
basis::LobattoLegendreBasis, number_bounds)
function create_cache(limiter::Type{SubcellLimiterIDP}, equations::AbstractEquations{2},
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,
Expand All @@ -26,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
Expand All @@ -50,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]
Expand Down