Skip to content

Commit 3b038fb

Browse files
DanielDoehringDanielDoehring
authored
Merge branch 'main' into 3d-subcell-limiting-first
2 parents 5f271e7 + 970db9f commit 3b038fb

File tree

14 files changed

+258
-200
lines changed

14 files changed

+258
-200
lines changed

‎examples/p4est_2d_dgsem/elixir_navierstokes_blast_reflective.jl‎

Lines changed: 1 addition & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -94,10 +94,6 @@ analysis_callback = AnalysisCallback(semi, interval = analysis_interval)
9494

9595
alive_callback = AliveCallback(analysis_interval = analysis_interval)
9696

97-
callbacks = CallbackSet(summary_callback,
98-
analysis_callback,
99-
alive_callback)
100-
10197
callbacks = CallbackSet(summary_callback,
10298
analysis_callback,
10399
alive_callback)
@@ -108,7 +104,5 @@ callbacks = CallbackSet(summary_callback,
108104
# https://docs.sciml.ai/DiffEqDocs/stable/solvers/ode_solve/#Low-Storage-Methods
109105
ode_alg = CKLLSRK65_4M_4R()
110106

111-
abstol = 1e-6
112-
reltol = 1e-4
113-
sol = solve(ode, ode_alg; abstol = abstol, reltol = reltol,
107+
sol = solve(ode, ode_alg; abstol = 1e-6, reltol = 1e-4,
114108
ode_default_options()..., callback = callbacks);

‎examples/tree_2d_dgsem/elixir_diffusion_steady_state_linear_map.jl‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -71,7 +71,7 @@ using Krylov
7171
# symmetric and positive definite with respect to the inner product induced by
7272
# the DGSEM quadrature, not the standard Euclidean inner product. Standard CG
7373
# would require multiplying the result of `A_map * u` (and `b`) by the mass matrix.
74-
u_ls, stats = gmres(A_map, b)
74+
u_ls, stats = gmres(A_map, b, atol = 1.0e-11, rtol = 1.0e-10)
7575

7676
###############################################################################
7777

‎ext/TrixiSparseConnectivityTracerExt.jl‎

Lines changed: 17 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,4 +10,21 @@ import SparseConnectivityTracer: AbstractTracer
1010
# we switch back to the Base implementation here which does not contain an if-clause.
1111
Trixi.sqrt(x::AbstractTracer) = Base.sqrt(x)
1212

13+
# if-clause free (i.e., non-optimized) implementations of some helper functions
14+
# that compute specialized mean values used in advanced flux functions
15+
16+
@inline function Trixi.ln_mean(x::AbstractTracer, y::AbstractTracer)
17+
return (y - x) / log(y / x)
18+
end
19+
20+
@inline function Trixi.inv_ln_mean(x::AbstractTracer, y::AbstractTracer)
21+
return log(y / x) / (y - x)
22+
end
23+
24+
@inline function Trixi.stolarsky_mean(x::AbstractTracer, y::AbstractTracer, gamma::Real)
25+
yg = exp((gamma - 1) * log(y)) # equivalent to y^(gamma - 1) but faster for non-integers
26+
xg = exp((gamma - 1) * log(x)) # equivalent to x^(gamma - 1) but faster for non-integers
27+
return (gamma - 1) * (yg * y - xg * x) / (gamma * (yg - xg))
28+
end
29+
1330
end

‎src/auxiliary/containers.jl‎

Lines changed: 28 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -350,6 +350,34 @@ function storage_type(C::Type{<:AbstractContainer})
350350
return storage_type(Adapt.unwrap_type(C))
351351
end
352352

353+
abstract type AbstractTreeElementContainer <: AbstractContainer end
354+
355+
# Return number of elements
356+
@inline nelements(elements::AbstractTreeElementContainer) = length(elements.cell_ids)
357+
# TODO: Taal performance, 1:nelements(elements) vs. Base.OneTo(nelements(elements))
358+
"""
359+
eachelement(elements::AbstractTreeElementContainer)
360+
361+
Return an iterator over the indices that specify the location in relevant data structures
362+
for the elements in `elements`.
363+
In particular, not the elements themselves are returned.
364+
"""
365+
@inline eachelement(elements::AbstractTreeElementContainer) = Base.OneTo(nelements(elements))
366+
@inline Base.real(elements::AbstractTreeElementContainer) = eltype(elements.node_coordinates)
367+
@inline nvariables(elements::AbstractTreeElementContainer) = size(elements.surface_flux_values,
368+
1)
369+
@inline nnodes(elements::AbstractTreeElementContainer) = size(elements.node_coordinates,
370+
2)
371+
@inline Base.eltype(elements::AbstractTreeElementContainer) = eltype(elements.surface_flux_values)
372+
373+
abstract type AbstractTreeBoundaryContainer <: AbstractContainer end
374+
375+
@inline nvariables(boundaries::AbstractTreeBoundaryContainer) = size(boundaries.u, 2)
376+
# Return number of boundaries
377+
@inline nboundaries(boundaries::AbstractTreeBoundaryContainer) = length(boundaries.orientations)
378+
# For 2D and 3D. 1D Hard-coded to 1
379+
@inline nnodes(boundaries::AbstractTreeBoundaryContainer) = size(boundaries.u, 3)
380+
353381
# backend handling
354382
"""
355383
trixi_backend(x)

‎src/auxiliary/math.jl‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -434,8 +434,8 @@ Given ε = 1.0e-4, we use the following algorithm.
434434
yg = y^(gamma - 1)
435435
xg = x^(gamma - 1)
436436
else
437-
yg = exp((gamma - 1) * log(y)) # equivalent to y^gamma but faster for non-integers
438-
xg = exp((gamma - 1) * log(x)) # equivalent to x^gamma but faster for non-integers
437+
yg = exp((gamma - 1) * log(y)) # equivalent to y^(gamma - 1) but faster for non-integers
438+
xg = exp((gamma - 1) * log(x)) # equivalent to x^(gamma - 1) but faster for non-integers
439439
end
440440
return (gamma - 1) * (yg * y - xg * x) / (gamma * (yg - xg))
441441
end

‎src/auxiliary/precompile.jl‎

Lines changed: 12 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -397,54 +397,54 @@ function _precompile_manual_()
397397
# 1D, serial
398398
@assert Base.precompile(Tuple{typeof(Trixi.init_boundaries), Array{Int, 1},
399399
TreeMesh{1, Trixi.SerialTree{1}, RealT},
400-
Trixi.ElementContainer1D{RealT, uEltype}})
400+
Trixi.TreeElementContainer1D{RealT, uEltype}})
401401
@assert Base.precompile(Tuple{typeof(Trixi.init_interfaces), Array{Int, 1},
402402
TreeMesh{1, Trixi.SerialTree{1}, RealT},
403-
Trixi.ElementContainer1D{RealT, uEltype}})
403+
Trixi.TreeElementContainer1D{RealT, uEltype}})
404404
@assert Base.precompile(Tuple{typeof(Trixi.save_mesh_file),
405405
TreeMesh{1, Trixi.SerialTree{1}, RealT}, String})
406406

407407
# 2D, serial
408408
@assert Base.precompile(Tuple{typeof(Trixi.init_boundaries), Array{Int, 1},
409409
TreeMesh{2, Trixi.SerialTree{2}, RealT},
410-
Trixi.ElementContainer2D{RealT, uEltype}})
410+
Trixi.TreeElementContainer2D{RealT, uEltype}})
411411
@assert Base.precompile(Tuple{typeof(Trixi.init_interfaces), Array{Int, 1},
412412
TreeMesh{2, Trixi.SerialTree{2}, RealT},
413-
Trixi.ElementContainer2D{RealT, uEltype}})
413+
Trixi.TreeElementContainer2D{RealT, uEltype}})
414414
@assert Base.precompile(Tuple{typeof(Trixi.init_mortars), Array{Int, 1},
415415
TreeMesh{2, Trixi.SerialTree{2}, RealT},
416-
Trixi.ElementContainer2D{RealT, uEltype},
416+
Trixi.TreeElementContainer2D{RealT, uEltype},
417417
mortar_type})
418418
@assert Base.precompile(Tuple{typeof(Trixi.save_mesh_file),
419419
TreeMesh{2, Trixi.SerialTree{2}, RealT}, String})
420420

421421
# 2D, parallel
422422
@assert Base.precompile(Tuple{typeof(Trixi.init_boundaries), Array{Int, 1},
423423
TreeMesh{2, Trixi.ParallelTree{2}, RealT},
424-
Trixi.ElementContainer2D{RealT, uEltype}})
424+
Trixi.TreeElementContainer2D{RealT, uEltype}})
425425
@assert Base.precompile(Tuple{typeof(Trixi.init_interfaces), Array{Int, 1},
426426
TreeMesh{2, Trixi.ParallelTree{2}, RealT},
427-
Trixi.ElementContainer2D{RealT, uEltype}})
427+
Trixi.TreeElementContainer2D{RealT, uEltype}})
428428
@assert Base.precompile(Tuple{typeof(Trixi.init_mortars), Array{Int, 1},
429429
TreeMesh{2, Trixi.ParallelTree{2}, RealT},
430-
Trixi.ElementContainer2D{RealT, uEltype},
430+
Trixi.TreeElementContainer2D{RealT, uEltype},
431431
mortar_type})
432432
@assert Base.precompile(Tuple{typeof(Trixi.init_mpi_interfaces), Array{Int, 1},
433433
TreeMesh{2, Trixi.ParallelTree{2}, RealT},
434-
Trixi.ElementContainer2D{RealT, uEltype}})
434+
Trixi.TreeElementContainer2D{RealT, uEltype}})
435435
@assert Base.precompile(Tuple{typeof(Trixi.save_mesh_file),
436436
TreeMesh{2, Trixi.ParallelTree{2}, RealT}, String})
437437

438438
# 3D, serial
439439
@assert Base.precompile(Tuple{typeof(Trixi.init_boundaries), Array{Int, 1},
440440
TreeMesh{3, Trixi.SerialTree{3}, RealT},
441-
Trixi.ElementContainer3D{RealT, uEltype}})
441+
Trixi.TreeElementContainer3D{RealT, uEltype}})
442442
@assert Base.precompile(Tuple{typeof(Trixi.init_interfaces), Array{Int, 1},
443443
TreeMesh{3, Trixi.SerialTree{3}, RealT},
444-
Trixi.ElementContainer3D{RealT, uEltype}})
444+
Trixi.TreeElementContainer3D{RealT, uEltype}})
445445
@assert Base.precompile(Tuple{typeof(Trixi.init_mortars), Array{Int, 1},
446446
TreeMesh{3, Trixi.SerialTree{3}, RealT},
447-
Trixi.ElementContainer3D{RealT, uEltype},
447+
Trixi.TreeElementContainer3D{RealT, uEltype},
448448
mortar_type})
449449
@assert Base.precompile(Tuple{typeof(Trixi.save_mesh_file),
450450
TreeMesh{3, Trixi.SerialTree{3}, RealT}, String})

‎src/solvers/dgsem_structured/containers.jl‎

Lines changed: 15 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -5,8 +5,8 @@
55
@muladd begin
66
#! format: noindent
77

8-
struct ElementContainer{NDIMS, RealT <: Real, uEltype <: Real,
9-
NDIMSP1, NDIMSP2, NDIMSP3}
8+
struct StructuredElementContainer{NDIMS, RealT <: Real, uEltype <: Real,
9+
NDIMSP1, NDIMSP2, NDIMSP3}
1010
# Physical coordinates at each node
1111
node_coordinates::Array{RealT, NDIMSP2} # [orientation, node_i, node_j, node_k, element]
1212
# ID of neighbor element in negative direction in orientation
@@ -49,22 +49,25 @@ function init_elements(mesh::Union{StructuredMesh{NDIMS, RealT},
4949
NDIMS - 1)..., NDIMS * 2,
5050
nelements)
5151

52-
elements = ElementContainer{NDIMS, RealT, uEltype, NDIMS + 1, NDIMS + 2, NDIMS + 3}(node_coordinates,
53-
left_neighbors,
54-
jacobian_matrix,
55-
contravariant_vectors,
56-
inverse_jacobian,
57-
surface_flux_values)
52+
elements = StructuredElementContainer{NDIMS, RealT, uEltype,
53+
NDIMS + 1, NDIMS + 2, NDIMS + 3}(node_coordinates,
54+
left_neighbors,
55+
jacobian_matrix,
56+
contravariant_vectors,
57+
inverse_jacobian,
58+
surface_flux_values)
5859

5960
init_elements!(elements, mesh, basis)
6061
return elements
6162
end
6263

63-
@inline nelements(elements::ElementContainer) = size(elements.left_neighbors, 2)
64-
@inline Base.ndims(::ElementContainer{NDIMS}) where {NDIMS} = NDIMS
64+
@inline nelements(elements::StructuredElementContainer) = size(elements.left_neighbors,
65+
2)
6566

66-
function Base.eltype(::ElementContainer{NDIMS, RealT, uEltype}) where {NDIMS, RealT,
67-
uEltype}
67+
function Base.eltype(::StructuredElementContainer{NDIMS, RealT, uEltype}) where {NDIMS,
68+
RealT,
69+
uEltype
70+
}
6871
return uEltype
6972
end
7073

0 commit comments

Comments
 (0)