diff --git a/src/auxiliary/precompile.jl b/src/auxiliary/precompile.jl index 4d5399b5ba3..962b2c377c9 100644 --- a/src/auxiliary/precompile.jl +++ b/src/auxiliary/precompile.jl @@ -244,10 +244,12 @@ function _precompile_manual_() end for TreeType in (SerialTree, ParallelTree), NDIMS in 1:3 @assert Base.precompile(Tuple{typeof(Trixi.initialize!), - TreeMesh{NDIMS, TreeType{NDIMS}}, Int, Tuple{}, + TreeMesh{NDIMS, TreeType{NDIMS}, Float64}, Int, + Tuple{}, Tuple{}}) @assert Base.precompile(Tuple{typeof(Trixi.save_mesh_file), - TreeMesh{NDIMS, TreeType{NDIMS}}, String, Int}) + TreeMesh{NDIMS, TreeType{NDIMS}, Float64}, String, + Int}) end # Constructors: linear advection @@ -396,58 +398,58 @@ function _precompile_manual_() # 1D, serial @assert Base.precompile(Tuple{typeof(Trixi.init_boundaries), Array{Int, 1}, - TreeMesh{1, Trixi.SerialTree{1}}, + TreeMesh{1, Trixi.SerialTree{1}, RealT}, Trixi.ElementContainer1D{RealT, uEltype}}) @assert Base.precompile(Tuple{typeof(Trixi.init_interfaces), Array{Int, 1}, - TreeMesh{1, Trixi.SerialTree{1}}, + TreeMesh{1, Trixi.SerialTree{1}, RealT}, Trixi.ElementContainer1D{RealT, uEltype}}) @assert Base.precompile(Tuple{typeof(Trixi.save_mesh_file), - TreeMesh{1, Trixi.SerialTree{1}}, String}) + TreeMesh{1, Trixi.SerialTree{1}, RealT}, String}) # 2D, serial @assert Base.precompile(Tuple{typeof(Trixi.init_boundaries), Array{Int, 1}, - TreeMesh{2, Trixi.SerialTree{2}}, + TreeMesh{2, Trixi.SerialTree{2}, RealT}, Trixi.ElementContainer2D{RealT, uEltype}}) @assert Base.precompile(Tuple{typeof(Trixi.init_interfaces), Array{Int, 1}, - TreeMesh{2, Trixi.SerialTree{2}}, + TreeMesh{2, Trixi.SerialTree{2}, RealT}, Trixi.ElementContainer2D{RealT, uEltype}}) @assert Base.precompile(Tuple{typeof(Trixi.init_mortars), Array{Int, 1}, - TreeMesh{2, Trixi.SerialTree{2}}, + TreeMesh{2, Trixi.SerialTree{2}, RealT}, Trixi.ElementContainer2D{RealT, uEltype}, mortar_type}) @assert Base.precompile(Tuple{typeof(Trixi.save_mesh_file), - TreeMesh{2, Trixi.SerialTree{2}}, String}) + TreeMesh{2, Trixi.SerialTree{2}, RealT}, String}) # 2D, parallel @assert Base.precompile(Tuple{typeof(Trixi.init_boundaries), Array{Int, 1}, - TreeMesh{2, Trixi.ParallelTree{2}}, + TreeMesh{2, Trixi.ParallelTree{2}, RealT}, Trixi.ElementContainer2D{RealT, uEltype}}) @assert Base.precompile(Tuple{typeof(Trixi.init_interfaces), Array{Int, 1}, - TreeMesh{2, Trixi.ParallelTree{2}}, + TreeMesh{2, Trixi.ParallelTree{2}, RealT}, Trixi.ElementContainer2D{RealT, uEltype}}) @assert Base.precompile(Tuple{typeof(Trixi.init_mortars), Array{Int, 1}, - TreeMesh{2, Trixi.ParallelTree{2}}, + TreeMesh{2, Trixi.ParallelTree{2}, RealT}, Trixi.ElementContainer2D{RealT, uEltype}, mortar_type}) @assert Base.precompile(Tuple{typeof(Trixi.init_mpi_interfaces), Array{Int, 1}, - TreeMesh{2, Trixi.ParallelTree{2}}, + TreeMesh{2, Trixi.ParallelTree{2}, RealT}, Trixi.ElementContainer2D{RealT, uEltype}}) @assert Base.precompile(Tuple{typeof(Trixi.save_mesh_file), - TreeMesh{2, Trixi.ParallelTree{2}}, String}) + TreeMesh{2, Trixi.ParallelTree{2}, RealT}, String}) # 3D, serial @assert Base.precompile(Tuple{typeof(Trixi.init_boundaries), Array{Int, 1}, - TreeMesh{3, Trixi.SerialTree{3}}, + TreeMesh{3, Trixi.SerialTree{3}, RealT}, Trixi.ElementContainer3D{RealT, uEltype}}) @assert Base.precompile(Tuple{typeof(Trixi.init_interfaces), Array{Int, 1}, - TreeMesh{3, Trixi.SerialTree{3}}, + TreeMesh{3, Trixi.SerialTree{3}, RealT}, Trixi.ElementContainer3D{RealT, uEltype}}) @assert Base.precompile(Tuple{typeof(Trixi.init_mortars), Array{Int, 1}, - TreeMesh{3, Trixi.SerialTree{3}}, + TreeMesh{3, Trixi.SerialTree{3}, RealT}, Trixi.ElementContainer3D{RealT, uEltype}, mortar_type}) @assert Base.precompile(Tuple{typeof(Trixi.save_mesh_file), - TreeMesh{3, Trixi.SerialTree{3}}, String}) + TreeMesh{3, Trixi.SerialTree{3}, RealT}, String}) end # various `show` methods @@ -456,16 +458,16 @@ function _precompile_manual_() for NDIMS in 1:3 # serial @assert Base.precompile(Tuple{typeof(show), Base.TTY, - TreeMesh{NDIMS, Trixi.SerialTree{NDIMS}}}) + TreeMesh{NDIMS, Trixi.SerialTree{NDIMS}, RealT}}) @assert Base.precompile(Tuple{typeof(show), IOContext{Base.TTY}, MIME"text/plain", - TreeMesh{NDIMS, Trixi.SerialTree{NDIMS}}}) + TreeMesh{NDIMS, Trixi.SerialTree{NDIMS}, RealT}}) # parallel @assert Base.precompile(Tuple{typeof(show), Base.TTY, - TreeMesh{NDIMS, Trixi.ParallelTree{NDIMS}}}) + TreeMesh{NDIMS, Trixi.ParallelTree{NDIMS}, RealT}}) @assert Base.precompile(Tuple{typeof(show), IOContext{Base.TTY}, MIME"text/plain", - TreeMesh{NDIMS, Trixi.ParallelTree{NDIMS}}}) + TreeMesh{NDIMS, Trixi.ParallelTree{NDIMS}, RealT}}) end # equations diff --git a/src/meshes/parallel_tree.jl b/src/meshes/parallel_tree.jl index 7175ed47743..52c84af6169 100644 --- a/src/meshes/parallel_tree.jl +++ b/src/meshes/parallel_tree.jl @@ -25,12 +25,12 @@ # function, which is required for implementing level-wise refinement in a sane # way. Also, depth-first ordering *might* not be guaranteed during # refinement/coarsening operations. -mutable struct ParallelTree{NDIMS} <: AbstractTree{NDIMS} +mutable struct ParallelTree{NDIMS, RealT <: Real} <: AbstractTree{NDIMS} parent_ids::Vector{Int} child_ids::Matrix{Int} neighbor_ids::Matrix{Int} levels::Vector{Int} - coordinates::Matrix{Float64} + coordinates::Matrix{RealT} original_cell_ids::Vector{Int} mpi_ranks::Vector{Int} @@ -38,11 +38,11 @@ mutable struct ParallelTree{NDIMS} <: AbstractTree{NDIMS} length::Int dummy::Int - center_level_0::SVector{NDIMS, Float64} - length_level_0::Float64 + center_level_0::SVector{NDIMS, RealT} + length_level_0::RealT periodicity::NTuple{NDIMS, Bool} - function ParallelTree{NDIMS}(capacity::Integer) where {NDIMS} + function ParallelTree{NDIMS, RealT}(capacity::Integer) where {NDIMS, RealT <: Real} # Verify that NDIMS is an integer @assert NDIMS isa Integer @@ -71,13 +71,14 @@ mutable struct ParallelTree{NDIMS} <: AbstractTree{NDIMS} end # Constructor for passing the dimension as an argument -ParallelTree(::Val{NDIMS}, args...) where {NDIMS} = ParallelTree{NDIMS}(args...) +ParallelTree(::Val{NDIMS}, args...) where {NDIMS, RealT <: Real} = ParallelTree{NDIMS}(args...) # Create and initialize tree -function ParallelTree{NDIMS}(capacity::Int, center::AbstractArray{Float64}, - length::Real, periodicity = true) where {NDIMS} +function ParallelTree{NDIMS}(capacity::Int, center::AbstractArray{RealT}, + length::Real, + periodicity = true) where {NDIMS, RealT <: Real} # Create instance - t = ParallelTree{NDIMS}(capacity) + t = ParallelTree{NDIMS, RealT}(capacity) # Initialize root cell init!(t, center, length, periodicity) @@ -86,13 +87,14 @@ function ParallelTree{NDIMS}(capacity::Int, center::AbstractArray{Float64}, end # Constructor accepting a single number as center (as opposed to an array) for 1D -function ParallelTree{1}(cap::Int, center::Real, len::Real, periodicity = true) - ParallelTree{1}(cap, [convert(Float64, center)], len, periodicity) +function ParallelTree{1}(cap::Int, center::RealT, len::RealT, + periodicity = true) where {RealT <: Real} + ParallelTree{1}(cap, [center], len, periodicity) end # Clear tree with deleting data structures, store center and length, and create root cell -function init!(t::ParallelTree, center::AbstractArray{Float64}, length::Real, - periodicity = true) +function init!(t::ParallelTree, center::AbstractArray{RealT}, length::RealT, + periodicity = true) where {RealT} clear!(t) # Set domain information @@ -222,12 +224,13 @@ function raw_copy!(target::ParallelTree, source::ParallelTree, first::Int, last: end # Reset data structures by recreating all internal storage containers and invalidating all elements -function reset_data_structures!(t::ParallelTree{NDIMS}) where {NDIMS} +function reset_data_structures!(t::ParallelTree{NDIMS, RealT}) where {NDIMS, + RealT <: Real} t.parent_ids = Vector{Int}(undef, t.capacity + 1) t.child_ids = Matrix{Int}(undef, 2^NDIMS, t.capacity + 1) t.neighbor_ids = Matrix{Int}(undef, 2 * NDIMS, t.capacity + 1) t.levels = Vector{Int}(undef, t.capacity + 1) - t.coordinates = Matrix{Float64}(undef, NDIMS, t.capacity + 1) + t.coordinates = Matrix{RealT}(undef, NDIMS, t.capacity + 1) t.original_cell_ids = Vector{Int}(undef, t.capacity + 1) t.mpi_ranks = Vector{Int}(undef, t.capacity + 1) diff --git a/src/meshes/serial_tree.jl b/src/meshes/serial_tree.jl index 143ac19f6ee..74201f85365 100644 --- a/src/meshes/serial_tree.jl +++ b/src/meshes/serial_tree.jl @@ -25,23 +25,23 @@ # function, which is required for implementing level-wise refinement in a sane # way. Also, depth-first ordering *might* not by guaranteed during # refinement/coarsening operations. -mutable struct SerialTree{NDIMS} <: AbstractTree{NDIMS} +mutable struct SerialTree{NDIMS, RealT <: Real} <: AbstractTree{NDIMS} parent_ids::Vector{Int} child_ids::Matrix{Int} neighbor_ids::Matrix{Int} levels::Vector{Int} - coordinates::Matrix{Float64} + coordinates::Matrix{RealT} original_cell_ids::Vector{Int} capacity::Int length::Int dummy::Int - center_level_0::SVector{NDIMS, Float64} - length_level_0::Float64 + center_level_0::SVector{NDIMS, RealT} + length_level_0::RealT periodicity::NTuple{NDIMS, Bool} - function SerialTree{NDIMS}(capacity::Integer) where {NDIMS} + function SerialTree{NDIMS, RealT}(capacity::Integer) where {NDIMS, RealT <: Real} # Verify that NDIMS is an integer @assert NDIMS isa Integer @@ -69,13 +69,15 @@ mutable struct SerialTree{NDIMS} <: AbstractTree{NDIMS} end # Constructor for passing the dimension as an argument -SerialTree(::Val{NDIMS}, args...) where {NDIMS} = SerialTree{NDIMS}(args...) +SerialTree(::Val{NDIMS}, args...) where {NDIMS, RealT <: Real} = SerialTree{NDIMS, + RealT}(args...) # Create and initialize tree -function SerialTree{NDIMS}(capacity::Int, center::AbstractArray{Float64}, - length::Real, periodicity = true) where {NDIMS} +function SerialTree{NDIMS}(capacity::Int, center::AbstractArray{RealT}, + length::RealT, + periodicity = true) where {NDIMS, RealT <: Real} # Create instance - t = SerialTree{NDIMS}(capacity) + t = SerialTree{NDIMS, RealT}(capacity) # Initialize root cell init!(t, center, length, periodicity) @@ -84,13 +86,14 @@ function SerialTree{NDIMS}(capacity::Int, center::AbstractArray{Float64}, end # Constructor accepting a single number as center (as opposed to an array) for 1D -function SerialTree{1}(cap::Int, center::Real, len::Real, periodicity = true) - SerialTree{1}(cap, [convert(Float64, center)], len, periodicity) +function SerialTree{1}(cap::Int, center::RealT, len::RealT, + periodicity = true) where {RealT <: Real} + SerialTree{1, RealT}(cap, [center], len, periodicity) end # Clear tree with deleting data structures, store center and length, and create root cell -function init!(t::SerialTree, center::AbstractArray{Float64}, length::Real, - periodicity = true) +function init!(t::SerialTree, center::AbstractArray{RealT}, length::RealT, + periodicity = true) where {RealT} clear!(t) # Set domain information @@ -203,12 +206,13 @@ function raw_copy!(target::SerialTree, source::SerialTree, first::Int, last::Int end # Reset data structures by recreating all internal storage containers and invalidating all elements -function reset_data_structures!(t::SerialTree{NDIMS}) where {NDIMS} +function reset_data_structures!(t::SerialTree{NDIMS, RealT}) where {NDIMS, RealT <: + Real} t.parent_ids = Vector{Int}(undef, t.capacity + 1) t.child_ids = Matrix{Int}(undef, 2^NDIMS, t.capacity + 1) t.neighbor_ids = Matrix{Int}(undef, 2 * NDIMS, t.capacity + 1) t.levels = Vector{Int}(undef, t.capacity + 1) - t.coordinates = Matrix{Float64}(undef, NDIMS, t.capacity + 1) + t.coordinates = Matrix{RealT}(undef, NDIMS, t.capacity + 1) t.original_cell_ids = Vector{Int}(undef, t.capacity + 1) invalidate!(t, 1, capacity(t) + 1) diff --git a/src/meshes/tree_mesh.jl b/src/meshes/tree_mesh.jl index 458a43203d2..051b5468c3a 100644 --- a/src/meshes/tree_mesh.jl +++ b/src/meshes/tree_mesh.jl @@ -24,16 +24,19 @@ get_name(mesh::AbstractMesh) = mesh |> typeof |> nameof |> string A Cartesian mesh based on trees of hypercubes to support adaptive mesh refinement. """ -mutable struct TreeMesh{NDIMS, TreeType <: AbstractTree{NDIMS}} <: AbstractMesh{NDIMS} +mutable struct TreeMesh{NDIMS, TreeType <: AbstractTree{NDIMS}, RealT <: Real} <: + AbstractMesh{NDIMS} tree::TreeType current_filename::String unsaved_changes::Bool first_cell_by_rank::OffsetVector{Int, Vector{Int}} n_cells_by_rank::OffsetVector{Int, Vector{Int}} - function TreeMesh{NDIMS, TreeType}(n_cells_max::Integer) where {NDIMS, - TreeType <: - AbstractTree{NDIMS}} + function TreeMesh{NDIMS, TreeType, RealT}(n_cells_max::Integer) where {NDIMS, + TreeType <: + AbstractTree{NDIMS}, + RealT <: + Real} # Create mesh m = new() m.tree = TreeType(n_cells_max) @@ -46,14 +49,13 @@ mutable struct TreeMesh{NDIMS, TreeType <: AbstractTree{NDIMS}} <: AbstractMesh{ end # TODO: Taal refactor, order of important arguments, use of n_cells_max? - # TODO: Taal refactor, allow other RealT for the mesh, not just Float64 - # TODO: Taal refactor, use NTuple instead of domain_center::AbstractArray{Float64} - function TreeMesh{NDIMS, TreeType}(n_cells_max::Integer, - domain_center::AbstractArray{Float64}, - domain_length, - periodicity = true) where {NDIMS, - TreeType <: - AbstractTree{NDIMS}} + function TreeMesh{NDIMS, TreeType, RealT}(n_cells_max::Integer, + domain_center::AbstractArray{RealT}, + domain_length, + periodicity = true) where {NDIMS, + TreeType <: + AbstractTree{NDIMS}, + RealT <: Real} @assert NDIMS isa Integer && NDIMS > 0 # Create mesh @@ -87,30 +89,32 @@ function TreeMesh(::Type{TreeType}, end # Constructor accepting a single number as center (as opposed to an array) for 1D -function TreeMesh{1, TreeType}(n::Int, center::Real, len::Real, - periodicity = true) where {TreeType <: AbstractTree{1}} - # TODO: Taal refactor, allow other RealT for the mesh, not just Float64 - return TreeMesh{1, TreeType}(n, SVector{1, Float64}(center), len, periodicity) +function TreeMesh{1, TreeType, RealT}(n::Int, center::RealT, len::RealT, + periodicity = true) where { + TreeType <: + AbstractTree{1}, + RealT <: Real} + return TreeMesh{1, TreeType, RealT}(n, SVector{1, RealT}(center), len, periodicity) end -function TreeMesh{NDIMS, TreeType}(n_cells_max::Integer, - domain_center::NTuple{NDIMS, Real}, - domain_length::Real, - periodicity = true) where {NDIMS, - TreeType <: - AbstractTree{NDIMS}} - # TODO: Taal refactor, allow other RealT for the mesh, not just Float64 - TreeMesh{NDIMS, TreeType}(n_cells_max, SVector{NDIMS, Float64}(domain_center), - convert(Float64, domain_length), periodicity) +function TreeMesh{NDIMS, TreeType, RealT}(n_cells_max::Integer, + domain_center::NTuple{NDIMS, RealT}, + domain_length::RealT, + periodicity = true) where {NDIMS, + TreeType <: + AbstractTree{NDIMS}, + RealT <: Real} + TreeMesh{NDIMS, TreeType, RealT}(n_cells_max, SVector{NDIMS, RealT}(domain_center), + domain_length, periodicity) end -function TreeMesh(coordinates_min::NTuple{NDIMS, Real}, - coordinates_max::NTuple{NDIMS, Real}; +function TreeMesh(coordinates_min::NTuple{NDIMS, RealT}, + coordinates_max::NTuple{NDIMS, RealT}; n_cells_max, periodicity = true, initial_refinement_level, refinement_patches = (), - coarsening_patches = ()) where {NDIMS} + coarsening_patches = ()) where {NDIMS, RealT} # check arguments if !(n_cells_max isa Integer && n_cells_max > 0) throw(ArgumentError("`n_cells_max` must be a positive integer (provided `n_cells_max = $n_cells_max`)")) @@ -139,10 +143,10 @@ function TreeMesh(coordinates_min::NTuple{NDIMS, Real}, end # Create mesh - mesh = @trixi_timeit timer() "creation" TreeMesh{NDIMS, TreeType}(n_cells_max, - domain_center, - domain_length, - periodicity) + mesh = @trixi_timeit timer() "creation" TreeMesh{NDIMS, TreeType, RealT}(n_cells_max, + domain_center, + domain_length, + periodicity) # Initialize mesh initialize!(mesh, initial_refinement_level, refinement_patches, coarsening_patches)