Skip to content

Commit

Permalink
Merge branch 'main' into sc/converters_coupling
Browse files Browse the repository at this point in the history
  • Loading branch information
SimonCan authored Oct 18, 2023
2 parents cdd984a + dfdaef5 commit 4668b26
Show file tree
Hide file tree
Showing 8 changed files with 74 additions and 6 deletions.
11 changes: 9 additions & 2 deletions src/callbacks_step/save_solution.jl
Original file line number Diff line number Diff line change
Expand Up @@ -222,21 +222,28 @@ end
end
end

node_variables = Dict{Symbol, Any}()
@trixi_timeit timer() "get node variables" get_node_variables!(node_variables,
semi)

@trixi_timeit timer() "save solution" save_solution_file(u_ode, t, dt, iter, semi,
solution_callback,
element_variables,
node_variables,
system = system)
end

@inline function save_solution_file(u_ode, t, dt, iter,
semi::AbstractSemidiscretization, solution_callback,
element_variables = Dict{Symbol, Any}();
element_variables = Dict{Symbol, Any}(),
node_variables = Dict{Symbol, Any}();
system = "")
mesh, equations, solver, cache = mesh_equations_solver_cache(semi)
u = wrap_array_native(u_ode, mesh, equations, solver, cache)
save_solution_file(u, t, dt, iter, mesh, equations, solver, cache,
solution_callback,
element_variables; system = system)
element_variables,
node_variables; system = system)
end

# TODO: Taal refactor, move save_mesh_file?
Expand Down
18 changes: 16 additions & 2 deletions src/callbacks_step/save_solution_dg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,9 @@ function save_solution_file(u, time, dt, timestep,
UnstructuredMesh2D, SerialP4estMesh,
SerialT8codeMesh},
equations, dg::DG, cache,
solution_callback, element_variables = Dict{Symbol, Any}();
solution_callback,
element_variables = Dict{Symbol, Any}(),
node_variables = Dict{Symbol, Any}();
system = "")
@unpack output_directory, solution_variables = solution_callback

Expand Down Expand Up @@ -73,6 +75,16 @@ function save_solution_file(u, time, dt, timestep,
var = file["element_variables_$v"]
attributes(var)["name"] = string(key)
end

# Store node variables
for (v, (key, node_variable)) in enumerate(node_variables)
# Add to file
file["node_variables_$v"] = node_variable

# Add variable name as attribute
var = file["node_variables_$v"]
attributes(var)["name"] = string(key)
end
end

return filename
Expand All @@ -81,7 +93,9 @@ end
function save_solution_file(u, time, dt, timestep,
mesh::Union{ParallelTreeMesh, ParallelP4estMesh}, equations,
dg::DG, cache,
solution_callback, element_variables = Dict{Symbol, Any}();
solution_callback,
element_variables = Dict{Symbol, Any}(),
node_variables = Dict{Symbol, Any}();
system = "")
@unpack output_directory, solution_variables = solution_callback

Expand Down
4 changes: 4 additions & 0 deletions src/semidiscretization/semidiscretization.jl
Original file line number Diff line number Diff line change
Expand Up @@ -335,6 +335,10 @@ function get_element_variables!(element_variables, u_ode,
get_element_variables!(element_variables, u, mesh_equations_solver_cache(semi)...)
end

function get_node_variables!(node_variables, semi::AbstractSemidiscretization)
get_node_variables!(node_variables, mesh_equations_solver_cache(semi)...)
end

# To implement AMR and use OrdinaryDiffEq.jl etc., we have to be a bit creative.
# Since the caches of the SciML ecosystem are immutable structs, we cannot simply
# change the underlying arrays therein. Hence, to support changing the number of
Expand Down
21 changes: 21 additions & 0 deletions src/solvers/dg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,11 @@ function get_element_variables!(element_variables, u, mesh, equations,
nothing
end

function get_node_variables!(node_variables, mesh, equations,
volume_integral::AbstractVolumeIntegral, dg, cache)
nothing
end

"""
VolumeIntegralStrongForm()
Expand Down Expand Up @@ -214,6 +219,18 @@ function Base.show(io::IO, mime::MIME"text/plain",
end
end

function get_node_variables!(node_variables, mesh, equations,
volume_integral::VolumeIntegralSubcellLimiting, dg, cache)
# While for the element-wise limiting with `VolumeIntegralShockCapturingHG` the indicator is
# called here to get up-to-date values for IO, this is not easily possible in this case
# because the calculation is very integrated into the method.
# See also https://github.com/trixi-framework/Trixi.jl/pull/1611#discussion_r1334553206.
# Therefore, the coefficients at `t=t^{n-1}` are saved. Thus, the coefficients of the first
# stored solution (initial condition) are not yet defined and were manually set to `NaN`.
get_node_variables!(node_variables, volume_integral.limiter, volume_integral,
equations)
end

# TODO: FD. Should this definition live in a different file because it is
# not strictly a DG method?
"""
Expand Down Expand Up @@ -403,6 +420,10 @@ function get_element_variables!(element_variables, u, mesh, equations, dg::DG, c
dg, cache)
end

function get_node_variables!(node_variables, mesh, equations, dg::DG, cache)
get_node_variables!(node_variables, mesh, equations, dg.volume_integral, dg, cache)
end

const MeshesDGSEM = Union{TreeMesh, StructuredMesh, UnstructuredMesh2D, P4estMesh,
T8codeMesh}

Expand Down
1 change: 1 addition & 0 deletions src/solvers/dgsem_tree/containers_2d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1372,6 +1372,7 @@ function Base.resize!(container::ContainerSubcellLimiterIDP2D, capacity)
(; _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)
resize!(_alpha1, (n_nodes + 1) * n_nodes * capacity)
container.alpha1 = unsafe_wrap(Array, pointer(_alpha1),
(n_nodes + 1, n_nodes, capacity))
Expand Down
7 changes: 7 additions & 0 deletions src/solvers/dgsem_tree/subcell_limiters.jl
Original file line number Diff line number Diff line change
Expand Up @@ -101,4 +101,11 @@ function Base.show(io::IO, ::MIME"text/plain", limiter::SubcellLimiterIDP)
summary_box(io, "SubcellLimiterIDP", setup)
end
end

function get_node_variables!(node_variables, limiter::SubcellLimiterIDP,
::VolumeIntegralSubcellLimiting, equations)
node_variables[:limiting_coefficient] = limiter.cache.subcell_limiter_coefficients.alpha

return nothing
end
end # @muladd
5 changes: 5 additions & 0 deletions src/time_integration/methods_SSP.jl
Original file line number Diff line number Diff line change
Expand Up @@ -31,6 +31,11 @@ struct SimpleSSPRK33{StageCallbacks} <: SimpleAlgorithmSSP
stage_callbacks::StageCallbacks

function SimpleSSPRK33(; stage_callbacks = ())
# Mathematically speaking, it is not necessary for the algorithm to split the factors
# into numerator and denominator. Otherwise, however, rounding errors of the order of
# the machine accuracy will occur, which will add up over time and thus endanger the
# conservation of the simulation.
# See also https://github.com/trixi-framework/Trixi.jl/pull/1640.
numerator_a = SVector(0.0, 3.0, 1.0) # a = numerator_a / denominator
numerator_b = SVector(1.0, 1.0, 2.0) # b = numerator_b / denominator
denominator = SVector(1.0, 4.0, 3.0)
Expand Down
13 changes: 11 additions & 2 deletions utils/trixi2txt.jl
Original file line number Diff line number Diff line change
Expand Up @@ -70,7 +70,7 @@ function trixi2txt(filename::AbstractString...;
center_level_0, length_level_0, leaf_cells, coordinates, levels = read_meshfile(meshfile)

# Read data
labels, data, n_elements, n_nodes, element_variables, time = read_datafile(filename)
labels, data, n_elements, n_nodes, element_variables, node_variables, time = read_datafile(filename)

# Check if dimensions match
if length(leaf_cells) != n_elements
Expand Down Expand Up @@ -263,7 +263,16 @@ function read_datafile(filename::String)
index += 1
end

return labels, data, n_elements, n_nodes, element_variables, time
# Extract node variable arrays
node_variables = Dict{String, Union{Vector{Float64}, Vector{Int}}}()
index = 1
while haskey(file, "node_variables_$index")
varname = read(attributes(file["node_variables_$index"])["name"])
node_variables[varname] = read(file["node_variables_$index"])
index += 1
end

return labels, data, n_elements, n_nodes, element_variables, node_variables, time
end
end

Expand Down

0 comments on commit 4668b26

Please sign in to comment.