Skip to content

Commit

Permalink
Subcell limiting: Revise order of bounds using a Dict (#1649)
Browse files Browse the repository at this point in the history 
* Implement order of bounds as Dict

* Implement `variable_bounds` directly as Dictionary

* Remove unnecessary line of code

* Use Symbols instead of Strings
  • Loading branch information
@bennibolm
bennibolm authored Oct 12, 2023
1 parent 1727b07 commit 7fd4503
Show file tree
Hide file tree
Showing 3 changed files with 32 additions and 35 deletions.
31 changes: 16 additions & 15 deletions src/solvers/dgsem_tree/containers_2d.jl
View file
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
View file
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
View file
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

0 comments on commit 7fd4503

Please sign in to comment.